part3_Lesson12-Hypothesis_Testing.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Put these at the top of every notebook, to get automatic reloading and inline plotting\n",
    "%reload_ext autoreload\n",
    "%autoreload 2\n",
    "%matplotlib inline\n",
    "from IPython.core.interactiveshell import InteractiveShell\n",
    "InteractiveShell.ast_node_interactivity = \"all\"\n",
    "def pv(variable):\n",
    "    if(len(variable)>1): print(variable, \":\", eval(variable))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# low, upper\n",
    "# plt.axvline(x=low, color='r', linewidth=2);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulating From the Null Hypothesis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "%matplotlib inline\n",
    "np.random.seed(42)\n",
    "\n",
    "full_data = pd.read_csv('coffee_dataset.csv')\n",
    "sample_data = full_data.sample(200)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. If you were interested in studying whether the average height for coffee drinkers is the same as for non-coffee drinkers, what would the null and alternative hypotheses be? Write them in the cell below, and use your answer to answer the first quiz question below.\n",
    "\n",
    "2. If you were interested in studying whether the average height for coffee drinkers is less than non-coffee drinkers, what would the null and alternative hypotheses be? Write them in the cell below, and use your answer to answer the second quiz question below.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "average height for coffee drinkers is the same as for non-coffee drinkers\n",
    "\n",
    "    H0 : height(coffee) = height(no coffee)\n",
    "    H1 : height(coffee) <> height(no coffee)\n",
    "    \n",
    "average height for coffee drinkers is less than non-coffee drinkers\n",
    "\n",
    "    H0 : height(coffee) >= height(no coffee)\n",
    "    H1 : height(coffee) < height(no coffee)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "3. For 10,000 iterations: bootstrap the sample data, calculate the mean height for coffee drinkers and non-coffee drinkers, and calculate the difference in means for each sample. You will want to have three arrays at the end of the iterations - one for each mean and one for the difference in means. Use the results of your sampling distribution, to answer the third quiz question below.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "np.mean(drink_means) : 68.1170203438573\n",
      "np.mean(no_drink_means) : 66.7802393564418\n",
      "np.mean(boot_diffs) : 1.3367809874155205\n",
      "np.std(drink_means) : 0.24199030603946156\n",
      "np.std(no_drink_means) : 0.40740156905046576\n",
      "np.std(boot_diffs) : 0.47461751257882373\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(array([   3.,   76.,  420., 1408., 2764., 2909., 1718.,  602.,   92.,\n",
       "           8.]),\n",
       " array([67.14733174, 67.33718939, 67.52704705, 67.71690471, 67.90676237,\n",
       "        68.09662002, 68.28647768, 68.47633534, 68.666193  , 68.85605065,\n",
       "        69.04590831]),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAFB5JREFUeJzt3X+w5XV93/HnS1DS+ItFFrsu2yza1Qan44q3iKVJjKT8cibgVBroVHfQzqYZmEKaznRJO4WJZQbbqBmnhoaELWASCdFYtmFT3BIyTmbCj8Ws/BS54kbW3cBGELROteC7f5zPLYfLuT/33nNWPs/HzJnzPe/v53u+7/Pdc+/rfn+cs6kqJEn9edmkG5AkTYYBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASerUkZNuYD7HHntsbdy4cdJtSNKPlHvuuedvqmrtQuMWDIAkPwZ8ETiqjf9sVV2e5ATgRuAY4EvAB6rqB0mOAm4A3gF8C/iFqtrbnusy4MPAc8C/qqpb51v3xo0b2b1790ItSpKGJPmrxYxbzCGg7wPvqaq3AZuBM5OcAnwU+ERVbQKeYvCLnXb/VFX9XeATbRxJTgTOB94KnAn8ZpIjFv+SJEkracEAqIHvtocvb7cC3gN8ttWvB85t0+e0x7T5pyVJq99YVd+vqq8D08DJK/IqJElLtqiTwEmOSLIHeALYBXwN+HZVPduG7APWt+n1wGMAbf7TwOuG6yOWGV7X1iS7k+w+ePDg0l+RJGlRFhUAVfVcVW0GjmfwV/tPjhrW7jPHvLnqs9d1TVVNVdXU2rULnsOQJC3Tki4DrapvA38GnAIcnWTmJPLxwP42vQ/YANDmvxZ4crg+YhlJ0pgtGABJ1iY5uk3/LeDngIeA24H3t2FbgJvb9I72mDb/T2vwv87sAM5PclS7gmgTcNdKvRBJ0tIs5nMA64Dr2xU7LwNuqqo/TvIgcGOS/wj8JXBtG38t8Okk0wz+8j8foKoeSHIT8CDwLHBRVT23si9HkrRYOZz/S8ipqanycwCStDRJ7qmqqYXG+VUQktSpw/qrIKTD2cZtt0xkvXuveu9E1quXHvcAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpU34SWD/SJvVpXOmlwD0ASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktSpBQMgyYYktyd5KMkDSS5p9SuSfDPJnnY7e2iZy5JMJ3k4yRlD9TNbbTrJttV5SZKkxVjM/wfwLPArVfWlJK8G7kmyq837RFX9+vDgJCcC5wNvBd4A/K8kb26zPwX8Y2AfcHeSHVX14Eq8EEnS0iwYAFV1ADjQpr+T5CFg/TyLnAPcWFXfB76eZBo4uc2brqpHAZLc2MYaAJI0AUs6B5BkI/B24M5WujjJvUm2J1nTauuBx4YW29dqc9UlSROw6ABI8irgc8ClVfUMcDXwJmAzgz2Ej80MHbF4zVOfvZ6tSXYn2X3w4MHFtidJWqJFBUCSlzP45f97VfVHAFX1eFU9V1U/BH6b5w/z7AM2DC1+PLB/nvoLVNU1VTVVVVNr165d6uuRJC3SYq4CCnAt8FBVfXyovm5o2PuA+9v0DuD8JEclOQHYBNwF3A1sSnJCklcwOFG8Y2VehiRpqRZzFdCpwAeA+5LsabVfBS5IspnBYZy9wC8CVNUDSW5icHL3WeCiqnoOIMnFwK3AEcD2qnpgBV+LJGkJFnMV0J8z+vj9znmWuRK4ckR953zLSZLGx08CS1KnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdWjAAkmxIcnuSh5I8kOSSVj8mya4kj7T7Na2eJJ9MMp3k3iQnDT3Xljb+kSRbVu9lSZIWspg9gGeBX6mqnwROAS5KciKwDbitqjYBt7XHAGcBm9ptK3A1DAIDuBx4J3AycPlMaEiSxm/BAKiqA1X1pTb9HeAhYD1wDnB9G3Y9cG6bPge4oQbuAI5Osg44A9hVVU9W1VPALuDMFX01kqRFW9I5gCQbgbcDdwKvr6oDMAgJ4Lg2bD3w2NBi+1ptrvrsdWxNsjvJ7oMHDy6lPUnSEiw6AJK8CvgccGlVPTPf0BG1mqf+wkLVNVU1VVVTa9euXWx7kqQlWlQAJHk5g1/+v1dVf9TKj7dDO7T7J1p9H7BhaPHjgf3z1CVJE7CYq4ACXAs8VFUfH5q1A5i5kmcLcPNQ/YPtaqBTgKfbIaJbgdOTrGknf09vNUnSBBy5iDGnAh8A7kuyp9V+FbgKuCnJh4FvAOe1eTuBs4Fp4HvAhQBV9WSSjwB3t3G/VlVPrsirkCQt2YIBUFV/zujj9wCnjRhfwEVzPNd2YPtSGpQkrQ4/CSxJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOrWY/w9A0mFk47ZbJrbuvVe9d2Lr1spzD0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASerUggGQZHuSJ5LcP1S7Isk3k+xpt7OH5l2WZDrJw0nOGKqf2WrTSbat/EuRJC3FYvYArgPOHFH/RFVtbredAElOBM4H3tqW+c0kRyQ5AvgUcBZwInBBGytJmpAFvw20qr6YZOMin+8c4Maq+j7w9STTwMlt3nRVPQqQ5MY29sEldyxJWhGHcg7g4iT3tkNEa1ptPfDY0Jh9rTZXXZI0IcsNgKuBNwGbgQPAx1o9I8bWPPUXSbI1ye4kuw8ePLjM9iRJC1lWAFTV41X1XFX9EPhtnj/Msw/YMDT0eGD/PPVRz31NVU1V1dTatWuX054kaRGWFQBJ1g09fB8wc4XQDuD8JEclOQHYBNwF3A1sSnJCklcwOFG8Y/ltS5IO1YIngZN8Bng3cGySfcDlwLuTbGZwGGcv8IsAVfVAkpsYnNx9Frioqp5rz3MxcCtwBLC9qh5Y8VcjSVq0xVwFdMGI8rXzjL8SuHJEfSewc0ndSZJWjZ8ElqROGQCS1KkFDwFJi7Fx2y2TbkHSErkHIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6tSCAZBke5Inktw/VDsmya4kj7T7Na2eJJ9MMp3k3iQnDS2zpY1/JMmW1Xk5kqTFWswewHXAmbNq24DbqmoTcFt7DHAWsKndtgJXwyAwgMuBdwInA5fPhIYkaTIWDICq+iLw5KzyOcD1bfp64Nyh+g01cAdwdJJ1wBnArqp6sqqeAnbx4lCRJI3Rcs8BvL6qDgC0++NafT3w2NC4fa02V12SNCErfRI4I2o1T/3FT5BsTbI7ye6DBw+uaHOSpOctNwAeb4d2aPdPtPo+YMPQuOOB/fPUX6SqrqmqqaqaWrt27TLbkyQtZLkBsAOYuZJnC3DzUP2D7WqgU4Cn2yGiW4HTk6xpJ39PbzVJ0oQcudCAJJ8B3g0cm2Qfg6t5rgJuSvJh4BvAeW34TuBsYBr4HnAhQFU9meQjwN1t3K9V1ewTy5KkMVowAKrqgjlmnTZibAEXzfE824HtS+pOkrRq/CSwJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ06ctINSPrRsXHbLRNZ796r3juR9b7UuQcgSZ0yACSpUwaAJHXqkAIgyd4k9yXZk2R3qx2TZFeSR9r9mlZPkk8mmU5yb5KTVuIFSJKWZyX2AH62qjZX1VR7vA24rao2Abe1xwBnAZvabStw9QqsW5K0TKtxCOgc4Po2fT1w7lD9hhq4Azg6ybpVWL8kaREONQAK+EKSe5JsbbXXV9UBgHZ/XKuvBx4bWnZfq71Akq1JdifZffDgwUNsT5I0l0P9HMCpVbU/yXHAriRfmWdsRtTqRYWqa4BrAKampl40X5K0Mg5pD6Cq9rf7J4DPAycDj88c2mn3T7Th+4ANQ4sfD+w/lPVLkpZv2QGQ5JVJXj0zDZwO3A/sALa0YVuAm9v0DuCD7WqgU4CnZw4VSZLG71AOAb0e+HySmef5/ar6n0nuBm5K8mHgG8B5bfxO4GxgGvgecOEhrFuSdIiWHQBV9SjwthH1bwGnjagXcNFy1ydJWll+GdxLzKS+rEvSjx6/CkKSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXqyEk3IEkL2bjtlomte+9V753YulebAbAKJvlmlaTF8hCQJHVq7AGQ5MwkDyeZTrJt3OuXJA2MNQCSHAF8CjgLOBG4IMmJ4+xBkjQw7j2Ak4Hpqnq0qn4A3AicM+YeJEmM/yTweuCxocf7gHeOuQdJWrRJXdQxjquPxh0AGVGrFwxItgJb28PvJnl41bua27HA30xw/Quxv0Njf4fucO/xR7a/fPSQnvcnFjNo3AGwD9gw9Ph4YP/wgKq6BrhmnE3NJcnuqpqadB9zsb9DY3+H7nDv0f7mN+5zAHcDm5KckOQVwPnAjjH3IElizHsAVfVskouBW4EjgO1V9cA4e5AkDYz9k8BVtRPYOe71LtNhcShqHvZ3aOzv0B3uPdrfPFJVC4+SJL3k+FUQktSprgIgydFJPpvkK0keSvKuJH+QZE+77U2yZ8RyG5Lc3pZ5IMklQ/OuSPLNoec4e9z9tWX3Jrmvjds9VD8mya4kj7T7NePuL8lbhsbsSfJMkkvbvNXefpuT3DGzXZKcPMeyW9o2eiTJlqH6O9p2nU7yySSjLmVe1f7amL9o7717k/zC0Lzrknx9aPttHnd/bdnnhnrYMVQ/Icmdbbv+Qbv4Y6z9JfnZWe+//5Pk3DZvtbff29q/3X1J/keS18yx7MivyFnJ7TdSVXVzA64H/kWbfgVw9Kz5HwP+w4jl1gEntelXA18FTmyPrwD+zST7a/P2AseOqP8nYFub3gZ8dBL9DY05Avhr4CfGsf2ALwBntdrZwJ+NWO4Y4NF2v6ZNr2nz7gLexeAzLH8y81xj7u/NwKY2/QbgwMy2B64D3j/J7dfmfXeO+k3A+W36vwK/NIn+Zv1bPwn8+Ji2393Az7Tah4CPzPEz8TXgjW25L/P875cV236jbt3sAbTk/WngWoCq+kFVfXtofoB/Cnxm9rJVdaCqvtSmvwM8xOBTzYdFfws4h8Ebk3Z/7oT7Ow34WlX91XL6WEZ/Bcz81fVaZn3upDkD2FVVT1bVU8Au4Mwk64DXVNVf1OAn8AZWfvst2F9VfbWqHmnT+4EngLXL6WM1+pvnOQO8B/hsK63G+2+p/b0f+JOq+t5y+lhGf28BvtiG7QL+yYjFR35Fzkpuv7l0EwAM0vUg8N+S/GWS30nyyqH5PwU8PvODNpckG4G3A3cOlS9uu+bbs/xDLIfaXwFfSHJPBp+mnvH6qjoAgyADjptQfzPO58UhsZrb71LgPyd5DPh14LIRy476ipL17bZvRH3c/f1/7RDHKxj8xTjjyrb9PpHkqAn192PtEMwdM4dXgNcB366qZ9vjiW8/Rr//VnP73Q/8fBtzHi/8IOyMud5/K7n9RuopAI4ETgKurqq3A/+bwSGRGRewwF+vSV4FfA64tKqeaeWrgTcBmxnsmn9sQv2dWlUnMfim1YuS/PQy+1it/mjHL38e+MOh8mpvv18CfrmqNgC/TPsLbXZrI2o1T33c/Q2aHOyRfBq4sKp+2MqXAX8P+AcMDm/82wn193dq8InWfwb8RpI3cXhuv7/P4HNIM1Z7+32Iwc/jPQwOH/9gVGsjaiv9/httJY8nHc434G8De4ce/xRwS5s+EngcOH6e5V/O4I3zr+cZsxG4fxL9zXquK2jH1YGHgXVteh3w8KT6Y3A46gvj3H7A0zx/uXOAZ0YsewHwW0OPf6vV1gFfmWvcuPpr814DfAk4b551vBv440n0N+u5rmNwqCUMvufmyFZ/F3DrpPoDLgGuGef2mzXmzcBdI5Z9wXZhEEqXreT2m+vWzR5AVf018FiSt7TSacCDbfrnGPyg7xu1bDsWdy3wUFV9fNa8dUMP38dgl2/c/b0yyatnpoHTh/rYAcxc1bIFuHnc/Q150V7CGLbffuBnWu09wKhDVLcCpydZ0w5Bnc7gB+0A8J0kp7T3wAdZ+e23YH9tz+nzwA1V9Yez5q1r92FwfHjs269tt6Pa9LHAqcCDNfitdTuDMIDVef8t5t93xpzvv9XafkmOa8//MuDfMziRO9vIr8hZye03X+Pd3BgcZtgN3Av8d56/0uM64F/OGvsGYGeb/kcMdr3uBfa029lt3qeB+9q8HbS/tsfc3xsZXDnwZeAB4N8NjXsdcBuDH4zbgGPG3V97/OPAt4DXzhq3qtuv/dvd07bNncA72tgp4HeGlv0QMN1uFw7Vpxj8Uvga8F9of22Osz/gnwP/d+i9twfY3Ob9adt+9wO/C7xqAv39w9bDl9v9h4ee840MrqSaZnDo76gJ/ftuBL4JvGzWc6729ruEwVWDXwWu4vm9ldk/H2e3MV/jhT+/K7b9Rt38JLAkdaqbQ0CSpBcyACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6tT/A3m8mFOtxQx/AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "(array([  16.,  154.,  649., 1747., 2760., 2647., 1453.,  450.,  110.,\n",
       "          14.]),\n",
       " array([65.31124228, 65.61188904, 65.9125358 , 66.21318257, 66.51382933,\n",
       "        66.81447609, 67.11512285, 67.41576961, 67.71641638, 68.01706314,\n",
       "        68.3177099 ]),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEF1JREFUeJzt3X+s3XV9x/HnayD+4Y9RQiFYmhVNXYZ/rLIOWZhOZStQ/qgucYE/tGEudQYWNW5JdckgGhK2+SMxcSRVOnE/RKY4G2mGlbgZ/wC5OAQqMu6ws9d29Dqcuplo0Pf+OJ+bHdpzf9/ec9vP85GcnHPe5/O55/PO97avfr/fc75NVSFJ6s8vjHsBkqTxMAAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTpz3AuYy7nnnlubNm0a9zIk6ZTy0EMPfa+q1s83bk0HwKZNm5iYmBj3MiTplJLkPxYyzkNAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUqTX9TWBpPpt23zO29z506zVje29pJbgHIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSp/wegLRE4/oOgt8/0EpxD0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE7NGwBJNib5cpLHkxxM8o5WvznJd5M83G7bh+a8J8lkkieSXDlUv6rVJpPsPjktSZIWYiFfBHsWeHdVfT3Ji4CHkhxor324qj4wPDjJxcC1wCuAlwBfSvLy9vJHgd8BpoAHk+yrqm+uRCOSpMWZNwCq6ihwtD3+UZLHgQ1zTNkB3FlVPwG+nWQSuLS9NllVTwEkubONNQAkaQwWdQ4gySbglcADrXRjkkeS7E2yrtU2AIeHpk212mz1499jV5KJJBPT09OLWZ4kaREWHABJXgh8FnhnVf0QuA14GbCFwR7CB2eGjphec9SfW6jaU1Vbq2rr+vXrF7o8SdIiLehicEmex+Av/7+rqrsBqurpodc/BnyhPZ0CNg5NvxA40h7PVpckrbKFfAoowO3A41X1oaH6BUPD3gg81h7vA65N8vwkFwGbga8BDwKbk1yU5CwGJ4r3rUwbkqTFWsgewOXAm4FHkzzcau8FrkuyhcFhnEPA2wCq6mCSuxic3H0WuKGqfgaQ5EbgXuAMYG9VHVzBXiRJi7CQTwF9ldHH7/fPMecW4JYR9f1zzZMkrR6/CSxJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROzRsASTYm+XKSx5McTPKOVj8nyYEkT7b7da2eJB9JMpnkkSSXDP2snW38k0l2nry2JEnzWcgewLPAu6vqV4DLgBuSXAzsBu6rqs3Afe05wNXA5nbbBdwGg8AAbgJeBVwK3DQTGpKk1TdvAFTV0ar6env8I+BxYAOwA7ijDbsDeEN7vAP4ZA3cD5yd5ALgSuBAVT1TVd8HDgBXrWg3kqQFW9Q5gCSbgFcCDwDnV9VRGIQEcF4btgE4PDRtqtVmq0uSxmDBAZDkhcBngXdW1Q/nGjqiVnPUj3+fXUkmkkxMT08vdHmSpEVaUAAkeR6Dv/z/rqrubuWn26Ed2v2xVp8CNg5NvxA4Mkf9OapqT1Vtraqt69evX0wvkqRFWMingALcDjxeVR8aemkfMPNJnp3A54fqb2mfBroM+EE7RHQvsC3Junbyd1urSZLG4MwFjLkceDPwaJKHW+29wK3AXUneCnwHeFN7bT+wHZgEfgxcD1BVzyR5P/BgG/e+qnpmRbqQJC3avAFQVV9l9PF7gCtGjC/ghll+1l5g72IWKEk6OfwmsCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnq1EL+PwBpXpt23zPuJUhaJPcAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKn5g2AJHuTHEvy2FDt5iTfTfJwu20feu09SSaTPJHkyqH6Va02mWT3yrciSVqMhewBfAK4akT9w1W1pd32AyS5GLgWeEWb81dJzkhyBvBR4GrgYuC6NlaSNCbz/ocwVfWVJJsW+PN2AHdW1U+AbyeZBC5tr01W1VMASe5sY7+56BVLklbEcs4B3JjkkXaIaF2rbQAOD42ZarXZ6idIsivJRJKJ6enpZSxPkjSXpQbAbcDLgC3AUeCDrZ4RY2uO+onFqj1VtbWqtq5fv36Jy5MkzWdJ/ydwVT098zjJx4AvtKdTwMahoRcCR9rj2eqSpDFY0h5AkguGnr4RmPmE0D7g2iTPT3IRsBn4GvAgsDnJRUnOYnCieN/Sly1JWq559wCSfAp4LXBukingJuC1SbYwOIxzCHgbQFUdTHIXg5O7zwI3VNXP2s+5EbgXOAPYW1UHV7wbqQObdt8ztvc+dOs1Y3tvrbyFfArouhHl2+cYfwtwy4j6fmD/olYnSTpp/CawJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOjVvACTZm+RYkseGauckOZDkyXa/rtWT5CNJJpM8kuSSoTk72/gnk+w8Oe1IkhZqIXsAnwCuOq62G7ivqjYD97XnAFcDm9ttF3AbDAIDuAl4FXApcNNMaEiSxmPeAKiqrwDPHFfeAdzRHt8BvGGo/skauB84O8kFwJXAgap6pqq+DxzgxFCRJK2ipZ4DOL+qjgK0+/NafQNweGjcVKvNVpckjclKnwTOiFrNUT/xByS7kkwkmZienl7RxUmS/t9SA+DpdmiHdn+s1aeAjUPjLgSOzFE/QVXtqaqtVbV1/fr1S1yeJGk+Sw2AfcDMJ3l2Ap8fqr+lfRroMuAH7RDRvcC2JOvayd9trSZJGpMz5xuQ5FPAa4Fzk0wx+DTPrcBdSd4KfAd4Uxu+H9gOTAI/Bq4HqKpnkrwfeLCNe19VHX9iWZK0iuYNgKq6bpaXrhgxtoAbZvk5e4G9i1qdJOmk8ZvAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ06c9wL0MratPuecS9B0inCPQBJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjq1rABIcijJo0keTjLRauckOZDkyXa/rtWT5CNJJpM8kuSSlWhAkrQ0K7EH8Lqq2lJVW9vz3cB9VbUZuK89B7ga2Nxuu4DbVuC9JUlLdDIOAe0A7miP7wDeMFT/ZA3cD5yd5IKT8P6SpAVYbgAU8MUkDyXZ1WrnV9VRgHZ/XqtvAA4PzZ1qtedIsivJRJKJ6enpZS5PkjSb5V4L6PKqOpLkPOBAkm/NMTYjanVCoWoPsAdg69atJ7wuSVoZy9oDqKoj7f4Y8DngUuDpmUM77f5YGz4FbByafiFwZDnvL0lauiUHQJIXJHnRzGNgG/AYsA/Y2YbtBD7fHu8D3tI+DXQZ8IOZQ0WSpNW3nENA5wOfSzLzc/6+qv4pyYPAXUneCnwHeFMbvx/YDkwCPwauX8Z7SxqDcV1u/NCt14zlfU93Sw6AqnoK+NUR9f8CrhhRL+CGpb6fJGll+U1gSeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1asn/Kbxmt2n3PeNegiTNyz0ASeqUASBJnTIAJKlTngOQtOaN87zaoVuvGdt7n2zuAUhSpwwASeqUASBJnVr1AEhyVZInkkwm2b3a7y9JGljVAEhyBvBR4GrgYuC6JBev5hokSQOrvQdwKTBZVU9V1U+BO4Edq7wGSRKr/zHQDcDhoedTwKtO1pt5SQZJyzWuv0dW4+Onqx0AGVGr5wxIdgG72tP/SfLESV/VwpwLfG/ci1gBp0sfYC9r0enSB4y5l/z5sqb/0kIGrXYATAEbh55fCBwZHlBVe4A9q7mohUgyUVVbx72O5Tpd+gB7WYtOlz7g9OplNqt9DuBBYHOSi5KcBVwL7FvlNUiSWOU9gKp6NsmNwL3AGcDeqjq4mmuQJA2s+rWAqmo/sH+133cFrLnDUkt0uvQB9rIWnS59wOnVy0ipqvlHSZJOO14KQpI61X0AJDk7yWeSfCvJ40l+I8nNSb6b5OF22z7L3ENJHm1jJlZ77SPWc0Ivrf5H7fIbB5P8xSxz19QlOpbZy5rZLrP8fn166HfrUJKHZ5m75rfJInpZM9ukrWdUL1uS3D+zxiSXzjJ3Z5In223naq99RVVV1zfgDuAP2uOzgLOBm4E/XsDcQ8C54+5hnl5eB3wJeH6rnzdi3hnAvwMvbfO+AVx8Kvay1rbLqD6Oe/2DwJ+dqttkIb2stW0yx+/XF4GrW2078M8j5p0DPNXu17XH68bdz1JvXe8BJHkx8BrgdoCq+mlV/fd4V7U0c/TyduDWqvpJqx8bMX1NXaJjmb2sGfP9fiUJ8HvAp0ZMP1W2yczrc/WypszRSwEvbsN+keO+o9RcCRyoqmeq6vvAAeCqk7/qk6PrAGDwr6tp4K+T/GuSjyd5QXvtxiSPJNmbZN0s8wv4YpKH2jeYx2m2Xl4OvDrJA0n+Jcmvj5g76hIdG07+kme1nF5g7WyXuX6/AF4NPF1VT46Ye6pskxlz9QJrZ5vA7L28E/jLJIeBDwDvGTF3rW2XZek9AM4ELgFuq6pXAv8L7AZuA14GbAGOMti1HeXyqrqEwdVNb0jympO/5FnN1suZDHZVLwP+BLir/Wtt2LyX6Fhly+kF1s52ma2PGdcx+7+YT5VtMmOuXmDtbBOYvZe3A++qqo3Au2h7CMdZa9tlWXoPgClgqqoeaM8/A1xSVU9X1c+q6ufAxxjsjp+gqo60+2PA52Ybt0pG9tLqd9fA14CfM7jGyfFz57xExypbTi9rabvM1gdJzgR+F/j0HHNPhW2ykF7W0jaB2XvZCdzdav/A6DWute2yLF0HQFX9J3A4yS+30hXAN5NcMDTsjcBjx89N8oIkL5p5DGwbNW61zNYL8I/A6wGSvJzBCa/jL3C1pi7RsZxe1tJ2maMPgN8GvlVVU7NMP1W2CczTy1raJjBnL0eA32q11wOjDmfdC2xLsq4dGt7WaqemcZ+FHveNwWGeCeARBn/BrAP+Bni01fYBF7SxLwH2t8cvZfDJjG8AB4E/XaO9nAX8LYM/cF8HXn98L+35duDfGHzy5JTtZa1tl1F9tPongD88buwpt00W0sta2yZz/H79JvBQW+cDwK+1sVuBjw/N/X1gst2uH3cvy7n5TWBJ6lTXh4AkqWcGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnfo/C3IbsDY9p20AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "(array([  21.,  124.,  595., 1645., 2667., 2709., 1544.,  555.,  123.,\n",
       "          17.]),\n",
       " array([-0.38962334, -0.04257622,  0.3044709 ,  0.65151802,  0.99856513,\n",
       "         1.34561225,  1.69265937,  2.03970649,  2.3867536 ,  2.73380072,\n",
       "         3.08084784]),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af7701d30>"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af827a1d0>"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEQNJREFUeJzt3X+s3XV9x/Hna4C4DDPKekEsxaKp2dBsiDfIZFnYUH6UxWI2E/hDqmGpZpBo4j9VE3EaMlz8kRANC45GSBzI/DG6Uccqshj/ANqyCtTKuGAjtQ2t4kBiwga+98f5Vg/l/jj39t5zTvt5PpKT8z3v7+ec8/5+7+153fP5fs9pqgpJUnt+a9QNSJJGwwCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNerYUTcwm+XLl9eqVatG3cbRYfv23vVb3jLaPiTw93GJbd++/adVNTHXuIzzV0FMTk7Wtm3bRt3G0SHpXY/xz1sN8fdxSSXZXlWTc41zCkiSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkho11l8FIY2zVRvuGsnz7r7+0pE8r44+vgOQpEYZAJLUKKeAdEQb1TSMdDTwHYAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUXMGQJKVSe5NsivJziQf7OqfSPKTJDu6y5q++3wkyVSSR5Nc1Fe/uKtNJdmwNJskSRrEIJ8EfgH4cFU9mORVwPYkW7p1n6+qz/QPTnImcDnwRuA1wLeTvKFb/UXgHcAeYGuSTVX1g8XYEEnS/MwZAFW1D9jXLf8iyS5gxSx3WQvcXlXPAz9KMgWc062bqqonAJLc3o01ACRpBOZ1DCDJKuDNwP1d6ZokDyXZmGRZV1sBPNl3tz1dbaa6JGkEBg6AJCcAXwc+VFXPAjcCrwfOovcO4bMHh05z95qlfujzrE+yLcm2AwcODNqeJGmeBgqAJMfRe/H/SlV9A6CqnqqqF6vqV8CX+M00zx5gZd/dTwP2zlJ/iaq6qaomq2pyYmJivtsjSRrQIGcBBbgZ2FVVn+urn9o37F3AI93yJuDyJMcnOQNYDTwAbAVWJzkjySvoHSjetDibIUmar0HOAjoPeA/wcJIdXe2jwBVJzqI3jbMbeD9AVe1Mcge9g7svAFdX1YsASa4B7gaOATZW1c5F3BZJ0jwMchbQ95h+/n7zLPe5Drhumvrm2e4nSRoePwksSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktSoOQMgycok9ybZlWRnkg929ZOSbEnyWHe9rKsnyQ1JppI8lOTsvsda141/LMm6pdssSdJcBnkH8ALw4ar6A+Bc4OokZwIbgHuqajVwT3cb4BJgdXdZD9wIvcAArgXeCpwDXHswNCRJwzdnAFTVvqp6sFv+BbALWAGsBW7pht0CXNYtrwVurZ77gBOTnApcBGypqqer6ufAFuDiRd0aSdLA5nUMIMkq4M3A/cApVbUPeiEBnNwNWwE82Xe3PV1tprokaQQGDoAkJwBfBz5UVc/ONnSaWs1SP/R51ifZlmTbgQMHBm1PkjRPAwVAkuPovfh/paq+0ZWf6qZ26K73d/U9wMq+u58G7J2l/hJVdVNVTVbV5MTExHy2RZI0D4OcBRTgZmBXVX2ub9Um4OCZPOuAO/vqV3ZnA50LPNNNEd0NXJhkWXfw98KuJkkagWMHGHMe8B7g4SQ7utpHgeuBO5JcBfwYeHe3bjOwBpgCfgm8D6Cqnk7yKWBrN+6TVfX0omyFJGne5gyAqvoe08/fA1wwzfgCrp7hsTYCG+fToCRpafhJYElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1KhjR92Ajg6rNtw16hYkzZPvACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKj5gyAJBuT7E/ySF/tE0l+kmRHd1nTt+4jSaaSPJrkor76xV1tKsmGxd8USdJ8DPI5gC8DXwBuPaT++ar6TH8hyZnA5cAbgdcA307yhm71F4F3AHuArUk2VdUPDqN3qUmj/MzF7usvHdlza/HNGQBV9d0kqwZ8vLXA7VX1PPCjJFPAOd26qap6AiDJ7d1YA0CSRuRwjgFck+ShbopoWVdbATzZN2ZPV5upLkkakYUGwI3A64GzgH3AZ7t6phlbs9RfJsn6JNuSbDtw4MAC25MkzWVBAVBVT1XVi1X1K+BL/GaaZw+wsm/oacDeWerTPfZNVTVZVZMTExMLaU+SNIAFBUCSU/tuvgs4eIbQJuDyJMcnOQNYDTwAbAVWJzkjySvoHSjetPC2JUmHa86DwEluA84HlifZA1wLnJ/kLHrTOLuB9wNU1c4kd9A7uPsCcHVVvdg9zjXA3cAxwMaq2rnoWyNJGtggZwFdMU355lnGXwdcN019M7B5Xt1JkpaMnwSWpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjZozAJJsTLI/ySN9tZOSbEnyWHe9rKsnyQ1JppI8lOTsvvus68Y/lmTd0myOJGlQg7wD+DJw8SG1DcA9VbUauKe7DXAJsLq7rAduhF5gANcCbwXOAa49GBqSpNGYMwCq6rvA04eU1wK3dMu3AJf11W+tnvuAE5OcClwEbKmqp6vq58AWXh4qkqQhWugxgFOqah9Ad31yV18BPNk3bk9Xm6n+MknWJ9mWZNuBAwcW2J4kaS6LfRA409RqlvrLi1U3VdVkVU1OTEwsanOSpN9YaAA81U3t0F3v7+p7gJV9404D9s5SlySNyEIDYBNw8EyedcCdffUru7OBzgWe6aaI7gYuTLKsO/h7YVeTJI3IsXMNSHIbcD6wPMkeemfzXA/ckeQq4MfAu7vhm4E1wBTwS+B9AFX1dJJPAVu7cZ+sqkMPLEuShmjOAKiqK2ZYdcE0Ywu4eobH2QhsnFd3kqQl4yeBJalRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRx466AS2uVRvumra+e471ktrjOwBJapQBIEmNcgpI0sAWawpx9zwfb/f1ly7K8+qlDusdQJLdSR5OsiPJtq52UpItSR7rrpd19SS5IclUkoeSnL0YGyBJWpjFmAL6s6o6q6omu9sbgHuqajVwT3cb4BJgdXdZD9y4CM8tSVqgpTgGsBa4pVu+Bbisr35r9dwHnJjk1CV4fknSAA43AAr4jyTbk6zvaqdU1T6A7vrkrr4CeLLvvnu6miRpBA73IPB5VbU3ycnAliQ/nGVspqnVywb1gmQ9wOmnn36Y7UmSZnJY7wCqam93vR/4JnAO8NTBqZ3uen83fA+wsu/upwF7p3nMm6pqsqomJyYmDqc9SdIsFhwASX4nyasOLgMXAo8Am4B13bB1wJ3d8ibgyu5soHOBZw5OFUmShu9wpoBOAb6Z5ODj/FNV/XuSrcAdSa4Cfgy8uxu/GVgDTAG/BN53GM8tSTpMCw6AqnoC+KNp6j8DLpimXsDVC30+SdLi8qsgJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWrUsaNu4Gi0asNdo25BOqqM8t/U7usvHdlzLzXfAUhSowwASWqUASBJjTIAJKlRBoAkNWroAZDk4iSPJplKsmHYzy9J6hlqACQ5BvgicAlwJnBFkjOH2YMkqWfYnwM4B5iqqicAktwOrAV+sBRP5vn4kjSzYQfACuDJvtt7gLcOuQdJGtio/pAcxgfQhh0AmaZWLxmQrAfWdzefS/LoNPdZDvx0kXtbSiPv99c7/tN/Mcjwkfc7T0dSv0dSr7BE/c7z93E+jpr9m08f1uO+dpBBww6APcDKvtunAXv7B1TVTcBNsz1Ikm1VNbn47S0N+11aR1K/R1KvYL9LbdT9DvssoK3A6iRnJHkFcDmwacg9SJIY8juAqnohyTXA3cAxwMaq2jnMHiRJPUP/NtCq2gxsPsyHmXWKaAzZ79I6kvo9knoF+11qI+03VTX3KEnSUcevgpCkRh0RAZDkpCRbkjzWXS+bYdyLSXZ0l6EfXJ7ray6SHJ/kq936+5OsGnaPfb3M1et7kxzo259/PYo++/rZmGR/kkdmWJ8kN3Tb81CSs4fd4yH9zNXv+Ume6du/Hx92j329rExyb5JdSXYm+eA0Y8Zm/w7Y7zjt31cmeSDJ97t+/3aaMaN5baiqsb8Afw9s6JY3AJ+eYdxzI+zxGOBx4HXAK4DvA2ceMuZvgH/oli8HvjrGvb4X+MKof/Z9/fwpcDbwyAzr1wDfoneK+bnA/WPe7/nAv416v3a9nAqc3S2/CvjvaX4fxmb/DtjvOO3fACd0y8cB9wPnHjJmJK8NR8Q7AHpfF3FLt3wLcNkIe5nJr7/moqr+Fzj4NRf9+rfja8AFSab7cNxSG6TXsVJV3wWenmXIWuDW6rkPODHJqcPp7uUG6HdsVNW+qnqwW/4FsIvep/b7jc3+HbDfsdHts+e6m8d1l0MPvo7kteFICYBTqmof9H74wMkzjHtlkm1J7ksy7JCY7msuDv2l/PWYqnoBeAb4vaF0N0Mfnel6BfjL7u3+15KsnGb9OBl0m8bJH3fTAt9K8sZRNwPQTT28md5fqf3Gcv/O0i+M0f5NckySHcB+YEtVzbh/h/naMDb/KXySbwOvnmbVx+bxMKdX1d4krwO+k+Thqnp8cTqc05xfczHgmGEYpI9/BW6rqueTfIDeXyd/vuSdLdy47NtBPQi8tqqeS7IG+Bdg9SgbSnIC8HXgQ1X17KGrp7nLSPfvHP2O1f6tqheBs5KcCHwzyZuqqv/40Ej279i8A6iqt1fVm6a53Ak8dfDtZne9f4bH2NtdPwH8J72/DIZlzq+56B+T5FjgdxnNNMEgX8nxs6p6vrv5JeAtQ+ptoQbZ/2Ojqp49OC1Qvc/GHJdk+aj6SXIcvRfTr1TVN6YZMlb7d65+x23/HlRV/0PvteniQ1aN5LVhbAJgDpuAdd3yOuDOQwckWZbk+G55OXAeS/Q10zMY5Gsu+rfjr4DvVHfUZ8jm7PWQ+d130ptnHWebgCu7s1XOBZ45OG04jpK8+uAcb5Jz6P1b/NmIeglwM7Crqj43w7Cx2b+D9Dtm+3ei+8ufJL8NvB344SHDRvPaMMyj4Qu90JsLuwd4rLs+qatPAv/YLb8NeJjeGS0PA1eNoM819M5IeBz4WFf7JPDObvmVwD8DU8ADwOtGuE/n6vXvgJ3d/rwX+P0R/w7cBuwD/o/eX0tXAR8APtCtD73/bOjx7uc/Oeb9XtO3f+8D3jbCXv+E3nTDQ8CO7rJmXPfvgP2O0/79Q+C/un4fAT7e1Uf+2uAngSWpUUfKFJAkaZEZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNer/AULcCPyLzkIkAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "lower, upper : (0.3933301883768994, 2.2627035091189343)\n",
      "upper-lower : 1.8693733207420349\n"
     ]
    }
   ],
   "source": [
    "# pv(\"sample_data.mean()\")\n",
    "# pv(\"sample_data.std()\")\n",
    "# pv(\"sample_data.var()\")\n",
    "\n",
    "boot_diffs = []\n",
    "drink_means = []\n",
    "no_drink_means = []\n",
    "for _ in range(10_000):\n",
    "    nb_sample = 200\n",
    "    boot_sample = sample_data.sample(nb_sample, replace=True)\n",
    "    \n",
    "    drink_mean = boot_sample.query(\"drinks_coffee == True\")['height'].mean()\n",
    "    drink_means.append(drink_mean)\n",
    "    no_drink_mean = boot_sample.query(\"drinks_coffee == False\")['height'].mean()\n",
    "    no_drink_means.append(no_drink_mean)\n",
    "    \n",
    "    boot_diffs.append(drink_mean - no_drink_mean)\n",
    "\n",
    "# pv(\"len(drink_means)\")\n",
    "# pv(\"len(no_drink_means)\")\n",
    "# pv(\"len(boot_diffs)\")\n",
    "\n",
    "\n",
    "pv(\"np.mean(drink_means)\")\n",
    "pv(\"np.mean(no_drink_means)\")\n",
    "pv(\"np.mean(boot_diffs)\")\n",
    "\n",
    "pv(\"np.std(drink_means)\")\n",
    "pv(\"np.std(no_drink_means)\")\n",
    "pv(\"np.std(boot_diffs)\")\n",
    "\n",
    "plt.hist(drink_means)\n",
    "plt.show()\n",
    "\n",
    "plt.hist(no_drink_means)\n",
    "plt.show()\n",
    "\n",
    "lower, upper = np.percentile(boot_diffs, 2.5), np.percentile(boot_diffs, 97.5) # np.percentile(mb, [2.5, 97.5])\n",
    "plt.hist(boot_diffs);\n",
    "plt.axvline(x=lower, color='r', linewidth=2);\n",
    "plt.axvline(x=upper, color='r', linewidth=2);\n",
    "plt.show()\n",
    "\n",
    "pv(\"lower, upper\")\n",
    "pv(\"upper-lower\")\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 26 - Drawing Conclusions - Calculating Errors\n",
    "\n",
    "Here are two datasets that represent two of the examples you have seen in this lesson.  \n",
    "\n",
    "One dataset is based on the parachute example, and the second is based on the judicial example.  Neither of these datasets is based on real people.\n",
    "\n",
    "Use the exercises below to assist in answering the quiz questions at the bottom of this page."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 66,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "jud_data = pd.read_csv('judicial_dataset_predictions.csv')\n",
    "par_data = pd.read_csv('parachute_dataset.csv')\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`1.` Above, you can see the actual and predicted columns for each of the datasets.  Using the **jud_data**, find the proportion of errors for the dataset, and furthermore, the percentage of errors of each type.  Use the results to answer the questions in quiz 1 below.  \n",
    "\n",
    "**Hint for quiz:** an error is any time the prediction doesn't match an actual value.  Additionally, there are Type I and Type II errors to think about.  We also know we can minimize one type of error by maximizing the other type of error.  If we predict all individuals as innocent, how many of the guilty are incorrectly labeled?  Similarly, if we predict all individuals as guilty, how many of the innocent are incorrectly labeled?\n",
    "\n",
    "Total percentage of errors\n",
    "\n",
    "Percentage of Type I Errors\n",
    "\n",
    "Percentage of Type II Errors\n",
    "\n",
    "If everyone was predicted to be guilty, the proportion of Type I Errors made.\n",
    "\n",
    "If everyone was predicted to be guilty, the proportion of Type II Errors made.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 67,
   "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>defendant_id</th>\n",
       "      <th>actual</th>\n",
       "      <th>predicted</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>22574</td>\n",
       "      <td>innocent</td>\n",
       "      <td>innocent</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>35637</td>\n",
       "      <td>innocent</td>\n",
       "      <td>innocent</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>39919</td>\n",
       "      <td>innocent</td>\n",
       "      <td>innocent</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>29610</td>\n",
       "      <td>guilty</td>\n",
       "      <td>guilty</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>38273</td>\n",
       "      <td>innocent</td>\n",
       "      <td>innocent</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   defendant_id    actual predicted\n",
       "0         22574  innocent  innocent\n",
       "1         35637  innocent  innocent\n",
       "2         39919  innocent  innocent\n",
       "3         29610    guilty    guilty\n",
       "4         38273  innocent  innocent"
      ]
     },
     "execution_count": 67,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "nb_errors/nb_cases : 0.042152958945489497\n",
      "errors_type_1.shape[0]/nb_cases : 0.001510366607167376\n",
      "errors_type_2.shape[0]/nb_cases : 0.04064259233832212\n",
      "errors_type_1.shape[0]/nb_cases : 0.45159961554304545\n"
     ]
    }
   ],
   "source": [
    "jud_data.head()\n",
    "nb_cases = jud_data.shape[0]\n",
    "errors = jud_data.query(\"actual != predicted\")\n",
    "nb_errors = errors.shape[0]\n",
    "pv(\"nb_errors/nb_cases\")\n",
    "errors_type_1 = jud_data.query(\"actual == 'innocent' and predicted == 'guilty'\")\n",
    "pv(\"errors_type_1.shape[0]/nb_cases\")\n",
    "errors_type_2 = jud_data.query(\"actual == 'guilty' and predicted == 'innocent'\")\n",
    "pv(\"errors_type_2.shape[0]/nb_cases\")\n",
    "\n",
    "errors_type_1 = jud_data.query(\"actual == 'innocent'\")\n",
    "pv(\"errors_type_1.shape[0]/nb_cases\")\n",
    "errors_type_2 = 0 # because no one is innocent !"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`2.` Using the **par_data**, find the proportion of errors for the dataset, and furthermore, the percentage of errors of each type.  Use the results to answer the questions in quiz 2 below.\n",
    "\n",
    "These should be very similar operations to those you performed in the previous question.\n",
    "Total percentage of errors\n",
    "\n",
    "Percentage of Type I Errors\n",
    "\n",
    "Percentage of Type II Errors\n",
    "\n",
    "If every parachute was predicted to not open, the proportion of Type I Errors made.\n",
    "\n",
    "If every parachute was predicted to not open, the proportion of Type II Errors made.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array(['opens', 'fails'], dtype=object)"
      ]
     },
     "execution_count": 68,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "nb_errors/nb_cases : 0.039972551037913875\n",
      "errors_type_1.shape[0]/nb_cases : 0.00017155601303825698\n",
      "errors_type_2.shape[0]/nb_cases : 0.03980099502487562\n",
      "errors_type_2.shape[0]/nb_cases : 0.9917653113741637\n"
     ]
    }
   ],
   "source": [
    "par_data.actual.unique()\n",
    "\n",
    "nb_cases = par_data.shape[0]\n",
    "errors = par_data.query(\"actual != predicted\")\n",
    "nb_errors = errors.shape[0]\n",
    "pv(\"nb_errors/nb_cases\")\n",
    "errors_type_1 = par_data.query(\"actual == 'fails' and predicted == 'opens'\")\n",
    "pv(\"errors_type_1.shape[0]/nb_cases\")\n",
    "errors_type_2 = par_data.query(\"actual == 'opens' and predicted == 'fails'\")\n",
    "pv(\"errors_type_2.shape[0]/nb_cases\")\n",
    "\n",
    "errors_type_2 = par_data.query(\"actual == 'opens'\")\n",
    "pv(\"errors_type_2.shape[0]/nb_cases\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Impact of Large Sample Sizes\n",
    "\n",
    "When we increase our sample size, even the smallest of differences may seem significant.  \n",
    "\n",
    "To illustrate this point, work through this notebook and the quiz questions that follow below.\n",
    "\n",
    "Start by reading in the libraries and data.\n",
    "\n",
    "`1.` In this case, imagine we are interested in testing if the mean height of all individuals in `full_data` is equal to 67.60 inches or **different** than 67.60.  First, use **quiz 1** below to state the null and alternative hypotheses for these cases in the cell below. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$H_0: \\mu = 67.60$$\n",
    "\n",
    "$$H_1: \\mu \\neq 67.60$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "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>user_id</th>\n",
       "      <th>age</th>\n",
       "      <th>drinks_coffee</th>\n",
       "      <th>height</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>4509</td>\n",
       "      <td>&lt;21</td>\n",
       "      <td>False</td>\n",
       "      <td>64.538179</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1864</td>\n",
       "      <td>&gt;=21</td>\n",
       "      <td>True</td>\n",
       "      <td>65.824249</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   user_id   age  drinks_coffee     height\n",
       "0     4509   <21          False  64.538179\n",
       "1     1864  >=21           True  65.824249"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "\n",
    "full_data = pd.read_csv('coffee_dataset.csv')\n",
    "\n",
    "full_data.head(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`2.` Find:\n",
    "* What is the population mean height?  \n",
    "* What is the standard deviation of the population heights?  Create a sample set of data using the code below.  \n",
    "* What is the sample mean height?  Simulate the sampling distribution for the mean of five values to see the shape and plot a histogram. \n",
    "* What is the standard deviation of the sampling distribution of the mean of five draws? \n",
    "\n",
    "Use **quiz 2** below to assure your answers are correct.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Simulating from the null\n",
    "1. Assume the null is true, we know what the sampling distribuntion would look like if we simulate from the closest value under the null to the alternative (67.60 ?)\n",
    "2. We simulate a normal distribution to get the standard deviation of the sampling distribution (for loop, df.sample(replace=True), append)\n",
    "3. We know it must follow a central distribution if the null is true (central limit theorem)\n",
    "4. We use NumPy to simulate draws from the normal using the hypothesized mean and the standard deviation of the sampling  distribution : np.random.normal(loc=null_hypothesis, scale=np.std(means_sampling_distribution), size = number of sample range)\n",
    "5. Each of the simulated draws represents a possible mean from the null hypothesis\n",
    "6. Where the sample mean falls in this distribution ? Below or above the null ?\n",
    "\n",
    "#### Calculate Confidence Interval\n",
    "Determine a confidence interval and degrees of freedom\n",
    "This is what we call alpha (α). The typical value of α is 0.05. This means that there is 95% confidence that the conclusion of this test will be valid. The degree of freedom can be calculated by the the following formula:\n",
    "alpha = 0.05\n",
    "degree_freedom = 1-alpha\n",
    "\n",
    "ordered = sort(statistics)\n",
    "lower = percentile(ordered, alpha/2)\n",
    "upper = percentile(ordered, (1-alpha) + alpha/2))\n",
    "\n",
    "=> Lower et upper sont les p-value ?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "np.mean(full_data.height) : 67.59748697307937\n",
      "np.std(full_data.height) : 3.1189087116821472\n"
     ]
    },
    {
     "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>user_id</th>\n",
       "      <th>age</th>\n",
       "      <th>drinks_coffee</th>\n",
       "      <th>height</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>2402</th>\n",
       "      <td>2874</td>\n",
       "      <td>&lt;21</td>\n",
       "      <td>True</td>\n",
       "      <td>64.357154</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2864</th>\n",
       "      <td>3670</td>\n",
       "      <td>&gt;=21</td>\n",
       "      <td>True</td>\n",
       "      <td>66.859636</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2167</th>\n",
       "      <td>7441</td>\n",
       "      <td>&lt;21</td>\n",
       "      <td>False</td>\n",
       "      <td>66.659561</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>507</th>\n",
       "      <td>2781</td>\n",
       "      <td>&gt;=21</td>\n",
       "      <td>True</td>\n",
       "      <td>70.166241</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1817</th>\n",
       "      <td>2875</td>\n",
       "      <td>&gt;=21</td>\n",
       "      <td>True</td>\n",
       "      <td>71.369120</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      user_id   age  drinks_coffee     height\n",
       "2402     2874   <21           True  64.357154\n",
       "2864     3670  >=21           True  66.859636\n",
       "2167     7441   <21          False  66.659561\n",
       "507      2781  >=21           True  70.166241\n",
       "1817     2875  >=21           True  71.369120"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "(5, 4)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "np.mean(sample1.height) : 67.88234252049084\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(array([   9.,   65.,  346., 1016., 2096., 2852., 2257., 1067.,  260.,\n",
       "          32.]),\n",
       " array([62.0966947 , 63.10054507, 64.10439543, 65.1082458 , 66.11209617,\n",
       "        67.11594654, 68.1197969 , 69.12364727, 70.12749764, 71.131348  ,\n",
       "        72.13519837]),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD9CAYAAAC1DKAUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEbVJREFUeJzt3X+s3Xddx/Hny5WhTmUFujm7agcWcPxhmdcxQ1RkunWbYTORpIuBBpdUcTNi/FU0cfiDZPhrkQSXVFYpBpnll2ukOOpAkT/2o5tj0I1l11HppXWrdgyVOFJ8+8f5VM66e3vPvT33Hm4/z0dycs73/f18z/fzyW3P634/3+/53lQVkqT+fNOkOyBJmgwDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpU/MGQJJvTnJPkk8n2Z/kd1r9giR3J3k0yV8nObPVn9uWp9v69UPv9ZZWfyTJ5Us1KEnS/EY5AngaeE1VfT+wEdiU5BLg7cDNVbUBeBK4rrW/Dniyqr4XuLm1I8mFwGbg5cAm4M+SnDHOwUiSRjdvANTAf7XF57RHAa8BPtDqO4Fr2uur2zJt/aVJ0uq3VdXTVfV5YBq4eCyjkCQt2EjnAJKckeQB4AlgL/AvwJeq6lhrMgOsba/XAgcB2vqngBcM12fZRpK0zFaN0qiqvgZsTHI28GHg+2Zr1p4zx7q56s+QZCuwFeCss876gZe97GWjdFGS1Nx3333/XlVr5ms3UgAcV1VfSvIPwCXA2UlWtd/yzwcOtWYzwDpgJskq4HnA0aH6ccPbDO9jO7AdYGpqqvbt27eQLkpS95L86yjtRrkKaE37zZ8k3wL8OPAw8Angp1uzLcDt7fXutkxb//Ea3HFuN7C5XSV0AbABuGe04UiSxm2UI4DzgJ3tip1vAnZV1d8meQi4LcnvA/8M3Nra3wr8ZZJpBr/5bwaoqv1JdgEPAceA69vUkiRpAvKNfDtop4AkaeGS3FdVU/O185vAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVML+iawpK9bv+0jE9nvgZuumsh+dfrxCECSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpU/MGQJJ1ST6R5OEk+5P8Uqu/NckXkzzQHlcObfOWJNNJHkly+VB9U6tNJ9m2NEOSJI1i1QhtjgG/UlX3J/l24L4ke9u6m6vqj4YbJ7kQ2Ay8HPgu4O+TvKStfifwE8AMcG+S3VX10DgGIklamHkDoKoOA4fb6/9M8jCw9iSbXA3cVlVPA59PMg1c3NZNV9VjAElua20NAEmagAWdA0iyHngFcHcr3ZDkwSQ7kqxutbXAwaHNZlptrvqJ+9iaZF+SfUeOHFlI9yRJCzByACT5NuCDwJur6svALcCLgY0MjhD++HjTWTavk9SfWajaXlVTVTW1Zs2aUbsnSVqgUc4BkOQ5DD7831tVHwKoqseH1v858LdtcQZYN7T5+cCh9nquuiRpmY1yFVCAW4GHq+pPhurnDTX7KeCz7fVuYHOS5ya5ANgA3APcC2xIckGSMxmcKN49nmFIkhZqlCOAVwGvBz6T5IFW+03g2iQbGUzjHAB+DqCq9ifZxeDk7jHg+qr6GkCSG4A7gDOAHVW1f4xjkSQtwChXAX2K2efv95xkm7cBb5ulvudk20mSlo/fBJakThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSp0b6ewCSvnGs3/aRie37wE1XTWzfGj+PACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE55LyCtaJO8L4600nkEIEmdMgAkqVPzBkCSdUk+keThJPuT/FKrPz/J3iSPtufVrZ4k70gyneTBJBcNvdeW1v7RJFuWbliSpPmMcgRwDPiVqvo+4BLg+iQXAtuAO6tqA3BnWwa4AtjQHluBW2AQGMCNwCuBi4Ebj4eGJGn5zRsAVXW4qu5vr/8TeBhYC1wN7GzNdgLXtNdXA++pgbuAs5OcB1wO7K2qo1X1JLAX2DTW0UiSRragcwBJ1gOvAO4Gzq2qwzAICeCc1mwtcHBos5lWm6suSZqAkQMgybcBHwTeXFVfPlnTWWp1kvqJ+9maZF+SfUeOHBm1e5KkBRopAJI8h8GH/3ur6kOt/Hib2qE9P9HqM8C6oc3PBw6dpP4MVbW9qqaqamrNmjULGYskaQFGuQoowK3Aw1X1J0OrdgPHr+TZAtw+VH9DuxroEuCpNkV0B3BZktXt5O9lrSZJmoBRvgn8KuD1wGeSPNBqvwncBOxKch3wBeB1bd0e4EpgGvgK8EaAqjqa5PeAe1u7362qo2MZhSRpweYNgKr6FLPP3wNcOkv7Aq6f4712ADsW0kFJ0tLwm8CS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqfmDYAkO5I8keSzQ7W3Jvlikgfa48qhdW9JMp3kkSSXD9U3tdp0km3jH4okaSFGOQJ4N7BplvrNVbWxPfYAJLkQ2Ay8vG3zZ0nOSHIG8E7gCuBC4NrWVpI0Iavma1BVn0yyfsT3uxq4raqeBj6fZBq4uK2brqrHAJLc1to+tOAeS5LG4lTOAdyQ5ME2RbS61dYCB4fazLTaXHVJ0oQsNgBuAV4MbAQOA3/c6pmlbZ2k/ixJtibZl2TfkSNHFtk9SdJ8FhUAVfV4VX2tqv4X+HO+Ps0zA6wbano+cOgk9dnee3tVTVXV1Jo1axbTPUnSCBYVAEnOG1r8KeD4FUK7gc1JnpvkAmADcA9wL7AhyQVJzmRwonj34rstSTpV854ETvI+4NXAC5PMADcCr06ykcE0zgHg5wCqan+SXQxO7h4Drq+qr7X3uQG4AzgD2FFV+8c+GknSyEa5CujaWcq3nqT924C3zVLfA+xZUO8kSUvGbwJLUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnq1KpJd0DSyrF+20cmst8DN101kf2e7uY9AkiyI8kTST47VHt+kr1JHm3Pq1s9Sd6RZDrJg0kuGtpmS2v/aJItSzMcSdKoRpkCejew6YTaNuDOqtoA3NmWAa4ANrTHVuAWGAQGcCPwSuBi4MbjoSFJmox5A6CqPgkcPaF8NbCzvd4JXDNUf08N3AWcneQ84HJgb1Udraongb08O1QkSctosecAzq2qwwBVdTjJOa2+Fjg41G6m1eaq6zQxqblhSYs37quAMkutTlJ/9hskW5PsS7LvyJEjY+2cJOnrFhsAj7epHdrzE60+A6wbanc+cOgk9Wepqu1VNVVVU2vWrFlk9yRJ81lsAOwGjl/JswW4faj+hnY10CXAU22q6A7gsiSr28nfy1pNkjQh854DSPI+4NXAC5PMMLia5yZgV5LrgC8Ar2vN9wBXAtPAV4A3AlTV0SS/B9zb2v1uVZ14YlmStIzmDYCqunaOVZfO0raA6+d4nx3AjgX1TpK0ZLwVhCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjp1SgGQ5ECSzyR5IMm+Vnt+kr1JHm3Pq1s9Sd6RZDrJg0kuGscAJEmLM44jgB+rqo1VNdWWtwF3VtUG4M62DHAFsKE9tgK3jGHfkqRFWoopoKuBne31TuCaofp7auAu4Owk5y3B/iVJIzjVACjgY0nuS7K11c6tqsMA7fmcVl8LHBzadqbVniHJ1iT7kuw7cuTIKXZPkjSXVae4/auq6lCSc4C9ST53kraZpVbPKlRtB7YDTE1NPWu9JGk8TukIoKoOtecngA8DFwOPH5/aac9PtOYzwLqhzc8HDp3K/iVJi7foAEhyVpJvP/4auAz4LLAb2NKabQFub693A29oVwNdAjx1fKpIkrT8TmUK6Fzgw0mOv89fVdXfJbkX2JXkOuALwOta+z3AlcA08BXgjaewb0nSKVp0AFTVY8D3z1L/D+DSWeoFXL/Y/UmSxstvAktSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROneqtIPQNZv22j0y6C5JWCI8AJKlTBoAkdcoAkKROGQCS1CkDQJI65VVAkr7hTfLqtgM3XTWxfS81jwAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI65d1Al4B/l1fSSrDsRwBJNiV5JMl0km3LvX9J0sCyHgEkOQN4J/ATwAxwb5LdVfXQcvZDkkY1qSP65fg7BMt9BHAxMF1Vj1XVV4HbgKuXuQ+SJJb/HMBa4ODQ8gzwyqXamXPxkjS35Q6AzFKrZzRItgJb2+J/JXlkgft4IfDvi+jbStXbeMEx96C38cIJY87bT+m9vmeURssdADPAuqHl84FDww2qajuwfbE7SLKvqqYWu/1K09t4wTH3oLfxwmTGvNznAO4FNiS5IMmZwGZg9zL3QZLEMh8BVNWxJDcAdwBnADuqav9y9kGSNLDsXwSrqj3AniXcxaKnj1ao3sYLjrkHvY0XJjDmVNX8rSRJpx3vBSRJnVrRAZDk7CQfSPK5JA8n+aEkf9iWH0zy4SRnT7qf4zTbmIfW/WqSSvLCSfZxnOYab5JfbLcU2Z/kDybdz3Ga49/1xiR3JXkgyb4kF0+6n+OS5KVtXMcfX07y5iTPT7I3yaPtefWk+zoOJxnvsn92regpoCQ7gX+qqne1q4q+lcG3jT/eTji/HaCqfmOS/Ryn2cZcVV9Ksg54F/Ay4Aeq6rS4hnqOn/ErgN8Crqqqp5OcU1VPTLSjYzTHmHcBN1fVR5NcCfx6Vb16kv1cCu12MV9k8AXR64GjVXVTu2/Y6tPp/zI8a7wvZZk/u1bsEUCS7wB+BLgVoKq+WlVfqqqPVdWx1uwuBt81OC3MNea2+mbg1znhi3Ur2UnG+ybgpqp6utVPpw//ucZcwHe0Zs/jhO/PnEYuBf6lqv6VwW1idrb6TuCaifVq6fz/eCfx2bViAwB4EXAE+Isk/5zkXUnOOqHNzwIfXf6uLZlZx5zktcAXq+rTE+7fuM31M34J8MNJ7k7yj0l+cLLdHKu5xvxm4A+THAT+CHjLJDu5hDYD72uvz62qwwDt+ZyJ9WrpDI932LJ8dq3kAFgFXATcUlWvAP4b+P/bSyf5LeAY8N7JdG9JzDbmtzKYDvntCfZrqcz1M14FrAYuAX4N2JVkttuMrERzjflNwC9X1Trgl2lHCKeTNt31WuD9k+7LcphrvMv52bWSA2AGmKmqu9vyBxj8xyHJFuAngZ+plXyS49nmGvMFwKeTHGBw2Hh/ku+cTBfHaq7xzgAfqoF7gP9lcB+V08FcY94CfKjV3s/gXNfp5grg/qp6vC0/nuQ8gPZ82kz1NSeOd9k/u1ZsAFTVvwEHk7y0lS4FHkqyCfgN4LVV9ZWJdXAJzDHm+6vqnKpaX1XrGXyAXNTarmhz/YyBvwFeA5DkJcCZnCY3DjvJmA8BP9pqrwEenUD3ltq1PHM6ZDeD4KM9377sPVpazxjvJD67VvpVQBsZXPlyJvAY8EYG9xt6LvAfrdldVfXzk+nh+M025qp6cmj9AWDqNLoKaLaf8X8DO4CNwFeBX62qj0+sk2M2x5hfDvwpgymi/wF+oarum1gnxyzJtzK4VfyLquqpVnsBg6ufvhv4AvC6qjo6uV6OzxzjnWaZP7tWdABIkhZvxU4BSZJOjQEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKn/g9WGU7sVwU8qgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "np.mean(means) : 67.59984472303822\n",
      "np.std(means) : 1.3873129885457822\n"
     ]
    }
   ],
   "source": [
    "np.random.seed(42)\n",
    "pv(\"np.mean(full_data.height)\")\n",
    "pv(\"np.std(full_data.height)\")\n",
    "nb_sample = 5\n",
    "sample1 = full_data.sample(nb_sample)\n",
    "sample1\n",
    "sample1.shape\n",
    "pv(\"np.mean(sample1.height)\")\n",
    "null_hypothesis = 67.60\n",
    "sampling_times = 10_000\n",
    "means = []\n",
    "for _ in range(sampling_times):\n",
    "    boot_sample = full_data.sample(nb_sample)\n",
    "    means.append(boot_sample.height.mean())\n",
    "plt.hist(means);\n",
    "plt.show();\n",
    "pv(\"np.mean(means)\")\n",
    "pv(\"np.std(means)\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`3.` Using the null and alternative hypotheses set up in question 1 and the results of your sampling distribution in question 2, simulate the mean values you would expect from the null hypothesis.  Use these simulated values to determine a p-value to make a decision about your null and alternative hypotheses.  Check your solution using **quiz 3** and **quiz 4** below.\n",
    "\n",
    "**Hint:** Use the numpy documentation [here](https://docs.scipy.org/doc/numpy-1.13.0/reference/generated/numpy.random.normal.html) to assist with your solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "alpha/2  * 100 : 2.5\n",
      "((1-alpha) + alpha/2)  * 100 : 97.5\n",
      "lower, upper : (64.93167771424777, 70.3297488018013)\n",
      "np.mean(means)-np.mean(null_hypothesis_values) : -0.025500615201551113\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af7550390>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af7417b38>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af6d1d668>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af77c9828>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af70594e0>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af7059a90>"
      ]
     },
     "execution_count": 49,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEeFJREFUeJzt3XuwXvVd7/H3Ryg9XqoE2SANwdBzQpX+YaBbitM52pYjV6fgzGFOGKeNtU68gGO9NtUZW63MxCtzOlNxsMSmTi0n1vaQOaRipPXSP7gkSGkDIlsam90gRNNStVOc0K9/PL+cPsC+X56Hze/9mnlmrfVdv7XW75dk9ifr8qydqkKS1J+vG3cHJEnjYQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOnXyuDswl9NPP702btw47m6sXQcODKavfvV4+7EGHDhygPNeNpj/+78f/Hn5x7YA/ht7QTpw4MA/V9XEfO3yQn4VxOTkZO3fv3/c3Vi7ksH0Bfx3/EKRXw2f+L7B/OtfP/jz8o9tAfw39oKU5EBVTc7XzktAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE7NGwBJ/kuSe5N8KsnBJL/a6ucmuSfJo0n+T5JTWv2lbXmqrd84tK93tPojSS5brUFJkua3kFdBPA28oar+LclLgE8m+Rjws8BNVXVbkt8H3grc3KZfqKr/lmQL8BvA/0pyPrAFeBXwcuAvkpxXVc+swrjUiY3b71iZHX394vd9aMdVK3NsaUzmPQOogX9riy9pnwLeAHy41XcB17T5q9sybf0lSdLqt1XV01X1WWAKuGhFRiFJWrQF3QNIclKSB4AngX3APwBfrKrjrck0sL7NrwcOA7T1TwHfOlyfYRtJ0ogtKACq6pmq2gyczeB/7d85U7M2zSzrZqs/S5JtSfYn2X/06NGFdE+StASLegqoqr4I/CVwMXBqkhP3EM4GjrT5aWADQFv/LcCx4foM2wwf45aqmqyqyYmJed9mKklaooU8BTSR5NQ2//XA/wAeBj4B/M/WbCtwe5vf05Zp6z9eg3dO7wG2tKeEzgU2Afeu1EAkSYuzkKeAzgJ2JTmJQWDsrqr/l+Qh4LYkvw78LXBra38r8EdJphj8z38LQFUdTLIbeAg4DlzvE0CSND7zBkBVPQhcMEP9MWZ4iqeqvgJcO8u+bgRuXHw3JUkrzW8CS1KnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKn5g2AJBuSfCLJw0kOJvnpVn9Xks8neaB9rhza5h1JppI8kuSyofrlrTaVZPvqDEmStBAnL6DNceDnqur+JC8DDiTZ19bdVFW/Pdw4yfnAFuBVwMuBv0hyXlv9XuD7gWngviR7quqhlRiIJGlx5g2AqnoceLzN/2uSh4H1c2xyNXBbVT0NfDbJFHBRWzdVVY8BJLmttTUAJGkMFnUPIMlG4ALgnla6IcmDSXYmWddq64HDQ5tNt9psdUnSGCzkEhAASb4J+FPgbVX1pSQ3A+8Gqk1/B/gRIDNsXswcNjXDcbYB2wDOOeechXZPGrmN2+8Yy3EP7bhqLMfVi8+CzgCSvITBD/8PVtVHAKrqiap6pqq+CvwBX7vMMw1sGNr8bODIHPVnqapbqmqyqiYnJiYWOx5J0gIt5CmgALcCD1fV7w7Vzxpq9oPAZ9r8HmBLkpcmORfYBNwL3AdsSnJuklMY3CjeszLDkCQt1kIuAb0WeBPw6SQPtNovAdcl2czgMs4h4McAqupgkt0Mbu4eB66vqmcAktwA3AmcBOysqoMrOBZJ0iIs5CmgTzLzdf29c2xzI3DjDPW9c20nSRodvwksSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ2aNwCSbEjyiSQPJzmY5Kdb/bQk+5I82qbrWj1J3pNkKsmDSS4c2tfW1v7RJFtXb1iSpPmcvIA2x4Gfq6r7k7wMOJBkH/DDwF1VtSPJdmA78HbgCmBT+7wGuBl4TZLTgHcCk0C1/eypqi+s9KA0ehu33zHuLkhapHnPAKrq8aq6v83/K/AwsB64GtjVmu0CrmnzVwMfqIG7gVOTnAVcBuyrqmPth/4+4PIVHY0kacEWdQ8gyUbgAuAe4MyqehwGIQGc0ZqtBw4PbTbdarPVJUljsOAASPJNwJ8Cb6uqL83VdIZazVF/7nG2JdmfZP/Ro0cX2j1J0iItKACSvITBD/8PVtVHWvmJdmmHNn2y1aeBDUObnw0cmaP+LFV1S1VNVtXkxMTEYsYiSVqEhTwFFOBW4OGq+t2hVXuAE0/ybAVuH6q/uT0NdDHwVLtEdCdwaZJ17YmhS1tNkjQGC3kK6LXAm4BPJ3mg1X4J2AHsTvJW4HPAtW3dXuBKYAr4MvAWgKo6luTdwH2t3a9V1bEVGYUkadHmDYCq+iQzX78HuGSG9gVcP8u+dgI7F9NBSdLq8JvAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktSpeQMgyc4kTyb5zFDtXUk+n+SB9rlyaN07kkwleSTJZUP1y1ttKsn2lR+KJGkxFnIG8H7g8hnqN1XV5vbZC5DkfGAL8Kq2ze8lOSnJScB7gSuA84HrWltJ0picPF+DqvrrJBsXuL+rgduq6mngs0mmgIvauqmqegwgyW2t7UOL7rEkaUUs5x7ADUkebJeI1rXaeuDwUJvpVputLkkak6UGwM3AfwU2A48Dv9PqmaFtzVF/niTbkuxPsv/o0aNL7J4kaT5LCoCqeqKqnqmqrwJ/wNcu80wDG4aang0cmaM+075vqarJqpqcmJhYSvckSQuwpABIctbQ4g8CJ54Q2gNsSfLSJOcCm4B7gfuATUnOTXIKgxvFe5bebUnScs17EzjJh4DXAacnmQbeCbwuyWYGl3EOAT8GUFUHk+xmcHP3OHB9VT3T9nMDcCdwErCzqg6u+GgkSQu2kKeArpuhfOsc7W8EbpyhvhfYu6jeSZJWjd8ElqROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkderkcXdA0uJs3H7H2I59aMdVYzu2Vt68ZwBJdiZ5MslnhmqnJdmX5NE2XdfqSfKeJFNJHkxy4dA2W1v7R5NsXZ3hSJIWaiGXgN4PXP6c2nbgrqraBNzVlgGuADa1zzbgZhgEBvBO4DXARcA7T4SGJGk85g2Aqvpr4NhzylcDu9r8LuCaofoHauBu4NQkZwGXAfuq6lhVfQHYx/NDRZI0Qku9CXxmVT0O0KZntPp64PBQu+lWm60uSRqTlX4KKDPUao7683eQbEuyP8n+o0ePrmjnJElfs9QAeKJd2qFNn2z1aWDDULuzgSNz1J+nqm6pqsmqmpyYmFhi9yRJ81lqAOwBTjzJsxW4faj+5vY00MXAU+0S0Z3ApUnWtZu/l7aaJGlM5v0eQJIPAa8DTk8yzeBpnh3A7iRvBT4HXNua7wWuBKaALwNvAaiqY0neDdzX2v1aVT33xrIkaYTmDYCqum6WVZfM0LaA62fZz05g56J6J0laNb4KQpI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpfyPYi8zwb4s6NENNkk7wDECSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1KllBUCSQ0k+neSBJPtb7bQk+5I82qbrWj1J3pNkKsmDSS5ciQFIkpZmJc4AXl9Vm6tqsi1vB+6qqk3AXW0Z4ApgU/tsA25egWNLkpZoNS4BXQ3savO7gGuG6h+ogbuBU5OctQrHlyQtwHIDoIA/T3IgybZWO7OqHgdo0zNafT1weGjb6VaTJI3Bcn8p/Gur6kiSM4B9Sf5ujraZoVbPazQIkm0A55xzzjK7J0mazbLOAKrqSJs+CXwUuAh44sSlnTZ9sjWfBjYMbX42cGSGfd5SVZNVNTkxMbGc7kmS5rDkAEjyjUledmIeuBT4DLAH2NqabQVub/N7gDe3p4EuBp46calIkjR6y7kEdCbw0SQn9vPHVfVnSe4Ddid5K/A54NrWfi9wJTAFfBl4yzKOLUlapiUHQFU9BnzXDPV/AS6ZoV7A9Us9niRpZflNYEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnlvsyOEkd2bj9jmctH5qlvtIO7bhqVfffK88AJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlN8DWAWr/Uy0JK0EzwAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnRp5ACS5PMkjSaaSbB/18SVJAyMNgCQnAe8FrgDOB65Lcv4o+yBJGhj1qyAuAqaq6jGAJLcBVwMPjbgfktaQcb5e5cX86yhHHQDrgcNDy9PAa1brYL6TR5Jml6oa3cGSa4HLqupH2/KbgIuq6qeG2mwDtrXFVwKPjKyDo3U68M/j7sQqe7GP0fGtfS/WMX57VU3M12jUZwDTwIah5bOBI8MNquoW4JZRdmockuyvqslx92M1vdjH6PjWvh7GOJdRPwV0H7ApyblJTgG2AHtG3AdJEiM+A6iq40luAO4ETgJ2VtXBUfZBkjQw8l8IU1V7gb2jPu4L0Iv+Mhcv/jE6vrWvhzHOaqQ3gSVJLxy+CkKSOmUAjEiSU5N8OMnfJXk4yfcMrfv5JJXk9HH2cTlmG1+Sn2qv/jiY5DfH3c+lmml8STYnuTvJA0n2J7lo3P1cqiSvbOM48flSkrclOS3JviSPtum6cfd1KeYY32+1v9MHk3w0yanj7usoeQloRJLsAv6mqt7XnoD6hqr6YpINwPuA7wBeXVVr8pnkmcYHXAD8MnBVVT2d5IyqenKsHV2iWca3G7ipqj6W5ErgF6vqdePs50por2z5PIMvaV4PHKuqHe3dXeuq6u1j7eAyPWd8rwQ+3h5Q+Q2AtT6+xfAMYASSfDPwvcCtAFX1H1X1xbb6JuAXgTWbxHOM7yeAHVX1dKuv1R/+s42vgG9uzb6F53ynZQ27BPiHqvpHBq9q2dXqu4BrxtarlfP/x1dVf15Vx1v9bgbfTeqGATAarwCOAn+Y5G+TvC/JNyZ5I/D5qvrUmPu3XDOODzgP+O9J7knyV0m+e7zdXLLZxvc24LeSHAZ+G3jHODu5grYAH2rzZ1bV4wBtesbYerVyhsc37EeAj424L2NlAIzGycCFwM1VdQHw78C7GFwe+ZUx9mulzDS+7a2+DrgY+AVgd5KMrZdLN9v4fgL4maraAPwM7QxhLWuXt94I/Mm4+7IaZhtfkl8GjgMfHEe/xsUAGI1pYLqq7mnLH2bwA+Vc4FNJDjE49bw/ybeNp4vLMtv4poGP1MC9wFcZvHtlrZltfFuBj7TanzB42+1adwVwf1U90ZafSHIWQJuuyct4Q547PpJsBX4A+KHq7KaoATACVfVPwOEkr2ylSxj8IzyjqjZW1UYGP2QubG3XlFnG9xDwf4E3ACQ5DziFNfjirTnGdwT4vlZ7A/DoGLq30q7j2ZdH9jAIOtr09pH3aGU9a3xJLgfeDryxqr48tl6NiU8BjUiSzQye9jkFeAx4S1V9YWj9IWByDT8F9LzxMbhUshPYDPwH8PNV9fGxdXIZZhnfq4D/zeAS0VeAn6yqA2Pr5DIl+QYGr2t/RVU91WrfyuBpp3OAzwHXVtWx8fVy6WYZ3xTwUuBfWrO7q+rHx9TFkTMAJKlTXgKSpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdeo/AUgEPCOA/APSAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "95.0% confident that the value is between 64.9 and 70.3\n",
      "pval : 1.0081\n",
      "p-value= 0.6226\n"
     ]
    }
   ],
   "source": [
    "\n",
    "null_hypothesis_values = np.random.normal(loc=null_hypothesis, scale=np.std(means), size=sampling_times)\n",
    "\n",
    "\n",
    "alpha = 0.05\n",
    "pv(\"alpha/2  * 100\")\n",
    "pv(\"((1-alpha) + alpha/2)  * 100\")\n",
    "lower = np.percentile(null_hypothesis_values, alpha/2  * 100)\n",
    "upper = np.percentile(null_hypothesis_values, ((1-alpha) + alpha/2)  * 100) # np.percentile(mb, [2.5, 97.5])\n",
    "\n",
    "\n",
    "pv(\"lower, upper\")\n",
    "pv(\"np.mean(means)-np.mean(null_hypothesis_values)\")\n",
    "\n",
    "_ = plt.hist(null_hypothesis_values);\n",
    "plt.axvline(x=lower, color='r', linewidth=2);\n",
    "plt.axvline(x=upper, color='r', linewidth=2);\n",
    "plt.axvline(x=sample1.height.mean(), color='b', linewidth=2);\n",
    "plt.axvline(x=null_hypothesis, color='r', linewidth=2);\n",
    "plt.axvline(x=np.mean(means), color='g', linewidth=2);\n",
    "plt.axvline(x=np.mean(null_hypothesis_values), color='y', linewidth=2);\n",
    "plt.show();\n",
    "\n",
    "print('%.1f%% confident that the value is between %.1f and %.1f' % ((1-alpha)*100, lower, upper))\n",
    "\n",
    "# p-value\n",
    "# for a two sided hypothesis, we want to look at anything \n",
    "# more extreme from the null in both directions\n",
    "obs_mean = sample1.height.mean()\n",
    "\n",
    "# probability of a statistic higher than observed\n",
    "prob_more_extreme_high = (null_hypothesis_values > obs_mean).mean()\n",
    "    \n",
    "# probability a statistic is more extreme lower\n",
    "prob_more_extreme_low = (null_hypothesis - (obs_mean - null_hypothesis) < null_hypothesis_values).mean()\n",
    "\n",
    "pval = prob_more_extreme_low + prob_more_extreme_high\n",
    "pv(\"pval\")\n",
    "\n",
    "# this is another way to compute the standard deviation of the sampling distribution theoretically  \n",
    "std_sampling_dist = full_data.height.std()/np.sqrt(nb_sample)  \n",
    "null_sims = np.random.normal(null_hypothesis, std_sampling_dist, sampling_times)  \n",
    "low_ext = (null_hypothesis - (sample1.height.mean() - null_hypothesis))  \n",
    "high_ext = sample1.height.mean()  \n",
    "print(\"p-value=\",(null_sims > high_ext).mean() + (null_sims < low_ext).mean())\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`4.` Now imagine if you received the same sample mean as you calculated from the sample in question 1 above, but that you actually retrieved it from a sample of 300.  What would the new standard deviation be for your sampling distribution for the mean of 300 values?  Additionally, what would your new p-value be for choosing between the null and alternative hypotheses you set up? Simulate the sampling distribution for the mean of five values to see the shape and plot a histogram.  Use your solutions here to answer the second to last quiz question below.\n",
    "\n",
    "**Hint:**  If you get stuck, notice you can use the solution from quiz regarding finding the p-value earlier to assist with obtaining this answer with just a few small changes.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "np.mean(full_data.height) : 67.59748697307937\n",
      "np.std(full_data.height) : 3.1189087116821472\n",
      "np.mean(sample1.height) : 67.687610765017\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEMZJREFUeJzt3X+s3XV9x/Hna0Xcpi6UcGHY1l00dRn+scIaxDgzJrNAWYL+QaSL0jGXugUSUZOluGQ4DRnLFBcTJaJ0wiYwNnU00Aw75mJMBlJcrdTCuGJHSyvUYcCNxAx974/zrRzb++PcX+fc6+f5SE7OOe/z+Zzv+9N72tf9fr/nnKaqkCS15+dG3YAkaTQMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjThh1A9M55ZRTanx8fNRtSNKy8tBDD32vqsZmGrekA2B8fJxdu3aNug1JWlaS/Ncg4zwEJEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjVrSnwSWZjK+9Z6RbXv/9RePbNvSQnAPQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSo2YMgCRrknw5yb4ke5O8p6t/MMmTSXZ3l419c65JMpHk0SQX9NUv7GoTSbYuzpIkSYMY5P8EfgF4f1V9PckrgIeS7Owe+1hVfaR/cJIzgcuA1wGvBP4lyWu7hz8BvAU4CDyYZHtVfWshFiJJmp0ZA6CqDgOHu9s/SLIPWDXNlEuAO6rqh8B3kkwA53SPTVTV4wBJ7ujGGgCSNAKzOgeQZBw4C3igK12VZE+SbUlWdrVVwIG+aQe72lR1SdIIDBwASV4OfB64uqqeA24EXgOso7eH8NGjQyeZXtPUj93OliS7kuw6cuTIoO1JkmZpoABI8hJ6//h/rqq+AFBVT1XVj6rqx8CnefEwz0FgTd/01cChaeo/papuqqr1VbV+bGxstuuRJA1oxnMASQLcDOyrqhv66qd35wcA3gY83N3eDtyW5AZ6J4HXAl+jtwewNskZwJP0ThT/3kItRBq28a33jGS7+6+/eCTb1c+eQd4F9EbgncA3k+zuah8ANiVZR+8wzn7g3QBVtTfJnfRO7r4AXFlVPwJIchVwL7AC2FZVexdwLZKkWRjkXUBfZfLj9zummXMdcN0k9R3TzZMkDY+fBJakRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjZgyAJGuSfDnJviR7k7ynq5+cZGeSx7rrlV09ST6eZCLJniRn9z3X5m78Y0k2L96yJEkzGWQP4AXg/VX1a8C5wJVJzgS2AvdV1Vrgvu4+wEXA2u6yBbgReoEBXAu8HjgHuPZoaEiShm/GAKiqw1X19e72D4B9wCrgEuCWbtgtwFu725cAt1bP/cBJSU4HLgB2VtUzVfV9YCdw4YKuRpI0sFmdA0gyDpwFPACcVlWHoRcSwKndsFXAgb5pB7vaVHVJ0ggMHABJXg58Hri6qp6bbugktZqmfux2tiTZlWTXkSNHBm1PkjRLAwVAkpfQ+8f/c1X1ha78VHdoh+766a5+EFjTN301cGia+k+pqpuqan1VrR8bG5vNWiRJszDIu4AC3Azsq6ob+h7aDhx9J89m4K6++uXdu4HOBZ7tDhHdC2xIsrI7+buhq0mSRuCEAca8EXgn8M0ku7vaB4DrgTuTvAt4Ari0e2wHsBGYAJ4HrgCoqmeSfBh4sBv3oap6ZkFWIUmatRkDoKq+yuTH7wHOn2R8AVdO8VzbgG2zaVCStDj8JLAkNcoAkKRGGQCS1CgDQJIaNci7gKQZjW+9Z9QtSJol9wAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJatSMAZBkW5KnkzzcV/tgkieT7O4uG/seuybJRJJHk1zQV7+wq00k2brwS5EkzcYgewCfBS6cpP6xqlrXXXYAJDkTuAx4XTfnk0lWJFkBfAK4CDgT2NSNlSSNyAkzDaiqryQZH/D5LgHuqKofAt9JMgGc0z02UVWPAyS5oxv7rVl3LElaEPM5B3BVkj3dIaKVXW0VcKBvzMGuNlX9OEm2JNmVZNeRI0fm0Z4kaTpzDYAbgdcA64DDwEe7eiYZW9PUjy9W3VRV66tq/djY2BzbkyTNZMZDQJOpqqeO3k7yaeDu7u5BYE3f0NXAoe72VHVJ0gjMaQ8gyel9d98GHH2H0HbgsiQvTXIGsBb4GvAgsDbJGUlOpHeiePvc25YkzdeMewBJbgfOA05JchC4FjgvyTp6h3H2A+8GqKq9Se6kd3L3BeDKqvpR9zxXAfcCK4BtVbV3wVcjSRrYIO8C2jRJ+eZpxl8HXDdJfQewY1bdSZIWjZ8ElqRGGQCS1CgDQJIaZQBIUqPm9DkASaMzvvWekW17//UXj2zbWnjuAUhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUqBkDIMm2JE8nebivdnKSnUke665XdvUk+XiSiSR7kpzdN2dzN/6xJJsXZzmSpEENsgfwWeDCY2pbgfuqai1wX3cf4CJgbXfZAtwIvcAArgVeD5wDXHs0NCRJozFjAFTVV4BnjilfAtzS3b4FeGtf/dbquR84KcnpwAXAzqp6pqq+D+zk+FCRJA3RXM8BnFZVhwG661O7+irgQN+4g11tqvpxkmxJsivJriNHjsyxPUnSTBb6JHAmqdU09eOLVTdV1fqqWj82NragzUmSXjTXAHiqO7RDd/10Vz8IrOkbtxo4NE1dkjQicw2A7cDRd/JsBu7qq1/evRvoXODZ7hDRvcCGJCu7k78bupokaUROmGlAktuB84BTkhyk926e64E7k7wLeAK4tBu+A9gITADPA1cAVNUzST4MPNiN+1BVHXtiWZI0RDMGQFVtmuKh8ycZW8CVUzzPNmDbrLqTJC0aPwksSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEbN+GVwWl7Gt94z6hYkLRPuAUhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRs0rAJLsT/LNJLuT7OpqJyfZmeSx7nplV0+SjyeZSLInydkLsQBJ0twsxB7Ab1fVuqpa393fCtxXVWuB+7r7ABcBa7vLFuDGBdi2JGmOFuMQ0CXALd3tW4C39tVvrZ77gZOSnL4I25ckDWC+AVDAl5I8lGRLVzutqg4DdNendvVVwIG+uQe7miRpBOb7X0K+saoOJTkV2JnkkWnGZpJaHTeoFyRbAF71qlfNsz1J0lTmtQdQVYe666eBLwLnAE8dPbTTXT/dDT8IrOmbvho4NMlz3lRV66tq/djY2HzakyRNY84BkORlSV5x9DawAXgY2A5s7oZtBu7qbm8HLu/eDXQu8OzRQ0WSpOGbzyGg04AvJjn6PLdV1T8neRC4M8m7gCeAS7vxO4CNwATwPHDFPLYtSZqnOQdAVT0O/Pok9f8Gzp+kXsCVc92eJGlh+UlgSWqUASBJjZrv20AlNWR86z0j2e7+6y8eyXZ/1rkHIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqP8/wAWwai+M12SZsM9AElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGuXnACQteaP8bM3+6y8e2bYXm3sAktSooQdAkguTPJpkIsnWYW9fktQz1ABIsgL4BHARcCawKcmZw+xBktQz7HMA5wATVfU4QJI7gEuAby3GxvxOHkma2rADYBVwoO/+QeD1Q+5BkgY2ql8kh3HyedgBkElq9VMDki3Alu7u/yR5dB7bOwX43jzmj9py7x+W/xqWe//gGpaCWfefv5zX9n5lkEHDDoCDwJq++6uBQ/0Dquom4KaF2FiSXVW1fiGeaxSWe/+w/New3PsH17AULNX+h/0uoAeBtUnOSHIicBmwfcg9SJIY8h5AVb2Q5CrgXmAFsK2q9g6zB0lSz9A/CVxVO4AdQ9rcghxKGqHl3j8s/zUs9/7BNSwFS7L/VNXMoyRJP3P8KghJatSyDIAkJyX5xySPJNmX5A1J/j7J7u6yP8nuKeYuia+imOsakqxJ8uVuzt4k71lO/ffNX5HkP5LcPcy+j+lhPq+j4+YOu/+p+pjFGt7bvYYeTnJ7kp9fIv2vS3J/1/+uJOdMMXdzkse6y+Zh997Xx5zW0I359+5nsCfJ24fefFUtuwtwC/CH3e0TgZOOefyjwJ9NMm8F8G3g1d28bwBnLrM1nA6c3d1+BfCfo1jDXPvve/x9wG3A3cvtdTTI3KW+BnofyvwO8Avd/TuB318K/QNfAi7qahuBf5tk3snA4931yu72yqXyMxhwDa8F1na3XwkcHvbraOh/WAvwh/1L3Qs3Uzweep82XjvJY28A7u27fw1wzXJawyRj7wLespz6p/f5j/uAN48qAOb5Opp27jJZw9FP5Z9M780gdwMblkL/9N4l+Pbu9ibgtknmbgI+1Xf/U8CmpfIzGGQNkzzXNwb5O7+Ql+V4COjVwBHgb7pDCJ9J8rK+x98EPFVVj00yd7Kvoli1eK1OaT5r+Ikk48BZwAOL1egU5tv/XwN/Avx4kfucznzWMNPcYZnzGqrqSeAjwBP0fvN8tqq+NIym+0zV/9XAXyU50PV4zSRzl/rf5UHW8BPdIaIT6R2hGJrlGAAnAGcDN1bVWcD/Av3H8jcBt08xd8avohiS+awBgCQvBz4PXF1Vzy1Wo1OYc/9Jfhd4uqoeWvQupzefn8FMc4dlPj+HlfS+iPEMeocfXpbkHYvb7nGm6v+PgfdW1RrgvcDNk8xd6n+XB1kDAElOB/4WuKKqhvtL0bB3mRZgl+uXgf19998E3NPdPgF4Clg9xdylcghozmvoxryE3i7m+5bhz+Av6P22th/4LvA88HfLbA1Tzl1Ga7gUuLnv/uXAJ5dC/8CzvPgW9QDPTTJ3qRwCmvMausd+Cfg6cOmwe69ahoeAquq7wIEkv9qVzufFr5P+HeCRqjo4xfQl8VUU81lDktD7bWJfVd2w6M1OYj79V9U1VbW6qsbp/fn/a1UN+zfP+a5hurlDM8+/C08A5yb5xe41dT6wb1EbPsY0/R8CfqurvRmY7DDcvcCGJCu7vZkNXW2o5rOG7t+gLwK3VtU/DKHd440idRYgddcBu4A9wD/Rnf0HPgv80TFjXwns6Lu/kd47Z74N/OlyWwPwm/R2dfcAu7vLxuXS/zH18xjtu4Dm8zqadO4yW8OfA48AD9M7BPHSpdB/9xp/iN5J0QeA3+jGrgc+0zf3D4CJ7nLFUnodDbIG4B3A//X9Pd4NrBtm734SWJIatewOAUmSFoYBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSo/4f+Zf3wtzjotkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "np.mean(means) : 67.5988283911597\n",
      "np.std(means) : 0.17076382206590648\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEjZJREFUeJzt3X+s3fV93/HnqyakXZIWM24yYju7tHWnkj/qsCtClVXLQgsGJpFKQ8VSG48xuZtAStJKk+mk0R+LRrUm2SIlqLR4cdomlDXNsMAq8WiqqtJCMKlDMA7jhnj4xhTckZJ2kbJB3/vjfJwczP1x7q9zfP15PqSj8/2+v5/P/b6/l2Ne9/s933NvqgpJUn++Z9INSJImwwCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdeq8STewmIsuuqimp6cn3YYkbSiPPvroX1bV1FLjzuoAmJ6e5vDhw5NuQ5I2lCT/a5RxXgKSpE4tGQBJvjfJF5J8KcnRJL/S6pckeTjJU0l+P8n5rf7atj7btk8Pfa3bWv3JJFev10FJkpY2yhnAt4F3VdWPATuAnUmuAH4d+HBVbQe+Adzcxt8MfKOqfhj4cBtHkkuBG4G3AjuBjyXZtJYHI0ka3ZIBUAN/01Zf0x4FvAv4g1bfD7y7LV/f1mnbr0ySVr+nqr5dVV8DZoHL1+QoJEnLNtJ7AEk2JTkCPA8cAr4K/FVVvdSGzAFb2vIW4ARA2/4i8HeH6/PMkSSN2UgBUFUvV9UOYCuDn9p/dL5h7TkLbFuo/gpJ9iQ5nOTwqVOnRmlPkrQCy7oLqKr+CvgT4ArggiSnbyPdCpxsy3PANoC2/QeAF4br88wZ3sddVTVTVTNTU0vexipJWqFR7gKaSnJBW/4+4CeBY8DngH/Whu0G7mvLB9o6bfsf1+DvTh4Abmx3CV0CbAe+sFYHIklanlE+CHYxsL/dsfM9wL1VdX+SJ4B7kvx74M+Bu9v4u4HfSTLL4Cf/GwGq6miSe4EngJeAW6rq5bU9HEnSqHI2/1H4mZmZ8pPAWsz03gcmtu/jd1w3sX1Li0nyaFXNLDXOTwJLUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdeq8STcgbVTTex+YyH6P33HdRParc49nAJLUKQNAkjplAEhSpwwASerUkgGQZFuSzyU5luRokve2+i8n+XqSI+1x7dCc25LMJnkyydVD9Z2tNptk7/ockiRpFKPcBfQS8ItV9cUkbwAeTXKobftwVf3G8OAklwI3Am8F3gz89yQ/0jZ/FPgpYA54JMmBqnpiLQ5EkrQ8SwZAVT0LPNuW/zrJMWDLIlOuB+6pqm8DX0syC1zets1W1dMASe5pYw0ASZqAZb0HkGQaeBvwcCvdmuSxJPuSbG61LcCJoWlzrbZQXZI0ASMHQJLXA58G3ldV3wTuBH4I2MHgDOGDp4fOM70WqZ+5nz1JDic5fOrUqVHbkyQt00gBkOQ1DP7n/3tV9YcAVfVcVb1cVX8L/BbfvcwzB2wbmr4VOLlI/RWq6q6qmqmqmampqeUejyRpRKPcBRTgbuBYVX1oqH7x0LCfBh5vyweAG5O8NsklwHbgC8AjwPYklyQ5n8EbxQfW5jAkScs1yl1A7wB+DvhykiOt9kvAriQ7GFzGOQ78PEBVHU1yL4M3d18CbqmqlwGS3Ao8CGwC9lXV0TU8FknSMoxyF9CfMf/1+4OLzPkA8IF56gcXmydJGh8/CSxJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTSwZAkm1JPpfkWJKjSd7b6hcmOZTkqfa8udWT5CNJZpM8luSyoa+1u41/Ksnu9TssSdJSRjkDeAn4xar6UeAK4JYklwJ7gYeqajvwUFsHuAbY3h57gDthEBjA7cDbgcuB20+HhiRp/JYMgKp6tqq+2Jb/GjgGbAGuB/a3YfuBd7fl64FP1MDngQuSXAxcDRyqqheq6hvAIWDnmh6NJGlky3oPIMk08DbgYeBNVfUsDEICeGMbtgU4MTRtrtUWqkuSJmDkAEjyeuDTwPuq6puLDZ2nVovUz9zPniSHkxw+derUqO1JkpZppABI8hoG//P/var6w1Z+rl3aoT0/3+pzwLah6VuBk4vUX6Gq7qqqmaqamZqaWs6xSJKWYZS7gALcDRyrqg8NbToAnL6TZzdw31D9Pe1uoCuAF9slogeBq5Jsbm/+XtVqkqQJOG+EMe8Afg74cpIjrfZLwB3AvUluBp4BbmjbDgLXArPAt4CbAKrqhSS/BjzSxv1qVb2wJkchSVq2JQOgqv6M+a/fA1w5z/gCblnga+0D9i2nQUnS+vCTwJLUKQNAkjplAEhSpwwASeqUASBJnRrlNlBpSdN7H5h0C5KWyTMASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4tGQBJ9iV5PsnjQ7VfTvL1JEfa49qhbbclmU3yZJKrh+o7W202yd61PxRJ0nKMcgbwcWDnPPUPV9WO9jgIkORS4EbgrW3Ox5JsSrIJ+ChwDXApsKuNlSRNyHlLDaiqP00yPeLXux64p6q+DXwtySxweds2W1VPAyS5p419YtkdS5LWxGreA7g1yWPtEtHmVtsCnBgaM9dqC9VfJcmeJIeTHD516tQq2pMkLWalAXAn8EPADuBZ4IOtnnnG1iL1Vxer7qqqmaqamZqaWmF7kqSlLHkJaD5V9dzp5SS/BdzfVueAbUNDtwIn2/JCdUnSBKzoDCDJxUOrPw2cvkPoAHBjktcmuQTYDnwBeATYnuSSJOczeKP4wMrbliSt1pJnAEk+BbwTuCjJHHA78M4kOxhcxjkO/DxAVR1Nci+DN3dfAm6pqpfb17kVeBDYBOyrqqNrfjSSpJGNchfQrnnKdy8y/gPAB+apHwQOLqs7SdK68ZPAktQpA0CSOmUASFKnDABJ6pQBIEmdWtEHwSRNzvTeBya27+N3XDexfWvteQYgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKnDABJ6pQBIEmdMgAkqVMGgCR1ygCQpE4ZAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROLRkASfYleT7J40O1C5McSvJUe97c6knykSSzSR5LctnQnN1t/FNJdq/P4UiSRjXKGcDHgZ1n1PYCD1XVduChtg5wDbC9PfYAd8IgMIDbgbcDlwO3nw4NSdJkLBkAVfWnwAtnlK8H9rfl/cC7h+qfqIHPAxckuRi4GjhUVS9U1TeAQ7w6VCRJY7TS9wDeVFXPArTnN7b6FuDE0Li5Vluo/ipJ9iQ5nOTwqVOnVtieJGkpa/0mcOap1SL1Vxer7qqqmaqamZqaWtPmJEnftdIAeK5d2qE9P9/qc8C2oXFbgZOL1CVJE7LSADgAnL6TZzdw31D9Pe1uoCuAF9slogeBq5Jsbm/+XtVqkqQJOW+pAUk+BbwTuCjJHIO7ee4A7k1yM/AMcEMbfhC4FpgFvgXcBFBVLyT5NeCRNu5Xq+rMN5YlSWO0ZABU1a4FNl05z9gCblng6+wD9i2rO0nSuvGTwJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOrXkn4TUxjK994FJtyBpg/AMQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkTq0qAJIcT/LlJEeSHG61C5McSvJUe97c6knykSSzSR5LctlaHIAkaWXW4gzgn1TVjqqaaet7gYeqajvwUFsHuAbY3h57gDvXYN+SpBVaj0tA1wP72/J+4N1D9U/UwOeBC5JcvA77lySNYLUBUMBnkzyaZE+rvamqngVoz29s9S3AiaG5c60mSZqA1f4uoHdU1ckkbwQOJfnKImMzT61eNWgQJHsA3vKWt6yyPUnSQlZ1BlBVJ9vz88BngMuB505f2mnPz7fhc8C2oelbgZPzfM27qmqmqmampqZW054kaRErDoAkr0vyhtPLwFXA48ABYHcbthu4ry0fAN7T7ga6Anjx9KUiSdL4reYS0JuAzyQ5/XU+WVV/lOQR4N4kNwPPADe08QeBa4FZ4FvATavYtyRplVYcAFX1NPBj89T/N3DlPPUCblnp/iRJa8tPAktSpwwASeqUASBJnTIAJKlT/lF4SSOb3vvARPZ7/I7rJrLfc51nAJLUKQNAkjplAEhSpwwASeqUASBJnTIAJKlTBoAkdcoAkKROGQCS1CkDQJI6ZQBIUqcMAEnqlAEgSZ0yACSpUwaAJHXKAJCkTvkHYdbBpP5ohiQth2cAktQpA0CSOmUASFKnDABJ6tTY3wROshP4z8Am4Ler6o5x9yBpY5nkjRXH77huYvteb2M9A0iyCfgocA1wKbAryaXj7EGSNDDuS0CXA7NV9XRV/V/gHuD6MfcgSWL8l4C2ACeG1ueAt6/XzrwfX5IWNu4AyDy1esWAZA+wp63+TZIn17Gfi4C/XMevv142at+wcXu37/E6a/rOry9r+NnS998fZdC4A2AO2Da0vhU4OTygqu4C7hpHM0kOV9XMOPa1ljZq37Bxe7fv8bLv8Rj3ewCPANuTXJLkfOBG4MCYe5AkMeYzgKp6KcmtwIMMbgPdV1VHx9mDJGlg7J8DqKqDwMFx73cBY7nUtA42at+wcXu37/Gy7zFIVS09SpJ0zvFXQUhSp87ZAEhyQZI/SPKVJMeS/HiS309ypD2OJzmywNydSZ5MMptk70boO8m2JJ9rc44mee9G6Hto/qYkf57k/o3S93xzN0jf72+vkceTfCrJ90647x1JPt/6Ppzk8gXm7k7yVHvsHlfPq+m7jfkf7fv9WJKfGWffS6qqc/IB7Af+ZVs+H7jgjO0fBP7dPPM2AV8FfrDN+xJw6Qbo+2Lgsrb8BuB/boS+h7b/AvBJ4P6N8DoZZe7Z2DeDD2N+Dfi+tn4v8M8n2TfwWeCaVrsW+JN55l0IPN2eN7flzRug7x8BtrflNwPPjvN1suRxTbqBdfqP9f3tRZ4FtofBJ5K3z7Ptx4EHh9ZvA2472/ueZ+x9wE9thL4ZfB7kIeBd4wyAVb5OFp17Fvd9+tP4FzK4CeR+4KpJ9s3grsCfacu7gE/OM3cX8JtD678J7Drb+57na31plH+/43qcq5eAfhA4BfyXdlnht5O8bmj7TwDPVdVT88yd79dVbFm/Vl9hNX1/R5Jp4G3Aw+vV6BlW2/d/Av4N8Lfr3OeZVtP3UnPX04r7rqqvA78BPMPgp9EXq+qz42iahft+H/Afk5xovd02z9yz8d/lKH1/R7tEdD6DKwxnhXM1AM4DLgPurKq3Af8HGL6Wvwv41AJzl/x1FetoNX0DkOT1wKeB91XVN9er0TOsuO8k/xR4vqoeXfcuX2013++l5q6n1Xy/NzP4BYyXMLgk8bokP7u+7X7HQn3/a+D9VbUNeD9w9zxzz8Z/l6P0DUCSi4HfAW6qqnH/oLOwSZ+CrNMp298Djg+t/wTwQFs+D3gO2LrA3EleAlpx323Maxiclv7CBvp+/wcGP80dB/4C+Bbwuxug7wXnnuV93wDcPbT+HuBjk+wbeJHv3pIe4JvzzJ3kJaAV9922fT/wReCGcfS7nMc5eQZQVX8BnEjyD1rpSuCJtvyTwFeqam6B6RP7dRWr6TtJGPwEcqyqPrTuzQ5ZTd9VdVtVba2qaQbf6z+uqrH8RLrKvhebu65W+fp+Brgiyd9pr5krgWPr2nCzSN8ngX/cau8C5rvk9iBwVZLN7SzmqlZbd6vpu/0/5DPAJ6rqv46h3eWZdAKtY2rvAA4DjwH/jXbHAPBx4F+dMfbNwMGh9WsZ3EXzVeDfboS+gX/E4JT4MeBIe1x7tvd9Rv2djP8uoNW8TuaduwH6/hXgK8DjDC5LvHaSfbfX7qMM3iB9GPiHbewMg78aeHruvwBm2+OmSb9ORukb+Fng/w39mzwC7Bhn74s9/CSwJHXqnLwEJElamgEgSZ0yACSpUwaAJHXKAJCkThkAktQpA0CSOmUASFKn/j9Zc4sw4661+wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "np.random.seed(42)\n",
    "pv(\"np.mean(full_data.height)\")\n",
    "pv(\"np.std(full_data.height)\")\n",
    "nb_sample = 300\n",
    "sample1 = full_data.sample(nb_sample)\n",
    "pv(\"np.mean(sample1.height)\")\n",
    "null_hypothesis = 67.60\n",
    "sampling_times = 10_000\n",
    "\n",
    "means = []\n",
    "for _ in range(sampling_times):\n",
    "    boot_sample = full_data.sample(nb_sample) # Why not with replace=True ?\n",
    "    means.append(boot_sample.height.mean())\n",
    "plt.hist(means);\n",
    "plt.show();\n",
    "pv(\"np.mean(means)\")\n",
    "pv(\"np.std(means)\")\n",
    "\n",
    "\n",
    "null_hypothesis_values = np.random.normal(loc=null_hypothesis, scale=np.std(means), size=sampling_times)\n",
    "\n",
    "plt.hist(null_hypothesis_values);\n",
    "plt.show();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "alpha/2  * 100 : 2.5\n",
      "((1-alpha) + alpha/2)  * 100 : 97.5\n",
      "lower, upper : (67.25222466486355, 67.95660471529004)\n",
      "np.mean(means)-np.mean(null_hypothesis_values) : -0.005220997801814065\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(array([  15.,   98.,  575., 1562., 2654., 2687., 1681.,  586.,  125.,\n",
       "          17.]),\n",
       " array([66.92853397, 67.06220868, 67.19588338, 67.32955809, 67.46323279,\n",
       "        67.5969075 , 67.73058221, 67.86425691, 67.99793162, 68.13160632,\n",
       "        68.26528103]),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af699da58>"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af69653c8>"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af68fc978>"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af68fcc50>"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.lines.Line2D at 0x7f5af690f1d0>"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD8CAYAAAB+UHOxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAEPFJREFUeJzt3X/sXXV9x/Hna0WcUwwlFIZts6KpZvjHCn5TMc7IZBYoS9A/iDRROqapW2ARNVmKS4Y/QkYyfywmSqzSCZvC2NTRQDPsGIsxGci3rlZqYXzFhtZWqMOAGYkZ+t4f9xQu7ff3j3u/Xz/PR3Jyz3nfz+ee97255fU95557SVUhSWrPbwy7AUnScBgAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEadNOwGJnP66afXmjVrht3Gr6/du3u3b3jDcPtYgnYf3s1rT+mtn3KKr9+s+R5cELt37/5pVa2YalwW809BjIyM1Ojo6LDb+PWV9G4X8XtgscrHwn1v7a1fcIGv36z5HlwQSXZX1chU4zwFJEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjVrUPwUhLQVrtt490P0duPHSge5Pv748ApCkRnkEoCVt0H99P+9lw9mtNJ88ApCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVFTBkCS1UnuS7I/yb4kH+jqH03y4yR7umVj35zrkowleSTJRX31i7vaWJKtC/OUJEnTMZ0vgj0HfLiqvpvkFGB3kl3dfZ+pqk/2D05yDnAF8HrgVcC/JXltd/fngLcDh4AHk+yoqh/MxxORJM3MlAFQVUeAI936z5PsB1ZOMuUy4Paq+gXwoyRjwPruvrGqegwgye3dWANAkoZgRp8BJFkDnAs80JWuSbI3yfYky7vaSuBg37RDXW2iuiRpCKYdAEleAXwNuLaqngFuAl4DrKN3hPCpY0PHmV6T1I/fz5Yko0lGjx49Ot32JEkzNK0ASPISev/x/0pVfR2gqp6oql9W1a+AL/LCaZ5DwOq+6auAw5PUX6SqtlXVSFWNrFixYqbPR5I0TdO5CijAzcD+qvp0X/2svmHvBB7q1ncAVyR5aZKzgbXAd4AHgbVJzk5yMr0PinfMz9OQJM3UdK4CejPwHuD7SfZ0tY8Am5Kso3ca5wDwfoCq2pfkDnof7j4HXF1VvwRIcg1wD7AM2F5V++bxuUiSZmA6VwF9m/HP3++cZM4NwA3j1HdONk+SNDh+E1iSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1asoASLI6yX1J9ifZl+QDXf20JLuSPNrdLu/qSfLZJGNJ9iY5r++xNnfjH02yeeGeliRpKtM5AngO+HBV/S5wPnB1knOArcC9VbUWuLfbBrgEWNstW4CboBcYwPXAG4H1wPXHQkOSNHhTBkBVHamq73brPwf2AyuBy4BbumG3AO/o1i8Dbq2e+4FTk5wFXATsqqqnqupnwC7g4nl9NpKkaZvRZwBJ1gDnAg8AZ1bVEeiFBHBGN2wlcLBv2qGuNlFdkjQE0w6AJK8AvgZcW1XPTDZ0nFpNUj9+P1uSjCYZPXr06HTbkyTN0LQCIMlL6P3H/ytV9fWu/ER3aofu9smufghY3Td9FXB4kvqLVNW2qhqpqpEVK1bM5LlIkmZgOlcBBbgZ2F9Vn+67awdw7EqezcCdffUru6uBzgee7k4R3QNsSLK8+/B3Q1eTJA3BSdMY82bgPcD3k+zpah8BbgTuSPJe4HHg8u6+ncBGYAx4FrgKoKqeSvIJ4MFu3Mer6ql5eRaSpBmbMgCq6tuMf/4e4MJxxhdw9QSPtR3YPpMGJUkLw28CS1KjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY2azv8UXtIismbr3UPb94EbLx3avjX/PAKQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcovgmleDPPLSZJmxyMASWrUlAGQZHuSJ5M81Ff7aJIfJ9nTLRv77rsuyViSR5Jc1Fe/uKuNJdk6/09FkjQT0zkC+DJw8Tj1z1TVum7ZCZDkHOAK4PXdnM8nWZZkGfA54BLgHGBTN1aSNCRTfgZQVd9Ksmaaj3cZcHtV/QL4UZIxYH1331hVPQaQ5PZu7A9m3LEkaV7M5TOAa5Ls7U4RLe9qK4GDfWMOdbWJ6idIsiXJaJLRo0ePzqE9SdJkZhsANwGvAdYBR4BPdfWMM7YmqZ9YrNpWVSNVNbJixYpZtidJmsqsLgOtqieOrSf5InBXt3kIWN03dBVwuFufqC5JGoJZHQEkOatv853AsSuEdgBXJHlpkrOBtcB3gAeBtUnOTnIyvQ+Kd8y+bUnSXE15BJDkNuAC4PQkh4DrgQuSrKN3GucA8H6AqtqX5A56H+4+B1xdVb/sHuca4B5gGbC9qvbN+7ORJE3bdK4C2jRO+eZJxt8A3DBOfSewc0bdSZIWjN8ElqRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElq1JQBkGR7kieTPNRXOy3JriSPdrfLu3qSfDbJWJK9Sc7rm7O5G/9oks0L83QkSdM1nSOALwMXH1fbCtxbVWuBe7ttgEuAtd2yBbgJeoEBXA+8EVgPXH8sNCRJwzFlAFTVt4CnjitfBtzSrd8CvKOvfmv13A+cmuQs4CJgV1U9VVU/A3ZxYqhIkgZotp8BnFlVRwC62zO6+krgYN+4Q11torokaUjm+0PgjFOrSeonPkCyJcloktGjR4/Oa3OSpBfMNgCe6E7t0N0+2dUPAav7xq0CDk9SP0FVbauqkaoaWbFixSzbkyRNZbYBsAM4diXPZuDOvvqV3dVA5wNPd6eI7gE2JFneffi7oatJkobkpKkGJLkNuAA4Pckhelfz3AjckeS9wOPA5d3wncBGYAx4FrgKoKqeSvIJ4MFu3Mer6vgPliVJAzRlAFTVpgnuunCcsQVcPcHjbAe2z6g7SdKC8ZvAktQoA0CSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRhkAktQoA0CSGmUASFKjThp2A5pfa7bePe2xB2YxR22b7/fKgWk+7oEbL53X/arHIwBJapQBIEmNmlMAJDmQ5PtJ9iQZ7WqnJdmV5NHudnlXT5LPJhlLsjfJefPxBCRJszMfRwB/UFXrqmqk294K3FtVa4F7u22AS4C13bIFuGke9i1JmqWFOAV0GXBLt34L8I6++q3Vcz9wapKzFmD/kqRpmGsAFPDNJLuTbOlqZ1bVEYDu9oyuvhI42Df3UFd7kSRbkowmGT169Ogc25MkTWSul4G+uaoOJzkD2JXk4UnGZpxanVCo2gZsAxgZGTnhfknS/JjTEUBVHe5unwS+AawHnjh2aqe7fbIbfghY3Td9FXB4LvuXJM3erAMgycuTnHJsHdgAPATsADZ3wzYDd3brO4Aru6uBzgeePnaqSJI0eHM5BXQm8I0kxx7nq1X1r0keBO5I8l7gceDybvxOYCMwBjwLXDWHfUuS5mjWAVBVjwG/N079f4ALx6kXcPVs9ydJml9+E1iSGmUASFKjDABJapQBIEmNMgAkqVEGgCQ1ygCQpEYZAJLUKANAkhplAEhSowwASWqUASBJjTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqMMAElqlAEgSY0yACSpUQaAJDXKAJCkRp007AYkaSprtt49tH0fuPHSoe17oXkEIEmN8ghgAQzzrxVJmi6PACSpUQaAJDXKAJCkRg08AJJcnOSRJGNJtg56/5KknoEGQJJlwOeAS4BzgE1JzhlkD5KknkFfBbQeGKuqxwCS3A5cBvxgwH1I0rQM66q+QXz/YNABsBI42Ld9CHjjQu3MyzElaWKpqsHtLLkcuKiq3tdtvwdYX1V/3jdmC7Cl23wd8MjAGjzR6cBPh7j/2ViKPYN9D9JS7BmWZt/D6vl3qmrFVIMGfQRwCFjdt70KONw/oKq2AdsG2dREkoxW1ciw+5iJpdgz2PcgLcWeYWn2vdh7HvRVQA8Ca5OcneRk4Apgx4B7kCQx4COAqnouyTXAPcAyYHtV7RtkD5KknoH/FlBV7QR2Dnq/s7QoTkXN0FLsGex7kJZiz7A0+17UPQ/0Q2BJ0uLhT0FIUqOaDYAkpyb55yQPJ9mf5E1J/jHJnm45kGTPBHOH8nMWs+05yeok93Vz9iX5wKB6nkvfffOXJfmvJHcthZ7Hm7tE+v5g9/54KMltSX5ziD2vS3J/1/NokvUTzN2c5NFu2TyIfufadzfmP7vXem+Sdw2y7xepqiYX4Bbgfd36ycCpx93/KeCvxpm3DPgh8Opu3veAcxZ5z2cB53XrpwD/Paie59J33/0fAr4K3LUUep5q7mLsm96XNH8EvKzbvgP442H1DHwTuKSrbQT+Y5x5pwGPdbfLu/Xlw3ytp9n3a4G13fqrgCODfI+8qJdh7HTYC/DK7s2eCe4PvW8srx3nvjcB9/RtXwdct5h7HmfsncDbF/tr3d2/CrgXeNugAmCO749J5y7ivo99S/80eheH3AVsGFbP9K4UfFe3vgn46jhzNwFf6Nv+ArBpmK/1dPoe57G+N51/twuxtHoK6NXAUeDvulMLX0ry8r773wI8UVWPjjN3vJ+zWLlwrT5vLj0/L8ka4FzggYVq9Dhz7ftvgb8AfrXAffabS89TzV1Is+67qn4MfBJ4nN5fpE9X1TeH2PO1wN8kOdj1dd04c4f1bxHm1vfzulNEJ9M7qzBwrQbAScB5wE1VdS7wv0D/ufxNwG0TzM04tUFcSjWXngFI8grga8C1VfXMQjV6nFn3neSPgCeraveCd/lic3mtp5q7kObyWi+n98OMZ9M7LfHyJO9e2HaBiXv+M+CDVbUa+CBw83htj1Mb1GWNc+kbgCRnAX8PXFVVg/wD5wXDOOwY9gL8NnCgb/stwN3d+knAE8CqCeYO6xTQrHvuxryE3uHph5bQa/3X9P6qOwD8BHgW+IdF3vOEcxd535cDN/dtXwl8flg9A0/zwmXqAZ4ZZ+4wTwHNuu/uvlcC3wUuH0S/Ey1NHgFU1U+Ag0le15Uu5IWfpP5D4OGqOjTB9KH8nMVcek4Sen+J7K+qTy90r/3m0ndVXVdVq6pqDb3X+d+rasH/Kp1jz5PNXVBzfF8/Dpyf5Le698uFwP4FbZhJez4MvLWrvQ0Y73TbPcCGJMu7I5gNXW3BzaXv7r8b3wBurap/GkC7Extm+gxzAdYBo8Be4F/orh4Avgz86XFjXwXs7NveSO9Kmh8Cf7nYewZ+n96h8V5gT7dsXOx9H1e/gMFeBTSX98e4c5dA3x8DHgYeondq4qXD6rl7z+6m9wHpA8AburEjwJf65v4JMNYtVw3qdZ5L38C7gf/r+7e4B1g3yN6PLX4TWJIa1eQpIEmSASBJzTIAJKlRBoAkNcoAkKRGGQCS1CgDQJIaZQBIUqP+H8sjFk1lCZbPAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "95.0% confident that the value is between 67.3 and 68.0\n",
      "p-value= 0.6278\n",
      "p-value= 0.6276\n"
     ]
    }
   ],
   "source": [
    "\n",
    "alpha = 0.05\n",
    "pv(\"alpha/2  * 100\")\n",
    "pv(\"((1-alpha) + alpha/2)  * 100\")\n",
    "lower = np.percentile(null_hypothesis_values, alpha/2  * 100)\n",
    "upper = np.percentile(null_hypothesis_values, ((1-alpha) + alpha/2)  * 100) # np.percentile(mb, [2.5, 97.5])\n",
    "\n",
    "\n",
    "pv(\"lower, upper\")\n",
    "pv(\"np.mean(means)-np.mean(null_hypothesis_values)\")\n",
    "\n",
    "plt.hist(null_hypothesis_values);\n",
    "plt.axvline(x=lower, color='r', linewidth=2);\n",
    "plt.axvline(x=upper, color='r', linewidth=2);\n",
    "\n",
    "plt.axvline(x=null_hypothesis, color='r', linewidth=2);\n",
    "plt.axvline(x=np.mean(means), color='g', linewidth=2);\n",
    "plt.axvline(x=np.mean(null_hypothesis_values), color='y', linewidth=2);\n",
    "plt.show();\n",
    "\n",
    "print('%.1f%% confident that the value is between %.1f and %.1f' % ((1-alpha)*100, lower, upper))\n",
    "\n",
    "# p-value\n",
    "low_ext = (null_hypothesis - (sample1.height.mean() - null_hypothesis))  \n",
    "high_ext = sample1.height.mean()\n",
    "print(\"p-value=\",(null_hypothesis_values > high_ext).mean() + (null_hypothesis_values < low_ext).mean())\n",
    "\n",
    "# this is another way to compute the standard deviation of the sampling distribution theoretically  \n",
    "std_sampling_dist = full_data.height.std()/np.sqrt(nb_sample)  \n",
    "null_sims = np.random.normal(null_hypothesis, std_sampling_dist, sampling_times)  \n",
    "low_ext = (null_hypothesis - (sample1.height.mean() - null_hypothesis))  \n",
    "high_ext = sample1.height.mean()  \n",
    "print(\"p-value=\",(null_sims > high_ext).mean() + (null_sims < low_ext).mean())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multiple Tests\n",
    "In this notebook, you will work with a similar dataset to the judicial dataset you were working with before.  However, instead of working with decisions already being provided, you are provided with a p-value associated with each individual.  \n",
    "\n",
    "Use the questions in the notebook and the dataset to answer the questions at the bottom of this page.\n",
    "\n",
    "Here is a glimpse of the data you will be working with:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "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>defendant_id</th>\n",
       "      <th>actual</th>\n",
       "      <th>pvalue</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>22574</td>\n",
       "      <td>innocent</td>\n",
       "      <td>0.294126</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>35637</td>\n",
       "      <td>innocent</td>\n",
       "      <td>0.417981</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>39919</td>\n",
       "      <td>innocent</td>\n",
       "      <td>0.177542</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>29610</td>\n",
       "      <td>guilty</td>\n",
       "      <td>0.015023</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>38273</td>\n",
       "      <td>innocent</td>\n",
       "      <td>0.075371</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   defendant_id    actual    pvalue\n",
       "0         22574  innocent  0.294126\n",
       "1         35637  innocent  0.417981\n",
       "2         39919  innocent  0.177542\n",
       "3         29610    guilty  0.015023\n",
       "4         38273  innocent  0.075371"
      ]
     },
     "execution_count": 64,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "df = pd.read_csv('judicial_dataset_pvalues.csv')\n",
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`1.` Remember back to the null and alternative hypotheses for this example. Use that information to determine the answer for **Quiz 1** and **Quiz 2** below.\n",
    "`2.` If we consider each individual as a single hypothesis test, find the conservative Bonferroni corrected alpha level we should use to maintain a 5% type I error rate.\n",
    "\n",
    "`3.` What is the proportion of type I errors made if the correction isn't used?  How about if it is used?  \n",
    "\n",
    "Use your answers to find the solution to **Quiz 3** below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Bonferroni : 6.86530275985171e-06\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "0.001510366607167376"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "alpha = 0.05\n",
    "nb_test=jud_data.shape[0]\n",
    "Bonferroni = alpha / nb_test\n",
    "pv('Bonferroni')\n",
    "\n",
    "# p-value : The probability of us observing the facts about an individual's case that are in favor of them being guilty, assuming they are truly innocent.\n",
    "# In order to find the number of type I errors made without the correction - we need to find all those that are actually innocent with p-values less than 0.05.\n",
    "df.query(\"actual == 'innocent' and pvalue < 0.05\").count()[0]/df.shape[0] # If not used\n",
    "df.query(\"actual == 'innocent' and pvalue < @Bonferroni\").count()[0]/df.shape[0] # If used"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`4.` Think about how hypothesis tests can be used, and why this example wouldn't exactly work in terms of being able to use hypothesis testing in this way.  Check your answer with **Quiz 4** below."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**This is looking at individuals, and that is more of the aim for machine learning techniques.  Hypothesis testing and confidence intervals are for population parameters. Therefore, they are not meant to tell us about individual cases, and we wouldn't obtain p-values for individuals in this way. We could get probabilities, but that isn't the same as the probabilities associated with the relationship to sampling distributions as you have seen in these lessons.**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.1"
  },
  "widgets": {
   "state": {},
   "version": "1.1.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}