inferential_statistics_2-Q.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Inferential statistics II - Bootstrapping"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the previous frequentist mini-projects, you did frequentist calculations to perform inference from a sample of data. Such inference relies on theory largely developed from the 19th-Century onwards that is subject to certain assumptions or theoretical limits. These are fine if those assumptions hold for the particular case you're working on, and what you want to do has a known theoretical distribution (for example the mean of a sampling distribution that we looked at in the previous mini-project.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this mini-project, you'll use the same medical charge data you used in the frequentist inference mini-project, but this time you'll make inferences about the population using bootstrapping (ie. simulating repeated re-runs of an experiment.) If frequentism is about using assumptions and theoretical results to calculate what we expect to happen were an experiment to be run again and again and again, then bootstrapping is about using computing power to essentially re-run the sample draw again and again and again to see what actually happens."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Prerequisites"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "While these exercises do not strictly depend on these concepts, we encourage you to complete the previous mini-projects before starting this one so that you can approach this assignment with a good understanding of frequentist concepts like:\n",
    "* the _z_-statistic\n",
    "* the _t_-statistic\n",
    "* the difference and relationship between the two\n",
    "* the Central Limit Theorem, its assumptions and consequences\n",
    "* how to estimate the population mean and standard deviation from a sample\n",
    "* the concept of a sampling distribution of a test statistic, particularly for the mean\n",
    "* how to combine these concepts to calculate confidence intervals and p-values\n",
    "* how those confidence intervals and p-values allow you to perform hypothesis (or A/B) tests"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To complete mini-project, it's important that you first complete the bootstrap resources listed in this subunit, as they contain valuable information about how to calculate bootstrap replicates of summary statistics. Having an basic understanding of what confidence intervals and p-values are will also be helpful (we touch on them in this mini-project, but please speak to your mentor or conduct individual research if you'd like to learn more.) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "from numpy.random import seed\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Medical charge data set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "med_charges = pd.read_csv('data/insurance2.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>age</th>\n",
       "      <th>sex</th>\n",
       "      <th>bmi</th>\n",
       "      <th>children</th>\n",
       "      <th>smoker</th>\n",
       "      <th>region</th>\n",
       "      <th>charges</th>\n",
       "      <th>insuranceclaim</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>19</td>\n",
       "      <td>0</td>\n",
       "      <td>27.900</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>3</td>\n",
       "      <td>16884.92400</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>18</td>\n",
       "      <td>1</td>\n",
       "      <td>33.770</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>1725.55230</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>28</td>\n",
       "      <td>1</td>\n",
       "      <td>33.000</td>\n",
       "      <td>3</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "      <td>4449.46200</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>33</td>\n",
       "      <td>1</td>\n",
       "      <td>22.705</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>21984.47061</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>32</td>\n",
       "      <td>1</td>\n",
       "      <td>28.880</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>3866.85520</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   age  sex     bmi  children  smoker  region      charges  insuranceclaim\n",
       "0   19    0  27.900         0       1       3  16884.92400               1\n",
       "1   18    1  33.770         1       0       2   1725.55230               1\n",
       "2   28    1  33.000         3       0       2   4449.46200               0\n",
       "3   33    1  22.705         0       0       1  21984.47061               0\n",
       "4   32    1  28.880         0       0       1   3866.85520               1"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "med_charges.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the previous assignment, you used the frequentist approach to estimate the lower limit for the 95% confidence interval on the mean hospital charge. This approach relies on statistical theory that has been developed over the years and is also limited to statistics for which theoretical results on the sampling distribution exist. These results are remarkably useful and applicable much of the time and under a surprisingly wide range of conditions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Having calculated the 95% lower confidence interval using frequentist theory in the previous exercise, you'll now use bootstrap inference to verify your calculations and check that you get consistent results without making the assumptions required before. After all, the distribution of charges really was very non-normal."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Q:__ Use bootstrap sampling to estimate the same 95% confidence interval lower limit as before."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "12723.822458375787"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.random.seed(47)\n",
    "N_rep = 10000\n",
    "mreplicate = np.empty ( N_rep )\n",
    "for i in range ( N_rep ):\n",
    "    samples = np.random.choice ( med_charges['charges'] , len ( med_charges['charges'] ) )\n",
    "    mreplicate[i]  = np.mean( samples )\n",
    "\n",
    "mean , std = np.mean ( mreplicate ) , np.std ( mreplicate )\n",
    "lowerbound = mean - 1.64 * std \n",
    "lowerbound \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you performed 10000 replicates immediately after setting the random seed to 47, you should get the value 12724 here, which compares very well with the value 12725 obtained using the _t_-distribution confidence interval previously. It is a most pleasant result to see the predictions of classical frequentist theory match with results that are now possible through the number-crunching ability of computers."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Remember, in the previous mini-projects, we saw that there are two ways of performing a _t_-test from a sample, depending on whether we can assume the groups have equal variance or not. We can actually easily test this using the bootstrap approach!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Q:__ Calculate the 95% confidence interval for the difference between the standard deviations of insurance and non-insurance claim charges (insured - non-insured). Calculate the differences over 10000 replicates. Plot the histogram of values and mark the locations of the percentiles. State the null and alternative hypothesis and comment on whether you would retain or reject the null hypothesis in this case and why."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "insured = med_charges.charges[med_charges.insuranceclaim == 1.0]  \n",
    "uninsured = med_charges.charges[med_charges.insuranceclaim == 0.0]\n",
    "insured_number = len(insured)\n",
    "uninsured_number = len(uninsured)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[6689.776011985517, 8508.067362244858]\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(47)\n",
    "std_diff_replicate = np.empty(N_rep)\n",
    "for i in range(N_rep):\n",
    "    insured_samples = np.random.choice(insured, insured_number)\n",
    "    uninsured_samples = np.random.choice(uninsured, uninsured_number)\n",
    "    std_diff_replicate[i] = np.std(insured_samples) - np.std(uninsured_samples)\n",
    "\n",
    "std_diff_mean, std_diff_std = np.mean(std_diff_replicate) , np.std(std_diff_replicate)\n",
    "conf_interval = [std_diff_mean - 1.96*std_diff_std , std_diff_mean + 1.96*std_diff_std]\n",
    "print(conf_interval)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dfZxcVZ3n8c+XhAQCykMIDCbB4BhQhtHIZBBRtDGIkAGja5CgYFA0LoIjahYIMkNWlFGMA7gomhUkoMuDUYeI+ACRCK4ECBAgEJAWdkhLhEgQjZFkkvz2j3OaFJ3qx6rqul33+3697qtunXvq1q9un/7VqXNvnVJEYGZm5bFdswMwM7PB5cRvZlYyTvxmZiXjxG9mVjJO/GZmJTO82QH0ZI899ogJEyY0O4zW9+ij6Xb//Zsbx1Dh49V3PlZNcc899/whIsZ0t73QiX/ChAksW7as2WG0vra2dLtkSTOjGDp8vPrOx6opJP1nT9s91GNmVjKF7vHbIDn33GZHYK3KbauQnPgNjjii2RFYq3LbKiQP9RgsX54Ws3pz2yok9/gNzjgj3foEnNWb21YhucdvZlYyTvxmZiXjxG9mVjJO/GZmJeOTuwYXXNDsCArt2XUbWL9x84v3x2zaAsCateur1h81Yhijdx45KLEVnttWITnxGxx6aLMjKLT1Gzdz2IW3vnj/2iefA2BGRVml2888nNGDEtkQ4LZVSB7qMfj1r9NiVm9uW4XkHr/BOeekW19rbfXmtlVI7vGbmZWME7+ZWck48ZuZlYwTv5lZyfjkrsHFFzc7AmtVbluF5MRvMGlSsyOwVuW2VUge6jG45Za0mNWb21Yhucdv8PnPp1v/WpLVm9tWIbnHb2ZWMk78ZmYl48RvZlYyvSZ+SVdIekbSioqyL0t6RNIDkn4oadeKbXMktUt6VNI7K8qPymXtks6u/0sxM7O+6MvJ3SuBS4GrKspuBuZExCZJXwLmAGdJOgCYAfwd8ArgFkn75cd8DXgH0AHcLWlRRDxcn5dhNfnmN5sdgbUqt61C6jXxR8RtkiZ0Kft5xd2lwPS8Pg24NiI2AE9IagcOztvaI+JxAEnX5rpO/EWw//7NjmBQdf1hld5siWhgNC2uZG1rqKjH5ZwfBq7L62NJbwSdOnIZwKou5W+stjNJs4BZAPvss08dwrNe/ehH6fbYY5sbxyDp+sMqvVkyu61xwbS6krWtoaKmxC/ps8Am4LudRVWqBdXPJVTtRkXEfGA+wOTJk93VGgxf+Uq69T+n1ZvbViENOPFLmgkcA0yJePGzcAcwvqLaOOCpvN5duZmZDaIBXc4p6SjgLOBdEVH5i9OLgBmSRkraF5gI3AXcDUyUtK+kEaQTwItqC93MzAai1x6/pGuANmAPSR3AeaSreEYCN0sCWBoR/z0iHpJ0Pemk7SbgtIjYnPdzOvAzYBhwRUQ81IDXY2ZmvejLVT0nVCm+vIf6XwC+UKX8JuCmfkVnNgQN306sWru+94rZqBHDGL3zyAZGZPZSnqTN4Oqrmx1BS9mwaQtt85b0uf7tZx7O6MaF01xuW4XkxG8wfnzvdcwGwm2rkDxXj8F116XFrN7ctgrJPX6Dyy5Lt8cf39w4rPW4bRWSe/xmZiXjxG9mVjJO/GZmJePEb2ZWMj65a7BwYbMjsFbltlVITvwGe+zR7AisVbltFZKHegyuvDItZvXmtlVITvzmf05rHLetQnLiNzMrGSd+M7OSceI3MysZJ34zs5Lx5ZwGN/n3caxB3LYKyYnfYNSoZkdgrcptq5A81GPw9a+nxaze3LYKyYnf4Prr02JWb25bheTEb2ZWMh7jtyHt2XUbWL9xc78esyWiQdGYDQ1O/Dakrd+4mcMuvLVfj1kyu60xwZgNEb0O9Ui6QtIzklZUlO0u6WZJj+Xb3XK5JH1VUrukByQdVPGYmbn+Y5JmNublmJlZb/oyxn8lcFSXsrOBxRExEVic7wMcDUzMyyzgMkhvFMB5wBuBg4HzOt8srACWLEmLWb25bRVSr4k/Im4D1nYpngYsyOsLgHdXlF8VyVJgV0l7A+8Ebo6ItRHxHHAz276ZmJnZIBjoVT17RcRqgHy7Zy4fC6yqqNeRy7or34akWZKWSVq2Zs2aAYZn/TJvXlrM6s1tq5DqfTmnqpRFD+XbFkbMj4jJETF5zJgxdQ3OunHjjWkxqze3rUIaaOJ/Og/hkG+fyeUdwPiKeuOAp3ooNzOzQTbQyzkXATOBL+bbGyrKT5d0LelE7vMRsVrSz4ALKk7oHgnMGXjYZq1j+HZi1dr1fa4/asQwRu88soERWavrNfFLugZoA/aQ1EG6OueLwPWSTgGeBI7L1W8CpgLtwHrgQwARsVbS+cDdud7nIqLrCWOzUtqwaQtt85b0uf7tZx7O6MaFYyXQa+KPiBO62TSlSt0ATutmP1cAV/QrOhscO+7Y7AisVbltFZK/uWvwk580OwJrVW5bheRJ2szMSsaJ3+D889NiVm9uW4XkxG+weHFazOrNbauQnPjNzErGid/MrGSc+M3MSsaXcxqM9teBrEHctgrJid/g+99vdgTWqty2CslDPWZmJePEbzBnTlrM6s1tq5A81GNwxx3NjsBaldtWIbnHb2ZWMk78ZmYl48RvZlYyHuM3GDeu2RFYq3LbKiQnfoPvfKfZEVg/DKmfanTbKiQnfrMhxj/VaLXyGL/BGWekxaze3LYKyT1+g+XLmx2BtSq3rUJyj9/MrGSc+M3MSqamxC/pU5IekrRC0jWSdpC0r6Q7JT0m6TpJI3Ldkfl+e94+oR4vwMzM+mfAY/ySxgL/DBwQEX+VdD0wA5gKXBQR10r6BnAKcFm+fS4iXi1pBvAl4PiaX4HVbr/9mh3Bi55dt4H1Gzf3uf6WiAZGYzUrUNuyrWo9uTsc2FHSfwGjgNXA24H35+0LgLmkxD8trwMsBC6VpAj/5zbd/PnNjuBF6zdu5rALb+1z/SWz2xoXjNWuQG3LthrwUE9E/A6YBzxJSvjPA/cAf4yITblaBzA2r48FVuXHbsr1t7m8WNIsScskLVuzZs1AwzMzs24MOPFL2o3Ui98XeAWwE3B0laqdPXr1sG1rQcT8iJgcEZPHjBkz0PCsP2bNSotZvbltFVItQz1HAE9ExBoAST8ADgV2lTQ89+rHAU/l+h3AeKBD0nBgF2BtDc9v9fKb3zQ7AmtVbluFVMtVPU8Ch0gaJUnAFOBh4FZgeq4zE7ghry/K98nbf+HxfTOzwVfLGP+dpJO09wIP5n3NB84CPi2pnTSGf3l+yOXA6Fz+aeDsGuI2M7MBqumqnog4DzivS/HjwMFV6r4AHFfL85mZWe08V4/BpEnNjsBaldtWITnxG1x8cbMjsFbltlVInqvHzKxknPgNTjwxLWb15rZVSB7qMejoaHYE1qrctgrJPX4zs5Jx4jczKxknfjOzkvEYv8Gb3tTsCKxVuW0VkhO/wb/9W7MjsFbltlVIHuoxMysZJ36D9743LWb15rZVSB7qMXj22WZHYK3KbauQ3OM3MysZJ34zs5Jx4jczKxmP8RtMmdLsCKxVuW0VkhO/wb/8S7MjsFbltlVIHuoxMysZJ36Do49Oi1m9uW0Vkod6DP7612ZHYK3KbauQ3OM3MyuZmhK/pF0lLZT0iKSVkt4kaXdJN0t6LN/ulutK0lcltUt6QNJB9XkJZmbWH7X2+C8BfhoRrwFeD6wEzgYWR8REYHG+D3A0MDEvs4DLanxuMzMbgAGP8Ut6OfBW4GSAiNgIbJQ0DWjL1RYAS4CzgGnAVRERwNL8aWHviFg94OitPo45ptkRWKty2yqkWk7uvgpYA3xb0uuBe4BPAnt1JvOIWC1pz1x/LLCq4vEducyJv9lmz252BNaq3LYKqZahnuHAQcBlEfEG4C9sHdapRlXKYptK0ixJyyQtW7NmTQ3hmZlZNbUk/g6gIyLuzPcXkt4Inpa0N0C+faai/viKx48Dnuq604iYHxGTI2LymDFjagjP+qytLS1m9ea2VUgDTvwR8XtglaT9c9EU4GFgETAzl80Ebsjri4AP5qt7DgGe9/i+mdngq/ULXJ8AvitpBPA48CHSm8n1kk4BngSOy3VvAqYC7cD6XNda3LPrNrB+4+Y+198S24z+mVmd1ZT4I2I5MLnKpm2m5MtX85xWy/PZ0LN+42YOu/DWPtdfMrutccGYGeBv7pqZlY7n6jF43/uaHYG1KretQnLiN/j4x5sdgbUqt61C8lCPwfr1aTGrN7etQnKP32Dq1HS7ZElTw7DGGL6dWLW2f8l31IhhjN55ZO1P7rZVSE78Zi1uw6YttM1b0q/H3H7m4YxuTDhWAB7qMTMrGSd+M7OSceI3MysZj/EbnHxysyOwVuW2VUhO/OZ/Tmsct61C8lCPwR/+kBazenPbKiT3+A2mT0+3vtba6s1tq5Dc4zczKxknfjOzknHiNzMrGSd+M7OS8cldg1NPbXYE1qrctgrJid/g+OObHYG1KretQvJQj8GqVWkxqze3rUJyj9/gpJPSra+1tnpz2yok9/jNzErGid/MrGRqTvyShkm6T9KN+f6+ku6U9Jik6ySNyOUj8/32vH1Crc9tZmb9V48e/yeBlRX3vwRcFBETgeeAU3L5KcBzEfFq4KJcz8zMBllNJ3cljQP+CfgC8GlJAt4OvD9XWQDMBS4DpuV1gIXApZIUEVFLDFYHn/lMsyOwVuW2VUi1XtVzMXAm8LJ8fzTwx4jYlO93AGPz+lhgFUBEbJL0fK7/kjlbJc0CZgHss88+NYZnfXLssc2OwFqV21YhDXioR9IxwDMRcU9lcZWq0YdtWwsi5kfE5IiYPGbMmIGGZ/3x6KNpMas3t61CqqXH/2bgXZKmAjsALyd9AthV0vDc6x8HPJXrdwDjgQ5Jw4FdgLU1PL/Vy8c+lm59rbXVm9tWIQ048UfEHGAOgKQ2YHZEfEDS94DpwLXATOCG/JBF+f4defsvPL4/9Dy7bgPrN27uc/0t/hObFU4jvrl7FnCtpM8D9wGX5/LLgasltZN6+jMa8NzWYOs3buawC2/tc/0ls9saF4yZDUhdEn9ELAGW5PXHgYOr1HkBOK4ez2dmZgPnb+6amZWMJ2kzOPfcZkdgrcptq5Cc+A2OOKLZEVirctsqJA/1GCxfnhazenPbKiT3+A3OOCPd+lprqze3rUJyj9/MrGSc+M3MSsaJ38ysZJz4zcxKxid3DS64oNkRWMEM306sWru+z/VHjRjG6J1HbrvBbauQnPgNDj202RFYwWzYtIW2eUv6XP/2Mw9ndLUNbluF5KEeg1//Oi1m9ea2VUju8Rucc0669bXWVm9uW4XkHr+ZWck48ZuZlYwTv5lZyTjxm5mVjE/uGlx8cbMjsFbltlVITvwGkyY1OwJrVW5bheShHoNbbkmLWb25bRWSe/wGn/98uvWvJVm9uW0VkhN/yT27bgM7bdoCwJo+zM2yJaLRIZlZgw048UsaD1wF/A2wBZgfEZdI2h24DpgA/D/gfRHxnCQBlwBTgfXAyRFxb23hW63Wb9zMY08+B8CMC2/ttf6S2W0NjsjMGq2WMf5NwGci4rXAIcBpkg4AzgYWR8REYHG+D3A0MDEvs4DLanhuMzMboAH3+CNiNbA6r/9Z0kpgLDANaMvVFgBLgLNy+VUREcBSSbtK2jvvx8yGsO6mcR7TzTBit9M426Coyxi/pAnAG4A7gb06k3lErJa0Z642FlhV8bCOXPaSxC9pFukTAfvss089wrNenPPO05sdgg1x3U3j/KrXfgCAx7sMI3Y7jbMNipoTv6Sdge8DZ0TEn9JQfvWqVcq2OVMYEfOB+QCTJ0/2mcRB8Pjocc0OwVqU21Yx1XQdv6TtSUn/uxHxg1z8tKS98/a9gWdyeQcwvuLh44Cnanl+q48p7Xcypf3OZodhLchtq5gGnPjzVTqXAysj4t8rNi0CZub1mcANFeUfVHII8LzH94vho3f9kI/e9cNmh2EtyG2rmGoZ6nkzcBLwoKTluewc4IvA9ZJOAZ4EjsvbbiJdytlOupzzQzU8t5mZDVAtV/X8iurj9gBTqtQP4LSBPp+ZmdWH5+oxMysZJ34zs5LxXD3Gp475TLNDsBbltlVMTvzG6pePaXYI1qLctorJQz3GMStv45iVtzU7DGtBblvF5B6/ceJ9NwFw42vf2uRIrNW4bRWTe/xmZiXjxG9mVjIe6mkxz67bwPqNm/tc37+oZVY+TvwtZv3GzRzWh1/S6uRf1DIrHyd+49R3z2l2CNai3LaKyYnfeG7ULs0OwVqU21YxOfEb0x+8BYCFf39EkyOxVtNd2+rupxq7459qrC8nfnPit4bprm1191ON3fFPNdaXL+c0MysZJ34zs5LxUE+B9feafPB1+WbWOyf+AuvvNfng6/KtNflkcH058RsnHze32SFYi6pX2/LJ4Ppy4jde2H6HZodgLcptq5h8ctc48d4fc+K9P252GNaC3LaKyYnfOOaR2znmkdubHYa1ILetYvJQzyDyzJlmVgSDnvglHQVcAgwDvhURXxzsGOplIIn8bV9e0uf6vkLHbGD6exUQlOtKoEFN/JKGAV8D3gF0AHdLWhQRDw9mHPXiKZDNiqm/VwFBua4EGuwe/8FAe0Q8DiDpWmAaUIjE76EYs/Lq76eEEcO2Y+PmLQ2r38hPIIpBTF6SpgNHRcRH8v2TgDdGxOkVdWYBs/Ld/YFHBy3A7u0B/KHZQfTBUIjTMdbPUIjTMdZPf+J8ZUSM6W7jYPf4VaXsJe88ETEfmD844fSNpGURMbnZcfRmKMTpGOtnKMTpGOunnnEO9uWcHcD4ivvjgKcGOQYzs1Ib7MR/NzBR0r6SRgAzgEWDHIOZWakN6lBPRGySdDrwM9LlnFdExEODGcMAFWroqQdDIU7HWD9DIU7HWD91i3NQT+6amVnzecoGM7OSceI3MyuZUid+SbtKWijpEUkrJb1J0u6Sbpb0WL7dLdeVpK9Kapf0gKSDKvYzM9d/TNLMQYhxrqTfSVqel6kV9efkGB+V9M6K8qNyWbuks+sY3/4VcSyX9CdJZxTwOHYXZ2GOZd73pyQ9JGmFpGsk7ZAvhrgzH5fr8oURSBqZ77fn7RN6i72BMV4p6YmK4zgp123W3/uTOb6HJJ2RywrVJnuIs/FtMiJKuwALgI/k9RHArsCFwNm57GzgS3l9KvAT0ncRDgHuzOW7A4/n293y+m4NjnEuMLtK3QOA+4GRwL7Ab0kn0Yfl9VflfdwPHNCA4zkM+D3wyqIdxx7iLMyxBMYCTwA75vvXAyfn2xm57BvAqXn948A38voM4LqeYm9wjFcC06vUH/S/N3AgsAIYRbqA5RZgYtHaZA9xNrxNlrbHL+nlwFuBywEiYmNE/JE0hcSCXG0B8O68Pg24KpKlwK6S9gbeCdwcEWsj4jngZuCoBsfYnWnAtRGxISKeANpJ02S8OFVGRGwEOqfKqLcpwG8j4j8p0HHsJc7uNOtYDgd2lDSclBBWA28HFubtXY9l5zFeCEyRpB5ib1SMPX0Xpxl/79cCSyNifURsAn4JvIfitcnu4uxO3dpkaRM/6d1xDfBtSfdJ+paknYC9ImI1QL7dM9cfC6yqeHxHLuuuvJExApyeP5Ze0fmRtUkxVpoBXJPXi3Qce4oTCnIsI+J3wDzgSVLCfx64B/hjTgxdn+/FWPL254HRgx1jRPw8b/5CPo4XSeqcZKYZf+8VwFsljZY0itSjH0/x2mR3cUKD22SZE/9w4CDgsoh4A/AX0se/7nQ33USv01DUoLsYLwP+FphE+uf7ShNjTE+cxp3fBXyvt6rdxNLwGKFqnIU5lvkffBrpY/wrgJ2Ao3t4vkLEKOlEYA7wGuAfSUMjZzUrxohYCXyJ1EP/KWnoY1MPD2lKm+whzoa3yTIn/g6gIyLuzPcXkpLs0/ljHvn2mYr61aabaOQ0FFVjjIinI2JzRGwB/jdbP8Y3I8ZORwP3RsTT+X6RjmO3cRbsWB4BPBERayLiv4AfAIeShh46v2xZ+XwvxpK37wKsbUaMEbE6D5VsAL5Nk9tkRFweEQdFxFtJx+QxCtgmq8U5GG2ytIk/In4PrJK0fy6aQpoeehHQefZ+JnBDXl8EfDBfAXAI6SPuatK3kI+UtFvuDR2ZyxoWY2fjzd5D+sjYGeOMfLXHvqQTRXcxOFNlnMBLh08Kcxx7irNgx/JJ4BBJo/JYfWebvBWYnut0PZadx3g68ItIZwG7i71RMa6sSKgijZ1XHsdB/3tL2jPf7gP8N9LfvHBtslqcg9Im+3MWutUW0kepZcADwH+QztyPBhaTegiLgd1zXZF+ROa3wIPA5Ir9fJh0oqUd+NAgxHh1juGB/Afeu6L+Z3OMjwJHV5RPBX6Tt322zjGOAp4FdqkoK9Rx7CHOoh3L/wk8kv/ZryZdwfGq/A/eThqiGpnr7pDvt+ftr+ot9gbG+It8HFcA3wF2bvL/ze2kN837gSkFbpPV4mx4m/SUDWZmJVPaoR4zs7Jy4jczKxknfjOzknHiNzMrGSd+M7OSceIvAEkh6eqK+8MlrZF0Y5PieU2eFfA+SX/bZduHJT2Yv06+QtI0SV/L9R+W9FdtnVVwurbO2ni/pN9IukpSr197l9TW+folvUt5xkFJY5RmorxP0mGSjlOatfTWxhyNxsnxP5SP1Y5Vtr8nt43XVJRNkLQir794jKo89mBJtynN2PiI0nQfoySdLOnSXuJ68XgP8HVtL+mLSjNarpB0l6Sj87Z1A92v1c+g/vSidesvwIGSdoyIvwLvAH7XxHjeDdwQEedVFkoaR7qO+KCIeF7SzsCYiLghb58A3BgRkyoecwzwPyJiYf5yzxnArZIOjDShVK8iYhFbv5AyBXgkImbm/f8U+HhE9CnxSxoWEZv7UncQfACYFxHf7mb7CcCvSF/ImdvXnUrai3R9/4yIuCMf9/cCL+vL47sc74E4H9gbODAiNuR43lbD/oDUIYqtcxZZDdzjL46fAP+U17t+u3Qnpcma7s493Wm5fIKk2yXdm5dDc3mbpCXaOo//d/M//0tImiRpae69/zB/Q3EqKTl/pEovek/gz8A6gIhYF2mWwD6J5CLSlMjbzEGjNKf4I5J+RfoWY2f5yZIuVZrj/UJgau4lnwe8BfiGpC9LGpZv786v6WMVx+NWSf+H9MUYJJ2Ye6LLJX1T0rBcvk7SF/InlKU5aSFpr3yM7s/LoT3tp8vrmpL/bg/mv+NISR8B3gf8q6TvVnnMzsCbgVNIib8/TgMWRMQdFcd9YWydTqPzOY6t+PR0S8VrffFTgdIntsvy8Xtc0tvya1gp6coqcY8CPgp8ItL0DUSaguD6ijrVjm93scyVNF/Sz4Gr8qeW6/Pf97r8mMm57pGS7sj/C9/LxxClTx8P58fM6+exbE31/iaalwF9e28d8DrSXDw7AMuBNlLvGeAC4MS8vivpG3o7kb6JukMunwgsy+ttpJkax5He3O8A3lLleR8A3pbXPwdcnNfnUn0+8GGkr6w/SZqP5dgu2ycAK7qUXUmXedqBi4GzupTtQJphcCLpm5TXV7z+k4FLu67n+0vI37QEZgHn5vWRpG8875uPx1+AffO21wI/ArbP978OfDCvR+frIr3JdO7vOuCMiuOwS0/7qfK69sv3r6rYzzbHpuJxJwKX5/Vfkz5lveQYV7aRLo/9ATCtm/1WHsvd2Pq72x8BvlKlzpWkaX47p3v+E/D3pHZ1DzCpy/5fB9zXQ1vv7vh2F8vc/Dyd8//PBr6Z1w8kTWo2GdgDuA3YKW87C/hX0oRxj1bse9dm/78XYfFQT0FExAN5qOQE4KYum48E3iVpdr6/A7APaSKmzp7wZmC/isfcFREdAJKWkxLGrzo3StqF9E/wy1y0gF5m1oyIzZKOIs3AOAW4SNI/RMTcfr3Y6rMJvoY0+ddjOb7vkBJ5fxwJvE5S57w2u5DeSDaSjkfnp5MpwD8Ad+cPQjuydcKujUDnuPk9pGE3SHPifxDScQCel3RSD/vptH9+Xb/J9xeQeuQX9/JaTqioc22+f28vj+mvccB1SnPDjCD9wEo1P4qIkPQg8HREdH5qeojUrpb34zm7O749xbIo0hAopE94lwBExApJD+TyQ0g/VPJ/899iBKnD8yfgBeBbkn5c8dyl5sRfLItIc523keYV6STgvRHxaGVlSXOBp4HXk3pgL1Rs3lCxvpk6/a0jdZvuAu6SdDOp5z+3n7t5A2mulG12X1t0iDTE8JKJtCS1kXr8lfUWRMScKvv4r/waoffj1tN+Kuv0i6TRpDeaAyUF6RNGSDqzj7t4iPSGdEMv9f4X8O8RsSgfo7nd1OtsS1t4abvawrbHpx3YR9LLIuLPVfbV3fHtKZauf7tqRPrRlBO22SAdTHqznwGcTjq2peYx/mK5AvhcZ4+qws+AT3SO00t6Qy7fBVgdafrWk0gJok8i4nngOUmH5aKTSL8A1C1Jr1DF75GSJpDr6Vesuj5ekv6ZdOLvp102PwLsq61XEW3zD9wHPwNOlbR9fr79tPWHayotBqZr68yIu0t6ZS/7XgycmusPU/p1tL7s5xFggqRX5/u9HmfSTJtXRcQrI2JCRIwn9YDf0svjOl0KzJT0xs6CfC7ib7rU24WtFxHMpA4iYj3pF+O+qq2/Dby30pz9PelrLL8inRtB0gGkYSeApcCbO49zPhewXx7n3yUibiKdu5pUZZ+l48RfIBHRERGXVNl0PrA98IDSpXzn5/Kvk/7Bl5KGef5S5bE9mQl8OX9cnkQa5+/J9sA8pROwy4HjgU/24Xm+LOl+0rmJfwQOjy5X9ETEC6ShnR8rndzt8xtKhW+RZjq8Nx+nb1Klxx4RDwPnAj/Pr/1m0ptRTz4JHJ6HO+4B/q4v+8mv60PA9/Jjt5B+N7cnJwA/7FL2feD9vTyu8zmfJvVu5yldzrkSOIw07FFpbo7rduAPfdl3H51L+uW4h/Pf4T/y/Z70NZavA2Py8T6LdJ7q+YhYQzo3cU3etpQ0fPgy4MZc9kvgUwN9Ua3Es3Oa2ZChdNXU9hHxQv50uJh04rxPlwZb4jF+MxtKRpG+B7I9aVz/VCf9/mKyDvAAAAAoSURBVHOP38ysZDzGb2ZWMk78ZmYl48RvZlYyTvxmZiXjxG9mVjL/H38KyPZmzlVTAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "_ = plt.hist(std_diff_replicate,  edgecolor = 'white', linewidth = 0.75, bins = 25)\n",
    "_ = plt.xlabel('Mean of STD difference of All Claim Charges')\n",
    "_ = plt.axvline(std_diff_mean, color='r')\n",
    "_ = plt.axvline(conf_interval[0], color='r', linestyle='--')\n",
    "_ = plt.axvline(conf_interval[1], color='r', linestyle='--')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Confidence interval and p-value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The confidence interval above is often a useful quantity to estimate. If we wish to limit our expected probability of making a Type I error (where we wrongly reject the null hypothesis, and it is, instead, true) to $\\alpha$, the associated confidence interval is our estimate of the interval within which we expect the true population value to be found $100\\times(1 - \\alpha)$% of the time we do this test. In the above we performed bootstrap replicates to estimate the interval and reject the null hypothesis if this interval did not contain zero. You will sometimes see such an interval reported in the output of statistical functions.\n",
    "\n",
    "The partner of the confidence interval is the p-value. The p-value and the confidence interval are linked through our choice of $\\alpha$. The p-value tells us how likely it is, under the null hypothesis, to get an outcome at least as extreme as what was observed. If this fails to reach the level of our _pre-specified_ $\\alpha$, we decide the null hypothesis is sufficiently unlikely to be true and thus reject it. To calculate this p-value via the bootstrap, we have to put ourselves in a position where we are simulating the null hypothesis being true and then calculate the fraction of times we observe a result at least as extreme as that actually observed.\n",
    "\n",
    "Remember how, previously, you used the _t_-test to calculate the p-value for the observed difference between the means of insured and non-insured medical cases. We're now going to repeat this, this time using the bootstrap approach.\n",
    "\n",
    "__Q:__ Perform a bootstrapped hypothesis test at the 5% significance level ($\\alpha = 0.05$) to calculate the p-value of the observed difference between insurance and non-insurance charges, state your null and alternative hypotheses and whether you retain or reject the null hypothesis for the given significance level."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "7602.506384231368\n"
     ]
    }
   ],
   "source": [
    "\n",
    "# Calculate the mean difference between the insured and uninsured\n",
    "mean_diff_observed = np.mean(insured) - np.mean(uninsured)\n",
    "print(mean_diff_observed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "# Calculate with and without permutation with bootstrap\n",
    "np.random.seed(47)\n",
    "uninsured_shifted = uninsured - np.mean(uninsured) + np.mean(insured)\n",
    "\n",
    "# Mean difference with permutation\n",
    "def permutation_sample(data1, data2):\n",
    "    data = np.random.permutation( np.concatenate((data1, data2)) )\n",
    "    perm_sample_1 = data[:len(data1)]\n",
    "    perm_sample_2 = data[len(data1):]\n",
    "    return perm_sample_1, perm_sample_2\n",
    "\n",
    "perm_mean_replicates = np.empty(N_rep)\n",
    "for i in range(N_rep):\n",
    "    perm_insured, perm_uninsured = permutation_sample(insured, uninsured_shifted)\n",
    "    perm_mean_replicates[i] = np.mean(perm_insured) - np.mean(perm_uninsured)\n",
    "    \n",
    "# Mean difference without permutation\n",
    "mean_diff = np.empty(N_rep)\n",
    "for i in range(N_rep):\n",
    "    mean_diff[i] = np.mean(insured_samples) - np.mean(np.random.choice(uninsured_shifted, len(uninsured_shifted)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0\n",
      "0.0\n"
     ]
    }
   ],
   "source": [
    "\n",
    "p_val_permutation = np.sum(perm_mean_replicates >= mean_diff_observed)/len(perm_mean_replicates)\n",
    "print(p_val_permutation)\n",
    "p_val = np.sum(mean_diff >= mean_diff_observed)/len(mean_diff)\n",
    "print(p_val)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The p values are very small and less than 0.05 so we reject the null hypothesis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Q:__ To put the above result in perspective, plot the histogram of your bootstrapped differences along with lines marking the locations of the observed difference. (Why would we plot more than one line, given that we only have one observed difference?)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAD4CAYAAAAHHSreAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAdfklEQVR4nO3dfZhXZb3v8ffniJhEyShYxMMBTsC51ONBG9G0kxQ+4CPmpqRLDR/aeEpLd1mitbfag5u2tTMv2+ZsRfBhi0aaaKaiOXapCYJOBho6YUcGMTBGzayM+p4/1j30i/kNA/N0r5n5vK7rvtZa97p/v993sWb4zn3f67eWIgIzM7NK/y13AGZmVj5ODmZm1oqTg5mZteLkYGZmrTg5mJlZKwNyB7AtQ4cOjTFjxuQOo89Y/dvVAEzcY2KVncU+JlbZZ2a9yooVK16JiGGdeY9SJ4cxY8awfPny3GH0GVPmTwGg/rT6KjuLfdRX2WdmPWrt2mI5alTHXi/p/3U2hlInBzOz/ujUU4tlzr/VnBys8OUv547AzJIy/Do6OVjhsMNyR2BmSRl+HX21khUaGopiZtmtWVOUnNpNDpLmSdogaeVW9Z+RtFrSKkn/VlF/oaTGtO/Iivppqa5R0pyuPQzrtPPOK4qZZXfGGUXJaXuGleYDVwE3tFRI+hAwHdg3Iv4kac9UvxcwE9gbeA/wgKQJ6WXfBQ4HmoAnJC2OiGe66kDMzPqKSy/NHcF2JIeI+KmkMVtVfwqYGxF/Sm02pPrpwMJU/4KkRmBy2tcYEWsAJC1MbZ0czMy2cuihuSPo+JzDBOD/SFoq6WFJB6T6EcDainZNqa6t+lYkzZa0XNLyjRs3djA8M7Pea/Xqv30vNZeOXq00AKgBDgIOAG6TNA5QlbZB9SRU9UESEVEH1AHU1tb6YRNm1u+cdVax7I3fc2gCbo/iSUHLJP0VGJrqK7/TNxJ4Ka23VW9lcNllXPfkdaytvwSAS6ZckjUcs/7ssstyR9DxYaUfAh8GSBPOA4FXgMXATEm7SBoLjAeWAU8A4yWNlTSQYtJ6cWeDty508MGs3aeD39U3sy518MFFyandnoOkW4ApwFBJTcDFwDxgXrq89S1gVupFrJJ0G8VE82bg7Ij4S3qfc4D7gJ2AeRGxqhuOxzrqsccYtXKtE4RZCaxMXxzYZ598MWzP1Uofb2PXKW20/zrw9Sr19wD37FB01nMuuoipr/6a+VecljsSs37vnHOKZW+cczAzs25y+eW5I3ByMDMrnQMOaL9Nd/O9lczMSqYMtzpzz6GPaa5vbnPf5lc3t9lmcNpnZvm13ObMcw6W3R/OuYx7+a/cYZgZcMUVuSNwcrDkL+/9X7zMu3OHYWbApEm5I/CcgyUDVtQzbkXmG8ibGQBPPFGUnNxz6McWNCzYsj7j6qv44ODNrHnfuIwRmRnAF75QLD3nYGZmW1x1Ve4InBzMzEon520zWnjOwcysZB57rCg5uedgZlYyF11ULD3nYNk9+LGP8fr/fD13GGYGXHNN7gicHCxp3nNPXh09MHcYZgZMnJg7As85WDJ25UomPJb5obVmBsDDDxclJ/ccDID31dezeflmnju4BH+ymPVzF19cLD3nYGZmW8yblzuC7RhWkjRP0ob0SNCt950vKSQNTduSdKWkRklPS9q/ou0sSc+nMqtrD8PMrO8YN64oOW3PnMN8YNrWlZJGAYcDL1ZUHwWMT2U2cHVquzvFs6cPBCYDF0uq6UzgZmZ91QMPFCWndpNDRPwU2FRl17eBLwJRUTcduCEKjwNDJA0HjgSWRMSmiGgGllAl4ZiZGXzta0XJqUNzDpKOB9ZFxM8lVe4aAayt2G5KdW3VV3vv2RS9DkaPHt2R8KwD7j35ZF7f299zMCuDG2/MHUEHkoOkQcCXgCOq7a5SF9uob10ZUQfUAdTW1lZtY13vjZoaXt+z2mkys542alTuCDr2PYf/AYwFfi7p18BI4ElJ76boEVQe1kjgpW3UW0lMeOop9v5Jq2sOzCyDe+8tSk47nBwi4hcRsWdEjImIMRT/8e8fES8Di4FPpKuWDgJei4j1wH3AEZJq0kT0EanOSmLfRx/lgMXLc4dhZsDcuUXJqd1hJUm3AFOAoZKagIsj4ro2mt8DHA00Am8CpwNExCZJXwVanm30lYioNsltJXFJ/SV/W59ySZvtzKzrLVyYO4LtSA4R8fF29o+pWA/g7DbazQNK8NUOM7Nye3cJHufueyuZmZXMXXcVJSffPsPMrGS+9a1iedxx+WJwcjAA7j7tNACGNAxuta+ZZgBqpvhL7WY9YdGi3BE4OVjyx8Gtk4KZ5TF0aO4IPOdgyV7LlrHXsmW5wzAz4Pbbi5KTew4GsCUxPDN5cuZIzOzKK4vliSfmi8HJwcysZO68M3cETg69TnN9c+4QzKyb7bZb7gg852BmVjq33lqUnNxzMDMrmauvLpYnnZQvBicHA+CHs2fnDsHMknvuyR2Bk4MlmwcOzB2CmSWDBuWOwMmh33n5jZdZ0LCgVf2+jzwCwNMf+EBPh2RmW7nppmJ5yin5YvCEtAEwoaGBCQ0NucMwM+Daa4uSk3sOZmYls2RJ7gicHMzMSmfnnXNH4GElM7PSmT+/KDm1mxwkzZO0QdLKirrLJf1S0tOS7pA0pGLfhZIaJa2WdGRF/bRU1yhpTtcfiplZ31CG5LA9w0rzgauAGyrqlgAXRsRmSd8ALgQukLQXMBPYG3gP8ICkCek13wUOB5qAJyQtjohnuuYwrLMWnXNO7hDMLKmvzx3BdvQcIuKnwKat6u6PiM1p83FgZFqfDiyMiD9FxAtAIzA5lcaIWBMRbwELU1szMyuhrphzOAP4cVofAayt2NeU6tqqb0XSbEnLJS3fuHFjF4Rn22P/hx5i/4ceyh2GmQH/+Z9FyalTyUHSl4DNwM0tVVWaxTbqW1dG1EVEbUTUDhs2rDPh2Q4Yt2oV41atyh2GmdHLb7wnaRZwLDA1Ilr+o28CRlU0Gwm8lNbbqjczswoPPJA7gg72HCRNAy4Ajo+INyt2LQZmStpF0lhgPLAMeAIYL2mspIEUk9aLOxe6mZl1l3Z7DpJuAaYAQyU1ARdTXJ20C7BEEsDjEfF/I2KVpNuAZyiGm86OiL+k9zkHuA/YCZgXER7DMDOr4j/+o1h++tP5Ymg3OUTEx6tUX7eN9l8Hvl6l/h6gBDeitWo2l+ErmWYGwF13FctSJwfrH3541lm5QzCz5Mc/br9Nd/PtM8zMrBUnBwNg8v33M/n++3OHYWbAd75TlJycHAyA0c89x+jnnssdhpkBDz5YlJw852BmVjKLS3Chv3sOZmbWipODmVnJfPObRcnJw0oGwB/e/vbcIZhZ8rOf5Y7AycGSH51+eu4QzCz5wQ9yR+BhJTMzq8LJwQA45O67OeTuu3OHYWbA3LlFycnDSgbA8F//OncIZpY0NOSOwMnBzKx0Fi7MHYGHlczMrAonBzOzkvnqV4uSk4eVDIDfDRmSOwQzS1avzh2Bk4Ml951ySu4QzCy56abcEWzHsJKkeZI2SFpZUbe7pCWSnk/LmlQvSVdKapT0tKT9K14zK7V/XtKs7jkcMzPrCtsz5zAfmLZV3RzgwYgYDzyYtgGOAsanMhu4GopkQvHs6QOBycDFLQnFyuHQO+7g0DvuyB2GmQH/8i9Fyand5BARPwU2bVU9HViQ1hcAJ1TU3xCFx4EhkoYDRwJLImJTRDQDS2idcCyjYevWMWzdutxhmBmwdm1RcuronMO7ImI9QESsl7Rnqh8BVB5SU6prq74VSbMpeh2MHj26g+GZmfVe11+fO4Kuv5RVVepiG/WtKyPqIqI2ImqHDRvWpcGZmdn26Why+E0aLiItN6T6JmBURbuRwEvbqDczs61ceGFRcupoclgMtFxxNAu4s6L+E+mqpYOA19Lw033AEZJq0kT0EanOSqJ52DCa3VMzK4Xf/rYoObU75yDpFmAKMFRSE8VVR3OB2ySdCbwIfDQ1vwc4GmgE3gROB4iITZK+CjyR2n0lIrae5LaMHjzppDb3LWhI1x6km4HNmlT9SuSaKb4Azawr1NXljmA7kkNEfLyNXVOrtA3g7DbeZx4wb4eiMzOzLPwN6X6g5S//l994uc02U2+9Fdh2D8LMesb55xfLnM+RdnIwAGo2bswdgpklf/hD7gicHEqpub45dwhmltF3v5s7At+y28zMqnByMDMrmfPOK0pOHlYyADaOqHo3EzPrp5wcDICHP/KR3CGYWXLFFbkj8LCSmZlV4eRgABx5000cWYbHT5kZZ59dlJw8rGQAvOPVV3OHYGbJrrvmjsDJwcysdHJ+M7qFh5XMzKwVJwczs5KZPbsoOXlYyQBYP2ZM7hDMLNljj9wRODlY8uixx+YOwcySf/3X3BF4WMnMzKpwcjAAjrn+eo65/vrcYZgZcPrpRcmpU8lB0j9JWiVppaRbJL1N0lhJSyU9L+lWSQNT213SdmPaP6YrDsC6xq6//z27/v73ucMwM2DUqKLk1OHkIGkE8FmgNiL2AXYCZgLfAL4dEeOBZuDM9JIzgeaIeC/w7dTOzMy28pWvFCWnzg4rDQB2lTQAGASsBz4MLEr7FwAnpPXpaZu0f6okdfLzzcysG3Q4OUTEOuCbwIsUSeE1YAXwakRsTs2agJZ7QY8A1qbXbk7tW12wJWm2pOWSlm/0oyvNrB865ZSi5NThS1kl1VD0BsYCrwLfB46q0jRaXrKNfX+riKgD6gBqa2tb7bfu8eKECblDMLNk4sTcEXTuew6HAS9ExEYASbcDBwNDJA1IvYORwEupfRMwCmhKw1C7AZs68fnWhZYdcUTuEMws+ed/zh1B5+YcXgQOkjQozR1MBZ4BHgJmpDazgDvT+uK0Tdr/k4hwz8DMrIQ6M+ewlGJi+UngF+m96oALgM9JaqSYU7guveQ6YI9U/zlgTifiti52wjXXcMI11+QOw8yAmTOLklOnbp8RERcDF29VvQaYXKXtH4GPdubzrPsM+POfc4dgZsmkSbkj8L2VzMxKZ04JxlV8+wwzM2vFycHMrGT+4R+KkpOHlQyANXvvnTsEM0ve//7cETg5WPLkhz6UOwQzS84/P3cEHlYyM7MqnBwMgBlXXcWMq67KHYaZAccfX5ScPKxkZlYyU6fmjsDJwcysdM49N3cEHlYyM7MqnBzMzErmqKOKkpOHlQyA58pwMxczA+C443JH4ORgydMf+EDuEMws+fSnc0fgYSVLBrz1FgPeeit3GGZWEk4OBsAJdXWcUFeXOwwzAw47rCg5eVjJzKxkTjopdwRODmZmpfOP/5g7gk4OK0kaImmRpF9KelbS+yXtLmmJpOfTsia1laQrJTVKelrS/l1zCGZm1tU623P4DnBvRMyQNBAYBFwEPBgRcyXNoXhW9AXAUcD4VA4Erk5L6yOa65vb3FczpaYHIzHr3aZMKZb19fli6HBykPRO4IPAaQAR8RbwlqTpwJTUbAFQT5EcpgM3REQAj6dex/CIWN/h6K3LPDO51WO/zSyT007LHUHneg7jgI3A9ZL+N7ACOBd4V8t/+BGxXtKeqf0IYG3F65tS3d8lB0mzgdkAo0eP7kR4tiOcHMzKo7cnhwHA/sBnImKppO9QDCG1RVXqolVFRB1QB1BbW9tqv7VvQcOCHX7N2954A4A/Dh68Q+8/a9KsHf4sM9u2P/+5WO68c74YOjMh3QQ0RcTStL2IIln8RtJwgLTcUNF+VMXrRwIvdeLzrQsdO38+x86fnzsMMwMOP7woOXU4OUTEy8BaSRNT1VTgGWAx0PLn5CzgzrS+GPhEumrpIOA1zzeYmbX2yU8WJafOXq30GeDmdKXSGuB0ioRzm6QzgReBj6a29wBHA43Am6ltv7Stq3rMzE45JXcEnUwOEdEA1FbZ1eo5RukqpbM783lmZv3Bm28Wy0GD8sXgb0ibmZXM0UcXy175PQfrW54+5JDcIZhZ8qlP5Y7AycGS5/bbL3cIZpaU4cZ7vmW3ATC4uZnBzZ4oNyuD114rSk5ODgbAtJtvZtrNN+cOw8yA6dOLkpOHlczMSuazn80dgZODmVnpnHhi7gg8rGRmVjqvvFKUnNxzMDMrmRkziqW/52DZrWh5uoiZZff5z+eOwMnBkhf22Sd3CGaWHHdc7gg852BJzYYN1GzY0H5DM+t2L79clJycHAyAqbfdxtTbbssdhpkBM2cWJScPK5mZlcycbT1Ts4c4OZiZlcy0abkj8LCSmVnprF1blJzcczAzK5lTTy2Wvfp7DpJ2ApYD6yLiWEljgYXA7sCTwKkR8ZakXYAbgPcBvwVOiohfd/bzrWsszf00czPb4stfzh1B1wwrnQs8W7H9DeDbETEeaAbOTPVnAs0R8V7g26mdlcTaiRNZO3Fi7jDMDDjssKLk1KnkIGkkcAxwbdoW8GFgUWqyADghrU9P26T9U1N7K4Fh69YxbN263GGYGbBmTVFy6mzP4Qrgi8Bf0/YewKsRsTltNwEj0voIYC1A2v9aam8lcOgdd3DoHXfkDsPMgDPOKEpOHZ5zkHQssCEiVkia0lJdpWlsx77K950NzAYYPXp0R8MzM+u1Lr00dwSdm5A+BDhe0tHA24B3UvQkhkgakHoHI4GXUvsmYBTQJGkAsBuwaes3jYg6oA6gtra2VfIwM+vrDj00dwSdSA4RcSFwIUDqOZwfESdL+j4wg+KKpVnAnekli9P2z9L+n0SE//PvJ5rrt/186popNT0UiVn5rV5dLHNeI9Id33O4AFgo6WvAU8B1qf464EZJjRQ9hsx3DjEzK6ezziqWvfp7DgARUQ/Up/U1wOQqbf4IfLQrPs+63qPHHJM7BDNLLrssdwT+hrQl68eOzR2CmSUHH5w7At9byZLhL7zA8BdeyB2GmQErVxYlJ/ccDIBDfvQjABadc84OvW5Bw4It67MmzerSmMz6q5Zfw14/52BmZl3n8stzR+Dk0GdU/gVvZr3bAQfkjsBzDmZmpdPQUJSc3HMwMyuZ884rlp5z6GPa+zZwGT38kY/kDsHMkiuuyB2Bk4MlG0eMaL+RmfWISZNyR+A5B0tGrV7NqJYbuphZVk88UZSc3HMwAA5csgTAT4MzK4EvfKFYes7BzMy2uOqq3BE4OZiZlc4+++SOwHMOZmal89hjRcnJPQczs5K56KJi6TkHy+7Bj30sdwhmllxzTe4InBwsad5zz9whmFlShosGPedgAIxduZKxuW8gb2YAPPxwUXLqcHKQNErSQ5KelbRK0rmpfndJSyQ9n5Y1qV6SrpTUKOlpSft31UFY572vvp735RzgNLMtLr64KDl1ZlhpM/D5iHhS0juAFZKWAKcBD0bEXElzgDnABcBRwPhUDgSuTkszM6swb17uCDrRc4iI9RHxZFr/HfAsMAKYDrQ8XGABcEJanw7cEIXHgSGShnc4cjOzPmrcuKLk1CVzDpLGAPsBS4F3RcR6KBII0DLTOQJYW/GyplS39XvNlrRc0vKNGzd2RXhmZr3KAw8UJadOX60kaTDwA+C8iHhdUptNq9RFq4qIOqAOoLa2ttV+M7O+7mtfK5aHHZYvhk4lB0k7UySGmyPi9lT9G0nDI2J9GjbakOqbgFEVLx8JvNSZz7euc+/JJ+cOwcySG2/MHUHnrlYScB3wbET8e8WuxcCstD4LuLOi/hPpqqWDgNdahp8svzdqanijpiZ3GGYGjBpVlJw603M4BDgV+IWklqedXgTMBW6TdCbwIvDRtO8e4GigEXgTOL0Tn23AgoYF7TfaThOeegqA5/bbr8vec0e09/S8milOXNZ/3HtvsZw2LV8MHU4OEfEI1ecRAKZWaR/A2R39POte+z76KJAvOZjZ38ydWyx7ZXLo73rjc6K729Y9mVmTZrXR0sy2ZeHC3BE4OZiZlc673507At9bycysdO66qyg5uedgZlYy3/pWsTzuuHwxODkYAHefdlruEMwsWbQodwRODpb8cfDg3CGYWTJ0aO4IPOdgyV7LlrHXsmW5wzAz4Pbbi5KTew4GsCUxPDN5cuZIzOzKK4vliSfmi8HJwcysZO68s/023c3JwcysZHbbLXcEnnMwMyudW28tSk7uOVivsK3blfimfNbXXH11sTzppHwxODkYAD+cPTt3CGaW3HNP7gicHNpUxhvrdeUture2eeDALn/Pynh9Ez6z7TdoUO4IPOdgyb6PPMK+jzySOwwzA266qSg5OTkYABMaGpjQ0NB+QzPrdtdeW5ScPKxUct05lGRm5bRkSe4IMiQHSdOA7wA7AddGxNyejgHKOafQl3Xn/IMfMWp9zc47546gh5ODpJ2A7wKHA03AE5IWR8Qz3fF5TgAGTh7W+8yfXyxz3iy5p3sOk4HGiFgDIGkhMB3oluTQW3joyMwqlSE5KCJ67sOkGcC0iPhk2j4VODAizqloMxtoueh+IrC6G0MaCrzSje9fRv3xmMHH3d/0x+OuPOb/HhHDOvNmPd1zUJW6v8tOEVEH1PVIMNLyiKjtic8qi/54zODjzh1HT+uPx93Vx9zTl7I2AaMqtkcCL/VwDGZm1o6eTg5PAOMljZU0EJgJLO7hGMzMrB09OqwUEZslnQPcR3Ep67yIWNWTMWylR4avSqY/HjP4uPub/njcXXrMPTohbWZmvYNvn2FmZq04OZiZWSt9NjlIulzSLyU9LekOSUMq9l0oqVHSaklHVtRPS3WNkuZU1I+VtFTS85JuTZPpvU5bx9cbSRol6SFJz0paJencVL+7pCXpXC2RVJPqJenKdOxPS9q/4r1mpfbPS+oV9xaXtJOkpyTdnbar/oxK2iVtN6b9Yyreo+rvQVlJGiJpUfq9flbS+/v6+Zb0T+nne6WkWyS9rcfOdUT0yQIcAQxI698AvpHW9wJ+DuwCjAV+RTE5vlNaHwcMTG32Sq+5DZiZ1r8HfCr38XXg36PN4+uNBRgO7J/W3wE8l87tvwFzUv2civN+NPBjiu/aHAQsTfW7A2vSsiat1+Q+vu04/s8B/wXcnbar/owCnwa+l9ZnArem9aq/B7mPq51jXgB8Mq0PBIb05fMNjABeAHatOMen9dS57rM9h4i4PyI2p83HKb5TAcXtOhZGxJ8i4gWgkeK2Hltu7RERbwELgemSBHwYWJRevwA4oaeOowtVPb7MMXVYRKyPiCfT+u+AZyl+maZTnCP4+3M1HbghCo8DQyQNB44ElkTEpohoBpYA03rwUHaYpJHAMcC1aXtbP6OV/x6LgKmpfVu/B6Uk6Z3AB4HrACLirYh4lb5/vgcAu0oaAAwC1tND57rPJoetnEHxVwQU/4GsrdjXlOraqt8DeLUi0bTU9zZtHV+vl7rP+wFLgXdFxHooEgiwZ2q2o+e9zK4Avgj8NW1v62d0y/Gl/a+l9r3tuMcBG4Hr03DatZLeTh8+3xGxDvgm8CJFUngNWEEPnetenRwkPZDG4rYu0yvafAnYDNzcUlXlraID9b1NXzmOvyNpMPAD4LyIeH1bTavU9brzK+lYYENErKisrtI02tnXq46b4i/o/YGrI2I/4PcUw0ht6fXHneZPplMMBb0HeDtwVJWm3XKue/XDfiLisG3tT5NNxwJTIw2+se1beFSrf4WiSzogZePeesuPPnfrEkk7UySGmyPi9lT9G0nDI2J9GkbYkOrbOv4mYMpW9fXdGXcnHQIcL+lo4G3AOyl6Em39jLYcd1MamtgN2ETv+3loApoiYmnaXkSRHPry+T4MeCEiNgJIuh04mB46172657AtKh4qdAFwfES8WbFrMTAzzeyPBcYDy2jj1h4pqTwEzEivnwXc2VPH0YX61K1L0ljqdcCzEfHvFbsWU5wj+PtztRj4RLqK5SDgtTQMcR9whKSa9JfaEamulCLiwogYGRFjKM7hTyLiZNr+Ga3895iR2gdt/x6UUkS8DKyVNDFVTaW41X9fPt8vAgdJGpR+3luOuWfOde4Z+e4qFJMua4GGVL5Xse9LFDP2q4GjKuqPprjq5VfAlyrqx6V/zEbg+8AuuY+vg/8mVY+vNxbgAxRd46crzvHRFGOsDwLPp+Xuqb0oHjT1K+AXQG3Fe52Rzm0jcHruY9uBf4Mp/O1qpao/oxS9i++n+mXAuIrXV/09KGsBJgHL0zn/IcXVRn36fAOXAr8EVgI3Ulxx1CPn2rfPMDOzVvrssJKZmXWck4OZmbXi5GBmZq04OZiZWStODmZm1oqTg5mZteLkYGZmrfx/ixXLENxcby8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "_ = plt.hist(mean_diff, alpha = 0.5, bins = 20, color = 'g')\n",
    "_ = plt.hist(perm_mean_replicates,  linewidth = 0.5, alpha = 0.25, bins = 20, color = 'm')\n",
    "_ = plt.axvline(np.mean(perm_mean_replicates),color='g')\n",
    "_ = plt.axvline(np.mean(mean_diff),color='r', linestyle='--')\n",
    "_ = plt.axvline(np.mean(mean_diff_observed),color='b', linestyle=':')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Q:__ Compare your p-value above with that obtained using the _t_-test function in the previous assignment. Do you think you would want to try to perform enough bootstrap replicates to observe a random difference as large as that we did observe?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "# The p value here with a bootstrap approach is 0.\n",
    "# The p value from the previous exercise with a frequentist (parametric) approach was 4.46.\n",
    "# A lower p value is better, so this approach (bootstrap) is preferred for this dataset."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Q:__ Consider the two variants of the _t_-test we performed in the previous assignment. Which one would you use now?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "# In the previous assignment we used two variants of the t-test: assuming a normal\n",
    "# variance, and assuming a not-normal variance.\n",
    "# For this data here, we should assume that there is a normal variance."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Q:__ If, instead of being asked whether the means of two groups were different, you were working with an ad-recommender team who wanted to release a new recommendation algorithm, and you were looking at click-through rate both for the current algorithm (call it A) and from trials of their new algorithm (call it B), would you perform a two-sided test as above? What would be your null and alternative hypotheses and what would be the real-world consequence of rejecting the null hypothesis?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__A:__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "# If instead we were analyzing an A/B test for an ad-recommeder system. \n",
    "# We would perform a permutation test to determine if the redesign improved the click-through rate.\n",
    "# The null hypothesis is 'The results from the algorithm are impervious to the change.'\n",
    "# A low p value implies that change in strategy led to improved performance.\n",
    "# Your alternative hypothesis is that Algo B improved results."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Learning outcomes\n",
    "\n",
    "You've previously applied frequentist methods to calculate confidence intervals, p-values, and perform hypothesis tests. Frequentist methods use theoretical results to calculate what we expect would happen if experiments were to be run again and again and again. Now you've seen how you can do the same things using the bootstrap approach, which does not rely on such theory, and attendant assumptions, but instead literally does run experiments again and again and again.\n",
    "\n",
    "In these exercises, you have:\n",
    "* calculated the same confidence interval lower limit as you did previously\n",
    "* tested the assumption that the variances of the two groups (insured vs. non-insured) were equal - something a bit harder to do using the frequentist method because of the nature of the sampling distribution for variance\n",
    "* calculated the p-value for the difference between the means of the two groups and compared with the result obtained using the previous frequentist approach\n",
    "\n",
    "You are now well equipped to apply the bootstrap approach to a wide variety of problems. Just think about what conditions you wish to recreate in your simulated reruns."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "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.7.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}