analysis-I.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 8.2. [20 points] Analysis Part I¶\n",
    "Find 95% lower (one-sided) and two-sided confidence intervals for the reduction in risk corresponding to the primary endpoint (data “Through day 29”), using method 3 and also using the cruder conservative approach via simultaneous Bonferroni confidence bounds for N⋅1 and N1⋅ described in the notes on causal inference. (For the Bonferroni approach to two-sided intervals, use Sterne’s method for the underlying hypergeometric confidence intervals. Feel free to re-use your own code from the previous problem set.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from cibin.permute_utils import hypergeom_conf_interval\n",
    "from cibin.utils import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Data from Regeneron Press"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:06.391453Z",
     "start_time": "2021-05-15T04:59:06.376341Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(11, 742, 59, 693)"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "n=753\n",
    "m=752\n",
    "N=n+m\n",
    "n01 = 59\n",
    "n11 = 11\n",
    "n00 = m-n01\n",
    "n10 = n-n11\n",
    "alpha=0.05\n",
    "n11, n10, n01, n00"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## One-sided confidence intervals"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Upper one-sided"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:07.138447Z",
     "start_time": "2021-05-15T04:59:07.122074Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(397, [-0.5, 0.2637873754152824])"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "N_dot_1 = hypergeom_conf_interval(n11*N/n, n11, N, 1-alpha, alternative='upper')\n",
    "N_1_dot = hypergeom_conf_interval(n01*N/m, n01, N, 1-alpha, alternative='lower')\n",
    "\n",
    "upper = (N_dot_1[1] - N_1_dot[0])/N\n",
    "\n",
    "ci = [-0.5, upper]\n",
    "N_dot_1[1] - N_1_dot[0], ci"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Lower one-sided"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:08.311354Z",
     "start_time": "2021-05-15T04:59:08.284404Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(395, [-0.26245847176079734, 0.5])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "N_dot_1 = hypergeom_conf_interval(n11*N/n, n11, N, 1-alpha, alternative='lower')\n",
    "N_1_dot = hypergeom_conf_interval(n01*N/m, n01, N, 1-alpha, alternative='upper')\n",
    "\n",
    "lower = (N_dot_1[0] - N_1_dot[1])/N\n",
    "\n",
    "ci = [lower, 0.5]\n",
    "N_1_dot[1] - N_dot_1[0], ci"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "start_time": "2021-05-15T04:57:45.231Z"
    }
   },
   "outputs": [],
   "source": [
    "tau_lower_oneside(n11, n10, n01, n00, alpha=0.05, nperm=1) # result around (-0.15, 0.42), span 858"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Two-sided"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:10.495344Z",
     "start_time": "2021-05-15T04:59:10.469153Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(515, 516, [-0.34219269102990035, 0.34285714285714286])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "N_dot_1 = hypergeom_conf_interval(n11*N/n, n11, N, 1-alpha/2, alternative='two-sided')\n",
    "N_1_dot = hypergeom_conf_interval(n01*N/m, n01, N, 1-alpha/2, alternative='two-sided')\n",
    "lower = (N_dot_1[0] - N_1_dot[1])/N\n",
    "upper = (N_dot_1[1] - N_1_dot[0])/N\n",
    "ci = [lower, upper]\n",
    "N_1_dot[1] - N_dot_1[0], N_dot_1[1] - N_1_dot[0], ci"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:38.469606Z",
     "start_time": "2021-05-15T04:59:17.076519Z"
    }
   },
   "outputs": [],
   "source": [
    "two_side = hypergeom_conf_interval(n, n01, N, 1-alpha/2, alternative='two-sided', method='Sterne')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:38.980600Z",
     "start_time": "2021-05-15T04:59:38.976706Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "49"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "two_side[1]-two_side[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:59:39.622738Z",
     "start_time": "2021-05-15T04:59:39.617553Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.06312292358803986, 0.09568106312292358)"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(two_side[0]/N, two_side[1]/N)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-05-15T04:40:26.581612Z",
     "start_time": "2021-05-15T04:40:22.638324Z"
    }
   },
   "outputs": [],
   "source": [
    "tau_twosided_ci(n11, n10, n01, n00, 0.05, exact=False, max_combinations=10**5, reps=100) # result around (-0.15, 0.42), span 858"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Discuss the differences between the two sets of confidence intervals."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "According to Sterne's method, one-sided interval is wider than the two-sided one. On the other hand, for Li and Ding's method, the two-sided interval is wider than the one-sided interval. \n",
    "\n",
    "Overall, the interval according to Sterne's method is wider as Sterne is more conservative."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Is it statistically legitimate to use one-sided confidence intervals? Why or why not?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I think it is statistically legitimate to use one-sided confidence intervals if we support the Regeneron study's statement. This is because one-sided confidence interval favor one bound,upper or lower. The study claims that subjects with their antibody cocktail are less likely to be infected. He wants to show that the antibody cocktail can make the infection rate reduce certain percentage with high confidence. Hence, one-sided confidence interval works.\n",
    "\n",
    "However, I would choose two-sided confidence intervals if we would like to coduct a study in general. Two-sided confidence interval will value both upper and lower bound. As a result, if we want to study the changes in infection rate of COVID-19, using two-sided confidence intervals can be more statistically legitimate."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "- Are the 2-sided confidence intervals preferable to the one-sided intervals? Why or why not?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In general, I think it is preferable to use a two-sided confidence interval as one-sided can accept more extreme outcomes. Nonetheless, it also depends on the assumptions for the data. In this situation, a one-sided confidence interval works better because we assume that the sample risk is less than the population risk. In a more general case, we may prefer a two-sided confidence interval as I mentioned above. We can be able to understand the increase or decrease and overall changes using two-sided interval."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "hw7",
   "language": "python",
   "name": "hw7"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": true
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}