{ "cells": [ { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# Creative Application 3\n", "by Alvaro Lavalle and Asal Bastani" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
1st Dose: Pfizer1st Dose: Astrazeneca
Systematic Side Effects38155116473
No Systematic Side Effects243948228807
" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Table Used for Cocalc Analysis\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
1st Dose: Pfizer1st Dose: Astrazeneca
Systematic Side Effects38116
No Systematic Side Effects244229
" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#1 Find the description of a study from sources such as CDC, FDA, published articles on PubMed, or even a news article with actual data." ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# Website with Data:\n", "# https://www.thelancet.com/action/showPdf?pii=S1473-3099%2821%2900224-3\n", "# This study analyzed whether the Pfizer-BioNTech and Oxford-AstraZeneca COVID-19 vaccines are both safe and effective during the Phase 3 Trials. The main aim of the study was to analyze the safety and effectiveness of each virus within the United Kingdom" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#2 Perform a suitable statistical analysis to test whether the vaccine causes some positive or negative outcome (either prevention of a level of infection or a side effect). Briefly explain your analysis and outcome of interest and provide the results of your analysis." ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#We will be utilizing a chi squared statistical analysis with a relative risk calculation as well. The outcome of interest is to understand whether the vaccines produced systematic side effects on individuals who received the 1st dose of their respective vaccine. We will explain the results of the analysis at the very end." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#importing libraries and packages\n", "import numpy as np\n", "import seaborn as sns\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "array([[ 38, 116],\n", " [244, 229]])" ] }, "execution_count": 3, "metadata": { }, "output_type": "execute_result" } ], "source": [ "#our original table of analysis, but we quickly want to note that we had to divide all numbers by 1000 due to the fact that Cocalc would not have been able to process the Chi Squared Statistic for such high numbers\n", "reality=np.array([[38,116],[244,229]])\n", "reality" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[[ 69.26315789 84.73684211]\n", " [212.73684211 260.26315789]]\n" ] } ], "source": [ "#creating the expected table\n", "pfizer=38+244\n", "astra=116+229\n", "side_effect=38+116\n", "no_side_effect=244+229\n", "total=38+244+116+229\n", "pfizer_side= (pfizer*side_effect)/total\n", "pfizer_no_side=(pfizer*no_side_effect)/total\n", "astra_side=(astra*side_effect)/total\n", "astra_no_side=(astra*no_side_effect)/total\n", "expectation=np.array([[pfizer_side,astra_side],[pfizer_no_side,astra_no_side]])\n", "print(expectation)" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "33.99525049096289\n" ] } ], "source": [ "#calculation of chi squared statistic\n", "chi_squared_expected=np.sum((reality-expectation)**2/(expectation))\n", "print(chi_squared_expected)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#second way to check our chi squared statistic by making a function to check\n", "def chisquared (observed,expected):\n", " chi_squared_expected=np.sum((observed-expected)**2/(expected))\n", " return (chi_squared_expected)" ] }, { "cell_type": "code", "execution_count": 20, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "33.99525049096289" ] }, "execution_count": 20, "metadata": { }, "output_type": "execute_result" } ], "source": [ "#this our observed chi square statistic\n", "chisquared(reality,expectation)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#making a list of outcomes as well as a place to store your results\n", "#then making a for loop to make simulated chi square tables and calculate the chi squared statistic each time and storing it to the list of zeroes\n", "trial=[\"S\"]*side_effect + [\"N\"]*no_side_effect\n", "hist_data=np.zeros(10000)\n", "for i in range (10000):\n", " sim_data=np.zeros([2,2])\n", " sim_pfizer=np.random.choice(trial, pfizer)\n", " sim_data[0,0]=np.sum(sim_pfizer=='S')\n", " sim_data[1,0]=np.sum(sim_pfizer=='N')\n", " sim_astra=np.random.choice(trial,astra)\n", " sim_data[0,1]=np.sum(sim_astra=='S')\n", " sim_data[1,1]=np.sum(sim_astra=='N')\n", " hist_data[i]=chisquared(sim_data,expectation)" ] }, { "cell_type": "code", "execution_count": 22, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 22, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsMAAALACAYAAAB2GH/sAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAABYlAAAWJQFJUiTwAAAt5klEQVR4nO3de5hldX3n+88XegSlTWtaEyXEeEkQT8wNOJrgjO1lxvFC1CTkjCcYnbRpxkRDVHCSx0u8nGhubWJaEj2Yg5iYczDBRxMyXpgRgUSNxtvQKIJGUJnGGBtoBYQW/Z0/9qpmU1R1V3Xtql27f6/X89Sze6/LXr9azwLevVh7rWqtBQAAenTYtAcAAADTIoYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbonhKaiqv6yqv5z2OAAAerdh2gPo1HHHH3/88Ul+YdoDAQCYpKo7/tza2mxyJSs7MwwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLc2THsArL7TTj8zu3bvSZIcvXlTzt6xfcojAgBYH8RwB3bt3pONW7aO/nzJOVMeDQDA+uEyCQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDo1oZpD4C1dfnOnTn51G373h+9eVPO3rF9iiMCAJgeMdyZve2wbNyydd/7C8964b44FsYAQG/EcOfG43jXJedMeTQAAGvLNcMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRrJmO4qp5SVRdW1bVV9c2q+kJV/XVV/dQiy59UVe+uquuH5S+rqhdU1eH72cbJVXVxVe2pqpuq6iNV9ezV+60AAFhrMxfDVfV7Sf4uyfFJ3pvkj5N8IsnTknywqp45b/mnJbk0yaOTvDPJWUnuluSPkpy3yDaen+SCJA9P8rYkb05ydJJzq2r75H8rAACmYcO0B7AcVXW/JGcm+ZckP9pa++rYvMcmuSjJqzMK2FTVd2UUst9O8pjW2seG6S8flj2lqp7RWjtv7HMemGR7kuuTnNhau2aY/uok/5TkjKp6R2vtw6v72wIAsNpm7czwD2Q05o+Mh3CStNY+kOQbSe47NvmU4f15cyE8LHtrkpcNb39l3ja2JjkiyVlzITysc0OS1w5vn7vi3wQAgKmbqTPDST6XZG+SR1TVfVprX5ubUVWPTnLPJO8aW/5xw+t7F/isS5PckuSkqjqitXbbEtZ5z7xl9quqPr7IrOOWsj4AAKtrps4Mt9auT/IbSb43yWeq6uyq+p2q+qskFyb570n+y9gqDx1er1rgs25PcnVGfyF48BLXuS7JzUmOqap7rPDXAQBgymbtzHBaa6+vqmuSnJNk29iszyc5d97lE5uG1z2LfNzc9Hstc52jhuVuOcBYT1ho+nDG+Pj9rQsAwOqbqTPDSVJV/zXJ+UnOTfKQjML0hCRfSPKXVfX70xsdAACzZKZiuKoek+T3kvxta+1FrbUvtNZuaa19IsnPJPlfGd3tYe6yh7mzu5vu8mF3nn7j2LSlrrPYmWMAAGbETMVwkpOH1w/Mn9FauyXJRzP6nX5imHzl8Hrs/OWrakOSByW5PaOzylnCOvfP6Ez0tcP2AACYYbMWw0cMr/ddZP7c9L3D60XD6xMXWPbRSe6R5ENjd5I40DpPmrcMAAAzbNZi+O+H19Oq6vvGZ1TVk5I8KsmtST40TD4/ydeSPKOqThxb9sgkvz28feO8bbwlyW1Jnj88gGNunXsnecnw9k0r/k0AAJi6WbubxPlJ/keSf5/kiqp6Z5KvJHlYRpdQVJLfbK3tTpLW2teratuw3sVVdV5GT5Z7aka3UDs/ydvHN9Bau7qqXpxkR5KPVdXbMzrTfEqSY5K8ztPnAAAODTMVw62171TVk5M8L8kzMvrS3D0yCtx3J9nRWrtw3jrvqqotSV6a5OeSHJnRbdheNCzfFtjOG4bbt52Z5FkZnUH/TJKXtdbeukq/HgAAa2ymYjhJWmvfSvL64Wep63wwyZOXuZ0LklywnHUAAJgts3bNMAAATIwYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBubZj2AFg/Lt+5Myefui1JcvTmTTl7x/YpjwgAYHWJYfbZ2w7Lxi1bkyS7LjlnyqMBAFh9LpMAAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAurVh2gNg/Tvt9DOza/eefe+P3rwpZ+/YPsURAQBMhhjmgHbt3pONW7be8f6Sc6Y4GgCAyXGZBAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0SwwDANAtMQwAQLfEMAAA3RLDAAB0a2ZjuKoeX1XvrKqvVNVtVbWrqt5XVU9eYNmTqurdVXV9VX2zqi6rqhdU1eH7+fyTq+riqtpTVTdV1Ueq6tmr+1sBALCWZjKGq+r3k/yPJCcm+dskr0vy35LcN8lj5i37tCSXJnl0kncmOSvJ3ZL8UZLzFvn85ye5IMnDk7wtyZuTHJ3k3KraPvFfCACAqdgw7QEsV1VtS/LiJG9Nclprbe+8+f9m7M/flVHIfjvJY1prHxumvzzJRUlOqapntNbOG1vngUm2J7k+yYmttWuG6a9O8k9Jzqiqd7TWPrxqvyQAAGtips4MV9URSV6T5EtZIISTpLX2rbG3p2R0tvi8uRAelrk1ycuGt78y7yO2JjkiyVlzITysc0OS1w5vn7uy3wQAgPVg1s4M/4eM4vb1Sb5TVU/J6FKGW5N8dIGztY8bXt+7wGddmuSWJCdV1RGttduWsM575i0DAMAMm7UY/t+H11uTfDKjEN6nqi5Nckpr7V+HSQ8dXq+a/0Gttdur6uokP5zkwUmuWMI611XVzUmOqap7tNZu2d9gq+rji8w6bn/rAQCwNmbqMokk3zO8vjhJS/LvktwzyY8muTCjL8n99djym4bXPYt83tz0ex3EOpsWmQ8AwIyYtTPDc/F+e5Knjl3Tu7OqfibJlUm2VNVPrYcvuLXWTlho+nDG+Pg1Hg4AAPPM2pnhG4fXT45/uS1JhksW3je8fcTweqCzuHPTbxybttR1FjtzfEi4fOfOnHzqtpx86rZ8+oorpz0cAIBVMWsxPFdlNy4y/4bh9e7zlj92/oJVtSHJgzI6y/yFBbax0Dr3T3JUkmsPdL3wrNvbDsvGLVuzccvW3Pat26c9HACAVTFrMfz+jK4V/t+qaqGxz32h7urh9aLh9YkLLPvoJPdI8qGxO0kcaJ0nzVsGAIAZNlMx3Fr7YkZPhntAkl8fn1dVT0jyHzM6azx3W7Tzk3wtyTOq6sSxZY9M8tvD2zfO28xbktyW5PnDAzjm1rl3kpcMb9+08t8GAIBpm7Uv0CXJ85L8RJI/HO4z/MmMLnd4ekZPmvvl1tqeJGmtfX14Yt35SS6uqvMyerLcUzO6hdr5Sd4+/uGttaur6sVJdiT5WFW9PcnejB7gcUyS162HL+cBALByMxfDrbVrq+qEJL+VUdQ+OsnXMzpj/DuttY/OW/5dVbUlyUuT/FySI5N8PsmLkuxorbUFtvGGqromyZlJnpXRGfTPJHlZa+2tq/W7AQCwtmYuhpNkeKjGrw0/S1n+g0mevMxtXJBRYAMAcIiaqWuGAQBgksQwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRrojFcVQ+oqu86wDL3rKoHTHK7AABwMCZ9ZvjqJL9+gGVOH5YDAICpmnQM1/ADAADr3jSuGb5fkpunsF0AALiTDSv9gKp61rxJP77AtCQ5PMkDkjwzyc6VbhcAAFZqxTGc5NwkbfhzS/K04We+ucsnbknyqglsFwAAVmQSMfxLw2slOSfJu5L8zQLLfTvJ7iQfbq3dOIHtAgDAiqw4hltrb537c1U9O8m7Wmt/vtLPBQCA1TaJM8P7tNYeO8nPAwCA1eQJdAAAdGuiZ4aTpKq2JHlxkkckuXcWDu7WWpv4tlkbl+/cmZNP3ZYkOXrzppy9Y/uURwQAcHAmGqRV9ZSMvkB3eJIvJbkyye2T3AbTt7cdlo1btiZJdl1yzpRHAwBw8CZ9dvaVSb6V5CmttQsn/NkAADBRk75m+OFJ3i6EAQCYBZOO4ZuSXD/hzwQAgFUx6Rh+f5KfmvBnAgDAqph0DP9GkodU1cuqqg64NAAATNGkv0D3iiSfTvKqJFur6lNJblxgudZae86Etw0AAMsy6Rj+z2N/fuDws5CWRAwDADBVk47hB0348wAAYNVMNIZba1+c5OcBAMBqmvQX6AAAYGZM+nHMD1jqsq21L01y2wAAsFyTvmb4moy+HHcgbRW2DQAAyzLpIP3zLBzD90ry40l+IMnFSVxbDADA1E36C3T/ebF5VXVYkpcneW6SZ09yuwAAcDDW7At0rbXvtNZeldGlFL+7VtsFAIDFTONuEh9K8oQpbBcAAO5kGjH83UmOmsJ2AQDgTtY0hqvq3yf5T0kuX8vtAgDAQiZ9n+GL9rOd708ydx/iV09yuwAAcDAmfWu1xywyvSW5Icn7kmxvrS0WzQAAsGYmfWs1j3cGAGBmiFcAALq1qo9Erqp7ZvT0uT2tta+v5rYAAGC5Jn5muKo2VNVvVtXnk9yY0UM2bqiqzw/TVzXAAQBgqSZ9N4m7JXlvki0ZfWnuy0muS3L/JA9M8pokT6yqJ7TW9k5y2wAAsFyTPkv7oozuKPF3Sc5orX1ubkZVPSTJ65L89LCcRzIfYk47/czs2r1n3/ujN2/K2Tu2T3FEAAD7N+kY/oWMHqjx9Nbad8ZntNb+uap+NsmnkpwaMXzI2bV7TzZu2XrH+0vOmeJoAAAObNLXDP9gkvfMD+E5w/T3JHnIhLcLAADLNukY3ptk4wGWOSrJtya8XQAAWLZJx/BlSU6pqvsuNLOq7pPklCT/c8LbBQCAZZt0DJ+V5L5JPlpVz6mqB1fV3avqQVX1S0k+Msw/a8LbBQCAZZv045j/qqp+PMlvJjl7gUUqye+31v5qktsFAICDMfEHYLTWXlJVf5vkOUl+IsmmJHuSfDLJOa21D096mwAAcDBW5WlwrbV/TPKPq/HZAAAwKSu+Zriq7lZVH62q91fVvznAcu+vqn/c33IAALBWJvEFumcmOSHJ61pri94ybXj88h8keURGD90AAICpmkQM/2ySL7TW3n2gBVtr703yuSQ/P4HtAgDAikwihn8iycXLWP7SJD8+ge0CAMCKTCKG75PkX5ax/L8k2TyB7QIAwIpMIoa/mQM/gnncxiS3TmC7AACwIpOI4S8nOXEZy5+Y5EsT2C4AAKzIJO4zfHGSX62qE1trH9vfglV1QpKTkrxhAttlHbh8586cfOq2JMmnr7gyj9wy5QEBACzDJM4Mn5WkJfnrqnrYYgtV1XFJ/jrJt5P86QS2yzqwtx2WjVu2ZuOWrbntW7dPezgAAMuy4jPDrbUrq+rVSV6Z5JNVdX6Si5JcOyzyfUken+TnkhyR5Ldaa1eudLsAALBSE3kcc2vt1VV1e5JXJPmFJP/nvEUqybeSvLS19juT2CYAAKzURGI4SVprr62qv0yyNcmjktx/mHVdkn9I8pbW2hcntT0AAFipicVwkgyx+4pJfiYAAKyWSXyBDgAAZpIYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuzXwMV9Uzq6oNP7+8yDInV9XFVbWnqm6qqo9U1bMP8LnPrqqPDsvvGdY/eXV+CwAApmGmY7iqvj/JWUlu2s8yz09yQZKHJ3lbkjcnOTrJuVW1fZF1tic5N8n9h+XfluRHklwwfB4AAIeAmY3hqqokb0myO8mbFlnmgUm2J7k+yYmttee11l6Y5EeT/HOSM6rqp+atc1KSM4b5P9pae2Fr7XlJThg+Z/vwuQAAzLiZjeEkpyd5XJJfSnLzIstsTXJEkrNaa9fMTWyt3ZDktcPb585bZ+79a4bl5ta5JsmfDJ/3SyscOwAA68BMxnBVPSzJ7yb549bapftZ9HHD63sXmPeeecusZB0AAGbQhmkPYLmqakOSv0jypSQvOcDiDx1er5o/o7V2XVXdnOSYqrpHa+2Wqjoqyfcluam1dt0Cn/e54fXYJY7144vMOm4p6wMAsLpmLoaT/FaSn0jyb1tr3zzAspuG1z2LzN+T5KhhuVuWuHyS3GtJIwUAYF2bqRiuqkdmdDb4da21D097PAfSWjthoenDGePj13g4AADMMzPXDA+XR/x5Rpc8vHyJq82dyd20yPz5Z4KXuvyNS9w+AADr2CydGd6YO67VvXV0Z7W7eHNVvTmjL9a9IMmVSe4zrHenM8lVdf+MLpG4trV2S5K01m6uqv+V5Puq6v4LXDf8Q8PrXa5B5q4u37kzJ5+6LUly9OZNOXvHgrd1BgCYmlmK4duS/D+LzDs+o+uI/yGjAJ4L34uSPCrJEzMvhpM8aWyZcRcl+cVhnbcscR0WsLcdlo1btiZJdl1yzpRHAwBwVzMTw8OX5RZ73PIrM4rht7bW/mxs1luS/Nckz6+qt8zda7iq7p077kQx/4Edb8oohl9aVe+au9fw8KCN52UU5fMjGQCAGTQzMXwwWmtXV9WLk+xI8rGqenuSvUlOSXJMFvgiXmvtQ1X1h0lelOSyqjo/yd2S/Kck353k18Yf4AEAwOw6pGM4SVprb6iqa5KcmeRZGX1p8DNJXtZae+si65xRVTszOhN8WpLvJPlEkj9orf3dmgwcAIBVd0jEcGvtlUleuZ/5FyS5YJmfeW6Sc1cwLAAA1rlDIoZZ/9xZAgBYj8Qwa8KdJQCA9WhmHroBAACTJoYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFsbpj0A+nba6Wdm1+49+94fvXlTzt6xfYojAgB6IoaZql2792Tjlq13vL/knCmOBgDojRhmzV2+c2dOPnVbkuTTV1yZR26Z8oAAgG6JYdbc3nbYvrPBt112xpRHAwD0zBfoAADolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuiWGAQDo1oZpDwDGXb5zZ04+dVuS5OjNm3L2ju1THhEAcCgTw6wre9th2bhla5Jk1yXnTHk0AMChTgyzbjlLDACsNjHMuuUsMQCw2nyBDgCAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOjWTMVwVW2uql+uqndW1eer6ptVtaeq/qGqnlNVC/4+VXVSVb27qq4f1rmsql5QVYfvZ1snV9XFw+ffVFUfqapnr95vBwDAWtsw7QEs088neWOS65J8IMmXknxvkp9N8mdJnlRVP99aa3MrVNXTkrwjya1J3p7k+iQ/neSPkjxq+Mw7qarnJ3lDkt1J3pZkb5JTkpxbVT/SWjtztX5BDuy008/Mrt179r0/evOmnL1j+xRHBADMqlmL4auSPDXJf2utfWduYlW9JMlHk/xcRmH8jmH6dyV5c5JvJ3lMa+1jw/SXJ7koySlV9YzW2nljn/XAJNsziuYTW2vXDNNfneSfkpxRVe9orX14dX9VFrNr955s3LL1jveXnDPF0QAAs2ymLpNorV3UWrtgPISH6V9J8qbh7WPGZp2S5L5JzpsL4WH5W5O8bHj7K/M2szXJEUnOmgvhYZ0bkrx2ePvclf0mLNflO3fm5FO35eRTt+XTV1w57eEAAIeIWTszvD/fGl5vH5v2uOH1vQssf2mSW5KcVFVHtNZuW8I675m3DGtkbzts39ng2y47Y8qjAQAOFYdEDFfVhiTPGt6OR+xDh9er5q/TWru9qq5O8sNJHpzkiiWsc11V3ZzkmKq6R2vtlgOM6+OLzDpuf+sBALA2Zuoyif343SQPT/Lu1tr7xqZvGl733HWVO02/10Gss2mR+QAAzIiZPzNcVacnOSPJZ5P84pSHcyettRMWmj6cMT5+jYcDAMA8M31meLgF2h8n+UySx7bWrp+3yIHO4s5Nv/Eg1lnszDEAADNiZmO4ql6Q0b2AL88ohL+ywGJztx04doH1NyR5UEZfuPvCEte5f5Kjklx7oOuFAQBY/2YyhqvqNzJ6aManMgrhry6y6EXD6xMXmPfoJPdI8qGxO0kcaJ0nzVsGAIAZNnMxPDww43eTfDzJ41trX9vP4ucn+VqSZ1TViWOfcWSS3x7evnHeOm9JcluS5w8P4Jhb595JXjK8fVMAAJh5M/UFuqp6dpJXZ/REub9PcnpVzV/smtbauUnSWvt6VW3LKIovrqrzMnqy3FMzuoXa+Rk9onmf1trVVfXiJDuSfKyq3p47Hsd8TJLXefocAMChYaZiOKNrfJPk8CQvWGSZS5KcO/emtfauqtqS5KUZPa75yCSfT/KiJDtaa23+B7TW3lBV1yQ5M6P7Fx+W0Zf0XtZae+skfhEAAKZvpmK4tfbKJK88iPU+mOTJy1zngiQXLHdbAADMjpm7ZhgAACZFDAMA0K2ZukwCFnL5zp05+dRtSZKjN2/K2Tu2T3lEAMCsEMPMvL3tsGzcsjVJsuuSc6Y8GgBglrhMAgCAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFsbpj0AWC2nnX5mdu3es+/90Zs35ewd26c4IgBgvRHDHLJ27d6TjVu23vH+knOmOBoAYD1ymQQAAN0SwwAAdEsMAwDQLTEMAEC3fIGOQ8rlO3fm5FO3JUk+fcWVeeSWKQ8IAFjXxDCHlL3tsH13kLjtsjOmPBoAYL0Tw3Rj/Kyxew4DAIkYpiPjZ43dcxgASHyBDgCAjolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOiWGAYAoFtiGACAbolhAAC6JYYBAOjWhmkPAKbh8p07c/Kp25IkR2/elLN3bJ/yiACAaRDDdGlvOywbt2xNkuy65JwpjwYAmBaXSQAA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C1PoKN7Hs0MAP0Sw3TPo5kBoF8ukwAAoFtiGACAbolhAAC65ZphGDP+Zbok+cJVn82Djz3uLn/2RTsAODSIYRgz/mW6JLnxsjP2vR//sy/aAcChwWUSAAB0SwwDANAtMQwAQLfEMAAA3fIFOjgI8+864e4SADCbxDAchPl3nXB3CQCYTS6TAACgW2IYAIBuiWEAALolhgEA6JYYBgCgW+4mARN22ulnZtfuPUnccg0A1jsxDBO2a/eefbddc8s1AFjfXCYBAEC3xDAAAN0SwwAAdEsMAwDQLV+ggwm4fOfOnHzqtiTJp6+4Mo/cMuUBAQBLIoZhAva2w/bdQeK2y86Y8mgAgKVymQQAAN0SwwAAdEsMAwDQLTEMAEC3fIEOVtH4XSaO3rwpZ+/YPuURAQDjxDCsovG7TOy65JwpjwYAmE8Mw5ScdvqZ2bV7TxJnjQFgWsQwTMmu3XucNQaAKfMFOgAAuiWGAQDolhgGAKBbYhgAgG6JYQAAuuVuErBGxh/AkSSfvuLKPHLLFAcEAIhhWCvjD+BIktsuO2OKowEAEpdJAADQMTEMAEC3XCYB68D49cRfuOqzefCxx+2bN/6oZo9wBoDJEsOwDoxfT3zjZWfc6driC8964b5Q/vQVV+aRz/29JB7hDACTIIZhnRsPZV+6A4DJcs0wAADdEsMAAHRLDAMA0C3XDMOMGr8Dxfw7S7jrBAAsjRiGGTX+xbrxO04kd77rxPg8YQwAdyaG4RCwv0c9j89zOzYAuDPXDAMA0C1nhoEkrjMGoE9iGDo1Hr/J2j3dTnQDsJ6IYejUrt17Fr3OeKl3qvjCVZ/Ng4897i7LzQ/t8Xnj23UNMwDTJoaBu1jqnSpuvOyMBcN2fmiLXgDWKzG8iKo6JsmrkzwxyeYk1yV5V5JXtdZumOLQ4KCNn/H99BVX5pFbDrzO/u5UsZrbnYb9ndEG4NAkhhdQVQ9J8qEk35Pkb5J8Nskjkvx6kidW1aNaa7unOEQ4KONhezBRuz/7C96lbnexSzAWuxxj0pzRBuiPGF7Yn2YUwqe31t4wN7Gq/jDJC5O8JslzpzQ2WJcOJrTHAzpZ/BKM8T/Pv2zjYOLYl/gAmCOG5xnOCj8hyTVJ/mTe7FckOS3JL1bVGa21m9d4eHBIOZhLMOavMx7HSz2DPH4GeHz99XwJBwCrQwzf1WOH1wtba98Zn9Fa+0ZVfTCjWP7JJO9f68EBdzYex4t9oW/B28htuev682N8f3fVWM+c+QZYumqtTXsM60pV/UGSM5Oc2Vp73QLzz0ryvCS/2lp74wE+6+OLzPqxu9/97oc/7GEPW/F4l+LzV38xh9/zPkmSb3z12tzze47ZN2/8/VL+POnlDoXxzOpnr7fxTPqzb/7ql3Pk3e+eJLn11tty1H2/b0Xj+fY3vpYffNAPJEm+9OVrs/f2bydJ7rbh8Dzg++9YZ3ze3ttuzd2OOPKAy82ft1Lj/8yPj3s9G98fyZ33yUr31Wp+NnBXn/jEHX8+/vi12N4n/t/W2qkHu74Ynqeqzk6yLcm21tqfLTD/NUlekuQlrbXfOcBnLRbDD09yU0aXYqy244bXz67Btg5F9t/K2H8rY/8dPPtuZey/lbH/Vma5+++zK4lhl0msotbaCdMew1yQr4exzCL7b2Xsv5Wx/w6efbcy9t/K2H8rs9b777C12MiMmbuwcNMi8+em37j6QwEAYDWJ4bu6cng9dpH5PzS8XrUGYwEAYBWJ4bv6wPD6hKq60/6pqnsmeVSSW5L841oPDACAyRLD87TW/jnJhUkemNFdI8a9KslRSf7CPYYBAGafL9At7Fczehzzjqp6fJIrkjwyo3sQX5XkpVMcGwAAE+LWaouoqu9P8uokT0yyOcl1Sd6Z5FWttRumOTYAACZDDAMA0C3XDAMA0C0xDABAt8QwAADdEsMAAHRLDAMA0C0xDABAt8TwIaqqjqmqc6pqV1XdVlXXVNXrq+re0x7bejfsq7bIz1emPb71oKpOqao3VNXfV9XXh33ztgOsc1JVvbuqrq+qb1bVZVX1gqo6fK3GvV4sZ/9V1QP3czy2qjpvrcc/TVW1uap+uareWVWfH46lPVX1D1X1nKpa8L9rjr+R5e4/x99dVdXvVdX7q+rLw/67vqo+WVWvqKrNi6zj+BssZ/+t1fHnPsOHoKp6SEZP0PueJH+T5LNJHpHRE/SuTPKo1tru6Y1wfauqa5LcK8nrF5h9U2tt+1qOZz2qqk8l+bEkNyW5NslxSf6ytfbMRZZ/WpJ3JLk1yduTXJ/kp5M8NMn5rbWfX4NhrxvL2X9V9cAkVyf5n0netcDHXd5aO3+1xrreVNVzk7wxowchfSDJl5J8b5KfTbIpo+Ps59vYf9wcf3dY7v5z/N1VVe1N8okkn0ny1SRHJfnJJCcm2ZXkJ1trXx5b3vE3Zjn7b82Ov9aan0PsJ8n7krQkvzZv+h8O09807TGu558k1yS5ZtrjWM8/Gf3F6oeSVJLHDMfV2xZZ9ruGf+HdluTEselHZvSXtpbkGdP+ndbx/nvgMP/caY97PfwkeVxGIXHYvOn3yyjsWpKfG5vu+FvZ/nP83XUfHrnI9NcM++pPx6Y5/la2/9bk+HOZxCFmOCv8hIyC7k/mzX5FkpuT/GJVHbXGQ+MQ0lr7QGvtc234t9UBnJLkvknOa619bOwzbk3ysuHtr6zCMNetZe4/xrTWLmqtXdBa+8686V9J8qbh7WPGZjn+xhzE/mOe4dhZyF8Nrz80Ns3xN88y99+a2LDWG2TVPXZ4vXCBf9l9o6o+mFEs/2SS96/14GbIEVX1zCQPyOgvEJclubS19u3pDmsmPW54fe8C8y5NckuSk6rqiNbabWs3rJlzdFX9lySbk+xO8uHW2mVTHtN6863h9faxaY6/pVto/81x/B3YTw+v4/vF8bd0C+2/Oat6/InhQ89Dh9erFpn/uYxi+NiI4f25X5K/mDft6qr6pdbaJdMY0Axb9Jhsrd1eVVcn+eEkD05yxVoObMb8h+Fnn6q6OMmzW2tfmsqI1pGq2pDkWcPb8fBw/C3BfvbfHMffPFV1ZpKNGV1rfWKSf5tRyP3u2GKOv0Uscf/NWdXjz2USh55Nw+ueRebPTb/X6g9lZr0lyeMzCuKjkvxIkv87o2uX3lNVPza9oc0kx+TK3JLk/0pyQpJ7Dz9bMvry02OSvN9lT0lG/wF9eJJ3t9beNzbd8bc0i+0/x9/izszo8sMXZBRy703yhNbav44t4/hb3FL235ocf2IY5mmtvWq4ru5fWmu3tNYub609N6MvIN49ySunO0J60lr7amvtt1prn2it3Tj8XJrR/+H5SJIfTPLL0x3ldFXV6UnOyOjOOb845eHMnP3tP8ff4lpr92utVUYnTn42o7O7n6yq46c7stmwlP23VsefGD70zP0tc9Mi8+em37j6QznkzH255NFTHcXscUyugtba7Un+bHjb7TFZVc9P8scZ3abpsa216+ct4vjbjyXsvwU5/u4wnDh5Z0aBtjnJn4/NdvwdwAH232LrTPT4E8OHniuH12MXmT/3Lc3FrilmcXP/66bX/yV4sBY9JofrFB+U0Rd2vrCWgzpEdH1MVtULkrwhyeUZhdxCD8Vx/C1iiftvf7o+/uZrrX0xo79U/HBV3WeY7PhbokX23/5M7PgTw4eeDwyvT1jgSUL3TPKojK7B+ce1Htgh4CeH1+7/pbVMFw2vT1xg3qOT3CPJh3yT+qB0e0xW1W8k+aMkn8oo5L66yKKOvwUsY//tT7fH334cPbzO3XnI8bc88/ff/kzs+BPDh5jW2j8nuTCjL3s9b97sV2X0N6i/aK3dvMZDmwlV9bCFLsYfnoJz1vB2v48d5i7OT/K1JM+oqhPnJlbVkUl+e3j7xmkMbBZU1fELPWK4qh6f5IXD266Oyap6eUZf+Pp4kse31r62n8Udf/MsZ/85/u6sqo6tqrtc8lBVh1XVazJ68uuHWms3DLMcf2OWu//W6vjzOOZD0AKPY74iySMzugfxVUlOah7HvKCqemVGXyS5NMkXk3wjyUOSPCWjJwa9O8nPtNb2TmuM60FVPT3J04e390vyHzP62/nfD9O+1lo7c97y52f0ONLzMnoc6VMzPI40yf/R0wMolrP/htsH/VBG/0xfO8z/0dxx/9KXt9bm/qN6yKuqZyc5N6MzR2/Iwt/Sv6a1du7YOk+P4y/J8vef4+/OhktLfifJP2T0mODdGT3OektGXwD7SkZ/wfjM2DpPj+MvyfL335odf0t5TJ2f2ftJ8v0Z3SLsuiR7Mwq71ye597THtp5/hn8g/7+MvlV9Y0Y3of/XJP89o3tw1rTHuB5+MrqjRtvPzzULrPOojP4ycUOSbybZmdHf7A+f9u+znvdfkuck+buMnip5U0aPdf1Skrcn+XfT/l3W4b5rSS52/E1m/zn+7rL/Hp7R/yX8VEZnfG/P6C8U/zTs2+9eZD3H30Hsv7U6/pwZBgCgW64ZBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW2IYAIBuiWEAALolhgEA6JYYBgCgW/8/VOlh6t0GUGkAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "execution_count": 22, "metadata": { "image/png": { "height": 352, "width": 353 }, "needs_background": "light" }, "output_type": "execute_result" } ], "source": [ "#plotting the results and showing our effect size in blue\n", "sns.displot(hist_data)\n", "plt.axvline(chi_squared_expected, color=\"blue\")" ] }, { "cell_type": "code", "execution_count": 24, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "0.0" ] }, "execution_count": 24, "metadata": { }, "output_type": "execute_result" } ], "source": [ "#calculation of our p value\n", "a=np.sum(hist_data>=chi_squared_expected)\n", "pval=a/10000\n", "pval" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#The results of our analysis showed that the two proportions are statistically different since we reject the null hypothesis based on our p value (0<0.05). Thus, we can say that the proportion of individuals who got systemic side effects and those who did not from the 1st dose of the Pfizer vaccine is statistically different from the proportion of individuals who systemic side effects and those who did not from the 1st dose of the Astrazeneca vaccine. Conclusion: prevalence of systemic side effects is not the same between individuals who received the 1st dose of Pfizer and those who received the 1st dose of Astrazeneca." ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#3 Find a bootstrap pivotal confidence interval for the relative risk of your outcome. What could be the practical significance of your findings?" ] }, { "cell_type": "code", "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "0.4007703595011005" ] }, "execution_count": 25, "metadata": { }, "output_type": "execute_result" } ], "source": [ "#relative risk of our observed data\n", "og_relative=(reality[0,0]/(reality[0,0]+reality[1,0]))/(reality[0,1]/(reality[0,1]+reality[1,1]))\n", "og_relative" ] }, { "cell_type": "code", "execution_count": 26, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "array([[ 38, 116],\n", " [244, 229]])" ] }, "execution_count": 26, "metadata": { }, "output_type": "execute_result" } ], "source": [ "#reminder of our observed table\n", "reality" ] }, { "cell_type": "code", "execution_count": 27, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#again making a list of outcomes, creating an empty table for the relative risk calculation, and creating an empty list to store our risk calculations\n", "#making a for loop to simulate relative risk calculations and storing it to the empty list\n", "pfizerlist=(38*['S'])+(244*['N'])\n", "astralist=(116*['S'])+(229*['N'])\n", "sim_chi=np.zeros([2,2])\n", "results=np.zeros(10000)\n", "for i in range (10000):\n", " sim_pfizer_rr=np.random.choice(pfizerlist,pfizer)\n", " sim_chi[0,0]=np.sum(sim_pfizer_rr=='S')\n", " sim_chi[1,0]=np.sum(sim_pfizer_rr=='N')\n", " sim_astra_rr=np.random.choice(astralist,astra)\n", " sim_chi[0,1]=np.sum(sim_astra_rr=='S')\n", " sim_chi[1,1]=np.sum(sim_astra_rr=='N')\n", " sim=(sim_chi[0,0]/(sim_chi[0,0]+sim_chi[1,0]))/(sim_chi[0,1]/(sim_chi[0,1]+sim_chi[1,1]))\n", " results[i]=sim" ] }, { "cell_type": "code", "execution_count": 28, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#storing the results and sorting them out to compute the upper and lower bounds of the confidence interval\n", "results.sort()\n", "upper=(2*og_relative)-results[49]\n", "lower=(2*og_relative)-results[9949]" ] }, { "cell_type": "code", "execution_count": 32, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "[0.2018327507085006, 0.4007703595011005, 0.5568598679383712]" ] }, "execution_count": 32, "metadata": { }, "output_type": "execute_result" } ], "source": [ "#this is the upper and lower bounds of the relative risk confidence interval\n", "[lower,og_relative,upper]" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "#Conclusion: Having systemic side effects from a first dose is not independent of whether the individual received either the Pfizer or Astrazeneca vaccine. The relative risk is 0.401 with a 99% CI of [0.20,0.56]. Thus, it is 0.40 times more likely (or 60% less likely) to have systemic side effects from the Pfizer first dose compared to the Astrazeneca first dose." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (system-wide)", "language": "python", "metadata": { "cocalc": { "description": "Python 3 programming language", "priority": 100, "url": "https://www.python.org/" } }, "name": "python3", "resource_dir": "/ext/jupyter/kernels/python3" }, "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" } }, "nbformat": 4, "nbformat_minor": 4 }