Other Projects and Exercises/Statistics - Human Temperature/sliderule_dsi_inferential_statistics_exercise_1.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# What is the True Normal Human Body Temperature? \n",
    "\n",
    "#### Background\n",
    "\n",
    "The mean normal body temperature was held to be 37$^{\\circ}$C or 98.6$^{\\circ}$F for more than 120 years since it was first conceptualized and reported by Carl Wunderlich in a famous 1868 book. But, is this value statistically correct?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3>Exercises</h3>\n",
    "\n",
    "<p>In this exercise, you will analyze a dataset of human body temperatures and employ the concepts of hypothesis testing, confidence intervals, and statistical significance.</p>\n",
    "\n",
    "<p>Answer the following questions <b>in this notebook below and submit to your Github account</b>.</p> \n",
    "\n",
    "<ol>\n",
    "<li>  Is the distribution of body temperatures normal? \n",
    "    <ul>\n",
    "    <li> Although this is not a requirement for CLT to hold (read CLT carefully), it gives us some peace of mind that the population may also be normally distributed if we assume that this sample is representative of the population.\n",
    "    </ul>\n",
    "<li>  Is the sample size large? Are the observations independent?\n",
    "    <ul>\n",
    "    <li> Remember that this is a condition for the CLT, and hence the statistical tests we are using, to apply.\n",
    "    </ul>\n",
    "<li>  Is the true population mean really 98.6 degrees F?\n",
    "    <ul>\n",
    "    <li> Would you use a one-sample or two-sample test? Why?\n",
    "    <li> In this situation, is it appropriate to use the $t$ or $z$ statistic? \n",
    "    <li> Now try using the other test. How is the result be different? Why?\n",
    "    </ul>\n",
    "<li>  Draw a small sample of size 10 from the data and repeat both tests. \n",
    "    <ul>\n",
    "    <li> Which one is the correct one to use? \n",
    "    <li> What do you notice? What does this tell you about the difference in application of the $t$ and $z$ statistic?\n",
    "    </ul>\n",
    "<li>  At what temperature should we consider someone's temperature to be \"abnormal\"?\n",
    "    <ul>\n",
    "    <li> Start by computing the margin of error and confidence interval.\n",
    "    </ul>\n",
    "<li>  Is there a significant difference between males and females in normal temperature?\n",
    "    <ul>\n",
    "    <li> What test did you use and why?\n",
    "    <li> Write a story with your conclusion in the context of the original problem.\n",
    "    </ul>\n",
    "</ol>\n",
    "\n",
    "You can include written notes in notebook cells using Markdown: \n",
    "   - In the control panel at the top, choose Cell > Cell Type > Markdown\n",
    "   - Markdown syntax: http://nestacms.com/docs/creating-content/markdown-cheat-sheet\n",
    "\n",
    "#### Resources\n",
    "\n",
    "+ Information and data sources: http://www.amstat.org/publications/jse/datasets/normtemp.txt, http://www.amstat.org/publications/jse/jse_data_archive.htm\n",
    "+ Markdown syntax: http://nestacms.com/docs/creating-content/markdown-cheat-sheet\n",
    "\n",
    "****"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import stats\n",
    "\n",
    "df = pd.read_csv('data/human_body_temperature.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>temperature</th>\n",
       "      <th>heart_rate</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>130.000000</td>\n",
       "      <td>130.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>98.249231</td>\n",
       "      <td>73.761538</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>0.733183</td>\n",
       "      <td>7.062077</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>96.300000</td>\n",
       "      <td>57.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>97.800000</td>\n",
       "      <td>69.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>98.300000</td>\n",
       "      <td>74.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>98.700000</td>\n",
       "      <td>79.000000</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>100.800000</td>\n",
       "      <td>89.000000</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "       temperature  heart_rate\n",
       "count   130.000000  130.000000\n",
       "mean     98.249231   73.761538\n",
       "std       0.733183    7.062077\n",
       "min      96.300000   57.000000\n",
       "25%      97.800000   69.000000\n",
       "50%      98.300000   74.000000\n",
       "75%      98.700000   79.000000\n",
       "max     100.800000   89.000000"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mean = np.mean(df[\"temperature\"])\n",
    "std = np.std(df[\"temperature\"])\n",
    "n = len(df[\"temperature\"].dropna())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "samples = np.random.normal(mean, std, size=10000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def ecdf(data):\n",
    "    \"\"\"Compute ECDF for a one-dimensional array of measurements.\"\"\"\n",
    "\n",
    "    # Number of data points: n\n",
    "    n = len(data)\n",
    "\n",
    "    # x-data for the ECDF: x\n",
    "    x = np.sort(data)\n",
    "\n",
    "    # y-data for the ECDF: y\n",
    "    y = np.arange(1, n+1) / n\n",
    "\n",
    "    return x, y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x, y = ecdf(df[\"temperature\"])\n",
    "x_theor, y_theor = ecdf(samples)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEKCAYAAADpfBXhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xl8XNV5//HPM6N932VZiy1b8iLb\nGNuyTYBgQwhLQu1ACAGalTSkScnyyy9tya8NTUmbkqVNm4YsTgItSYAAbYlDTBwIEMABYxmvsi1b\nlmVttnbJ2peZ5/fHjKmwJWskj3RHo+f9es3Lc+8czXyvJT06c+bcc0VVMcYYE15cTgcwxhgTfFbc\njTEmDFlxN8aYMGTF3RhjwpAVd2OMCUNW3I0xJgwFVNxF5AYRqRCRShG5d5THvyMie/23oyLSEfyo\nxhhjAiXjzXMXETdwFHg3UAfsAu5Q1UNjtP8ssEpV7wpyVmOMMQEKpOe+DqhU1SpVHQQeBzZfoP0d\nwGPBCGeMMWZyIgJokwvUjtiuA9aP1lBE5gGFwAtjPH43cDdAfHz8miVLlkworDHGzHa7d+9uUdXM\n8doFUtxllH1jjeXcDjylqp7RHlTVLcAWgNLSUi0rKwvg5Y0xxpwlIicDaRfIsEwdkD9iOw9oGKPt\n7diQjDHGOC6Q4r4LKBaRQhGJwlfAt57bSEQWA6nAa8GNaIwxZqLGLe6qOgzcA2wHDgNPqGq5iNwv\nIptGNL0DeFxtmUljjHFcIGPuqOo2YNs5++47Z/urwYtljDHmYtgZqsYYE4asuBtjTBiy4m6MMWHI\nirsxxoQhK+7GGBOGApotY0w4eXRnzUU/x53rC4KQxJipYz13Y4wJQ1bcjTEmDFlxN8aYMGTF3Rhj\nwpAVd2OMCUNW3I0xJgxZcTfGmDBkxd0YY8KQFXdjjAlDVtyNMSYMWXE3xpgwZMXdGGPCkBV3Y4wJ\nQ1bcjTEmDFlxN8aYMGTF3RhjwlBAxV1EbhCRChGpFJF7x2hzm4gcEpFyEXk0uDGNMcZMxLhXYhIR\nN/Ag8G6gDtglIltV9dCINsXAl4ErVLVdRLKmKrAxxpjxBdJzXwdUqmqVqg4CjwObz2nzSeBBVW0H\nUNWm4MY0xhgzEYEU91ygdsR2nX/fSIuARSKyQ0ReF5EbRnsiEblbRMpEpKy5uXlyiY0xxowrkOIu\no+zTc7YjgGJgI3AH8BMRSTnvi1S3qGqpqpZmZmZONKsxxpgABVLc64D8Edt5QMMobX6lqkOqegKo\nwFfsjTHGOCCQ4r4LKBaRQhGJAm4Htp7T5mngagARycA3TFMVzKDGGGMCN25xV9Vh4B5gO3AYeEJV\ny0XkfhHZ5G+2HWgVkUPAi8BfqmrrVIU2xhhzYeNOhQRQ1W3AtnP23TfivgJf9N+MMcY4zM5QNcaY\nMBRQz92Yi/XozpqgPM+d6wuC8jzGhDvruRtjTBiy4m6MMWHIirsxxoQhK+7GGBOGrLgbY0wYsuJu\njDFhyIq7McaEISvuxhgThqy4G2NMGLLibowxYciKuzHGhCFbW8bMGo1n+tlb28Erx5rpHfQwMOwl\nPtpNckwk+WlxZCfFOB3RmKCx4m7Clqqyt7aD3+w/xW/LT1PX3vfWYy6BqAgX/UPet/ZlJUazqiCV\nyxemE+m2N7VmZrPibsLOsMfLtoOn+cFLxzl86gxRbhfvLM7grisKWZmfwr7aDqIjXIgIwx4vnX1D\nHGvqZl9dB9vLT7P7ZDs3r8qlMCPe6UMxZtKsuJuw8trxVv726QMcb+5hYWY8D9yyghtX5JAcG/lW\nm4rTXW/dj3C7SE+IJj0hmssWpHOssYun99bz41eqeM/yOVxZbBdyNzOTFXcTFrr6h/j7Xx/iqd11\n5KfF8sMPrea6kjm4XDKh5ynOTuRz7yrmv3bXse3gaUSEK4oypii1MVPHiruZ8apbevizR8o40dLD\nZzYu5LPXFBMb5Z7080VHuPng2gJ0Vw2/OXAKt0u4bEF6EBMbM/WsuJsZ7fWqVj71s92IwM8+sY7L\nFwanl+12CR9cm88vXq/hmf0N5KfGkZsaG5TnNmY62JQAM2PtPtnGxx/eRWZiNFv/4sqgFfazIlwu\nPlCaR0J0BE/srmXI4x3/i4wJEVbczYx0sL6Tjz28iznJMTz2ycsoSI+bkteJi4rg/avzaO4a4Hfl\np6fkNYyZCgEVdxG5QUQqRKRSRO4d5fGPiUiziOz13/4s+FGN8TnV2cfHHn6DxOgIfv5n68lMjJ7S\n1yvOTuSyBWnsON5KbVvvlL6WMcEybnEXETfwIHAjUALcISIlozT9pape6r/9JMg5jQF8c9g/99ge\negc9PPKJdeSmTM84+PXL5hAfHcH2Q9Z7NzNDID33dUClqlap6iDwOLB5amMZM7rvPH+UXdXtfP3m\nFRRlJU7b60ZHuLl6cSZVzT1UNnVP2+saM1mBFPdcoHbEdp1/37neLyL7ReQpEckf7YlE5G4RKROR\nsubm5knENbPZ8eZuvv/ScT5Yms/7Vo32Izi11s1PIzk2kt8dOo2qTvvrGzMRgRT30c4COfcn+9fA\nfFW9BHge+M/RnkhVt6hqqaqWZmbamX8mcEMeL0/vqWdeWhxf3bTMkQwRbhfvWpJFXXsf28sbHclg\nTKACKe51wMieeB7QMLKBqraq6oB/88fAmuDEM8bnlWPNtPYMcv/m5Rd1gtLFWlWQSnp8FD/8w3HH\nMhgTiECK+y6gWEQKRSQKuB3YOrKBiOSM2NwEHA5eRDPbtfUM8lJFMytyk7lqkbPv+Nwu4R0L09lb\n28G+2g5HsxhzIeMWd1UdBu4BtuMr2k+oarmI3C8im/zNPici5SKyD/gc8LGpCmxmn2f2N+ByCe9Z\nkTN+42mwuiCV+Cg3j7x20ukoxowpoOUHVHUbsO2cffeNuP9l4MvBjWYM1LT2cOR0F9eVZL9tZUcn\nxUS6uWV1Hr8sq+Vv3ruUtPgopyMZcx47Q9WEtOcONxIf5eYdC0Nr4a6PvGMeg8NefrmrdvzGxjjA\nirsJWVUt3Rxv7mHDokyiI5z7EHU0xdmJXL4wnZ+/fhKP16ZFmtBjxd2EJFXl+UNNJEZHsD5El9u9\nY10B9R197KxqdTqKMeex4m5C0omWHqpbe9iwODNkr2d67dJsEqIjeHpvvdNRjDlPaP7WmFnv1coW\n4qLcrJ2f5nSUMcVGubl+2RyePXCa/iGP03GMeRsr7ibktHQNUHG6i/WF6SHbaz/r5lW5dA0M88KR\nJqejGPM2of2bY2alHcdbcLmEyxaEbq/9rHcsTCcrMZqn99jQjAktVtxNSOkdHObNmnYuzUshMSY0\n5rVfiNsl/MnKubxY0URH76DTcYx5ixV3E1J2Vbcz5FGuKAruJfOm0s2rchnyKM8etLXeTeiw4m5C\nhleVN060siAjnjnJMU7HCdiyuUkUpMWx3S7DZ0KIFXcTMo43ddPeO8TawtAfax9JRLiuJJs/VrbS\n1T/kdBxjACvuJoS8Ud1GXJSbZTlJTkeZsOuXz2HQ4+WlCrsIjQkNVtxNSOjqH+LwqTOsLkglIsSn\nP45mdUEqGQlRNjRjQsbM+y0yYenNmg68SkiftHQhbpdw7dJsXqpoZmDYTmgyzrPibhznVWVXdRuF\nGfFkJkY7HWfSrl82h+6BYf543NaaMc6z4m4cd7K1l7aeQUrnpTod5aJcXpROfJSb39n1VU0IsOJu\nHLenpp0ot4tlc5OdjnJRoiPcbFycxfOHG1G1ZYCNs6y4G0cNebwcqO9keW4SUREz/8fx6iVZNHcN\nUN5wxukoZpab+b9NZkY7dOoMA8NeVhXM7CGZszb4L+D9UoUtJGacZcXdOGpPTTvJsZEUZsQ7HSUo\nMhOjuSQv2ea7G8dZcTeOOdM/xLHGblblp+AScTpO0GxcnMWbNe22kJhxlBV345j9tR0ocGlBitNR\ngurqxZl4FV4+1uJ0FDOLBVTcReQGEakQkUoRufcC7W4VERWR0uBFNOFqX10nuSmxZCXOnEXCAnFJ\nXgpp8VG8ZBfwMA4at7iLiBt4ELgRKAHuEJGSUdolAp8DdgY7pAk/Ld0D1Hf0cUnezJ7+OBq3S9iw\nKJOXjjbj9dqUSOOMQHru64BKVa1S1UHgcWDzKO2+BnwT6A9iPhOm9tV1IPh6ueFo4+JM2noG2V/f\n6XQUM0sFUtxzgdoR23X+fW8RkVVAvqo+E8RsJkypKvtrO5mXHk9ybOhfbWkyzl5s5NVjNmvGOCMi\ngDajTWN4672miLiA7wAfG/eJRO4G7gYoKCgILKEJO6fP9NPcPcDlRekT/tpHd9ZMQaLgy0iIpiQn\niVeOtXDPNcVOxzGzUCA99zogf8R2HtAwYjsRWA68JCLVwGXA1tE+VFXVLapaqqqlmZmZk09tZrR9\ntZ24BJbP8OUGxvPO4gzerGmnZ2DY6ShmFgqkuO8CikWkUESigNuBrWcfVNVOVc1Q1fmqOh94Hdik\nqmVTktjMaKrK/voOirISiI8O5I3jzHVlcQZDHuWNE21ORzGz0LjFXVWHgXuA7cBh4AlVLReR+0Vk\n01QHNOGlpq2Xjt4hVobpB6kjrZ2fRnSEi1dsvrtxQEBdJ1XdBmw7Z999Y7TdePGxTLjaV9dJhEso\nmYGX0puomEg36wrTeLXSPlQ108/OUDXTxuNVDtR3smROItGRbqfjTIsrizI42thN4xmbIWymV3gP\nepqQUtXSTc/AcFjMbQ901k5Xv+/D1G9tr2D1OStf3rneZoyZqWM9dzNt9td2Eh3hYvGcRKejTJs5\nyTHER7mpbOp2OoqZZay4m2kx7PFSfqqTZXOTiHTPnh87lwgLsxI43tRtV2cy02r2/JYZRx1t7KJ/\nyBsWQzITVZyVQNfAMI1nBpyOYmYRK+5mWuyt6yQ+OoKFmQlOR5l2RVm+YajKpi6Hk5jZxIq7mXJd\n/UMcOXWGS3KTcbvC56IcgUqOjSQzIZrKZht3N9PHiruZctvLGxn2KivzZ9+QzFlFWQmcaOlh2ON1\nOoqZJay4myn3q731pMZFkp8a63QUxxRlJTDkUU629TodxcwSVtzNlGruGmBHZQsr81OQMLpO6kQt\nyIjHJdiUSDNtrLibKfWb/Q14lVmxlsyFREe6yU+Ls+Jupo0VdzOlfrWvgaU5SWQnhdd1UiejKCuB\nho4+em0JYDMNrLibKVPT2suemg42XzrX6SghoTgzAQWOt/Q4HcXMAlbczZT51d56AP5kpRV3gNzU\nOKIjXDY0Y6aFFXczJVSVp/fWs64wjdyU2TtLZiS3S1iQmWAnM5lpYcXdTIlDp85wvLnHhmTOUZSV\nQHvvEG09g05HMWHOiruZElv3NhDhEt6zPMfpKCGlyL/8wjHrvZspZsXdBJ3Xq2zd18CGRZmkxkc5\nHSekZCREkRwbyXEbdzdTzIq7CbrXqlo51dnP5lW5TkcJOSJCUWYCx5t78HhtCWAzday4m6B7ancd\niTERXFeS7XSUkFSUlUDfkIfyhk6no5gwZsXdBFVX/xDPHjzFppVziZkl10mdqIVZvnH3V461OJzE\nhDMr7iaoth04Rf+Ql1vX5DkdJWQlREeQkxzDjkor7mbqBFTcReQGEakQkUoRuXeUx/9cRA6IyF4R\neVVESoIf1cwET+2uY2FmPJfO4uV9A7EwM4Gy6nb6Bj1ORzFhatziLiJu4EHgRqAEuGOU4v2oqq5Q\n1UuBbwL/EvSkJuRVt/Swq7qdW9fkz+oVIANRlJXAoMfLG9VtTkcxYSqQnvs6oFJVq1R1EHgc2Dyy\ngaqeGbEZD9g0gFnoibJaXAI32yyZcc1PjyfK7bKhGTNlIgJokwvUjtiuA9af20hE/gL4IhAFXBOU\ndGbGGPJ4eaKsjmuWZDEn2VaAHE9UhIs181LtQ1UzZQLpuY/2/vq8nrmqPqiqC4G/Bv521CcSuVtE\nykSkrLm5eWJJTUj7/eFGWroHuHN9gdNRZowrizM4fOoMLd0DTkcxYSiQ4l4H5I/YzgMaLtD+ceB9\noz2gqltUtVRVSzMzMwNPaULeL3bWMDc5hg2LspyOMmNcWZQBYEMzZkoEUtx3AcUiUigiUcDtwNaR\nDUSkeMTme4FjwYtoQl1tWy+vHGvhg2sLcLvsg9RALc9NJjk20oq7mRLjjrmr6rCI3ANsB9zAQ6pa\nLiL3A2WquhW4R0SuBYaAduCjUxnahJbHd9XgErhtrc1tnwi3S7iiKJ2Xj7agqjbDyARVIB+ooqrb\ngG3n7LtvxP3PBzmXmSEGhj38clct1yzJIifZ1m2fqA2LMtl24DQVjV0smZPkdBwTRuwMVXNRfrP/\nFC3dg3z08vlOR5mRzn5G8VKFTTAwwWXF3UyaqvLwjmqKshLe+nDQTMyc5BiWzEnkpYomp6OYMGPF\n3UzamzXtHKjv5GOXz7fx4ouwcXEWZdXtdPUPOR3FhBEr7mbSHt5RTWJMBLestjNSL8bGxZkMe5Ud\nla1ORzFhxIq7mZRTnX08e/A0t6/NJy4qoM/lzRjWzEslMTqCPxy1oRkTPFbczaT89JUTAHzkHfOd\nDRIGIt0urijK4KWKZlRtWSYTHFbczYR19g7x2Bs13HRJDvlpcU7HCQsbF2dyqrOfo412bVUTHFbc\nzYT97PVqegY9fOqqhU5HCRsbFvuW47BZMyZYrLibCekf8vDwjmo2LMqkZK6ddBMsOcmx/imRNt/d\nBIcVdzMhT+6uo7VnkD/fYL32YNuwOJOyk210Dww7HcWEASvuJmCDw15+8GIlqwpSuGxBmtNxws7G\nRVkMedQWEjNBYcXdBOzJ3bU0dPbzhWsX2UlLU6B0fioJ0RE2NGOCwiYom3E9urOGYY+Xf37uKPmp\nsdS19fLozhqnY4Ud35TIdP5Q0WSrRJqLZj13E5DdNe109g1x7dJsKzpTaOPiLBo6+znWZFMizcWx\n4m7GNeTx8lJFMwVpcRRlJTgdJ6xttCmRJkisuJtx7axqtV77NDk7JfL5w1bczcWx4m4u6Ez/EC9W\nNFOclWC99mny7pJsyqrbaOsZdDqKmcGsuJsL+tEfjtM35OH6ZXOcjjJrvLskG6/Ci0es924mz4q7\nGVPTmX5++uoJLslLZm6KXUJvuqzITWZOUgzPHWp0OoqZwWwqpBnTN7dX4PEq716a7XSUsHSh6aTz\n0uN44UgT//nHaiLdF+6D3bm+INjRTBiwnrsZ1d7aDp7aXcddVxSSnhDtdJxZZ2lOEoMeL8ebbUqk\nmRwr7uY8Xq/y1a3lZCREc881RU7HmZUWZMQTHeHi8Kkup6OYGcqKuznP/+ypZ29tB399w2ISYyKd\njjMrRbhdFGcncuTUGbx2AQ8zCQEVdxG5QUQqRKRSRO4d5fEvisghEdkvIr8XkXnBj2qmQ/fAMA/8\n9ggr81N4/+o8p+PMastykugaGKa2rdfpKGYGGre4i4gbeBC4ESgB7hCRknOa7QFKVfUS4Cngm8EO\naqbH916opLlrgK/+SQkul52w5KTFcxJxu4SD9Z1ORzEzUCA993VApapWqeog8DiweWQDVX1RVc92\nL14HrMs3A1W39PDQqye4ZXUuqwpSnY4z68VEuinOSuBgwxm7tqqZsECKey5QO2K7zr9vLJ8Anh3t\nARG5W0TKRKSsudmWNQ01//CbQ0S6hXtvWOJ0FOO3PDeZzr4h6tr7nI5iZphAivto781H7UaIyIeA\nUuBboz2uqltUtVRVSzMzMwNPaabcC0caef5wE/dcU0xWUozTcYzf0jlJuAQONtjQjJmYQIp7HZA/\nYjsPaDi3kYhcC/wNsElVB4ITz0yH3sFhvvJ0OcVZCXziykKn45gRYqPcFGUlUG5DM2aCAinuu4Bi\nESkUkSjgdmDryAYisgr4Eb7CbgtizDDfee4o9R19fP2WFURF2OzYULNsbjJtPYOc6ux3OoqZQcb9\nTVbVYeAeYDtwGHhCVctF5H4R2eRv9i0gAXhSRPaKyNYxns6EmIP1nTy0o5o71uWzdr5dFzUUleT4\nhmb219nQjAlcQGvLqOo2YNs5++4bcf/aIOcy08DjVf7f/xwgNS6Se29Y6nQcM4b46AiKshLYX9/B\ndcuycdma+iYA9h58FvvZa9Xsr+vkKzeVkBxnZ6KGspV5KXT0DtkJTSZgtipkmBtr5cHOviG+8/xR\nirMS6O4ftgteh7iSnCQiXMK+ug7mpcc7HcfMANZzn6V+va8BVWXzpbl26bwZIDrSzdKcJA7UdeLx\n2qwZMz4r7rPQgfpODp06wzVLskmLj3I6jgnQyrxkegY9tgywCYgV91mme2CYrXvryU2J5cqiDKfj\nmAlYlJ1ITKSLfbUdTkcxM4AV91nm1/sa6B/y8v41ebhtYbAZJcLtYvncZMobzjAw7HE6jglxVtxn\nkYP1nRyo7+TqJVnMsSUGZqQ181IZ9Hg5WH/G6SgmxFlxnyV6B4b51b4G5ibHsGGRreszUxWkxZGR\nEMXuk+1ORzEhzor7LPHMgVP0DQ7bcMwMJyKsLkilurWH1m5bwsmMzYr7LHD41Bn21nZw9eIscpJj\nnY5jLtKqglQEeLPGPlg1Y7PiHua6B4b5nz31zEmKYcNiG44JB8mxkRRlJbCnpt2ur2rGZMU9jKkq\n//1mHf1DHm4rzSfCZd/ucLFmXiodfUNUNtmcdzM6+20PY4++UcOR011cv2wOc5Jtdkw4KclJIj46\ngp1VrU5HMSHKinuYOt7czdeeOURRZgLvWJjudBwTZBFuF6XzUjlyuov6DrsEnzmfFfcwNDjs5QuP\n7yUm0s2ta/Jsidgwtc6//v7jb9iib+Z8VtzD0L8+f5QD9Z08cMsKkmJtKd9wlRofxaLsRB57o5bB\nYa/TcUyIseIeZnZWtfKDPxznttI8blie43QcM8XWL0ijpXuA3x067XQUE2KsuIeRpq5+PvvYHual\nxfF3f7LM6ThmGizKTiQ/LZaHXj3hdBQTYqy4h4lhj5fPPrqHM/1D/OBDa4iPtuuwzAYuET5xRSFv\n1nSw+2Sb03FMCLHiHia+ub2CnSfa+PrNK1iak+R0HDONblubT3JsJFternI6igkhVtzDwG8PnmLL\ny1V86LICblmd53QcM83ioiL48GXz+N2hRqrsQh7Gz4r7DHe8uZsvPbmflfkpfOWmEqfjGId89PL5\nRLpc/NTG3o1fQMVdRG4QkQoRqRSRe0d5/CoReVNEhkXk1uDHNKPpHRzm0z/fTaRb+MGfriY6wu10\nJOOQzMRoblmdy5O762g80+90HBMCxi3uIuIGHgRuBEqAO0Tk3C5iDfAx4NFgBzSj83iVzz++l8qm\nbr57xyrmpthqj7PdZzYW4fUq33+x0ukoJgQE0nNfB1SqapWqDgKPA5tHNlDValXdD9iZFNPkgWcP\n89yhRr5yUwnvLLbVHg0UpMfxgdI8Hnuj1pYkMAEV91ygdsR2nX/fhInI3SJSJiJlzc3Nk3kKA/xi\n50l+/MoJPvqOeXz8ikKn45gQcs81xQB87wXrvc92gRT30RYmmdQi0qq6RVVLVbU0M9N6m5Px24On\n+crTB7l6caZ9gGrOk5sSy+3r8nmyrJaa1l6n4xgHBVLc64D8Edt5QMPUxDEX8trxVj73+B5W5qfw\n4J+uJsJtk53M+f7i6iIi3S4e+O1hp6MYBwVyGuMuoFhECoF64HbgzilNZQB4dOf/rvZX29bLQztO\nkBIbyXtX5PD0Hvv7akaXnRTDZzYu5J+fO8ofK1u4vCjD6UjGAeN2/VR1GLgH2A4cBp5Q1XIRuV9E\nNgGIyFoRqQM+APxIRMqnMvRsU9fey8N/PEF8dAQfv6KQuChbWsBc2CevWkB+Wixf/XU5wx6b5zAb\nBfS+XlW3qeoiVV2oqv/o33efqm7139+lqnmqGq+q6apqq1YFSV27r8ceG+nmz64sJNmW8DUBiIl0\n87fvLeFoYzc/e/2k03GMA2zQNoQdb+7mJ6/6C/s7F5ASF+V0JDODXFeSzVWLMvn29gpq2+zD1dnG\ninuI2l5+mv/4YzUpsZHcfdVCUq2wmwkSEf7plhW4RPi/T+7D653UJDczQ1lxDzGqykOvnuDTP9/N\n3OQY7r5qgQ3FmEnLTYnl7zYt440TbTy0w9admU3sk7kQMuzxcv8zh3jktZNcV5LN5QsziIqwv7/m\n4rx/dS7by0/zze0VrC9MZ0VestORzDSwyhEimrsG+NBPd/LIayf51FUL+OGH1lhhN0EhIjxwywoy\nE6L55CNlNHXZwmKzgVWPEFBW3cZN//4Ke2s7+JfbVvLl9yzF5RrtxGBjJic9IZotH1lDZ98Qn/rZ\nbvqHPE5HMlPMiruDBoe9fGv7EW770WvERLr5709fYRfbMFNm2dxkvvPBleyp6eALj+9lyOa/hzUr\n7g451tjFzd/fwYMvHuf9q/N45rNXUjLXLo9nptYNy3O476YSflt+mnsefZPBYSvw4co+UJ1mXq/y\n8B+r+cZvj5AQHcGPPryG65fNcTqWmUXuurIQEfj7Xx/iM7/Yzb/dvsouqB6G7Ds6jQ6fOsNXnj5I\n2cl2rl2axT/dcgmZidFOxzKz0MevKCTCJfzd1nJu/v4Otny4lPkZ8U7HMkFkwzLToHtgmH945hA3\n/furVLX08O0PrOTHHym1wm4c9eF3zOeRu9bT1DXApu+9ytN76lG1E53ChTj1zSwtLdWysjJHXnuq\nnV3NUVU52HCG3+xvoKt/mLXz07huWbYt/GWC6s71BRf19bVtvdz549epbe9jUXYCm1bmkhY/8TOi\nLzaHCYyI7FbV0vHaWZWZIvXtfWw7eIoTLT3MTY7hT9fPIz8tzulYxpwnPy2OT21YyOtVrfzuUCPf\nee4opfNTuXpxFkl2dvSMZcU9yBo6+niyrJY9tR3ER7nZtHIua+en4bZ56yaEuUS4fGEGy+Ym82JF\nE7uq2yg72c6K3GTWF6ZRkBaHiP0MzyRW3IOks2+IH79cxY9fqcLjVTYsymTDokxiIt1ORzMmYMmx\nkbzv0lyuKs7k1coW9tS0s7e2gzlJMawrTGNVfgrR9jM9I1hxv0i9g8M8vKOaLS9X0dk3xOZL57Io\nO9FWcTQzWlp8FJtWzuX6Zdnsr+1k54lWtu5r4NmDp1gyJ4mVeSksyk6wSz2GMCvuk9Q/5OEXO2v4\nwUuVtHQPcs2SLL747kUsz01Ioj+XAAALkklEQVR+2+XxjJnJoiPcrC1Mo3R+KnXtfbxZ086B+k4O\n1HcSE+li2dxkVualsCDTplGGGivuE3Smf4jHdtbw0I4TNJ4Z4PKF6fzow4tZMy/V6WjGTBkRIT8t\njvy0OG66ZC7Hm7vZV9vBgfpOdp9sJzE6gsqmbjZdOpdV+Sk2Ph8CbCpkgBrP9PPQjhM8+noNXQPD\nXFGUzl9sLBr14sPWczezxZDHy5HTXeyv6+BYUzeDw17y02K5cXkO71qSxZp5qRMaugnW7044T8u0\nqZBB4PUqr1e18uTuOp7Z34DHq7xnRQ6fumqhrYltDBDpdrEiN5kVucnctDKH35U38ut9DTy84wRb\nXq4iJS6Sqxdn8a6lWVy+MGNS8+fN5FhxH0VNay9PvVnHf+2uo76jj8SYCO5YV8AnrixkXrqNLRoz\nmqSYSG5dk8eta/Lo6h/ilWMtPH+okRcqmvifPfUAFGclsKoghZKcJBbPSaIgPY7sxOiAevcerzIw\n5KF/2Ev/kIf+YQ8DQ2fve32PDfkWQmvuGiAuyk1MlJu4SDcZidHMTY4hJyWWhFmyjs7sOMpxeLzK\noYYzvFjRxPOHG9lf14kIXFmUwV/dsJjrl82xKY3GTEBiTCTvWZHDe1bkMOzxsq+ug50n2thZ1cbz\nh5t4oqzurbZul5AaF0lSTCR9Qx68qnjV985ZgaFhL/3DHoY84w8hu0VA4OVjzRfIFkFuSiwFaXEU\nZsQzPyOe+enxzM+IIzsxJmyupRBQcReRG4B/A9zAT1T1gXMejwYeAdYArcAHVbU6uFGDw+NVTnX2\ncbSxiwN1Z9hX18Gu6ja6+ocRgUvzU/irGxZz86pccpJjnY5rzIwX4XaxZl4aa+al8ZmNvmU5mroG\nONrYRV17H3XtvbT3DtHZN8Txpm5cIrhcgkt8J1dFuoWYCDfRkS5iIt1j3o+OcBHpfwdwW2ke/cNe\negeH6R3w0NQ1wKnOPho6+v3/9lHV0sNLFc0MjljXPibS5Sv06fHMy4ijMN1X/Asz4slKjJ5RHxSP\nW9xFxA08CLwbqAN2ichWVT00otkngHZVLRKR24FvAB+cisAerzI47GXQ42XI42Vw2PfvkMfLwLDv\n1tk3RGfvEB29g3T0DdHSPUBNWx81rT3Ud/S91QMQgQUZ8dx0yVwuW5DG5QszbDEvY6aYiJCdFEN2\nUsx5jwXrA9UIt4sEt8s3BJPImCteerxKQ0cfJ1t7OdHaQ3WL73a0qYvfH2l827uF2Eg389LjyEmO\nIS0+mvSEKNLifbeE6AhiI93ERLqJjXITF+V+azvSffaPleAWweWCCJdrys9aD6Tnvg6oVNUqABF5\nHNgMjCzum4Gv+u8/BXxPRESnYCrOlper+MZvj0zoa5JjIylIi2PZ3GRuXJFDQVocCzMTKJmbNGvG\n34wx53O7/neK55XFb5/5drbwn2jp4WRrDydaeqlu7aHxTD9HTnfR2jM46Yud/MP7lvOhy+YF4xDG\nFEhlywVqR2zXAevHaqOqwyLSCaQDLSMbicjdwN3+zW4RqZhM6MnYP/EvyeCc/DOUHUdoCdvj+FOH\ngoxmAlkc+X58+Bvw4cl/eUB/FQIp7qO9dzi3Rx5IG1R1C7AlgNd0nIiUBTKXNNTZcYQWO47QEi7H\nMZpAzi6oA/JHbOcBDWO1EZEIIBloC0ZAY4wxExdIcd8FFItIoYhEAbcDW89psxX4qP/+rcALUzHe\nbowxJjDjDsv4x9DvAbbjmwr5kKqWi8j9QJmqbgV+CvxMRCrx9dhvn8rQ02RGDB8FwI4jtNhxhJZw\nOY7zOLa2jDHGmKljizEbY0wYsuJujDFhyIq7n4h8XkQOiki5iHxhxP7PikiFf/83ncw4ntGOQUR+\nKSJ7/bdqEdnrdM7xjHEcl4rI6/7jKBORdU7nHM8Yx7FSRF4TkQMi8msRSXI657lE5CERaRKRgyP2\npYnIcyJyzP9vqn+/iMh3RaRSRPaLyGrnkr/dBI9jif/7MiAiX3IudRCp6qy/AcuBg0Acvg+ZnweK\ngav996P97bKczjrRYzinzT8D9zmddZLfi98BN/rbvAd4yemskzyOXcAGf5u7gK85nXWU7FcBq4GD\nI/Z9E7jXf/9e4BsjvhfP4jvX5TJgp9P5J3kcWcBa4B+BLzmdPRg367n7LAVeV9VeVR0G/gDcDHwa\neEBVBwBUtcnBjOMZ6xgAXw8LuA14zKF8gRrrOBQ428tN5vxzLULNWMexGHjZ3+Y54P0O5RuTqr7M\n+eepbAb+03//P4H3jdj/iPq8DqSISM70JL2wiRyHqjap6i5gaPoSTi0r7j4HgatEJF1E4vD1RvKB\nRcA7RWSniPxBRNY6mvLCxjqGs94JNKrqMUfSBW6s4/gC8C0RqQW+DXzZwYyBGOs4DgKb/G0+wNu/\nR6EsW1VPAfj/zfLvH215ktxpzjYRYx1H2LFVswBVPSwi38DXk+oG9gHD+P5/UvG93VwLPCEiC9T/\nPi6UXOAYzrqD0O+1X+g4Pg38H1X9LxG5Dd+5Fdc6l/TCLnAcdwHfFZH78J38N+hcyqAIaOkRM/2s\n5+6nqj9V1dWqehW+t3LH8PVC/tv/lvMNwItvoaGQNMYxnF0S4hbgl07mC9QYx/FR4L/9TZ7Et1pp\nSBvtOFT1iKpep6pr8P2xPe5syoA1nh1u8f97dogykOVJQslYxxF2rLj7iUiW/98CfIXwMeBp4Br/\n/kVAFCG8ot8YxwC+Hu4RVa0b62tDyRjH0QBs8De5Bv8frlA22nGM2OcC/hb4oXMJJ2TkEiMfBX41\nYv9H/LNmLgM6zw57hKixjiP8OP2JbqjcgFfwrVG/D3iXf18U8HN846RvAtc4nXOix+Df/x/Anzud\n7yK/F1cCu/37dgJrnM45yeP4PHDUf3sA/1nioXTD98f0FL4PF+vwXYwnHfg9vj+qvwfS/G0F38V8\njgMHgFKn80/yOOb425wBOvz3k5w+hou52fIDxhgThmxYxhhjwpAVd2OMCUNW3I0xJgxZcTfGmDBk\nxd0YY8KQFXcTkkQkRUQ+43SO8YjIF/zLCxgTUqy4m1CVAjhe3P0n51zo9+QL+FZ+nMhz2rIfZspZ\ncTeh6gFgoX/99m+JyF+KyC7/muF/DyAi80XkiIj8xL9u+i9E5FoR2eFfr3udv91XReRnIvKCf/8n\nz77IBZ73sIh8H9/Ja/ki8gP/OvLlI9p9DpgLvCgiL/r3dY947ltF5D/89/9DRP7F3+4bIhLvX298\nl4jsEZHN0/B/amYR60GYUHUvsFxVLxWR64Bb8a0nI8BWEbkKqAGK8K2ueDe+tdLvxHc26ybg//G/\nS9Negm8BuHhgj4j8Bt+a68VjPO9i4OOq+hkAEfkbVW0TETfwexG5RFW/KyJfBK5W1UCWpVgEXKuq\nHhH5OvCCqt4lIinAGyLyvKr2XMT/mTFvseJuZoLr/Lc9/u0EfEW5BjihqgcARKQc+L2qqogcAOaP\neI5fqWof0OfvPa/D90dgrOc9qb71yc+6TUTuxvc7kwOUAPsneBxPqqpnxDFtGnHVnxigADg8wec0\nZlRW3M1MIMA/qeqP3rZTZD4wMGKXd8S2l7f/fJ+7zoaO87w9I7YLgS8Ba1W13T/UEjNG1pGvc26b\nkb1yAd6vqhVjPI8xF8XG3E2o6gIS/fe3A3eJSAKAiOSeXV1xAjaLSIyIpAMb8Q3hBPq8SfgKc6eI\nZAM3jpETfEvKLvV/CHszY9sOfFZExP/aqyZ4PMZckPXcTUhS1Vb/B6MH8V2j81HgNX8t7AY+BHgu\n8BTnegP4Db6hj6+pagPQICJLx3teVd0nInuAcqAK2DHi4S3AsyJySlWvxvdZwTP4rk50EN9Qz2i+\nBvwrsN9f4KuBmyZwPMZckK0KacKeiHwV6FbVbzudxZjpYsMyxhgThqznbowxYch67sYYE4asuBtj\nTBiy4m6MMWHIirsxxoQhK+7GGBOG/j8NdY71Eyzn7AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x2b28e7a6a90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.distplot(df[\"temperature\"])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEFCAYAAADqujDUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xl8U1XawPFfliZpmy50YacFWiqL\nQGWRHUXEBXVcUAEVl+EVZ7TgCIqICgwvIKCD+zivjiKio+iIiLuCCAKCLBZkXy2l0I3uaZM2uef9\noxJpaUvRpkmT5/vRT3PvTe59TkPvk3ty7nN0SimFEEKIgKP3dgBCCCG8QxKAEEIEKEkAQggRoCQB\nCCFEgJIEIIQQAcro7QDqKyenGIBmzULIzy/1cjR/nLTDd/hDG0Da4Wt8pR2xsWG1bmtyVwBGo8Hb\nITQIaYfv8Ic2gLTD1zSFdjS5BCCEEKJhSAIQQogAJQlACCEClCQAIYQIUJIAhBAiQEkCEEKIAOXR\nBLBjxw7GjRt31vpvv/2WUaNGMXr0aN5//31PhiCEEA3G7nRwtPAYdqfD26E0CI/dCPbaa6+xcuVK\ngoODq6yvqKjgqaee4r///S/BwcGMHTuWYcOGERsb66lQhBB+TClFuVPDUe7CXuGivMKFy6VwujRc\nWrWfLoVTq3xueUXleqUULk2hKYX260+X9utjDfe6cq2cPUErsesLMWsRdLZfi56gX2NwR8Pph2ZT\nEA5HhXv5zMr7px+qaivO2I37sdGg49qB7Wkba23oX53nEkBcXBwvvvgiU6dOrbL+8OHDxMXFERER\nAUDv3r3ZunUrV199dZ37a9YsxH1jRV13tjUl0g7f4Q9tgKbfDlt5GRsP7cFlCyUnr5zi0gpK7RWU\n2p0Ul5ZjL3dR5nBidzixlzspc7iwlztpjFlNdKEFWLoVAuDQF7L56GGULdLzBwYGJbflIg+8tx5L\nAFdeeSXHjx8/a31JSQlhYb81JDQ0lJKSknPu7/Qt1bGxYe6yEE2ZtMN3+EMboGm1w1Hh4lhOAfuy\nj+EoCiE3v4KMU0WcavENumAbWlkojt0DQKt6itIBZpMBs8mAJciANSKoyrLZZMAUZMCo12M06DAY\nfv2p12E06DEa9Bj0OkxBesxBBvT6ym16nQ79mT+rr9frcGrlvHHwELmOXGLMMUy47TLMBnNlUIDu\n1we6X5ejo63k5dmqxq/7rR2nF3RnNu6MZd2v2w16HcFm4+9+b+v6UNDotYCsVis222+/FJvNViUh\nCCH8S0GJg8MZRaRlFZOZV8qxrGJyCosxdfsB/emT/f4BmMNK0QdXnhv0wTZuvboFic3aYzEbCDEb\nCbEYMQcZ3CfGxhfKY1EPctKWRavQFliM5jqfHRMZjKpwNlJsv0+jJ4CEhATS0tIoKCggJCSErVu3\nMn78+MYOQwhxnuxOxzlPfpqmOJ5Twv70Ag5nFHI4o5BTRVW/MLWYDMS315F1xsl+0h0duCC2HU9v\nO0RWaTYtQppzaZfO5zzJNjaL0UyHiDhvh9FgGi0BfPLJJ5SWljJ69GimTZvG+PHjUUoxatQoWrRo\n0VhhCCHqofrJ3u50sHDri+6T89Q+E7EYzWhKcTijkL2/5HPweAGHTxRhL3e592MNDqJnQjQd20TQ\nsVU4raJDiAwzU+4qZ+HWHe79dW4eh8VoZmqfidhNxVjKw3zu5O+PdE1lUvjT/V9NqZ+zLtIO3+EP\nbYCGa0dNJ/uTtiye2faS+zkP93qAgpwQ/vvdYU6e+q3kccuoEBLbRpDUNpJObSNo3iy41i6b2q4o\n5P1o+Dhq02TmAxBCNI6TtiyySrMByCrNdp+kW4Q0J6s0mwhDFG9/coKjGWXodToGdGtBr6TmJLWL\nICzEVO/j+Ft3SlMkCUCIAFf9k3iUJRKDzoBLuTDoDERZIrEYzYxqfScf/JDKsTRAK+OiTjHcNLQj\nbTwwPl00DkkAQgSwmrp78uwFuFRlP75LuTicm8WmH9PYuj8HMNIrKZZbLk2gRVSId4MXf5gkACEC\n2Lm6e0KI5NX30yh36OnQKpwxwxPp1LZxbn4SnicJQIgAUt/untvi7mHxmh85mWHAajYz7ppEBnRr\niV7vrTH4whMkAQgRIOrT3ZNXls+6raf4eP0vaMrM0J6tGDs8CbPJ9+e3FedPEoAQAeJc3T2xlhg+\n+y6P7fvyiQ43M+7KC+iREOPlqIUnSQIQIkDU1t0zpdcDfLxlJ9//WMIxez4JbcJJuakHEaH1H9Ip\nmiZJAEIEiExbdpXunkxbNkYVzIv/3cOB4zaswSZuvbIjl/RsLX39AUISgBB+6ly1e4pLy5m/fDsZ\nuTYu6hTDnVd1lk/9AUYSgBB+yO50sGDL82SX5dI8OIZH+z5IfHg7YoNjyCnLpVlQFO9+kkN+oYvL\ne7dlzOWd0HutyqbwFpkTWAg/lFaUTnZZLgDZZbmkFaVjMZqZ1vdBro29nVNb+1JQ6OLGIR0YKyf/\ngCVXAEIEkB0HCvjv56cwGgw8cFM3eiXJVKyBTBKAEH7A7nRw8FQuFmdlGeWWoc3Ro0dDQ4+elqHN\n2fDzSd74bC8Ws4HJtyaT0CbC22ELL5MEIEQTV9sNXhoaABoa3+w4wOffFhFqMTJlTDLtW4Z7OWrh\nC+Q7ACGauLpu8AKIMEbx+Xd5hIeamDxaTv7iN5IAhGjiTt/gBVS5wWtqn4lcE3M7OT/2xmQwMfnW\nnnRoJSd/8RtJAEI0cWfV87EXALD/l2I+/OIUBoJIuak7cS1qnxlKBCb5DkCIJq6mEg8H0gt4ZcUu\njAY9j9x2EQmt5QtfcTZJAEI0cdVLPBzMOc6SD3NxaYqUm7rLyV/USrqAhPAz7317CJvdyZjhneiZ\nKNU8Re0kAQjRxNidDo4WHsPudAC4x/wDoHTkZRu5ql8cl/Vq48UoRVMgXUBCNCHnGvOPTtGpo5lb\nLk1AJ+UdxDnIFYAQTUhtY/6bBUUDoHdYmXBFPzn5i3qRKwAhfFh95vA16oLQHRyM3ZbFwzcMIdpq\n9XbYoomQBCCEj6rXHL72Ar7YmkVGtoMr+vWkSzsp7ibqT7qAhPBR5yrx0CKkOYWnTHy9JZ2WUSGM\n/1M3b4YrmiC5AhDCR9U2h+/UPhM5acsi0hjNnDdTMeh13HtdV0IsQdiK7d4OWzQhcgUghI+qaQ5f\nAIvRTIeIOD75Pp38YgfXDmwvNX7E7yIJQIgmaF9aPmtTT9AmJpSR/eO9HY5ooiQBCOGjzrzB6/Sk\nLgA2ewWvfboHnU7H3Vd3Jsgof8bi9/HYvxxN05gxYwajR49m3LhxpKWlVdn++uuvc9NNNzFq1Ci+\n+eYbT4UhRJOVacuuMqlLpi0bpRRLv9pPfrGD6we3l1m9xB/isS+BV61aRXl5OcuWLSM1NZX58+fz\nyiuvAFBUVMTSpUv5+uuvKSsr44YbbmDEiBGeCkUIv7FxVyY/7s0msU0EIwdI14/4Yzx2BbBt2zaG\nDBkCQHJyMrt27XJvCw4OpnXr1pSVlVFWViZ3LQpRg+pdQOW2EJZ8uR+zycC913XFoJeuH/HHeOwK\noKSkBOsZdyQaDAacTidGY+UhW7VqxTXXXIPL5eK+++475/6aNQvBaKyc9Sg21j8mtpB2+A5fbEN2\n1okqXUDvb0jF6bLw+D396dqpeY2v8cV2/B7SjsbhsQRgtVqx2WzuZU3T3Cf/devWkZ2dzerVqwEY\nP348vXr1okePHrXuLz+/FKj8hebkFHsq7EYj7fAdvtqG/ILSKssnTtkY0C2e+JiQGuP11XacL2lH\nw8dRG49dQ/bq1Yt169YBkJqaSlJSkntbREQEFosFk8mE2WwmLCyMoqIiT4UiRJNUpcyzBoZyKzcO\n7ejdoIRf8dgVwIgRI9iwYQNjxoxBKcW8efNYvHgxcXFxDB8+nI0bN3Lrrbei1+vp1asXgwYN8lQo\nQjRJZ44CQg/J3YKJiQj2blDCr3gsAej1embPnl1lXUJCgvvxpEmTmDRpkqcOL4Tf6X9hC2+HIPyM\nDCMQwkdUn+mrRUgsqMoRcjp0JETLDF+iYUkxOCF8QE2ln1ftPAg6BYBCkWcvIMIsNX9Ew5ErACF8\nQPXSz8cKT/D9lhKUPRSoLP3cKlS6gETDkisAIXxA9dLP23fZKC7WuDp0DH2Sg90zggnRkCQBCOED\nqpd+/m7vPiKtLbm2fwLBZvkzFZ4hXUBC+KAKl+LGIR3l5C88ShKAED7gzJu+lAYtg5szqHsrL0cl\n/J0kACF8wJk3fen0MKivFb1eiiQKz5IEIIQPkjr/ojFIAhDCB5QWWlDuqg+/zf4lhCdJAhDCy1ya\nxvsbU9G5675p7gnghfAkSQBCeNn3O0+SXVDm7TBEAJIEIIQX2cudrPj+KEHOsBongBfCkyQBCOEF\npwu/rdr+C0W2cvolR1SZ/SvPXuDlCEUgkLtMhGhkdqeDBVueJ7ssFxyhWCyDub5Pd9J/3uAuBid1\nf0RjkAQgRCNLK0qvPPkDmG0k9wgiKtTK1D4TOWnLkro/otFIAhDCyy7uUvlp32I00yEizsvRiEAi\n3wEI0chahjZ3T/QCOjpESckH4R2SAIRoZLszjrknegElY/6F10gCEKIRKaVYtT3d22EIAUgCEKJR\n7T6axy+/KHcXkIz5F94kCUCIRuJ0abz37SH0wTZ3F5CUfRDeJAlAiEay8/ApTuTa6N4h2tuhCAFI\nAhCi0azedhyAET26SNkH4RMkAQjRCDJyStiblk/nuEgMIbYqZR+kC0h4iyQAIRrB8nVHALiir9zo\nJXyH3AkshIedPGXjp4O5dGgVTs/EaIrKTejRo6FJF5DwKrkCEMLD1mzPAOCqfnHodLoq8/9KF5Dw\nJrkCEMKD8orsfLczjYhYO107hnk7HCGqkCsAITzo29Q0DJ03UN7he57Z9iJ2p4P48HbEBscAEBsc\nQ3x4Oy9HKQKVXAEI4SGFJQ5W7f0ZQ6dSALLLckkrSueCqESm9X1QSj8Lr5MEIISHfLjuCBVOhaGG\nbVL6WfgCjyUATdOYNWsW+/fvx2QyMWfOHOLj493b165dy8svvwxA165dmTlzJjqdrrbdCdEk2J0O\nTtqyCCWSTbuziImMoVRG/Agf5bEEsGrVKsrLy1m2bBmpqanMnz+fV155BYCSkhKefvpp3nrrLaKi\nonjttdfIz88nKirKU+EI4XF2p4OFW18kqzSbECJxqj706RHFd7aqc/1GmMO9HKkQlTyWALZt28aQ\nIUMASE5OZteuXe5tP/30E0lJSSxYsID09HRuueWWc578mzULwWisvJiOjfWP0RTSDt/REG04eCqX\nrNLKIZ2lFBAe7WDssIs5uHYdGcWZtAlrSY/4BCxBlj98rNr4w3sB0o7G4rEEUFJSgtVqdS8bDAac\nTidGo5H8/Hw2b97MihUrCAkJ4fbbbyc5OZkOHTrUur/8/Mov0mJjw8jJKfZU2I1G2uE7GqoNOocJ\ng86AS7lQmo7LuiXiLNWYfNED7i98iwsqKKaiAaI+mz+8FyDt8EQctfHYMFCr1YrNZnMva5qG0ViZ\nbyIjI+nevTuxsbGEhobSp08f9u7d66lQhGgUmbZsXMoFgE6v6JhQecV6+gtfGe0jfI3HEkCvXr1Y\nt24dAKmpqSQlJbm3XXjhhRw4cIC8vDycTic7duwgMTHRU6EI4RVBBrnNRvg2j3UBjRgxgg0bNjBm\nzBiUUsybN4/FixcTFxfH8OHDmTJlCv/zP/8DwFVXXVUlQQjRFIUZolAa6PRS5lk0DR5LAHq9ntmz\nZ1dZl5CQ4H58zTXXcM0113jq8EI0us9/2oXu1w/9p2v8yIgf4cvkGlWIBlBSVsH2/TneDkOI8yIJ\nQIgG8M2WdOzFwYBM9i6aDkkAQvxBjnIX324/TkikA5DJ3kXTIQlAiD/o+50nsNmd9EqM8XYoQpwX\nSQBC/AFOl8ZXP6ZjMuq5MrmrTPYumpQ6E0BWVlZjxSFEk7RlXzaniuwM6dGaEi1fZvoSTUqdCeAv\nf/mL+/Ebb7zh8WCEaErsFXZW/vQTeoOLKy+WSV1E01NnAlBKuR9/8sknHg9GiKbC7nQw+4fnKGq9\nhvCLNmO16mWmL9Hk1Hkj2Jn1+c9MBkIEurSidAqdeQA49EUy05dokup9J7BM1iLEb9Kzaq7yKDN9\niaakzgRw8OBBhg8fDlR+IXz6sVIKnU7H6tWrPR+hED5GKcW6zcWotlL3RzRtdSaAr776qrHiEKLJ\n+OlgLhnFuVjOqPsjM32JpqjOBNCmTRsADhw4wJEjR7BYLCQkJNCunXy5JQLXmp8yUA4LegxouDDo\nDERZIr0dlhDnrc4EcOrUKSZNmsTBgweJj49Hp9Nx9OhRkpOTWbRoEWFhvj3dmRAN7VBGIbuP5hHX\nEXKonPzFpVxS+VM0SXUOA/3HP/5B79692bBhAx988AHvv/8+GzZsoHPnzsydO7exYhTC6+xOB0cL\nj/HZpkMADLywlZcjEuKPq/MK4KeffuKLL76oss5kMjF58mSuv/56jwYmhK+wOx0s2PI82WW5aJYQ\n4ltfydBOXdhYGENOWa6M+RdNVp0JwGyueRyzTqdDr5cyQiIwpBWlk12WC4A+uJRePYOwBFlkzL9o\n8up9I9j5bBPCn3VqGwHImH/R9NX7PoAzKaXIyZHZj0RgqBzjr6Oy1r+OVtYWXo5IiIZxzvsACgoK\ncLlcREdHA7B582YSExPdy0L4u6N5Jzg90QsoGfEj/EadHfmFhYXcd999FBQU0KZNG9q0aUNaWhoP\nPPAAxcU13wovhL/5ca+UdRb+qc4EsGDBAv7xj38wdOhQ97qHHnqIefPmMX/+fI8HJ4S3FdrK2b6z\nlF/L/EvZB+FX6kwARUVF9OvX76z1Q4YMIT8/32NBCeErPl5/lIqgYvdfikz0IvxJnQnA6XSiadpZ\n6zVNo6KiwmNBCeELSsoq2PjzSSJCTd4ORQiPqDMB9O3bl5deeums9f/85z+58MILPRaUEL5gzU8Z\nlDs1hl6QJHP9Cr9U5yigyZMnM2HCBFasWEHnzp0xm83s2bOHqKgoXnnllcaKUYhG53RpfLv9OMFm\nAwkd9Xyzu+pcvzIKSPiDOhOA1WrlnXfeYdOmTezduxe9Xs/tt99Onz59Gis+Ibxi6/5sCkvKGdGn\nHaYguelR+Kdzzgim0+kYMGAAAwYMaIx4hPA6pRSrth4H4LLebbCEONGjR0OTLiDhV6SgjxDVHEgv\n4MiJIpITY2jRLIRMWzYaVbuAhPAHkgCEqGblDwfRhRYwop+UfBb+TRKAEGc4mpnHEevnWLpt4v3j\ni7E7HcSHtyM2OAZASj8Lv3LO7wB+L03TmDVrFvv378dkMjFnzhzi4+PPes6ECRMYPnw4Y8eO9VQo\nQtTbBz9uQ28tBSC7LJe0onQuiEqU0s/CL3nsCmDVqlWUl5ezbNkypkyZUmPpiOeee47CwkJPhSDE\neUnLLGZ/ekGN206XfpaTv/AnHrsC2LZtG0OGDAEgOTmZXbt2Vdn+5ZdfotPpqtQZEqKx2Z0ODp7K\nxeIM4/NNaaiyUHToUCgZ8SP8nscSQElJCVar1b1sMBhwOp0YjUYOHDjAp59+ygsvvMDLL79cr/01\naxaC0WgAIDbWPyajl3Z4l73CzrxvFpFRnEnL0OakH06mVWs9Bb+WftbQUMHlxEY3nfY11feiOmlH\n4/BYArBardhsNveypmkYjZWHW7FiBVlZWdx1111kZGQQFBREmzZt6rwayM+v7JeNjQ0jJ6fpl6KW\ndnjf0cJjZBRnApBpy8ZlKmJIp5786NpFVmk2LUKaYylvOu1ryu/FmaQdDR9HbTyWAHr16sWaNWsY\nOXIkqampJCUlubdNnTrV/fjFF18kJiZGuoJEo4uyRGLQGXApF2g6jK4QhnSP41LjRPnCVwQEjyWA\nESNGsGHDBsaMGYNSinnz5rF48WLi4uJqnGZSiMaWacuuPPkD6BW9e4QQagkCkLl+RUDwWALQ6/XM\nnj27yrqEhISznjdx4kRPhSDEeenftaW3QxCiUcmNYCJgtQxtjo7Thd50xDeTBCACiyQAEbAyS7JQ\n1SZ7FyKQSAIQASstq8TbIQjhVZIARMBK/dmOksneRQCTBCACUlpmMftzj6OTyd5FAJMEIALSF5vT\nvB2CEF4nCUAEnOz8Urbsy6ZNaAuZ7F0ENEkAIuB8+kMaSkGvnsEy05cIaJIAREA5VWjnh12ZtIoO\noXNcpLfDEcKrPHYnsBC+6Ksfj+HSFCP7x9PKGiqTvYuAJlcAImAUlZazbscJosLN9Ovagjx7QZUu\noDx7zZPBCOGvJAGIgPHB2v1UmPMY3qclRoPeXQ0UwKAzEGWRLiERWCQBiIBw/FQB27TlWLptYou2\nHLvTQZ69wF0N1KVccgUgAo4kABEQPt2+C11w5QRFWaU57nr/LUIq+/1bhDSnVWgLb4YoRKOTL4GF\n3ytzONm514auqw70yt3dYzGamdpnInZTMZbyMJn8RQQcuQIQfm/jrkwchiLQV1b+dCmXe8y/xWim\nU3QHOfmLgCQJQPi1UnsFK74/QpBe/qkLUZ38VQi/tmL9UWx2J1f07CJlH4SoRhKA8FvZ+aWs2Z5B\nTISFzklBUvZBiGokAQi/ZHc6eH3ND7io4OZLEzDoded+kRABRkYBCb9jdzr43x+epSAyD2tPKz0S\nh4IuktjgGHLKcokNjiE+vJ23wxTC6yQBCL/zS+ExCiryAHAFlZBWfJwLohKZ1vdB9/h/GfUjhCQA\n4YcOHC+scb3FaKZDRFwjRyOE75LvAIRfcZS72PBjscz1K0Q9SAIQfuWLzWnklufKXL9C1IMkAOE3\n8orsfPnjMULM0rMpRH1IAhB+4+P1Rymv0Liu94Vy05cQ9SAJQPiFzLxSNvxcOdVjXBxy05cQ9SAJ\nQDR5Sin+880BNKW4cUhHdDq56UuI+pDOUtHkfb8rnd3ZR+jaoR29L4ilqNwic/0KUQ+SAESTlltc\nzLL0N7F0K6HYfAiHq0eNc/1GmMO9HKkQvke6gESTtnTtNrCUAJDryJWZvoQ4Dx67AtA0jVmzZrF/\n/35MJhNz5swhPj7evf3NN9/ks88+A+CSSy4hJSXFU6EIP7V+50l+3ldK8EU60J0905eUfRCibh67\nAli1ahXl5eUsW7aMKVOmMH/+fPe29PR0Vq5cyXvvvceyZctYv349+/bt81Qowg/lFdn5z6oDBIeX\nga7mmb46RMTJyV+IOnjsCmDbtm0MGTIEgOTkZHbt2uXe1rJlS/79739jMBgAcDqdmM11/6E2axaC\n0Vj5/NjYMA9F3bikHefPXmHnWOEJ3lmTgb3cxU0DO/JF9mb39sjIkN8Vj7wXvkXa0Tg8lgBKSkqw\nWq3uZYPBgNPpxGg0EhQURFRUFEopFi5cSNeuXenQoUOd+8vPLwUqf6E5OcWeCrvRSDvOn93pYMGW\n58kuy0UzhnBB+5FcfkEXthb/VuY5Qos+73jkvfAt0o6Gj6M2HksAVqsVm83mXtY0DaPxt8M5HA6m\nT59OaGgoM2fO9FQYwo+kFaWTXZYLgD64lEu7h2IJskiZZyF+J499B9CrVy/WrVsHQGpqKklJSe5t\nSinuv/9+LrjgAmbPnu3uChKiTkpVWYy0Vp7spb9fiN/HY1cAI0aMYMOGDYwZMwalFPPmzWPx4sXE\nxcWhaRo//vgj5eXlfP/99wBMnjyZiy66yFPhCD9w4JALpYFOLzV+hGgIHksAer2e2bNnV1mXkJDg\nfvzzzz976tDCD+1Ly2fl1l2YulQun67xIzd4CfH7yZ3AwmfZnQ5O2rIIN0Txyse70Mm/ViEalPxJ\nCZ905ogfsxZOsf1ibhzak23qiEzsLkQDkQQgfNKZI34c+iJati1n5MUJjNBkxI8QDUVqAQmfYHc6\nOFp4DLvTUeP2kQPiMej1MuJHiAYkVwDC6+xOBwu3vkhWaTYtQpoztc9EWlhagSMUzDYijVH0aZd0\n7h0JIc6LJADhdSdtWWSVVtbwySrN5mRJJt+sK6Fs3wAGX2zltoG95RO/EB4gCUB4XZQlEoPOgEu5\nMOgM7D/sYNOeLDq2juLOIb0wGqSnUghPkL8s4XWZtmxcygVUVvRcvnkn1uAg/vKnbnLyF8KD5K9L\n+Byl4P4bLiQmMtjboQjh1yQBCK9rGdoc/a//FJUGf+p7IZ3jm3k5KiH8nyQA4XWZJVnuOXx1ekhK\nkOKAQjQGSQCi0Z055l/TFN9sO171CTqddwITIsDIKCDhcadr+pyenP30mP/mwbFEZ13O9gMVhPaw\noplKpMSDEI1IrgCER52+yeuZbS+xcOuLpBWlu8f8Z5fl8FP6URJbRzN7yGQe7p3CtL4Pyph/P7R9\n+1auuupSsrIy3eteeeVFPv/8E48e9+TJE0yYcPdZ6y+9tD8pKRNISZnAvffexb///S80rbIbcvr0\nR2rdX2ZmJuvXrztr/aZNG/n44+W1Hq82a9euITc3h1OncnnmmfnnfkEDkysA4VHVb/ICaBYUTX7F\nKbSyULq3iuf+65MxBRlohrWuXYkG8v63h9iyL7tB99m3c3NuvSyxzucYjUHMmzeb5557GZ2Xu/nC\nwyN46aVXgcoJqp5+eh7Ll7/PzTePYd68p2t93fbtW0hL+4XBg4dWWd+//0CgMuGcjw8+eJf27acT\nH9+ehx+edp6t+OMkAQiPqn6T18kMHSd/uAhdcAkjundl1FVJMtY/QPTu3QdNUyxf/j6jRo2usu3d\nd99m9eqvMRgMDBjQj7vv/guvv/5/7Nq1k7KyMqZNe5K5c2fRokULTp48yfDhV3D06GEOHNjPwIGD\nue++B/jpp20sXvwaAHa7nSee+DtBQUHnjEun0zFmzB089dRsbr55DH/605WsXPkVy5d/wBdffIpe\nr6dHj2T+8pcU3n77Tex2O9279+C9994hMrIZxcXFjBhxBenp6dxwwygKCvJ59NGHKCkpom/fAdx9\n9/8wd+4shg+/gv79B7Jp00ZWr/6aYcMu59ChA8yZM4Mnn/xf5syZyauvvsmWLZt49dVXMJvNhIdH\n8NhjMzh4cD/vvPMWQUFGTp48wWWXjeCuu8b/4fdEEoDwqOo3eb29bjsmY3Om3NSPxDYRXo4uMN16\nWeI5P617ysMPT+Pee+/i4ov6nDIzAAAXsElEQVQHuNcdPnyIb7/9hn/96w0MBgOzZ09nw4bKmQLj\n4zvwt789zMmTJzh5MoNnn30Zh8POLbdcz4oVn2M2W7j55uu4774HOHr0CDNm/C8xMbG89dYbrFmz\niiuuuLpecUVFRVNYWFBl3eeff8Lf/vYIF17YnY8++i9KKe644+5frwAu4b333mHEiKu45JJhVbqy\nysrKePLJ/yUurjm33jqGQYOGVj8cAAMHDiYxMYlHHpnuTlRKKRYunMc///lvYmOb8/7777JkyesM\nHDiYrKyTvPnmu1RUVHDDDVc1SAKQj16iQZ2rqqfFZGTK6GQ5+QeoiIhIJk2awrx5s1Cqss89Le0X\nunXrjtFoRKfT0adPH44ePQxAXFy8+7WtWrXBarVitYYRFRVFeHgEZrPZ3Z0UGxvLc889zdy5s9i+\nfStOp7PecWVmniQ2tkWVddOnz+Djjz8kJWUCmZkna3zdmfGdlpjYCavVisFgoEuXbqSnH6uyXVWb\n2/pMBQUFhISEEhtbOd1pcvJFHD16BICOHRMxGo0EBwdjNlvq3ba6SAIQDeb0JC7PbHuJBVuex+50\nYFUx6Mt/7dt3hPLI9ZfKyT/ADR48lHbt4vn8808BiI9vz549u3A6nSil2LJlC+3aVZ5Y9frfvis4\n1/cGCxbMYfr0mTz++CxiYmLrHY+mabz77lIuv/yKKutXrlzBww8/xksvvcrBg/v5+ecd6HQ6d+Kq\njO/sU2ha2i+UlpbidDrZs2cXHTp0xGQycepU5fwWBw7sq/L6018+A0RGRlJaaiM3t/K5qanbadcu\n7tf217tJ9SZdQKLBnDmJS3ZZLqv37OKL1TZKK/rTo5uJ24b0pnl4mJejFL7gwQensG3bFgASEhK5\n7LLL+etfx6OUol+/vgwdeimHDh04r31eeeVIJky4m7CwMJo1iyY3N6fW5xYVFZKSMgG9Xo/T6aRv\n335ce+31VZ6TkJDIvffeSWRkM2JjY+na9UJCQ0N56603SErqXOu+w8LCmTnzMUpKirj00uF06NCR\n6667gaeems3XX3/pPqEDXHhhD+bMmcnUqY8DlUlu6tTHefzxR9DrdYSFhTN9+iyOHDl0Xr+L+tKp\nuq5HfEhOTjEAsbFh7sdNmT+2Y3/eIV5IfdW9zb63L0Flsdw+IonBPVp5feRHbfzxvWjKpB0NH0dt\n5ApA/G52p4ODp3KxOMOwGM3umj4aGkqDGHMsE2/pQ9vmMrxTCF8kCUD8LtVn8Xq41wN8unVvlZo+\n91wfT9tYOfkL4avkS2BRL9VH91S/wWveB2tZvaEQHJUn/BYhsbRv1tpr8Qohzk2uAMRZzqzdYzGa\na5yzt1VoC6JM0eSVV97ReyJDz4Aubbhh8EOUaPnu1wohfJckAFFFTSf76p/2v9m1h9QdTjKyKu/o\nbRXagrG3dKFbhygAYpGRPkI0BZIAAlz1T/tnTdD+67YYcwy5jlx0DisffZ2DThnpFh/L9ZcOJKFF\nqM+O8BFC1E4SQAA7feNWdlkuzYNjeLTvg7QKbUGLkOZklWYTY4kh64Selbv2kX40GV1wCUHOcIYl\nt+PKi+NoHhnsM0PdRMOr/uHgj9i+fSszZjxG+/Yd3OsiI5sxZ86C373PmTMf44knZlep93Pw4H7W\nr1/HPffcW+/9fP75J6Sl/cJf/zqx1nidTie33DKW4cNH1LqfzMxMDh06cFahOF8mCSCAVb9x61B+\nGgZbDImlIynLSCP9pIFXtcqbcRLbRjO4+4Vc3KU5FpP8s/F3NXUF/tEk0Lt3H/7+96caKEJq3Fen\nThfQqdMFDbL/M+MtLS0lJWUCcXFxte6/tkqhvkz+kv3A7/2k5qjQqiw/9/4OtOJoAEzGYLq3b0bn\nuEjaNbfSrUOUdPMEkJq6AjtExJ3jVb9PSsoEEhOTOHr0MMHBwQwc2J9vv/2OkpISFi16ifXr1/L9\n92spLbVRUFDAPff8D5deOpybb76Od975L8888xSFhYUUFRUyduw4vv32a/7+96f49NMVfPTRh2ia\ni8GDL2H8+Pv48MNlrF27BqfTidVqZe7c2ks/nykkJITrr7+JNWtW07FjIk8/PY/s7CwKCwvp338g\nf/7zhCqVQkNDrbz99htUVLjclUlrqhvkbZIAmpj6jNCpngSUUhSUlJORU8Kx7BJOnrJxIreUo1l5\nmLuFoA8uBXsoiVFxxCU1o0dCNEntIjEFydy8gerMrsAWIc3ds7n9Edu2bSUlZYJ7eeDAwdx2250A\ndO3ajb/97WEmT56IxWLhuef+yZw5M0lN3Q5AWVkpzz77MgUF+dx7710MHnxJlX337t2H0aNvZ/v2\nrQDk5+fx9ttLWLLkXYKCTLz00rPYbCUUFhby3HP/RK/XM3lyCnv37q53/FFRURw4sI/s7Cy6devO\ntGlP4nA4uOmmkdx771+rVApdvvwDnn76afT6EHdl0oao3tnQPJYANE1j1qxZ7N+/H5PJxJw5c4iP\n/y0Dvv/++7z33nsYjUb++te/MmzYME+F4jfqM0Jn46GDBDtjKChxkJVfxi+ZRWTll+Eod1XZlw5I\nahdDfPBoWrZ2cXHHBEJNwV5olfBFFqPZ/e+roYb01tUFdLq2TliYlcTExF8fh1NeXnnfSXJyL/R6\nPVFR0YSFhVNQULV0c/VP1xkZGXTokOCumjlp0hQAgoKCmDXrcYKDg8nOzj7PiqGZxMY2Jzw8nL17\nd7N9+1ZCQ0MpL68467mxsbHMnTsXvT6InJxsunfvWe/jNCaPJYBVq1ZRXl7OsmXLSE1NZf78+bzy\nyisA5OTksHTpUj788EMcDge33XYbgwYNwmQyeSQWTSmUUpyuenT6sVKg+O0xKDT122uo/O+M56ua\nt4F7+6+b4IzHv63/LYa80gry80urrNc0RbnThb288n9HhYvyChflFRqOChd5rkyyjL+d7BcuX4fT\nFgotrGApQSsL5e2VJ0D7bbYng15Hq+hQWjQLpk1sKG1jrbRrYSXSasYsn/BFHSxGs8e6fao7V/fi\n/v2VFTTz8k5hs9lo1qxZtddXvae1TZu2HDv2C+Xl5ZhMJp54YiqjRo1m3brveO21JdjtdsaPv6Pe\n8ZWW2vjkk4+YM2cBn3/+KVZrGFOnPs7x4+msXPkRSqkqlUIXLJjD6tWrKStTzJkzs97HaWweSwDb\ntm1jyJAhACQnJ7Nr1y73tp07d3LRRRdhMpkwmUzExcWxb98+evTo0eBxfL/zBEu+2F950m7q9E7M\n3ULRB9vQykL55ReFxVhBhPMyQiLLiLXE0nJAOJFWMxGhJqIjLLSOCZUZt4RPqN4FBPCPf7xQr9fm\n5Z3iwQf/SklJCVOmPIrBUPeHl2bNmnH77XeRkjIBnU7HoEFD6NKlG8HBwYwfPw6TKYjo6Jg6K4ae\njtdgMOByuRg//j7i4trjdLqYNWs6O3emYrFYaNu2Hbm5OSQkJLorhV555UhuvfVWQkJCz1mZ1Js8\nVg308ccf54orruCSSyr76i699FJWrVqF0Wjk448/5sCBAzzySOXky1OnTuWGG25g4MCBte7P6XRh\nNJ7/J9Y9R0/x7tf70TSFTlf5SUMH6PS//tTpKtfz609dTevOWNZTy+t0VfdffV3159b0HMBsMhBi\nCSLYbMRiMmA2GTAHVf60mIygc5LryKZDVFvCg0Oq1EsXwh8tX76cI0eO8PDDD3s7FL/jsSsAq9WK\nzWZzL2uahtForHGbzWYjLKzuu0fz80uB8y+xGms1Memm7ucTeqP4/ePnjbS1tKaiVONUaUmDx3W+\n/OE+AH9oA/hvO4qL7ZSWlje5tvnK+1FXOWiP9Q306tWLdevWAZCamkpSUpJ7W48ePdi2bRsOh4Pi\n4mIOHz5cZbsQQpw2cuR1VW7SEg3HY1cAI0aMYMOGDYwZMwalFPPmzWPx4sXExcUxfPhwxo0bx223\n3YZSioceegizWQqHCSFEY5IZwbxE2uE7/KENIO3wNb7SDq90AQkhhPBtkgCEECJASQIQQogAJQlA\nCCEClCQAIYQIUE1mFJAQQoiGJVcAQggRoCQBCCFEgJIEIIQQAUoSgBBCBChJAEIIEaAkAQghRICS\nBCCEEAHKY+WgG0J5eTmPPfYY6enpWK1WZsyYwYEDB1i4cCGtWrUCYOLEiVx88cVejrRuNbVDp9Mx\nc+ZMKioqMJlMLFq06Kx5Tn1JTW148skn3duPHDnCjTfe6POzNtXUjhMnTvDMM89gNBoZMGAADz30\nkLfDPKea2nH8+HGeeeYZgoODGTJkCPfff7+3w6zTjh07eOaZZ1i6dClpaWlMmzYNnU5Hp06dmDlz\nJnq9npdeeonvvvsOo9HI9OnTPTJt7B9RnzYApKWl8cADD/Dpp596OeJqlA9bunSpeuKJJ5RSSh0+\nfFj9+c9/VosWLVJffvmllyM7PzW1Y9y4ceqnn35SSin15Zdfqu3bt3szxHOqqQ2nHTt2TN14442q\npKTEW+HVW03tuP7669XBgweVpmlqzJgxat++fV6O8tyqt+Puu+9Wl1xyiTp27JhSSqkpU6aoLVu2\neDPEOr366qvq2muvVbfccotSSqn77rtPbdq0SSml1JNPPqm+/vprtWvXLjVu3DilaZrKyMhQN910\nkzdDPkt92qCUUh999JG68cYb1cCBA70Wa218ugvo0KFDDB06FICOHTty+PBhdu/ezYcffshtt93G\n/PnzcTqdXo7y3Kq3Y/fu3eTl5bFmzRrGjRtHamqqz32yqa6m9+K0uXPn8sgjjxAaGuqt8OqtpnZ0\n6dKFgoICKioqcDgc55xw3BdUb8f27dsJDw+nXbt2QOWMfNu3b/dmiHWKi4vjxRdfdC/v3r3bfSU/\ndOhQNm7cyLZt2xg8eDA6nY7WrVvjcrnIy8vzVshnqU8bACIiInj77be9EuO5+HQC6NKlC2vWrEEp\nRWpqKllZWQwcOJAnn3ySd955h9LSUt577z1vh3lO1duRn5/PwYMHGTBgAG+99RaFhYV89NFH3g6z\nTjW9Fy6Xi3379mGz2RgwYIC3Q6yXmtrRqVMn/vKXvzBy5EhatWpFx44dvR3mOVVvR3l5OXa7ncOH\nD+NyuVi3bh2lpaXeDrNWV155pXuOcAClFDqdDoDQ0FCKi4spKSnBarW6n3N6va+oTxsAhg0bRkhI\niFdiPBefTgCjRo3CarVy5513smbNGrp168bNN99Mu3bt0Ol0DB8+nD179ng7zHOq3o4LL7yQ0NBQ\n+vfvj06nY9iwYezatcvbYdappvfCYDCwcuVKbrnlFm+HV2/V2xEXF8drr73GZ599xqpVq4iPj+eN\nN97wdpjnVNP7sXDhQmbNmsWkSZPo0KGDT3+nVN3pvnIAm81GeHg4VqsVm81WZX1YWO2zW3lbTW3w\ndT6dAH7++Wd69+7N0qVLufzyy2nbti1/+tOfyMzMBOCHH36gW7duXo7y3Kq3Iy4ujvbt27N161YA\ntmzZQqdOnbwcZd2qt+F0V8OmTZsYMmSIl6Orv+rtSExMJCQkxP0JrXnz5hQVFXk5ynOr6f1Yt24d\n//d//8dLL73EsWPHGDhwoLfDrLeuXbuyefNmANatW0efPn3o1asX69evR9M0Tpw4gaZpREVFeTnS\n2tXUBl/n06OA4uPjef7553njjTcICwtj7ty5HDx4kJSUFCwWCwkJCdx6663eDvOcampHfn4+f//7\n33G5XLRt29bnR8/U1AaAnJycJvVJs6Z27Ny5kz//+c+YzWbCwsKYP3++t8M8p5rasXbtWsaOHYvF\nYuG6667z+Q8VZ3r00Ud58sknWbRoER07duTKK6/EYDDQp08fRo8ejaZpzJgxw9th1qmmNvg6KQct\nhBAByqe7gIQQQniOJAAhhAhQkgCEECJASQIQQogAJQlACCEClCQA4TWbN29m3LhxjfL64uJiHnjg\ngd99rMawZs0aFi9e7JF9Z2Zm8thjjwGwdu1ahgwZwpQpU5g/f36TuJlSeIZP3wcgREMpLCxk7969\n3g6jTp68G3zevHk8+OCDAHz55ZekpKQwevRo8vLymDRpks/WqhGeJQlAeFV+fj7jx48nOzubHj16\nMHPmTEwmE2vWrOG5555D0zTatWvH7NmziYmJYf369Tz11FOYzWY6dOgAVJbaveuuu/j222/R6/Vs\n3ryZ1157jX//+9/u48yZM4fs7GweeOABXn75ZVasWMGSJUvQNI1u3boxc+ZMzGYzgwYNYvjw4ezc\nuZOYmBhGjRrF0qVLyczMZP78+Vx88cWMGzeOzp07s3XrVhwOB9OnT2fw4MHk5uYyY8YMMjMz0el0\nTJkyhYEDB/Liiy+SmprKyZMnueOOO0hMTOTZZ5/FbrdTVFTEY489Rvv27d11rVq3bs2JEyeAynLn\nAJdddhlvvfUWP/74Ix999BEFBQUMGzaMO++8s8ZjnunYsWNkZ2eTkJDABx98wOrVq/nhhx/Q6/Xc\ncsstREVFsWnTJvr3798Yb7nwJV6rQyoC3qZNm1TPnj3V0aNHlaZp6sEHH1Rvvvmmys3NVYMHD1bp\n6elKKaVee+01NXHiROVwONSgQYPUoUOHlFJKTZ8+Xd1xxx1KKaVuv/12tXHjRqWUUtOmTVOfffZZ\nlWOlp6erYcOGKaWUOnDggBo7dqyy2+1KKaWeeeYZ9fLLLyullEpKSlLffPONUkqpO+64Q02ePFkp\npdTy5cvV/fff714/bdo0pZRSe/bsUYMGDVIOh0P97W9/U6tWrVJKKZWVlaWGDx+uiouL1QsvvOCO\nUymlJk6c6G7Dxo0b1bXXXquUUuqFF15QL7zwwlmPlVJq2LBhKj09XX344YdqxIgRqqKiQimlaj3m\nmRYvXqwWLFjgXn700UfVhx9+6F5esmSJmjt3bp3vlfBPcgUgvKpPnz60b98egOuuu47ly5cTFxdH\njx49aNu2LQCjR4/m1VdfZf/+/TRv3pyEhAQAbrzxRp5//nmgsjjaypUrSU5OZtOmTcyaNavWY27e\nvJm0tDR3GZGKigq6du3q3n66zHKbNm3o3bs3UPmp/MwaQadf26VLF2JjY9m/fz8bN27kyJEjvPDC\nCwA4nU7S09MBqpT7fvrpp1mzZg1ffvklO3bsqFLwrD66du3qrkJZ2zG7dOnifn5aWpr7aqkmrVu3\nZsOGDecVg/APkgCEV1Uvp2s0GtE0rcpzlFI4nU50Oh3qjMolZ9btv+qqq3j22Wf56quvGDp0KGaz\nudZjulwurr76ap544gmgsnKjy+VybzeZTDUe40xnrtc0zR33kiVLiIyMBCA7O5vo6GhWrVqFxWJx\nP/+2226jX79+9OvXjwEDBtRYB0qn01X5PVRUVLgfn7mv2o5ZfV9n/p6rMxqN7jLGIrDIKCDhVdu2\nbXNXelyxYgUDBw6kZ8+e7Nixg+PHjwOwbNky+vXrxwUXXEBubi779u0D4LPPPnPvJzg4mKFDh7Jo\n0SJuuumms45jNBrdkwf169ePb775hlOnTqGUYtasWSxZsuS84v7888+ByqqcRUVFJCUl0b9/f/7z\nn/8AlRO2XHfddZSVlVV5XUFBAb/88gsPPvggQ4cOZfXq1e7kYzAY3DE2a9aMQ4cOAbBz505ycnJq\njKM+x4yLiyMjI6PWtmRkZBAfH39e7Rf+Qa4AhFclJiYyffp0cnJy6N+/PzfffDMGg4HZs2eTkpJC\nRUUFrVu3Zu7cuQQFBbFo0SIeeeQRjEZjlW4bgGuuuYbt27fTs2fPs44THR1N69atGTduHEuXLiUl\nJYW77roLTdPo0qULEyZMOK+409PTufHGGwF49tlnMRgMPPHEE8yYMYPrrrsOgIULF1aZ0AQgMjKS\nm2++mWuuuQaj0Uj//v2x2+2UlpbSt29fHn30UWJiYrj22mv56quvGDlyJN26dTurrafV55jDhg2r\ns9rs5s2bueOOO86r/cI/SDVQ4RdcLhfPPvss0dHR3HPPPR491rhx40hJSaFfv34ePU5DSklJYdKk\nSSQlJVVZf+rUKVJSUnj33Xe9FJnwJukCEn5h1KhR7N69m7Fjx3o7FJ/02GOP8frrr5+1/l//+hfT\np0/3QkTCF8gVgBBCBCi5AhBCiAAlCUAIIQKUJAAhhAhQkgCEECJASQIQQogA9f9gZ1msbTfm0wAA\nAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x224f9e30710>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.set()\n",
    "_ = plt.plot(x_theor, y_theor)\n",
    "_ = plt.plot(x, y, marker=\".\", linestyle=\"none\")\n",
    "_ = plt.xlabel(\"body temperature (f)\")\n",
    "_ = plt.ylabel(\"CDF\")\n",
    "_ = plt.legend(('Normal Distribution', 'Empirical Data'), loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "NormaltestResult(statistic=2.7038014333192031, pvalue=0.2587479863488254)"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stats.normaltest(df[\"temperature\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Is the distribution of body temperatures normal? \n",
    "\n",
    "The normal CDF overlaid with the data CDF shows promise that the data may be normally distributed. It's very close to the curve, but there are some spots with visible variance. A null hypothesis asserts that the temperature data is distributed normally. SciPy's handy normaltest function, produces a p-value of 0.2587, far too high to reject the null hypothesis that the data conforms to a normal distribution. The temperature data therefore is normally distributed.\n",
    "\n",
    "## 2. Is the sample size large? Are the observations independent? \n",
    "\n",
    "As far as the Central Limit Theorem is concerned, the sample size of 130 observations is sufficiently large, as it is above the commonly espoused minimum of 30 samples. The observations can be declared indepedent as well, because one person's temperature does not affect another's. The data is naturally unordered.\n",
    "\n",
    "## 3. Is the true population mean really 98.6 degrees F?\n",
    "\n",
    "Since the sample data not being compared to another group, a one-sample test can be used. The data is assumed to follow the normal distribution. That condition, plus the large sample size (> 30 observations) allow the preferable Z-Test.\n",
    "\n",
    "Null Hypothesis: The true population mean is equal to the conventional mean 98.6 F\n",
    "\n",
    "Alternative hypothesis: The true population mean is not equal to the conventional mean 98.6 F\n",
    "\n",
    "For this test, a 95% confidence level will be used"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The 95% confidence interval is (98.123682111456645, 98.37477942700491)\n"
     ]
    }
   ],
   "source": [
    "# z_critical value for normal distribution given a two-tailed error\n",
    "# the z-critical value for 95% confidence is commonly represented as 1.96 as per textbook lookup tables\n",
    "z_critical = stats.norm.ppf(q = 0.975)\n",
    "\n",
    "# sigma represents the population standard deviation for which we use the sample standard deviation\n",
    "sigma = std\n",
    "\n",
    "# margin of error for the sample mean is computed as such\n",
    "moe = z_critical * sigma / np.sqrt(n)\n",
    "\n",
    "# show the 95% confidence interval\n",
    "ci = (mean - moe, mean + moe)\n",
    "\n",
    "print(\"The 95% confidence interval is\", ci)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the given population mean of 98.6 is not found within the 95% confidence interval, the null hypothesis is rejected.\n",
    "\n",
    "A T-Test could be arguably more relevant, as it doesn't assume the population standard deviation is the same as the sample's. The same hypothesis and confidence requirement will be used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Ttest_1sampResult(statistic=-5.4548232923645195, pvalue=2.4106320415561276e-07)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stats.ttest_1samp(df[\"temperature\"], 98.6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "t_test_statistic: -5.45482329236\n",
      "pvalue: 2.41063204156e-07\n",
      "test_stat2: -1.97852449146\n",
      "test_stat2: 1.97852449146\n",
      "pvalue2: 2.41063204156e-07\n"
     ]
    }
   ],
   "source": [
    "# this t value shows how much the sample differs from the null hypothesis\n",
    "t_test_statistic = stats.ttest_1samp(df[\"temperature\"], 98.6).statistic\n",
    "pvalue           = stats.ttest_1samp(df[\"temperature\"], 98.6).pvalue\n",
    "\n",
    "print(\"t_test_statistic:\", t_test_statistic)\n",
    "print(\"pvalue:\", pvalue)\n",
    "\n",
    "# next, check the quantiles of the t-distribution of the with 95% confidence\n",
    "test_stat1 = stats.t.ppf(.025, n-1)\n",
    "print(\"test_stat2:\", test_stat1)\n",
    "\n",
    "test_stat2 = stats.t.ppf(.975, n-1)\n",
    "print(\"test_stat2:\", test_stat2)\n",
    "\n",
    "# now check the probability of the sample's occurrance given that the null hypothesis is true\n",
    "pvalue2 = stats.t.cdf(t_test_statistic, n-1) * 2\n",
    "\n",
    "print(\"pvalue2:\", pvalue2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A p-value of 2.41e-07 means that observing data as extreme as the sample would happen about 0.0000241% of the time if the null hypothesis were true. In this case, the p-value is much lower than the significance level of 0.05, so the null hypothesis is rejected, as intuitively expected.\n",
    "\n",
    "A T-distribution is similar to a normal distribution, but varies by sample size. As sample size increases, the curve narrows, approaching the normal curve. With a sample of 130 observations, the curve is still wider, producing a larger critical value and wider confidence interval, shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(98.12249318821641, 98.375968350245145)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stats.t.interval(alpha=.95, df=n-1, loc=mean, scale=std / np.sqrt(n))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Draw a small sample of size 10 from the data and repeat both tests.\n",
    "\n",
    "Small Sample"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 98.9  97.9  97.7  97.9  97.9  97.6  99.   98.8  97.8  99. ]\n"
     ]
    }
   ],
   "source": [
    "# drawing a small sample of size 10\n",
    "small_n = 10\n",
    "small = np.random.choice(df[\"temperature\"], small_n)\n",
    "small_mean = np.mean(small)\n",
    "small_std = np.std(small)\n",
    "print(small)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Z-Test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The 95% confidence interval is (97.902417089155946, 98.597582910844054)\n"
     ]
    }
   ],
   "source": [
    "# z_critical value for normal distribution given a two-tailed error\n",
    "# the z-critical value for 95% confidence is commonly represented as 1.96 as per textbook lookup tables\n",
    "small_z_critical = stats.norm.ppf(q = 0.975)\n",
    "\n",
    "# sigma represents the population standard deviation for which we use the sample standard deviation\n",
    "small_sigma = small_std\n",
    "\n",
    "# margin of error for the sample mean is computed as such\n",
    "small_moe = small_z_critical * small_sigma / np.sqrt(small_n)\n",
    "\n",
    "# show the 95% confidence interval\n",
    "small_ci = (small_mean - small_moe, small_mean + small_moe)\n",
    "\n",
    "print(\"The 95% confidence interval is\", small_ci)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "T-Test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "t_test_statistic: -5.45482329236\n",
      "pvalue: 2.41063204156e-07\n",
      "test_stat2: -2.26215716274\n",
      "test_stat2: 2.26215716274\n",
      "pvalue2: 0.000403129999163\n"
     ]
    }
   ],
   "source": [
    "# this t value shows how much the sample differs from the null hypothesis\n",
    "small_t_test_statistic = stats.ttest_1samp(df[\"temperature\"], 98.6).statistic\n",
    "small_pvalue           = stats.ttest_1samp(df[\"temperature\"], 98.6).pvalue\n",
    "\n",
    "print(\"t_test_statistic:\", small_t_test_statistic)\n",
    "print(\"pvalue:\", small_pvalue)\n",
    "\n",
    "# next, check the quantiles of the t-distribution of the with 95% confidence\n",
    "small_test_stat1 = stats.t.ppf(.025, small_n-1)\n",
    "print(\"test_stat2:\", small_test_stat1)\n",
    "\n",
    "small_test_stat2 = stats.t.ppf(.975, small_n-1)\n",
    "print(\"test_stat2:\", small_test_stat2)\n",
    "\n",
    "# now check the probability of the sample's occurrance given that the null hypothesis is true\n",
    "small_pvalue2 = stats.t.cdf(small_t_test_statistic, small_n-1) * 2\n",
    "\n",
    "print(\"pvalue2:\", small_pvalue2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(97.84882570413825, 98.65117429586175)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "stats.t.interval(alpha=.95, df=small_n-1, loc=small_mean, scale=small_std / np.sqrt(small_n))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that the sample size is reduced to below the magical, statistical threshold of 30, the T-Test should be more appropriate. Both tests now produce increased margins of error and often include the given population mean of 98.6 degrees, preventing the consistent rejection of the null hypothesis at 95% confidence in the sample. The confidence intervals above for Z and T tests reject and retain the null hypothesis respectively.\n",
    "\n",
    "This illustrates the importance of repeating experiments and the value of drawing more samples."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. At what temperature should we consider someone's temperature to be \"abnormal\"?\n",
    "\n",
    "As of 2018, it has become difficult to label someone as abnormal. In this exercise, a prediction interval is used to make a defensible attempt at doing so."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(98.123682111456645, 98.37477942700491)\n",
      "The 95% confidence interval would suggest 115 of 130 people are abnormal\n"
     ]
    }
   ],
   "source": [
    "print(ci)\n",
    "low, high = ci\n",
    "abnormals = (df[\"temperature\"] < low).sum() + (df[\"temperature\"] > high).sum()\n",
    "print(\"The 95% confidence interval would suggest\", str(abnormals), \"of 130 people are abnormal\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Disregarding the assignment's suggestion to use a confidence interval, a prediction interval is used to find the range that encompasses 95% of observations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(96.798652852736836, 99.699808685724719)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# prediction interval with t-distribution\n",
    "t_critical = stats.t.ppf(0.975, n-1)\n",
    "mean - t_critical * np.sqrt(1 + 1/n) * std, mean + t_critical * np.sqrt(1 + 1/n) * std"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(96.81775582674824, 99.680705711713316)"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# prediction interval with normal distribution\n",
    "stats.norm.interval(0.95, loc=mean, scale=std)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(97.04499999999999, 99.3)"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# good old fashioned quantile definition\n",
    "df[\"temperature\"].quantile(.05), df[\"temperature\"].quantile(.95)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Is there a significant difference between males and females in normal temperature?\n",
    "\n",
    " - What test did you use and why?\n",
    " - Write a story with your conclusion in the context of the original problem. \n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Null Hypothesis: The true population mean is equal to the conventional mean 98.6 F\n",
    "\n",
    "Alternative hypothesis: The true population mean is not equal to the conventional mean 98.6 F"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style>\n",
       "    .dataframe thead tr:only-child th {\n",
       "        text-align: right;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: left;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>temperature</th>\n",
       "      <th>gender</th>\n",
       "      <th>heart_rate</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>99.3</td>\n",
       "      <td>F</td>\n",
       "      <td>68.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>98.4</td>\n",
       "      <td>F</td>\n",
       "      <td>81.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>97.8</td>\n",
       "      <td>M</td>\n",
       "      <td>73.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>99.2</td>\n",
       "      <td>F</td>\n",
       "      <td>66.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>98.0</td>\n",
       "      <td>F</td>\n",
       "      <td>73.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   temperature gender  heart_rate\n",
       "0         99.3      F        68.0\n",
       "1         98.4      F        81.0\n",
       "2         97.8      M        73.0\n",
       "3         99.2      F        66.0\n",
       "4         98.0      F        73.0"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "mtemp = df[\"temperature\"][df[\"gender\"] == \"M\"]\n",
    "ftemp = df[\"temperature\"][df[\"gender\"] == \"F\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "count    65.000000\n",
       "mean     98.104615\n",
       "std       0.698756\n",
       "min      96.300000\n",
       "25%      97.600000\n",
       "50%      98.100000\n",
       "75%      98.600000\n",
       "max      99.500000\n",
       "Name: temperature, dtype: float64"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mtemp.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "count     65.000000\n",
       "mean      98.393846\n",
       "std        0.743488\n",
       "min       96.400000\n",
       "25%       98.000000\n",
       "50%       98.400000\n",
       "75%       98.800000\n",
       "max      100.800000\n",
       "Name: temperature, dtype: float64"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ftemp.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "xm, ym = ecdf(mtemp)\n",
    "xf, yf = ecdf(ftemp)\n",
    "x_theor, y_theor = ecdf(samples)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAEFCAYAAADnpECTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3XtYVOXeN/DvMCMozmSoQzvdgYhA\nWlsRLDyE+yIkO+izEw9gSe60N9uFVrI76NNWX980Ommp6d5ZmbnbSQfzst3BxHxkq0mPGJqpRGMh\nWQoWKgsUGOZ+/0AGFswMwzBrzen7ua6ua5jVrHVzO8Nv7t+679+tEUIIEBERXRbk6QYQEZF3YWAg\nIiIZBgYiIpJhYCAiIhkGBiIiktF5ugHOqqys9nQTPCYsLBRVVbWeboZXYF/IsT/k2B8tmvvCaDR0\n+rUcMfgAnU7r6SZ4DfaFHPtDjv3Roit9wcBAREQyDAxERCTDwEBERDIMDEREJMPAQEREMgwMREQk\no2hgOHToELKysto9/8UXX2Dy5MnIyMjAu+++q2QTiEhhkgQUFQVBkjzdEj8lAbqiIEDF/lVsgdv6\n9euxbds29OjRQ/Z8Q0MDnnnmGbz//vvo0aMHpk+fjpSUFBiNRqWaQkQKkSRg/PhQlJZqERPTiO3b\na6HXe7pVfkQCwsaHQleqhTmmEVXbawEV+lexwBAREYHVq1fj8ccflz1vMpkQERGBXr16AQASExNx\n4MAB3HbbbQ7PFxYWGtCLV1xZveiv2BdynuyPEyeA0tKmx6WlWlRUGBAV5bHmAPCz98cJAJf7V1eq\nhbHCAHSif13tC8UCw/jx4/HTTz+1e16SJBgMLY3t2bMnJCfGoIG8zN1oNAR0SZDW2BdyavWHJAEl\nJUGIi7PIRgTh4UBMTMuIITy8FpWVijfHTiMBY4UBleHVqnyrVkU4EBbTasQQXgs42b/N7w1XgoPq\ntZL0ej1qamqsP9fU1MgCBRF5F0fpIr0e2L691mbQULeRTSkXlDb9IVUr5aI4PVC1vRa6kiCY4yyq\n/U6qz0qKjo5GWVkZzp07h/r6ehw4cADDhw9XuxlE5KSSkiCUljalcUtLtSgpkf/Z0OuBxEQPBgUA\nupIg6C63UVeqha7EjyZc6gFzonpBAVBxxPDRRx+htrYWGRkZePLJJzF79mwIITB58mRcddVVajWD\niOywly6Ki7MgJqbROmKIi7N4sJGw+e3ZHGeBOabRmnIxe7KNfkAjhBCeboQzAjmvzLx6C/aFnLv6\no6PZRfaChqo6mqHjj/cYuqAr9xj8aLxFRK7yi3SRHkASGBTcgIGBKIDYW4zWnC4CYDtd5IFFVm01\np4sAMF2kMJ/ZwY2Iusbl2UUeWmTVjodm6AQijhiIAoSr6SKvmvHjgRk6gYiBgShAdJgusoMpnMDD\nVBKRH7I1i6jDxWh2poIyhRN4GBiI/ExH9xISE2184+/oPkJzCocCAlNJRH6mo3sJtnjVfQTyOP7r\nE/koSQIKC9G5qad2pp3yPgK1xlQSkQ9qSRc1VTd1auqpo3QR7yNQKxwxEPkgV6aeOrNymFNBCWBg\nIPJJcXEWREc3pX6io52bemqOs8B8+TXmaKaLyD4GBiIikmFgIPJBJSVBMJma0kImUydmHl1+jc7E\nmUdkH98ZRD6ow1SSjdlHnHlEzuKsJCJ/Y2/2EWcekZM4YiDyQY5SSQ5nH3HmETmBgYHIi7myfwJn\nH1FXMZVE5KWc2T+hosKA8PBaj+6sRv6HIwYiL+XMIrakJNjeP4Gzj6gL+I4h8rCubLeJQrDuEbkd\nU0lEHtTV7TZRCoTFhLLuEbkVRwxEHqTYdpucfURdwMBA5EGu1DwCOPOIlMXAQEREMgwMRB7kSs0j\ngDOPSFl8NxF5kCs1jwDOPCJlcVYSkbdyYsc1Y4UBVeG1vMlMbsURA5EHuVzzCGgKBklgUCC3Y2Ag\nUomthWxxcRYMi27EjQCGRduoecR0EXkAU0lEKrC3kE0P4CsAwQDqAZxv/SIuVCMP4YiBSAX2FrLp\nSoIQfDmVFGxrdhEXqpEHKBYYLBYLFi1ahIyMDGRlZaGsrEx2/PXXX0d6ejomT56MHTt2KNUMIq9g\nb/YRF6qRN1IslZSfn4/6+nrk5eWhuLgYubm5WLduHQDgwoUL2LRpEz7//HNcvHgRd955J9LS0pRq\nChERdYJigaGoqAjJyckAgPj4eBw5csR6rEePHujXrx8uXryIixcvQqPRdHi+sLBQ6HRapZrr9YxG\ng6eb4DV8sS9OnABMpqbHJpMWFRUGREUBOAHg8vM6kxbGCgMQ1blz+2J/KIn90cLVvlAsMEiSBH2r\nyl9arRZmsxk6XdMlr776atxxxx1obGzEnDlzOjxfVVWtUk31ekajAZWV1Z5uhlfw9r6QJNishhoe\nDsTEtNx8Dg+vRWUlgHAgLDoUOpMW5ujGpjUJlc5fz9v7Q23sjxbNfeFKcFAsMOj1etTU1Fh/tlgs\n1qBQUFCAiooK7Ny5EwAwe/ZsJCQkYOjQoUo1h0hxLpfQJvIyit18TkhIQEFBAQCguLgYsbGx1mO9\nevVC9+7dERwcjJCQEBgMBly4cEGpphCpwpUS2qx5RN5IsRFDWloa9u7di8zMTAghsHz5cmzYsAER\nERFITU3Fvn37MG3aNAQFBSEhIQFjxoxRqilEqmjeca15xOBMCe3mRWzNZS84K4m8gUYIITzdCGcE\nct6QedMW3t4XZ84A+fk6jBtnxlVXOfsiICRfh7pxZsDZ11zm7f2hNvZHC6+8x0AUaCQJSE+3fY/B\n/ouAsHQ7hfKIPIQJTSI36egegy0dFsoj8gC+C4ncxJVtOrnymbwRAwMREckwMBC5iSvbdHK6Knkj\nvguJ3ISpJPIXDAxERCTDwEDkJkwlkb/gu5DITZhKIn/BwEBERDIMDERuwlQS+Qu+C4ncpMNUkgTo\nioIAqeWp5iJ6AFhEj7wGayURqUECwsbbqImkB6q210JXEtQUFFgnibwARwxEbuIoleSwJpIeMCcy\nKJD3YGAg6iRJAoqKgiBJ8ueb92MA0G4/Bs4+Il/CVBJRJ3D7TgoEHDEQdYIr23cCnH1EvoXvTqJO\ncGURG8BUEvkWBgYiIpJhYCDqBFcWsQFMJZFv4buTyAZXZh45wlQS+RLOSiJqgzOPKNBxxEDUhqsz\njxxhKol8Cd+dRG24OvPIEaaSyJcwMBARkQwDA1Ebrs48coSpJPIlfHcStcFUEgU6BgYiIpJhYCBq\ng6kkCnR8dxK10aVUko1d2gDu1Ea+hQvciNzF3i5tAHdqI5/CEQNRG12qh2RvlzaAO7WRz1AsMFgs\nFixatAgZGRnIyspCWVmZ7Pju3bsxbdo0TJs2DUuWLIEQQqmmEHVKl+ohMV1EfkCxVFJ+fj7q6+uR\nl5eH4uJi5ObmYt26dQAASZLw/PPP46233kLv3r2xfv16VFVVoXfv3ko1h8gmSUK7ukcu10PSA1Vb\nahGSr0PdODNHBuSzFAsMRUVFSE5OBgDEx8fjyJEj1mNff/01YmNj8eyzz6K8vBxTp05lUCDVdVQs\nLzGxk9/4JSAs3c49BiIfolhgkCQJ+lZftbRaLcxmM3Q6HaqqqlBYWIitW7ciNDQUd999N+Lj4xEV\nFWX3fGFhodDptEo11+sZjQZPN8FruKsvTpwASkubHpeWalFRYYCDt6ATJwRw+Xy6Ui2MFQagK+dz\nEt8bcuyPFq72hWKBQa/Xo6amxvqzxWKBTtd0uSuvvBJ/+MMfYDQaAQAjRozAsWPHHAaGqqpapZrq\n9YxGAyorqz3dDK/gzr4IDweio0NhMmkRHd2I8PBaVFZ25YRAWHQodCYtzNGNqAqvBbpyPifwvSHH\n/mjR3BeuBAfFbj4nJCSgoKAAAFBcXIzY2Fjrseuvvx7fffcdfvvtN5jNZhw6dAiDBg1SqilERNQJ\nio0Y0tLSsHfvXmRmZkIIgeXLl2PDhg2IiIhAamoqcnJycN999wEAbr31VlngIFKDrWmpnb6v0Iqt\n1c3mLpyPyFM0wkfmiQby8JDD4xau9IWtmUfNz9u7+ewSCQhLa5VK2qH8zWe+N+TYHy26kkriymfy\na9ymk6jzuPKZ/JoS23Taw0J55C/4ziW/psTeCiyUR/6OqSSizmChPAoADkcMZ86cUasdRIpw994K\nLJRHgcDhp+SBBx6wPn7jjTcUbwyRu7mcSmK6iAKYw8DQeibrRx99pHhjiLzC5XRR2G09ETY+VB4c\nLqeLqj6tYS0k8lsOA4NGo7E+9pHlDkQyrqSSmC6iQOd0wrV1kCDyFa6kksxxFpgvv8YczXQRBR6H\ns5JKS0uRmpoKoOlGdPNjIQQ0Gg127typfAuJiEhVDgPD9u3b1WoHkSJcqYfEmkcU6Bymkvr374/+\n/fujpqYG33zzDUpLS2GxWKzPE3kLSQKKioIgtZlF5Mo2nUwlUaBzOGL49ddfMW/ePJSWliIyMhIa\njQY//PAD4uPjsWLFChgM3BCDPI/1kIjcy+GI4cUXX0RiYiL27t2L9957D++++y727t2La6+9FsuW\nLVOrjUQOubseEmseUaBz+I7/+uuvMX/+fHTr1s36XHBwMObPn4+jR48q3jgiZ7i7HhJTSRToHAaG\nkJAQm89rNBoEBfFbFBGRP3J6gVtnjhGpSZF6SEwlUQBzeh1Da0IIVHZp13SizpMk4MQJIDwcsvsF\nzTOPmm8+uyuV1LwTG1NJFGg6XMdw7tw5NDY2ok+fPgCAwsJCDBo0yPozkRpaZh4BMTGhnHlEpCCH\nY+Tz589jzpw5OHfunHXtQllZGR566CFUV3NfVVIPd2IjUo/Dd/yzzz6LF198EWPHjrU+9+ijj2L5\n8uXIzc1VvHFEzRTZiQ2wWV6bpbUp0DkMDBcuXEBSUlK755OTk1FVVaVYo4hUYa+8NktrU4BzGBjM\nZjMslvbfliwWCxoaGhRrFFFb7p55BHRQXpultSmAOfx03XDDDVizZk2759euXYvrr79esUYRtaVE\nKokL2Yhsczgraf78+bj//vuxdetWXHvttQgJCcHRo0fRu3dvrFu3Tq02EhGRihwGBr1ej7fffhv7\n9+/HsWPHEBQUhLvvvhsjRoxQq31EAFwrn90Rltcmss1hYACaVjiPGjUKo0aNUqM9RDY1p5JMJq3b\nU0lcyEYkxwnaREQkw8BAPkGxWUlcyEbUDj8J5BNc2YmtI1zIRmRbh/cYiLyBXg9s2VKLwkIDkpJq\nO1f6QmoaHZjj2qxLuLyQzeYxogDGwEA+QZKA9HTbRfQcv7BpdbOuVAtzTGP7lczNC9mIyEqxVJLF\nYsGiRYuQkZGBrKwslJWV2fx/7rvvPrzzzjtKNYP8REdF9OxxuLqZiGxS7FOSn5+P+vp65OXlIScn\nx2bRvZdeegnnz59XqgnkR1xd+czVzUSdp1hgKCoqQnJyMgAgPj4eR44ckR3/7LPPoNFoZJVbiYjI\n8xS7xyBJEvStksBarRZmsxk6nQ7fffcd/v3vf2PVqlV45ZVXnDpfWFgodDqtUs31ekajwdNN8KgT\nJwCTqemxyaRFRYUBUVHOvBDA5dfpTFoYKwyAM6/zIYH+3miL/dHC1b5QLDDo9XrU1NRYf7ZYLNDp\nmi63detWnDlzBjNnzsSpU6fQrVs39O/f3+HooaqqVqmmej2j0YDKysDZGEmS0G43tvDwppvOzdt3\nhofXwqndZcOBsOhQ6+rmqvBawI92pQ2090ZH2B8tmvvCleCgWGBISEjArl27cPvtt6O4uBixsbHW\nY48//rj18erVq9G3b1+mlAhA6y08mwJA8+yj5u07KyoMCA/v5HRVIuoUxe4xpKWlITg4GJmZmXjm\nmWewYMECbNiwATt37lTqkuQHHM0+0uuBpCR0KihwdTNR5yk2YggKCsLSpUtlz0VHR7f7/+bOnatU\nE8gHubtYHgvlEXUevz4REZEMAwN5FXcXy2Mqiajz+Ckhj5AkoKgoCJIkf97dxfK4wI2o81griVRn\nb+YR0DL7qO10VSJSD0cMpLqO6h7p9UBionuCAlNJRJ3HTwmpztW6R65gKomo8xgYiIhIhoGBVKfE\nNp32MJVE1Hn8lJDqFEklSYCuKAhoM8uJ23cSdR5nJZHvc7RLG7fvJOo0jhhIdYosYnO0S1vz9p0M\nCkROYWAg1bmcSpIAFILpIiKFMZVEvuFyugilQFhMKNNFRAriiIFU50oqiekiIvUwMJBi3FkPiQvV\niNTDVBIpgvWQiHwXRwykCHfXQ+JCNSL18NNFXcLy2UT+h6kkchnTRUT+iSMGchnLZxP5J366yGUs\nn03knxgYiIhIhoGBXMby2UT+iZ8ucpm7Zx45wnpIROrhrCRymV4PbNlSi/x8HcaNMys780gPVG2p\nhbHQgKqkWpa+IFIQAwO5TJKA9HTb01XdfzEgLN1OET0iciumkshlHU1XdacOi+gRkdvw00UdUmt1\nc8sF22/TyXsMROphKokcUn11s71tOi/vuWCsMKAqnGkkIiVxxEAOqbm6GeggZaQHkAQGBSKFMTCQ\nQ2qubga4wpnIGzAwEBGRjGKBwWKxYNGiRcjIyEBWVhbKyspkx998801MnToVU6dOxZo1a5RqBnWR\nmqubAa5wJvIGin3q8vPzUV9fj7y8POTk5CA3N9d6rLy8HNu2bcPmzZuRl5eHPXv24Pjx40o1hZxk\na/aRmqubAaaSiLyBYrOSioqKkJycDACIj4/HkSNHrMd+97vf4bXXXoNW2/TN0Gw2IyQkRKmmkBPs\nzT7ivgpEgUexwCBJEvSt/opotVqYzWbodDp069YNvXv3hhACzz33HIYMGYKoqCiH5wsLC4VOp1Wq\nuV7PaDQoev4TJ4DS0qbHpaVaVFQY0PxPYjQCHfzzuLEhAExND3UmLYwVBqDNtZXuC1/D/pBjf7Rw\ntS8UCwx6vR41NTXWny0WC3S6lsvV1dVh4cKF6NmzJxYvXtzh+aqqahVppy8wGg2orKxW9Brh4UB0\ndChMJi2ioxsRHl6LykoFLyg13U8wx1nk00/Dm0peWNcxhNcCrdqhRl/4EvaHHPujRXNfuBIcFLvH\nkJCQgIKCAgBAcXExYmNjrceEEHjwwQcRFxeHpUuXWlNKFCAuL2ILu60nwsaHylY4Ny9kq/q0hvWQ\niDxEsRFDWloa9u7di8zMTAghsHz5cmzYsAERERGwWCz46quvUF9fj//85z8AgPnz52P48OFKNYc6\nYGv2UWKiMjd+bS1iM7e+lh7yn4lIVYoFhqCgICxdulT2XHR0tPXxN998o9SlyQXNC9maU0lqbNOp\nM2k584jIC3GSOBERyTAwEABu00lELfiJJADq1kTiIjYi78bAQEREMgwMBICpJCJqwU9kgFF1NzYb\nO7EB3I2NyNtxB7cAoupubPZ2YgOsi9hsrnwmIo/jiCGAqLkbm8Od2ICWRWwMCkReh4HBD6maLrKD\nM4+IfBdTSX5G1XQREfkljhj8jJrpIkc484jId/HT6qO8IV3U1BDOPCLfYe9z44qDBw/gpptGYOfO\nz2XPz5yZiWXLlth8zSeffIR161Z3/eIKYyrJB3lNuogzj8iHOPrcuCoycgDy87cjNfUWAIDJ9D0u\nXrzohtZ6FgODD7KVLmpdIrs5XaQ0ls8mX9LR58YVgwbFoLz8JKqrq2EwGLB9+ye45ZbbcObMaXzw\nQR52794Fs9kMvV6PZcuel732/fc3Y8eO7dBoNEhNvQVTp2Z2qS3uxFSSD1KsrpGdtJA9TBeRL1Eq\nzTp2bAoKCnZBCIFjx77F9dcPhcViwfnz5/HSS2uxdu1rMJvNOHbsW+trfvjhBHbu3IG1a1/D2rWv\n4T//+R+cPPmjW9rjDhwxUBNHaSF7mC4iH6JUmjUt7Va8+GIu+vXrj2HDmjYbCwoKQrdu3bBkyX+j\nR48eqKiogNlstr7mxAkTzpw5jYcf/gsAoLq6Gj/99BMiIga4p1FdxMDgg5TYba3DtJA9TBeRD1Ei\nzdq//+9x8eJFvP/+ZsyZk42ffz6FmpoaFBT8D9av34hLly5h9uwZstdERERiwICBePHFVdBoNMjL\nexsDBw5ya7u6goHBBzUPiZtvorljSNycFmoeMTAtROS81NQ0bN/+CSIiIvHzz6eg1WrRo0cPzJ6d\nheDgbujTpy/Onq20/v8xMbEYMeIGPPjgbNTXN2Dw4OtgNBo9+BvIaYQQwtONcEZlZbWnm+AxRqOh\n3e9/5gyQn6/DuHFmXHVVmxdIcC294+rrVGSrLwIZ+0OO/dGiuS+MRkOnX8sRgw+SJCA93c60O1fu\nFTRjWoiIwFlJPsnR6uYOi9cREXWAfzW8nCQBhYWQrdR0NO2OxeuIqKuYSvJiLSs1gZiYUGvKiMXw\niEhJHDF4MUcpI3vF8Fi8joi6in811OJgVbE7C+JxNTIRdRVTSWpwMFPImYJ4FRUGhIc7WfCLq5GJ\n7HPjlOxffvkZM2dOR2xsnPW5xMQbcO+9/6eLjWyRnX0/HntsISIjB7jtnM5gYFCBo1XFzhTEi4oC\nKivbn9cuTjslaq8rU7ntGDAgCmvWvOqe9nkRBgY3qjkD/JKvw9XjzOjZatGZo1XFSqxiJqL2XC77\n0kl///saHDp0EBaLQEbG3bj55nHIzr4fgwbF4ocfTOjRoweGDh2Or776EpIkYcWKNdBqg5Cb+zQk\nqRrnz5/DxImTMGnSFOs5JUlCbu5SnD9/HgDwyCOPITpauRIaDAxuUnMGaEzQY1SDBqZuAjUHpZbg\n4CC9wxlGROpQouzLjz/+gOzs+60//9d/TcIvv5zCunVvoK6uDnPm3IsbbkgCAAwZch0eeeSvmD9/\nLrp3746XXlqLp59ejOLig7jqqt9h3Lhb8Mc/3oyzZyuRnX2/LDC89dYbSEy8EZMmTUF5+UksX/5/\nsW7d611uvz0MDG7yS74Ooxo0AIDoBg2+zNdh0N0t1RQdpXfU2j+BKKApcP+tbSrp7bc3oqTkuDVY\nmM1mnD79CwAgNvZaAIDBoMeAAVGXH1+B+vo69OnTB++++y/s3r0LoaE9ZZVYAeDEie9x8OAB625x\n1dXKlv3w+8AgSbD/bdyNN6KuHmeGqZtA9OURw9XjzB2/iIjUpfD9t8jIARg+fASeeOK/YbFY8Oab\nr6F///4AAI1GY/d177yzCddfPxSTJk3BwYMH8OWXe9qd95ZbhuCWW25FVdVv+OijrYr9DoCfBwaH\nW/m5+UZUz6uAmoMSvrRxj4GIAsOYMWPx9ddFePDB+3DxYi3Gjk1BaGhPp173wgvP4PPPP0WvXr2g\n1WpRX19vPX7PPbOQm/v/sG3bFtTW1mDWrPsdnK3rFKuuarFYsGTJEpSUlCA4OBhPP/00IiMjrcff\nffddbN68GTqdDn/5y1+QkpLi8HyuVEwsKgrCbbe1/KN8+mmNNWWjKwpCWKtjVZ/WeO1MHlaMbMG+\nkGN/yLE/WnSluqpiC9zy8/NRX1+PvLw85OTkIDc313qssrISmzZtwubNm/H6669jxYoVsujoLh3W\nFOJCMCKidhRLJRUVFSE5ORkAEB8fjyNHjliPHT58GMOHD0dwcDCCg4MRERGB48ePY+jQoW5tg8MZ\nP1wIRkRkk2KBQZIk6Fv9JdZqtTCbzdDpdJAkCQZDy/CmZ8+ekNrWg2gjLCwUOp220+0wGpsWiNk+\nCMDeMS/jynDQX7Ev5NgfcuyPFq72hWKBQa/Xo6amxvqzxWKBTqezeaympkYWKGypqqpVpqE+gHnT\nFuwLOfaHHPujhVfeY0hISEBBQQEAoLi4GLGxsdZjQ4cORVFREerq6lBdXQ2TySQ7TkREnqPYiCEt\nLQ179+5FZmYmhBBYvnw5NmzYgIiICKSmpiIrKwt33XUXhBB49NFHERISolRTiIioExSbrupugTw8\n5PC4BftCjv0hx/5o4ZWpJCIi8k0MDEREJOMzqSQiIlIHRwxERCTDwEBERDIMDEREJMPAQEREMgwM\nREQkw8BAREQyDAxERCTDwOBl6uvrkZOTg2nTpmHWrFn48ccfkZWVZf1vzJgxeOGFFzzdTFXY6ot9\n+/YhPT0d06ZNw8qVKz3dRFXZ6o89e/bgzjvvxPTp07F27VpPN1EVhw4dQlZWFgCgrKwM06dPx113\n3YXFixfDYmnacGvNmjWYMmUKMjMzcfjwYU82V3HO9EfzsQkTJjh3UkFeZdOmTeKpp54SQghhMpnE\nrFmzrMdOnjwpJk2aJCRJ8lTzVGWrL/70pz+J0tJSYbFYRGZmpjh+/LiHW6metv3x5z//Wfzxj38U\nJ0+eFEIIkZOTI/73f//Xk01U3KuvviomTJggpk6dKoQQYs6cOWL//v1CCCH+9re/ic8//1wcOXJE\nZGVlCYvFIk6dOiXS09M92WRFOdMfQgjx4YcfikmTJonRo0c7dV6OGLzM999/j7FjxwIABg4cCJPJ\nZD22bNkyPPbYY+jZs+PNxf2Brb4YPHgwzp07h4aGBtTV1UGr7fzmTb6qbX8cPHgQV1xxBa655hoA\nTaXuDx486MkmKi4iIgKrV6+2/vztt9/ixhtvBACMHTsW+/btQ1FREW666SZoNBr069cPjY2N+O23\n3zzVZEU50x8A0KtXL/zzn/90+rwMDF5m8ODB2LVrF4QQKC4uxpkzZ9DY2Ijjx4+jpqYGo0aN8nQT\nVWOrL2JiYvDAAw/g9ttvx9VXX42BAwd6upmqadsf9fX1uHTpEkwmExobG1FQUIDaWv/e0Gr8+PHW\nDb8AQAgBjUYDoGknyOrq6na7RzY/74+c6Q8ASElJQWhoqNPnZWDwMpMnT4Zer8c999yDXbt24brr\nroNWq8W2bdswdepUTzdPVW37IiIiAuvXr8fHH3+M/Px8REZG4o033vB0M1Vj673x3HPPYcmSJZg3\nbx6ioqIQFhbm6WaqKiio5U9YTU0NrrjiCpd2iPQXtvrDpfO4q0HkHt988w0SExOxadMmjBs3zpom\n2L9/P5KTkz3cOnW17YtBgwYhNDTU+s0nPDwcFy5c8HAr1WPrvVFQUIB//OMfWLNmDU6ePInRo0d7\nupmqGjJkCAoLCwEABQUFGDFiBBISErBnzx5YLBb8/PPPsFgs6N27t4dbqg5b/eEKxXZwI9dERkbi\n5ZdfxhtvvAGDwYBly5YBACoUmfShAAAHNklEQVQrKwPu26Ctvjh8+DBmzZqFkJAQGAwG5ObmerqZ\nqrHVH7t378b06dPRvXt3TJw4ETExMZ5upqqeeOIJ/O1vf8OKFSswcOBAjB8/HlqtFiNGjEBGRgYs\nFgsWLVrk6WaqxlZ/uIJlt4mISIapJCIikmFgICIiGQYGIiKSYWAgIiIZBgYiIpJhYCCPKSwstBb/\nUvr11dXVeOihh1y+lhp27dqFDRs2KHLu06dPY8GCBQCA3bt3Izk5GTk5OcjNzcXRo0cVuSb5Lq5j\noIBw/vx5HDt2zNPNcOjIkSOKnXv58uV4+OGHAQCfffYZsrOzkZGRgd9++w3z5s3rVB0d8n8MDORR\nVVVVmD17NioqKjB06FAsXrwYwcHB2LVrF1566SVYLBZcc801WLp0Kfr27Ys9e/bgmWeeQUhICKKi\nogA0lROeOXMmvvjiCwQFBaGwsBDr16/Ha6+9Zr3O008/jYqKCjz00EN45ZVXsHXrVmzcuBEWiwXX\nXXcdFi9ejJCQEIwZMwapqak4fPgw+vbti8mTJ2PTpk04ffo0cnNzceONNyIrKwvXXnstDhw4gLq6\nOixcuBA33XQTzp49i0WLFuH06dPQaDTIycnB6NGjsXr1ahQXF+OXX37BjBkzMGjQIKxcuRKXLl3C\nhQsXsGDBAgwYMACbN28GAPTr1w8///wzAGDu3LkAgJtvvhlvvfUWvvrqK3z44Yc4d+4cUlJScM89\n99i8ZmsnT55ERUUFoqOj8d5772Hnzp348ssvERQUhKlTp6J3797Yv38/Ro4cqcY/OfkCN1eBJXLa\n/v37xbBhw8QPP/wgLBaLePjhh8Wbb74pzp49K2666SZRXl4uhBBi/fr1Yu7cuaKurk6MGTNGfP/9\n90IIIRYuXChmzJghhBDi7rvvFvv27RNCCPHkk0+Kjz/+WHat8vJykZKSIoQQ4rvvvhPTp08Xly5d\nEkII8cILL4hXXnlFCCFEbGys2LFjhxBCiBkzZoj58+cLIYTYsmWLePDBB63PP/nkk0IIIY4ePSrG\njBkj6urqxCOPPCLy8/OFEEKcOXNGpKamiurqarFq1SprO4UQYu7cudbfYd++fWLChAlCCCFWrVol\nVq1a1e6xEEKkpKSI8vJy8cEHH4i0tDTR0NAghBB2r9nahg0bxLPPPmv9+YknnhAffPCB9eeNGzeK\nZcuWOfy3osDCEQN51IgRIzBgwAAAwMSJE7FlyxZERERg6NCh+P3vfw8AyMjIwKuvvoqSkhKEh4cj\nOjoaADBp0iS8/PLLAJoKzG3btg3x8fHYv38/lixZYveahYWFKCsrw7Rp0wAADQ0NGDJkiPV4c2nr\n/v37IzExEUDTt/jWdZmaXzt48GAYjUaUlJRg3759OHHiBFatWgUAMJvNKC8vBwAMHTrU+trnn38e\nu3btwmeffYZDhw7JCr45Y8iQIdaKmvauOXjwYOv/X1ZWZh1d2dKvXz/s3bu3U20g/8bAQB7VtmSw\nTqeT7TrV/LzZbIZGo4FoVcGl9V4Mt956K1auXInt27dj7NixCAkJsXvNxsZG3HbbbXjqqacANFWh\nbGxstB4PDg62eY3WWj9vsVis7d64cSOuvPJKAEBFRQX69OmD/Px8dO/e3fr/33XXXUhKSkJSUhJG\njRqFv/71r+3Or9FoZP3Q0NBgfdz6XPau2fZcrfu5LZ1OZy3VTARwVhJ5WFFRkbUC5tatWzF69GgM\nGzYMhw4dwk8//QQAyMvLQ1JSEuLi4nD27FkcP34cAPDxxx9bz9OjRw+MHTsWK1asQHp6ervr6HQ6\nmM1mAEBSUhJ27NiBX3/9FUIILFmyBBs3buxUuz/55BMATRVPL1y4gNjYWIwcORL/+te/ADRtqjNx\n4kRcvHhR9rpz587hxx9/xMMPP4yxY8di586d1qCk1WqtbQwLC8P3338PADh8+DAqKytttsOZa0ZE\nRODUqVN2f5dTp04hMjKyU78/+TeOGMijBg0ahIULF6KyshIjR47ElClToNVqsXTpUmRnZ6OhoQH9\n+vXDsmXL0K1bN6xYsQKPPfYYdDqdLP0DAHfccQcOHjyIYcOGtbtOnz590K9fP2RlZWHTpk3Izs7G\nzJkzYbFYMHjwYNx///2dand5eTkmTZoEAFi5ciW0Wi2eeuopLFq0CBMnTgQAPPfcc7INYwDgyiuv\nxJQpU3DHHXdAp9Nh5MiRuHTpEmpra3HDDTfgiSeeQN++fTFhwgRs374dt99+O6677rp2v2szZ66Z\nkpJic1TSrLCwEDNmzOjU70/+jdVVyS80NjZi5cqV6NOnD+69915Fr5WVlYXs7GwkJSUpeh13ys7O\nxrx58xAbGyt7/tdff0V2djbeeecdD7WMvBFTSeQXJk+ejG+//RbTp0/3dFO80oIFC/D666+3e/7v\nf/87Fi5c6IEWkTfjiIGIiGQ4YiAiIhkGBiIikmFgICIiGQYGIiKSYWAgIiKZ/w8IFkqTQu4e8AAA\nAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x224fa1c6828>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.set()\n",
    "_ = plt.plot(xm, ym, marker=\".\", linestyle=\"none\", color=\"blue\")\n",
    "_ = plt.plot(xf, yf, marker=\".\", linestyle=\"none\", color=\"magenta\")\n",
    "_ = plt.xlabel(\"body temperature (f)\")\n",
    "_ = plt.ylabel(\"CDF\")\n",
    "_ = plt.legend(('Male', 'Female'), loc='lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting the CDFs implies that there is significant difference betewen them. A permutation test will determine the significance and difference. A permutation test simulates what the data would look like if the temperatures for the two genders were identically distributed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def permutation_sample(data1, data2):\n",
    "    \n",
    "    # first construct a single array with all the observations in it\n",
    "    temps = np.concatenate((mtemp, ftemp))\n",
    "\n",
    "    # scramble the entries in the array\n",
    "    temps_perm = np.random.permutation(temps)\n",
    "\n",
    "    # create permutation samples\n",
    "    perm_m = temps_perm[:len(mtemp)]\n",
    "    perm_f = temps_perm[len(mtemp):]\n",
    "    \n",
    "    return perm_m, perm_f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEFCAYAAADqujDUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzsvXl4HNWd7/3pqlKr1WqtRtiWZdmW\njZEXbHlhMcZmMWAMARIYlixMMknum9zMJDMJSd6EuZBM7oWQhJtnnknycJ9kEpLhvSGQTJKBJARj\nDDiYxXiRjWwLY8mbLFuWrbW61Wp1Vb1/nGpVb5Kx49Zi/T7PU8/vnD461UfVVWf51jm/43Mcx0EQ\nBEGYcGijXQBBEARhdJAGQBAEYYIiDYAgCMIERRoAQRCECYo0AIIgCBMUY7QL8H5pb+8d7SLknLKy\nIJ2dkdEuxphErs3wyPUZmol+bSoqioZMkxHAGMIw9NEuwphFrs3wyPUZGrk2Q5PTBmDnzp3cd999\nGZ9v3LiRO++8k3vuuYdnnnkml0UQBGE4hlkFNOwKITv7x6YJW97SMM2/qlRnTP8wAoHZMnLlGBVM\n6PqlBmdxzXMmAf3kJz/h2WefpaCgIOXzgYEBvv3tb/Ob3/yGgoICPvzhD3PttddSUVGRq6IIgpCN\nOKoB8JFRE8TjqgHw+cBIryWiqAZAAwLex6YJ118fpLlZp6bGYsOGCKFQDsvv0rIdBqKQF4Cqpalp\n7/wCgjZENLjk47kvy4hjQl5dkIt6dLoeshioj8AZXPOcNQDV1dX84Ac/4Ktf/WrK501NTVRXV1NS\nUgLAsmXL2Lp1K+vWrRv2fGVlwQkxlBtOr5voyLUZnjO6Pg6qAUhgoBoCVMUfT0ozDNUQAKrijybl\nCzCoI+zbB83NKtzcrNPaWsTKlWfwD5wF/b1QHGCwISoOQL57GcwWVfkDBO0iCvohVJXb8ow07X+C\n0h4VLu3Rad9URMUZNHQ5awDWrl1LS0vm2Ms0TYqKvBu1sLAQ832MFyfCS5yKiqIJ8bL7bJBrMzxn\ndX3O8QigshJqarwRQGVlhPb2M/xHzoKeqDcC6IniNVD5qucftIuIaL305UPfCJRnRFkNXcVBSnt0\nuootWJ15zYfrGIz4LKBQKEQ4HB6Mh8PhlAZBEISz4DRa/mAPPhkDsIAsA2vDANsCLUuaGYdtb2os\nu8JOURtCIdiwIUL9Jo261faIyD+gZJ/2vVAxLzPtkjtB3wTW6pEpy4gTgoH6CO89q1F6m31G8g+M\nwiyg2bNnc+jQIbq6uojFYmzdupUlS5aMdDEE4fwhnnSkJ8W9I4Mo0E+qpJPIFwWrX9lkTBNWXxXk\nro8UsvqqYMbL3pMvwLQ2m5MvnNV/clZs/QU0bVA2BRMKVwYp/4SyZ/OSdFwQgtKPnHnlDyPYADz3\n3HM8/fTT5OXl8bWvfY1PfepT3Hvvvdx5551Mnjx5pIohCOcXDl7vPzmM6vknZvIkhwEl4dhZwoBj\nqyM9DLDtLY2WVjUsaGnV2faWV4VEj4Pj6tFOj4rnmp4WiLvfGe9R8QTxNzWCx1RZg8d04m/KrPd0\nfOPFHfRE0H9F5x4auTbDEIeKSUW0n+rNquUnpJxsWn48CkbSS9QE3adg6xsay1fYlEzyPk+MAFpa\ndaoqLTa9ljrT5+B/wkAL5FXBzDvP5T85NDt+DIFGjWitzZL/JynBHQEEj+lEplqEN5/ZDJnzhTH1\nDkAQhBEkjtL5HTKe9q4urwEoneJ9bppw7XVBWo7pVE212LTZq+RDIfjuzRF27dZYtCBT5y+tUPlD\nIzWr24RrfhLEf0gnNsOi+yNJlXwIwpsjaI1FhGsnZuV/OmRMJAjjmYTs4yNTAnLlG58vU8qJx9SB\nlhR22famRosrnbQc09mWJJ2074aBXphXbTPQq+KD5+yEuAmBCmXjnTn4f9Pw7dLwH1Jl9R/S8e1K\nq9JCEFiHVP5DIA2AIIxnfOpwbC88mKSpwxrwwgkMP/TH4fVNGv1xFU+w7AqbqqkWAFVTLZZd4bUc\nFQugrAhKmjXKilR88JxlYITAiitrlOXkP07BWWQTm6HKGpth4SwaYomykBWRgARhnJM80ydd54+a\nEO8HYwBCaat2b/tgkCPHdaZPsXj19VSZ55t/E6HhHY2Fl6TJPCbc8GyQwAmd6CGL3gdSpRVjEVjt\nYIyUBBSC7pcj+HZpqvKXnv4ZISMAQRjHDEo7WWQea8Dr/SfCCba/rnHkuJJOjhzX2f66VxWc2AcD\nJlw8y2bAVPHB73tNI3BC5Quc0HFe8/JZcbf3X+aFR4QQOFdK5X82SAMgCOOYZGknXebR89SRHgZY\neqXN9ClKOpk+xWLplV7LceFcKCpX4aJyFR/8vqtsoheqfNELLXxXefl0Qx3pYWHsIj+RIIxzjABY\nljudM50A7G7QWLQ8VRsPheDV1yO8+B8aN/xt5myeG/8eTv5K44J70zT1EPS+GaHjZY28azN73aFS\n1fOXyn98ICMAQRjntDTA0T3KJmOasObqIB+8u5A1V2eu2t30MMR32mx6OO2EJhSuCDLvgUIKV2Su\noO0yoXuxTdcQK2ul8h8/SAMgCOOY/jDE3CmcsZiKJ2io1zhwROn1B47oNNR7j/vxXRDuUOFwh4on\nsDZrBNvcFbRtOtZmL1/ylNH06aPC+EMaAEEYJzhZZjjmF0Lcgi1vK5tf6KUtrLOZNV3p9bOmWyys\n804wZRGUFEHJXo2SIhVPoK+0iUxW+SKTLfSVXj7Drw7H9sLC+EUGa4IwDohH3UVdWqrWb5rwsc8E\nOdgCM6uCbNyUOp3zWx+JsHuvxoJ5mdM5b3suSFGnTm+bRfSbqStoT/xnhPjrGsaVNoVZdP6E+2Vh\nfCMjAEEY4wznnO2dbRoHW5Rcc7BF551t3iN96gBYfVA708bqU/EEkRc0ijpVvqJOncgLSdM5+yFu\nAKtt4oaKp5clMQrINioRxg/SAAjCGGGoyjQxvdOOZ071vGSZTXWlkmuqKy0uWeadZNIsCJUoD5mh\nEhVPEFxrY5aqfGapRXBt0nTOfND97mwev4qnlyU9LIxPRAIShDHAUBJPcroVc3vfaenx/lSbjNYD\nmqlsOonzZJ0+GkK5lSjMTDICXlmF8Y38hIIwygwn8YDq+Sd6/4lwgq2varSeUlJO6ymdra96j/TJ\nRjA7oHiysicbk75zi0bAXQkcOK7jbMlc0avnD72iVyr/8wP5GQVhlBl02tafXVbRDHX09XrhBMuv\ntqmcpKScykkWy6/2Wo8LaqEkBIX7NUpCKj74nZfZRF3pKFpp4bssc0Wv48iK3vMd+WkFYQwQPen5\n5g9VZaa3NCv/PHmdMGdpapoTT7WDmHDdb4KE2nXM/RZ9X0md6XPytxHstzS0y23y02b6BEIQHwAj\nD+E8RkYAgjDKWP1K3/dpylppWn5fL1hR1fO3oiqe4K3nNY51KynnWLfOW897j3T/Jo1Qu0oLtev0\nb/LS7BhYBjgrbSxDxRMkto9MjALGx56BwtkgDYAgjDKJWTeQOesGoKAIdPdFrR5Q8QSXr7OZ5s7m\nmVZqcfk6T8rJX21jVrgzfSos8ld7aZpfHelhUBvI+HyZYeH8QxoAQRgDhKrACGaXfwCmzIWWLo0p\nc9PyheAv2yP8328om7LYKwR9GyO0frmPvo2ZWyIGpoC/Qtl0DMM7hPMXaQAEYQxwcBMc365sOqYJ\n16wK8uFPFHLNqkynbm89CL5dyqZmhKJ1QSofK6BoXaZTt6ip/AdFh3DqJj3/8x9pAARhlIl2guVW\nwpap4snUb9M4fFRp+YeP6tQnrfZt266meIKybdu9fM52jYC7SjjQouNsT3oHYKkjPSxMLKQBEIRR\nJlAG/TpseUujX1fxZOqW2VRPc1f7TrOoS1rtO3kplBZBYJ+yk5NmCPmW2kSr3KmeVRa+pUnvAHR1\nOLYXFiYeovAJwihjmvDRLwY5elJn2q8t/rIlVcsPheBnj0TYvUNjwZJMp25rNwQp6oDeY0GiZupU\nz+M/j+DUa/jqMp26GQb4HJnnP5GREYAgjDJbXtE4elJ1wY+e1NnySupj2dcB/gFYstDGP6DiCSLr\nNYo6XKduHTqR9UkrevuUUzdruevUrc/Ll1hxnBgFiFO3iYk0AIIwylx2Tepq3suuSa2NC8rBKIJY\nj7IF5V5a8EYbs8yd6llmEbwxaUVvgZpW6tju9NICL584dRNAJCBBGBMk77KVDcPV6Y1sWr2WZpMI\nVEE8DIY4dROyID+9IIwyf3lW40SvqtlP9Or85dnUxzLaqdxA5JcomzxLKPaqRsh1Bhc6pRNLcgaX\nkHb0gqFlHqn8Jzby8wvCCJLNs+aq22wmFykZZ3KRxarbUmvqQBnk61DQoJGfNkvIf7WNeYErAV1g\n4U9yBpeQdhxHZB4hOyIBCcIIYXa5bpYNta1iMomGIVsDgQkLHwiSf0Snf7pFz8rUmT5N341QsLOI\nvsURKtNm+mCg/PrLNE8hC9InEIQRINmvfrqP/Y1PaZzsUzX0yT6djU+lPZb1GvlHVHr+ER3qkxy+\nnYJoDPQVyvaf8rIlHLklRgHi1E1IRxoAQRgBkv3qp/vYv+7DNhcGlYxzYdDiug+nifV1Nv3TVXr/\ndAvqkhy+TQJ/sQr7i1U8gTh1E06HSECCMFIYsHunxqJlqRV8KARvvRXhpf9usObxeMaCLULQ82qE\ngTc08lbYGU7dqq6HkAXFWWQewwDbBk26ekIW5LYQhBHANGHN1UE+eFcha65Oc+hmQtGlQT71WgFF\nl2Y6bQM4fhCOTbI5fjAzLXpSyTvRk5lp8ThYlrKCkE7OGgDbtnnooYe45557uO+++zh06FBK+k9/\n+lPuuOMO7rzzTl588cVcFUMQxgQN9RoHXB3/wBGdhiQdv+0XGpOiKm1SVKftF6mPZSysvHaCsrGw\nl2YPqCM9DKm6v7wDELKRswZgw4YNxGIxnn76ae6//34effTRwbSenh6efPJJfvWrX/Gzn/2MRx55\nJFfFEIQxwcI6mxmuY7YZVRYLk3T8yR+36QiotI6AxeSPp0pE/kLw+8EaUNaftKhLy1OHY3vhBPIO\nQDgdOXsHsG3bNlatWgVAXV0dDQ0Ng2kFBQVUVlbS19dHX18fPrkzhQmANZBqU3DSbBqhKaB3Q0FJ\nZpoRco8sC70Md1tHecSEbOSsATBNk1CS20Jd14nH4xjuFkNTp07llltuwbIsPvOZz5z2fGVlQYys\n6+DPLyoqik7/RxOU8XxtGrZCS5sKt7TpHGws4tp1Kn70lzDN3Qe4vF/n6LNFTPtHL681AP444K4d\nCJWC7vb0HVttJg9wwaQijIAs+MrGeL53cknOGoBQKEQ47ImVtm0PVv6bNm3ixIkTvPTSSwB86lOf\nYunSpSxatGjI83V2RnJV1DFDRUUR7e29p//DCch4uTamCe82alxcm+q2eWYtVE8NcviYTvVUi5m1\nEdrbVZr/Nuj8WpCyqE5nwMJ/m5c2eN4uiPeDkQ+xtKc2HlWV/8lTvRiB3P5/45Hxcu/kiuEav5w1\nAEuXLuXll1/m5ptvpr6+nrlzvc1MS0pKCAQC+P1+fD4fRUVF9PT05KoogjAimCasvTHIe/t1Lppj\n8cJ6z69/KATf+9sI9Ts16han+fQHtVo32aYRCEE8TzUA6RgB7xCEMyFnDcANN9zA5s2buffee3Ec\nh0ceeYQnnniC6upq1qxZw+uvv87dd9+NpmksXbqUlStX5qoogjAivNuo8d5+JVO+t1/n3UaNZcuV\nMN9xAOw+WDTXxu5T8fJZKt/JpzTmuSuBy/p09j6lccF/8wT9QadueV44XeYR2Uc4G3yOMz4mh02E\nIdxEH6oOx3i4NqYJN1wfpKlZZ3aNxYsbUnf2evMH0LIDqpbAFZ9PzgjGgiBlfTqdBRbx3ZGMxV7x\nqFfxZ+vpj4frM1pM9GszKhKQIExEhpvM03UY+o8rm45mpNp0xHe/kAvkdhKEc0Rjo0Zzs5Jympt1\nGhu9x+tYPUQ6oGCqssfqvXzdz2qUuPsBlPTqdD+b/bGUyl8418gtJQjniNpamzlz1IKuOXMsams9\nHX9qHQTdrRyD5SqeoOQ2m+5ila+72KLkNtmgVxgZpAEQhHNEKAS/+12Ex77bx+9+F8mY6XPT92Da\nUmVTM0Jsc4R9H+0jtjlT/xeEXCENgCCcI0wTbr81yJe/WsDtt6Y5fAP+8+Nw6M/KpmaEwDVB5v7f\nAgLXZHcGJwi5QBoAQThH7GnQaD7ovgM4qLOnwXu8Dm8Gq0OFrQ4VT9C3UaPI3de36JRO30Z5LIWR\nQe40QThHzF9oM2uG0vJnzbCYv9DT8qtXgl4OnQeVrU5a9lJwnY1Z7u7rW25RcJ28AxBGBmkABOEc\nMhBNtSkYENfJOvk64cUz2ZunIOQaaQAE4Ryx8y2NljYl5bS06ex8K0kCegNiHWr1b6xDxRPE39QI\nuvmCbTrxN+WxFEYGudME4Ryx+HKb6VOUlDN9isXiy5MkoBXgL1eLwPzlKp7AuMImMlnli0y2MK4Q\nCUgYGaQBEIRziBVLtckE/WD4lE3ndCuBBSEXSAMgCOeItzdqtHYoKae1Q+ftpNk8rVugpwOKZyvb\nusXL59RrBI6qfIGjOk69PJbCyCB3miCcIy69zqZqkpJyqiZZXJo0m6fyMih2VwIXl6t4Al+dTXS6\nyhedbuGrEwlIGBmkARCEsyGLt7dQCJ5/JcKXP9bH869krgRe9z1YXquxLstK4N5XI5z4dZjeV2Ul\nsDByiOIoCGdKHNUA+Eh5gkwT1q4M0tqj88tnLV6rT2oETMhbEuSKbp2uX1oM7Eit6ONxcBbY+OLy\nUAojh4wABOFMcEj1+Zw0Etj8nEZrj/sOoEdn83Pe49X9R43SbpVW2q3T/UcvzYmrIz0sCLlGGgBB\nOBN8YEZg61YNM0LKFo4rb7WpdL16VhZbrLzV0/JLbrHpLnE9fpZYlNzipfkMdaSHBSHXyK0mCGeA\nacKNNwXZ36QzZ7bF+hdT9/39p8UR3tutcdGCzH1/NT3VJmOUqp6/VP7CSCIjAEE4Axr3auxvUjX4\n/iadxr3eI3TgVfCbsGCGjd9U8QR9GzWK3CmiRR3ZHb5J5S+MNNIACMJQZJmNWTvPpmaWknJqZlnU\nzvP+aNbV4CuFE43KzrrayycO34SxiPQ5BCEbUVQDoAFpm7Bb8VSbjBZTD5WWZSWwLy/VCsJoIyMA\nQUjHxuv9J4eBhl0ah44oKefQEZ2GXUkS0AYY6IbS2coe2ODlG3hLo9B1+FbYpjPwljx6wugjd6Eg\npKOB2QdbtmiYfaQ8JQsXpfr8X7goSQK6HvJKoLNR2VnXe/nyLreJuI7iIlMs8i4XCUgYfUQCEoQ0\nTBOuXxek+aBOzUyLDRtTZ/o8/r8i7HtHY+4lmTN9fA7YtrLpGP5UKwijjYwABCGN4bZ2jHRCvgWX\nzLfJt1Q8wYEXob8byucpe+BFL01r0PAfVuf0H9bRGuTRE0YfuQsFIY35C21qZrozfWambu0YLAPd\n7fXrIRVPMOsGyC9R4fwSFU9gL7SJueeMzbSwF4oEJIw+IgEJQhqhEPzXcxGe+43BrX8Tz5B5aler\nvX3LZmbmvfv38Mo/wzUPp58Uup+LkPe8wcC6uDh8E8YEMgIQhDRME25ZG+SBfynglrVBTDM1veVN\n6G5WNp2n74IjbyubelIo+1CQ0q8WUPahIJiZeQVhpJEGQBDSqN+icbhV6fWHW3Xqt3iPSX83xCMq\nHI+oeIJDr0DMjce6VTyB3qhh7FfnNPbr6I3y6Amjj9yFwsQmixRfd5nN9Knu3r5TLeou8/4ovwSM\nIMROKpvQ/AFmXAP+EuhpVXbGNV6aVWsTr1HnjNdYWLXyDkAYfaQBECYu0aQjjeFW+9IOWqey6RRW\ng69QWUEY60gDIExMhlntu3WzRmu769e/XWfr5iQJ6ARYPZAXUrb/hJfv2A4Id8Ck+coe2+Gl6Y0a\nRrMrATWLBCSMDeQuFCYmmjrsuBdOsHylTeWFrl//Cy2Wr0ySgC4EvRj6OpTNv9DLN3UJlBRBaZNG\nSZGKJ7BqbeJzXAlojkhAwthApoEKE5ZoFOwYaDYE0hy+ObFUm0y3Dv0GWDpUJSeY8MEXgxSd0ult\nt4iaSds+hqBzfQS9UVOVv0wDFcYAMgIQJiR23Ov9D4Zd3npB41iXkmuOdem89YL3mITbwe6BvDJl\nw0nvAfo2ahSdcn3+n8ri8z8E1nKp/IWxQ84aANu2eeihh7jnnnu47777OHToUEr6q6++yt13383d\nd9/NN7/5TRwni/MUQfgrGeq20gzo64ft2zT6+lU8weVrbaaVKblmWpnF5Ws9uaawArRi6OtStrDC\nyyc+/4XxRs4koA0bNhCLxXj66aepr6/n0Ucf5fHHHwfANE2+973v8R//8R+Ul5fzk5/8hM7OTsrL\ny3NVHGECEo+rBsDnAyPtTjdNuPXOIPubdebUWKzfkOrw7X+ujbBjq8aS5ZkO3yyg38nYJgDwGhJN\nxFVhHJCz23Tbtm2sWrUKgLq6OhoaGgbTduzYwdy5c/nOd77DkSNHuOuuu05b+ZeVBTGMLJupnmdU\nVBSNdhHGLGdybRxHNQAJDEM1BAn274P9zW64Wed4axErVqr4sR1gDMCli4EBiLd4L3TNdigACmao\neAEQckcB0bch4M4KCp7Q0fYWEbj5TP/Ls0funaGRa5OdnDUApmkSSuo66bpOPB7HMAw6Ozt56623\n+P3vf08wGOSjH/0odXV1zJo1a8jzdXZGclXUMUNFRRHt7b2jXYwxydlcm+FGAFMqYU6NNwKYUhmh\n3dXzjSq1yKurA0rLVbw9WesH4j1gFKtwXyJtHthTggSP60SmWITnRejNslYgF8i9MzQT/doM1/jl\nrAEIhUKEw+HBuG3bGO5TWFpayiWXXEJFheo6LV++nL179w7bAAjCmWIYaiGXnuUuD4Xgd/8V4U/P\nGdx8a6bDt+v/N2z5EVz295l5F90MHfugfG76SSH8YoTepwy0D4vDN2Hsk7OXwEuXLmXTpk0A1NfX\nM3eu97QsXLiQffv20dHRQTweZ+fOncyZMydXRREmKFET+sPKpmOacOstQf7fBwq49ZZMh29//Ds4\nvF7ZdI5vAnO/sqknheJbgkx+pIDiW8ThmzD2ydkI4IYbbmDz5s3ce++9OI7DI488whNPPEF1dTVr\n1qzh/vvv59Of/jQAN910U0oDIQh/LbaljuSwlvQK6Z16jYPuBi0HD+u8U6+x4io1a+fI62qhFyh7\n5HWYfqWKxzoh5lbsMVPF/Yk9Aeo18t1z5h/WoV6Dq2QmkDB28TnjZP7lRNDwJrpWORxnc21OtcHu\ndzQWXGIzaXJqmmnCddcEOXhYZ2a1xcZXIiky0B//Do5vhylL4ZYnUvMe3wSRVghWwpTVySeF4tVB\n8lt0+qssejZFRkwGkntnaCb6tRmVdwCCMJqYJtxya/Z9fUG9A/g/j0Zo3KVRuyhzqqcD+GxlM9BR\n4mmWSWlaXqoVhLGMrAQWzkuG29cXINIFfhsWLbTx2yqeoOV15cytsEbZlte9tFgXxMJglCobS8qn\n7dXIO6C+M++AjrZXHi9hbDPsHbphw4aRKocgnFPmL7SZNUO9BJg1I3VfX4BgKRiFYJ5SNljqpVVd\nCYXl0NekbNWVXpq/FAIaBPdqBDQVT2DPsxmYrb5zYLaFPU/0f2FsM2wD8KMf/Wgw/KUvfSnnhRGE\nc4ljpdp0OjrgVLey6eTFId+nbAomzP3nIPMeLmTuP6fN9AnBqT9HOPFfYU79eeT0f0E4W4ZtAJLf\nDx84cCDnhRGEc8Xueo2DLe4snxad3fWpt3rXCejvgWCRsl1Jfv1bX4NIBwSqlW19zUvTGtJkniRp\nyXHADkJ8uY0dHNoPkSCMFYZtAHzJa+cFYRyxoM5mZpXq+s+sslhQlyrHlF4I+cXQdUzZ0iS//pVX\nQXExlB3UKC5W8QT2QpuBWa7MM8vCTpKWfD7P3URyWBDGKsPOAhoYGODYsWPYtj0YTh4VVFZW5ryA\ngnC2JHwBxbNt6wj4DoPRBr6CtAQTbtkapKhbp3drpl//tj9F0Bo07IU2+Wkyj2F47icEYawzbAMQ\niUT42Mc+Nljpf/SjHx1M8/l8vPTSS7ktnSCcJfVva7QcV1JNy3Gd+rc1rrrW6613HlT+fErLle08\nCGUzVZr5nEZFt+vXv1vn5HMaoQ+rvLYNVgFYl3pxLW0cLZW/MF4YtgHYuHHjSJVDEM4pdZfaVFda\nHG7Vqa60qLs0VQIqm6mcuSWcuiUqf4DQrTa9D1lqBFBiEbrVy6tp6khU/OmVvyCMJ067EKyrq4s/\n/OEPNDc3k5+fz0UXXcS6desoKEgfNwvC2CEUgj++EOE3vzD4m49nOnsDWPH30PhnqL0pPTNEX4tw\n8nsGoa9kOnULBLL3/AVhvDHsLbx7925uuukmNm3aREFBAbqu8/zzz3PDDTfw7rvvjlQZBeGMMU1Y\ntybIvzxWwLo1mc7eAF58GPa+qGxqZgheE2TWfxQQvCa7Uzep/IXzgWFHAN///vf5zne+w9VXX53y\n+caNG/nOd77Dz372s5wWThDOlm2vabSccN8BnNDZ9prG1Td5Uk5bI5ju/H+zQ8Un16r4wMsahe7e\nvoWndLpf1si7VRZ1Cecfw/Zjjh8/nlH5A1x33XV0ZFs9IwgjjGnC21u0jB7+sqtsqi5U0zWrLrRY\nluaVc3IthMqhr03ZROUPkHetTdjd2zdcbpF3rVT+wvnJsCMAv98/ZJqsERBGG9OEG64P0tSsM7vG\n4sW0fX0fuSdCfb1GXV2mszeA0iCE/cqm4zNSrSCcj7yvdQDZPEYPDAzkrFCC8H7Yu0ejqVlJNU3N\nOnv3aFx6meqtt+9T/vrnz7GJmSpekbTlRGcTRDqhrErZziYom63SrDc0gq58FDyh0/OGhn6DjAKE\n848zWgeQjIwAhBEli1uFefNtZs2yOHBAZ9Ysi3nzvUq6Yi4Ul8PhLVB9WWrlD6qyD5ZB236YPMer\n/AH0FTZ9ky0K2nT6JlvoK6TyF85PZB2AMPaJJx1pd2yib5LN747ZAv0dymbDXw7+kLLp6HqqFYTz\nkdNOZnvjjTdoamoajD/55JMGHcVFAAAgAElEQVS88cYbOS2UIAzi4PX+k8PA7t0aB12f/wcP6uze\n7d3OR7cqX/7Flcoe3Zp62p4WiPZASZWyPUmNhG+7hr9VndffquPbLnM+hfOTYe/sP/3pTzz00EP0\n9fUNfjZp0iQefPBBXnjhhZwXThDwuUd6GFiwwKamRs3WqamxWLDAk2qmLVe+/EHZactTT1tcBYFi\nFQ4Uq3gCZ6lNzHUkF6uycJaKBCScnwy7J/Add9zB448/zuTJqRuqtrS08IUvfIHf/va3OS9ggomw\np+dE37t0KEwTjrcWMaWyN2M2j2nCznqNxUPM9Hn9f8OV9w997rYdMHlJtgTQf29gfTAOk7OkjzHk\n3hmaiX5thtsT+LT7AaRX/gBVVVXYtvSKhNxjmnDjDUGuvErZ9Pn+0S6oqbSJdmXm/e09cGi9stnY\n/yfofVfZ1C+FktuClD9YQMlt2VcCC8L5wGkbgHA4nPG5aZoyDVQYERr3auxvUnr8/iadxqR9dgf6\nYSDmhmMqnuDQq9DvNgr9XSqeTF870ONGety4i9ag4Xc3ffGnbfoiCOcTw97Zt99+O1/84hc5cuTI\n4GfHjx/ny1/+MjfdlO5BSxDOPbXzbGrcDVhqZlnUJu2zm5cPeX5V+ef5VTzBjKshv1TN788vVfFk\nCipA80Ngr4bmV/EE9kKbmPudsbRNXwThfGLYaaB/93d/R2dnJ7feeit5eXn4/X76+vr42Mc+xj/8\nwz+MVBmFCY4VT7XvFycOzoCyGZhQ97+C5Lfo9FdZ9NyYuunLid9G8DVoOFk2fRGE84VhG4Bf/vKX\nfOlLX+Kzn/0se/bsIRgMMnv2bPLz83n44Yf553/+55EqpzBBadilceiIkmMOHdFp2KVxxZWqR56Q\ngBKjgIF+bxRwcCMMdMOki5Q9uBFmXpd04u0a+e6ewfktOmzXYLW7yUscrHxgmRfXxCWEcB4yrAT0\n61//GoBgMMjDDz/M/Pnzyc9XT9jWrVuHyyoI54SFi2xmVCs5Zka1xcJFmRJQfzhTApp5HeSVQNd+\nZVMqf4ClNv3T1Hn7p1mQNNVTM9Rh215YEM5Hhr21k2eIDjNbVBByim2l2mSiJ6C/G4wBoCo1rciA\niK5sNnxaqk0hgOoeDe0PURDGPe97eoP4/hFGg107NI4cVVLNkaM6u3Z4t2x/L8Siqucfi6p4gpZX\nIdIFRXOUbUmbBeTbpeF3pSX/ER3fLu+8tu31/hNhQTgfGbYBkEpfGG0WLbGpdqWa6mkWi5Z4tXF+\nEfgDMNClbH7SepeqqyFYCv1Hla1KmwXkLLKJudJSrNrCWZS57296WBDON4aVgN577z3WrFkDQFtb\n22DYcRza29uHyyoI54x4PNUmY5hghJVNxymHPh1Kszh7Ox2y768wERi2ARB/P8Jos+1NjdY2JdW0\ntulse1Pj6utVbz3WqXz+5xUqG+sEf5nK17oF+jqgtFbZ1i1QeZl3Xq1Bw3/YlYAOq8Ve9hWpWo9U\n/sL5zrANwLRp00aqHIKQlWVX2FRNtWg5plM11WJZUiXtL1PunGOm69a5zMtXeRkUlKvKv6A8tfIH\nd7HXTAv/QZ3YTFnsJUxMhnUGN5aYCM6cJrrTqqFoa4Pnniri1g/3ksU1FV0NULowS0YTDv+bRvUX\nbG+RV1q6b6eGs3iI9HGE3DtDM9GvzVk7gxOE0cY0Ye01QR54RNl0Z3CN/x8c36RsakYoWB5k2b8W\nUrA8u0O3OBBbaHOGC4wF4bwhZw2Abds89NBD3HPPPdx3330cOnQo6998+tOf5qmnnspVMYRxzpsb\nNVpPue8ATum8udG7ZSNHSXHoFjnq5et7SSPUofKFOnT6Xkq91R1LHelhQZhI5KwB2LBhA7FYjKef\nfpr777+fRx99NONv/vVf/5Xu7u5cFUEYb2SR4a+4zmZKuaqdp5RbXHGd90fBaUAx9J9UNpj0yqpg\njY3p5jPLLQrWpJ7cp6sDJyksCBOMnC1y37ZtG6tWrQKgrq6OhoaGlPQ///nP+Hw+Vq9enasiCOOJ\nKKoB0FCrcJMYbhooBeDkK5uOpqfadIyAchTnE1cPwgQlZ7e+aZqEkrZo0nWdeDyOYRjs27ePP/zh\nD/zbv/0bP/rRj97X+crKghjG+d9NG+6FzXmLjWoAEiTcMACv/xlOujLPyR6dXa8Vcft9Kh45DoU2\nFNaocxRaEJyi0qJbIOAuVQm262h7igjckvQdDqSI/wYp202ORybkvfM+kWuTnZw1AKFQKGUzGdu2\nMQz1db///e9pa2vj4x//OEePHiUvL49p06YNOxro7Izkqqhjhok8W8E8BXsbNOYttAlN8j5fdBVc\nWBLkRLfOhSUWi66KMLgGUYewBrEW8FdBWIdwIm0+2BcGCZ7QiVxoEZ4foTd97WIc1RD4yOGTMDJM\n5HvndEz0azNc45ez237p0qW8/PLL3HzzzdTX1zN37tzBtK9+9auD4R/84AdccMEFIgVNYEwTbrg5\nSNMBndmzLF58KTK4v28oBH+/PMK+d4uYe3EkY9/fUDmYXcqmo+upNgMDrwEQhAlIzl4C33DDDfj9\nfu69916+/e1v8/Wvf50nnniCl156KVdfKYxT9jZoNLlbMDYd0NmbtAVj02vg64H5s5Rtes3L198B\n8V4IVCjb35F00m0a+cdcf//HdNg2xK0ulb8wgcnZCEDTNL71rW+lfDZ79uyMv/v85z+fqyII44R5\nC21mzbA4cEhn1gyLeUmrcmdfBa+WQ9sumLxIxRPkl4NRBP3HIX+Kig+yzKa/0iK/Vae/0hrc3EUQ\nBI9xrnwK5wuWlWpTMEELk3UxV6gcfDEozCIBDevvXxAEWQksjD4NOzUOu9szHm7Radjp3Zb7N0C8\nGy6oUXb/Bi9f3IR4P+SXKRtPaiC0nRp+95z+Fh1tp9zqgpCOPBXCqLNwsc1M1zf/zGqLhYs9uWbO\n9WCUqLBRouIJjBAY7jaQRr6KJ7AXp/r7txeLBCQI6YgEJIwYpgnvNmpcXGunzOYJhWDjKxG2bNa4\nbKWdMdPnk7+DHd+DJV/JPGfpYug/BfmT0hJC0P3HCHnPGQzcGh/3zt4EIRfICEAYEUwT1t4YZN3N\nhay9MdOpW9NOKPLZNO3MzPvL22H3C8pmnLcFot3KpiZA2W1BSh8ooOy27M7gBGGiIw2AMCK826jx\n3n6lyb+3X+fdRu/WM7sg5q4ZjIVVPEHzeoi57qJi3SqewOoHK+aGYyqeQG/QMNyppcYBHb1BbnVB\nSEeeCmFEuLjWZnaN0uRn11hcXOtp8qFS8BdC1yllQ6VevpobwV8CZrOyNTd6aXo+6H6wB5TV8700\na6FNfIb6vvgMC0s2fBGEDOQdgDBiOHaqTaanA7q6IdCRmTbJD7aubDoJh25GIDNtcJGXLPYShKzI\nCEAYERp3azQfVJJM80Gdxt3erdd5AsweCBUp23nCy9fyCkS6oXiusi2veGn2gNv7z/fCCfQ9Gob7\nfcZBHX2P3OqCkI48FcKIULvApmamkmRqZlrULvCGAWUXQqgYuo4pW3ahl6/qGgiWQOyIslXXeGla\nHugxyKvX0GMqnsCabxN3Jad4jYU1XyQgQUhHJCBhxLDiqTYZ8yhEToKZZbpmaSn0B5RNzQRT7g1i\nNOnEZ1t0vhjxpnuGoHNDBH2Ppip/mQYqCBnICEAYERp2ahxyV+YeSlvt235IST/F5cq2J+0eeuJN\nMDtg6kXKnnjTS9P3ahhNrszTpKPvTbudQ2BdJpW/IAyFNADCiLBwsc0sd2XurLTVvhUzlPQDylbM\n8PJdeIXn6jlUruIJrHk28dmuzDPbwponMo8gnAk+x3Gc0S7E+2EibOhw3mxcMYSP/bY2+K9nDG6/\nO87kyZnp76yHS27M/BwT+n9WRP4nezN78yZoezTsCS7znDf3Tg6Y6NdmVDaEESYoQ+yyZZpwy9og\nh1t1fvxTi1deS93c5Q/fht4OOLQNPvD1pPOZULAsSEUnmD8K0rctklLRxwPgLLHx+eRmFoQzRSQg\n4dzhuEd6GKjfonG41fX42apTv8W79Y6/qyp/UPb4u16+2HqNUKfKF+rUia338jmOOtLDgiC8P6QB\nEM6ObJWtD8wIbNumYUZIkYHqLrOprlR6fXWlRd1lnl4/5WIoKodom7JTLvby+W+0CZeqfOFSC/+N\nXj6fTx3pYUEQ3h8yahbOnGFknrXrgry3X+eiORYvrE/d2/eL90TY06Axf2Gmx88pQejPVzYdzUi1\nyRiG6vlL5S8IZ46MAIQzYxiZZziHb8f2Q8yEOTNtYqaKJ+hshmgXTK5StrPZS7Nf0yg4qc5ZcFLH\nfi3zlpXKXxDODmkAhDPDlXm2vp0p81xca3PRHCXXXDQn1eHb1Dlqnv+po8pOnePlK6uBQCmE25Qt\nq/HStKts+irUOfsqLLSrZKqnIJwrRAISzgjThBvXBtnfrDOnxmL9hlSZ55mnIjTu1aidlynz9Ecg\nFlU2Hb8GZZqy6QwnAQmCcPbICEA4Ixr3aOxvVpLM/madxiQna/EByPfD4sU2+X4VT9CyF8IdUDFV\n2Za9Xpp5AALdUDpFWfNA0hdu1cg/pr4v/5gOW+WWFYRzhTxNwhlRO9+mZpbr1G2WRW2SkzUjTx39\nUS+coGoeFJbDiWZlq+Z5aaFZEC2BSLuyoVlJX7jcpt+dPdRfacFykYAE4VwhDYBw5jhpNgnThHCE\njC0fAfw9UBBRNh2tELR8ZdPxaalWEIRzgzxSwhkxnF///ijEYqrnH4upeIKj9cqZW3m1skfrvbTw\nESjogcJKZcNHvDRtp4bfdSLnb9HRdsotKwjnCnmahDOidoHNHNfP/pyaVL/++QHwu7t2+f0qnmBa\nXapTt2l1XlrhdOhzncH1Fat4AnuxTcx1IhertrAXiwQkCOcKcQY3hhhLTqsce2jJpa0N/vxHg5tu\nye7UreMElF+Y+TnA24/Dpf89e1rXbihdkO0LoezlIjqv7YUs3yeMrXtnrDHRr81wzuBkBCBkEI96\nRzqmCbffGuQrXyvg9luDGVp/w9twsFHZdP78KTi4Xtl0Tr4NfaeUTf1CKLktiPEFZcnybkEQhLND\nGgAhBcdO3bw9fQP3PQ2p7wD2NHi3ULgHYn0qHOtT8QStb6Q6fGt9w0sb6IEBN99An4on0Bo0/Afc\ndwAHdLQGuWUF4VwhT9MEJr1yByX7RPpg23aNSF+mDDR/YerevvMXeicpLAZ/AfScULaw2MtXuUI5\neou1Klu5wkvLKwY/EGjQ8LvxBPZCm5g77TQ2y8JeKO8ABOFcIWsrJyjxqKfzG0kva00T1t0aZP9+\nnTlzLNavT/XbHwrBDx6J8O5ujYsXZK72DR+FvnYIG8ClqWkXXACRfGVTMGH2l4PkHdQZmGnRtTF1\nb98Tf4pQ1lJEZ1WE/Am86YsgnGtkBDABGU7maWzU2O86dNu/X6exMfUWCXdBng0L59nk2SqeoLsN\noj0QLFS2u81L69gJ8Q6ovFjZjp1emvaORp4rK+Ud1NHe8b7TtsEqAN+VytoyABCEc4Y0ABMQn6YO\nK+6FE9TW2tS40zxraixqa1Nr3MJSMAqhq03ZwlIvrWQyBIoh3KlsSdKMnfLFYJRDuFXZ8sVemn2J\nzYA71XOg2sK+xPtOTVNHelgQhL8eeZwmKNGod2QQT7NptB2GtuPKplNYAKUFyqaTXw4lQWXTcbRU\nm0wg4B2CIJw7pAGYgFhxr/efCCdobNBoPuzO8jms05g26+ZUK/T2QH6xsqdavbTwSbB7IFSmbPik\nl2bug4IOKK1R1tznpWm7NfyuBOQ/qKPtzrwtpecvCOeenD1Wtm3z0EMPcc8993Dfffdx6NChlPSf\n//zn3HXXXdx111388Ic/zFUxhCzoBkT7YesWjWi/iieoXWgzx53lM2emRW3arJtJlVBUDD0tyk6q\n9NIKLwCtWL0E1opVPEFoLvSVQ6xH2dBcL81eYBNzvzM208JeIEK/IIwEOZsFtGHDBmKxGE8//TT1\n9fU8+uijPP744wAcOXKEZ599ll//+tf4fD4+8pGPcP3111NbW5ur4ghJmCbc/sEgzYd0amZYbHg5\n1af/Uz+PsLdBY16WrRsB+k5Cf6+y6RjtkN+pbDoV8yEah8D8zLRknV8QhJEhZw3Atm3bWLVqFQB1\ndXU0NDQMpk2ZMoV///d/R9fVsD8ej5Ofnz/s+crKghiGnqvijhmGW7Z9rti3B5rdAVnzIZ3Ww0Ws\nvFrFB6IQnwzT3Be4pUWQl6S9tx0A24LJM5W1TZjsum+OHIVC2/XlY0NhDILTVFrsBPQ5UHIJ4ECJ\nA37XXYTTBD53G0ijWeeCtiJ8yS6hXUbi2oxn5PoMjVyb7OSsATBNk1BS91HXdeLxOIZhkJeXR3l5\nOY7j8N3vfpf58+cza1aWJz6Jzs4s20idZ4yUz5LKaqiZ4Y0AKqsjtCf12Lt6YSAGeX4weoGkImkh\n0HSIdEOwRMUH8/ohrCn9XyuGsB/CiTQfdPuUBOQvhgIfkEibDGVzghj7deJzLDonR7w0l4nuz+V0\nyPUZmol+bYZr/HLWAIRCIcLh8GDctm0Mw/u6/v5+HnjgAQoLC/nGN76Rq2JMeBwnc9P0UAj+6w8R\n/vAbgw/8TTxD5rlgJnQch/Ip2c/5gc/Dvr/A3FWZabUfg1PbYNKyzLQp10J8o4ZxbZrGH4LO9RH0\nRg2r1vYWgQmCkFNyprguXbqUTZs2AVBfX8/cud5bP8dx+NznPsfFF1/Mt771rUEpSDi3xOPekYxp\nwq1rg3z9Xwq4dW2mQ7cdm+HdXcpm4+1fQHu9suns/0/o3KVs6pdCydVBpn66kJKrszh1C4G1XCp/\nQRhJcjYCuOGGG9i8eTP33nsvjuPwyCOP8MQTT1BdXY1t22zZsoVYLMZf/vIXAL70pS+xZMmSXBVn\nwuE46kgOJ0YC77ytcbBVNboHW3XeeVtjhdsr7+2CqDtwi4ZVvChpsVd3i5JxQNnuFiipUvG+40DC\nkVuPihe4owjfDg3/EXeq5xEd3w4NZ5XM9hGE0SRnDYCmaXzrW99K+Wz27NmD4XfeeSdXXz3hyCbz\n+HzqcCzw6anpl1xqUz3V4vAxneqpFpdc6lXERaUQKISTR+CC6amVP6jK3l8MHe9C+cVe5Q9uZV8M\nThv4JnuVP4CzxCY2zcJ/VCc2zcJZIpW/IIw24gxunBOPew2AkfZrGnG3AXDI+KXjA6k2GcsEK6Js\nNowYFMSVTWfKxWAGITQ9M21QcJSpnoIwJpBHcRyTTeYZxFaHz+eFE2x7TaP1pJJjWk/qbHvNuw16\nTkJ/DxQVKduTNte/97Ca5VM4WdneJHcQVhisGBRMVtby5gDg25UmAe2SW08QRht5CscJ2TbuHJR5\n4l54EE0dAzEvnGDZVTbTKtTK22kVFsuu8lqH4guUm4f+sLLFaa6bi6rVFM/oCWWLqr00vRB0Pzgx\nZfXCpPIvSt3b11kkEpAgjDYiAY0DhpN5MIEBIA9I0+tPdqmFXXkBuCBtSmfC/4+VxeHblBLw90B5\nSfbyFIVAK4DCLDN2AhdAPA+MbHn1NCsIwqgiI4AxznAyjxP3ev+JcIKBmNf7Hwy7vLFe43inqoWP\nd+q8sd67Dfq7IRaB0CRl+7tTy9PXCk4PBKcq25fkDM6Jub3/Ai+cQNubtrXjXrn1BGG0kadwjDMo\n8ziZMo/PUEd/rxdOkOeHWBx27dSIxVU8wYobbSrLlRxTWW6x4kZPjskvAX8QYr3K5qf15AsqwVcM\ndljZgiRncD4/aHEw6jW0uIonsOfZxGa7EtBsC3ueSECCMNqIBDTOOd4M/aZywFa11PvcNOHDHw2y\nv0lnzmyL9S+mOnz7xtURGnZoLFySxeGbH8h3bRYumA+9hVA0Iy3BhEkfCZLXrDNQY9G1IXVrx44/\nR/Dt0XDm2+iy4EsQRh0ZAYxxkhdxpUtA/b0Qi6qefyyq4gka92rsb3K3dmzSaUySXI6+DfTAwtk2\n9LjxtHP6izPPCRA3Id6v5vjH+1U8gbZHI6/Z3dqxWUfb432n44AdVKt97WD2l9qCIIws0gCMcZJl\nn3QJKL8I/K6nTn9AxRPUzrOZ40ouc2Zb1CZJLtMuhYC7K1egXMXfzzkBjBAYruNWI1/FE9jzbQbc\n7SQHaizs+d53Dvd/CIIwOvgcZ3z0xSaCN7/hvBbG41lmAKGknq0bNZZflynlmCbsqdeYX5fdr/+2\nn8Cy/5a9LOF2KKwYuqz9pyB/UpaENjCeN4ivi8PkzORsq5bfDxPdo+PpkOszNBP92oyKN1Dh3GGa\nYFmg66RU5KYJ16wMDrp0eGVzJLWiPw4Xl9twHJiTes71XwGzA07tgxu/l5rWdRDiMRgIQ+nMLOVp\nUQu9BvogVJWcACW3B/E368T+j0V38jsAF+n5C8LYQSSgMY5lqSM9DLDzLY3Dx5TmfviYzs63vJ/T\n6lOVNLircvu8fMfrVeUPyh6v99LiUVX5g7LxtE3jrf608/Z7adoeDb/7DsCf9g5AEISxhzyhYwg7\ny8xIXYdoFLa+rhGNqniCxZcrp24A1VMtFl/unUAvUKtxI13uqtwCL9+UOgiVQ3+XslPqvDQjAIYf\n7AFljaTdwAD0fOVjyL9Lw4ir+GD559vE3HcAsbR3AIIgjD1EAhojRKPeEUiqdE0Tbrs5yIEjOrOm\nW7z0aup0zn/7SITd72gsuCRT529phfhJ6I5DbZoEVHEFDPig4vLMsgRCoDngzzZV04SqjwcxmnTi\nsy06X0yb6vlCBN9uDWeBTPUUhLGOjADGALbt9f6TwwANOzQOuE7UDhzRadjh/WTdTaCHYVGNjR5W\n8QSRTnBM0APKRjq9tN7jYPXA5HnK9h5PKsuA2/sPeOFk9L0ahju91GjS0fdmmep5qUz1FITxgDQA\nI0g2iQdA09Rh2144wcIlNjOmKVllxjSLhUl+9Etmg78M2g8rW+Jtt0CwDHwhpdH7QiqeoGgK+P0Q\nfFfD71fxwbLkqQM7KZyENc8m7k4vjc+2sObJVE9BGK+IBDRCRKNeBR8InP7vk4n3p9pkOrqgL6xs\nOlU10BuCogvTEkxY+f0ggSM60ekWvbenztYJhMDOAy2fTEJw8g8RtD0adpYVvYZx9lM9BUEYWWQE\nMAIMJ/Ekf5YYBSSn79ikcdT13X/0pM6OTd5PdmK3msXjn6Lsid1evsRsnWBp5mwdZ4dGwJWVAkd0\nnCRZKbF3QGIUQFpZHRvsAogvs7ELVDwdqfwFYXwgDcAIkCzrpEs8p0tfstpmuuu7f3qFxZLVXo17\n4QI1iweUvXCBl0/PV7N/wJ0FlNSb9y2xiU5X54xOt/Alb8+YvHdA2j4CAD5NHelhQRDGHyIBnWPi\nA2DkZX4eCLg+doaQfwIBYCBTHgqF4E8bIzz9rwb3/FM8Y6bPmq9A+680Ku7N7IqHqmAgAnnB9ATo\nfTVC55saxhV2xmItAq6b6SHuDiOgev5S+QvC+EYagHNI10mvAShN20mr67haWGX4oXRKZt6D26Eo\nAL1RmJnm1XPdqiAtXTo//0+LTduTVvuaUHhVkIpWncgPLcKvpWr5URNsS+3tG0h3E3EcrEob/TiE\n0qaIxqNuBR/PXAeQQCp/QRj/yGN8jogPpG60nrzZejyWtro2bTP1qOuBE5SNJrkteft5jZYupde3\ndOm8/bz3k8W3aARbVVqwVSe+xUuzLXWkh2H4VcKO7en6yWFBEM4/pAE4CwayzMYx8tQR7fXCg2l+\ndfR1e+FkAkUQd2D7dmUDSb6bLl1nU1Wqau+qUotL13k1snGZTcRdCRyZamFc5qVpujoc2wsnSKwS\ntq3MVcKi8QvCxEEkoDPkZIvaXjHPDxdUpaZ1HYBID8S7MyWgk3uVB83+k5kSkGnCJ/8xSPNhqKkO\nsuGVtNW+n4qw7XWNZVdmrvY19FSbkmaAzwE9y69sTFELxIwsq3VF4xeEiYE84mfAQL+3t+5ALHUk\nEO2BWJ/q+cf6VDxBXwcM9ILmV7avI/W8e+o1mg+rGrz5sM6eeu9n6WwCXxiWL7bxhVU8gW+nhr/F\ndb7WouPbmbQq1/Z6/+lSTkIS0gsy5aHBc8udIQjnPfKYnwF5+arn33FU2bykqZWBYggA/q0aATee\noKAcBvJg9x6NgTwVT2Z+nU1NtaqFa6ot5td5tXXZbAiUqQYlUKbiCZzFNjE3X6zawlmctCp3GCln\nUBJyMuUhQRAmDiIBnSG71kPPSWh7D677ZFKCCbWfD1J4TCc81SKy2ZuRY5rwmQeDNB3WmV1t8eKq\nVL/9oRD8+HsRjh4sYtrMSIbMM+liyItD8cVphQlB9wsRfPUaTl3mdM7hpJxAAOw4aHIHCMKERUYA\nZ8CpIxDuUS9Owz0qnmDgLY1C1zd/4TGdgSTf/HvrNZpciafpsM7e+tTL3tcN+Q4sWaxsX7eXFm0H\neqB4urLR9qSMNlAAzgobCshYtQtDSDkOg73/RFgQhInHed8AnDp27s41aTrkA90bDPLdeIK8y23C\n7oyc8FSLvCTf/PPqbGa7Us3saot5dak1dUEJ6EE1JVMPqniCQAVQjKrci914gsRKXYesq3aHxOce\n6WFBECYU57UA8NIzYPZAqBjW3P3Xny/cBpc9GuIi28d7Wx3CHzYpTOx7G4LOpyKc2KKRd5lNfprE\n86sfR9jzjsb8LH77AaZUQYEDfVkq45mrIXocAlkWkGGgKvAz1fENVMMhlb8gTFjO2xHAqWOq8gdl\nz8VIYN/PDS6yVY15ke1j38+99tOOgmWAdqWNZah4goEI+A2oW2LjN1Q8mYTjtoKSTMdtWOoIVHjh\nQRLyjcbZSTlS+QvChOa8bQAmTVU9f1B20tS//pxzPxHnPU3Vsu9pDnM/ER9M0wJqmicoqyW5UMgL\nqllD4M4eSvPNM5zjNnS83n1yGETKEQThr2LcNwCmCX/4pYZpZqatuRlWBjXW3HxuvqtwMuTtNHn+\ns33k7UySf1wCMyG/Sj3duCkAAAwtSURBVNl0LqiFSbOUzUaoyjsyE4FCMp22gZJydM5zMU8QhFww\nrhsA04TL64J88p8KubwumNoImBBcEWTePxUSXBGELA3E2TCwB2YujzOwJzMtHlXvauPRLGkm+Gxl\nh0LPtgELQBwl78SHSLOGSBMEQRiGnDUAtm3z0EMPcc8993Dfffdx6NChlPRnnnmGO+64g7vvvpuX\nX375rL7jlWc12nuUJtLeo/PKs96/M7BZo7DNnZbZpjOw+a//V+MdEHcdtcV7VTzBcE7UHEsd6eH3\nRbK2n67zD5cmCIJwGnLWAGzYsIFYLMbTTz/N/fffz6OPPjqY1t7ezpNPPsmvfvUrfvrTn/L973+f\nWCw2zNmyc81tNhXFqjatKLa45jav1s1baROe7E7LnGyRt/Kvd2tplIPhOmozilQ8wXArb326OtLD\n74vhdH55ByAIwl9BzpTjbdu2sWrVKgDq6upoaGgYTNu1axdLlizB7/fj9/uprq6msbGRRYsWndF3\nhELwVn2EV57VuOa2tOmVIYi8EaF7s6Yq/2z6+VlQeq3q+RvlmWnDrbw1Qqrnf0aV/2Bmhp6yKdM5\nBUE4S3LWAJimSSipRtZ1nXg8jmEYmKZJUZHn87iwsBAz21vcJMrKghhZXF5WVMCsfxwiUwUw62xK\nfxoqTv8nZ33qiqLT/9EERa7N8Mj1GRq5NtnJmQQUCoUIh8ODcdu2MQwja1o4HE5pELKRrfIXBEEQ\nzp6cNQBLly5l06ZNANTX1zN37tzBtEWLFrFt2zb6+/vp7e2lqakpJV0QBEHIPT7HcXIyd8S2bb75\nzW+yb98+HMfhkUceYdOmTVRXV7NmzRqeeeYZnn76aRzH4TOf+Qxr167NRTEEQRCEIchZAyAIgiCM\nbcb1QjBBEATh7JEGQBAEYYIiDYAgCMIERRoAQRCECYr4kBwFYrEYX//61zly5AihUIiHHnqIBx98\ncDC9ubmZD33oQ3z5y18exVKOHtmuT2trK4899hiGYbBixQq++MUvjnYxR41s16elpYXHHnuMgoIC\nVq1axec+97nRLuaIs3PnTh577DGefPJJDh06xNe+9jV8Ph8XXXQR3/jGN9A0jR/+8Ie88sorGIbB\nAw88cMbeB847HGHEefLJJ53/8T/+h+M4jtPU1OR88pOfHEw7fPiw86EPfcgxTXO0ijfqZLs+t99+\nu/Pee+85tm079957r9PY2DjKpRw90q/PJz7xCefqq692Dh8+7DiO49x///3O22+/PZpFHHF+/OMf\nOx/4wAecu+66y3Ecx/nMZz7jvPnmm47jOM6DDz7orF+/3mloaHDuu+8+x7Zt5+jRo84dd9wxmkUe\nE4gENArs37+f1atXA1BTU0NTU9Ng2sMPP8xXvvIVCgsLR6t4o0626zNv3jy6uroYGBigv78fXZ+4\nK8PTr8/27dspLi5m+nS1SfXSpUvZvn37aBZxxKmuruYHP/jBYHz37t1cdtllAKxevZrXX3+dbdu2\ncdVVV+Hz+aisrMSyLDo6OoY65YRAGoBRYN68ebz88ss4jkN9fT1tbW1YlkVjYyPhcJgVK1aMdhFH\nlWzX56KLLuKzn/0sN998M1OnTqWmpma0izlqpF+fWCxGNBqlqakJy7LYtGkTkUjk9Cc6j1i7du2g\nqxkAx3Hw+ZSHxMLCQnp7ezP8kyU+n8hIAzAK3HnnnYRCIf72b/+Wl19+mQULFqDrOs8++yx33XXX\naBdv1Em/PtXV1fzkJz/hj3/8Ixs2bGDGjBn87Gc/G+1ijhrZ7p/vfve7fPOb3+QLX/gCs2bNoqys\nbLSLOapomle1hcNhiouLz8oH2fmONACjwDvvvMOyZct48v9v795CouraAI7/zSOGYTlZWSaYdjBJ\nslHDQ9JoRNoBUvCAXhgVZpKXWlSgWCGSNpigBhZpollKzhSROZYYOUGSdlFBkIZWEJIOHWys8bsI\nh8as7H39vvFrnt+Vm3n2Wg9L2WuvPdtn1dQQGxtrXrp3dXWZS2jbssnj4+fnh6urK66u3zZT9vT0\nxGAwWDlL65nq76ejo4PKykrOnj3Ly5cvCQ8Pt3aaVhUQEIBerwego6MDpVJJcHAwnZ2dmEwmXr16\nhclkYsGCKeq62xB5C8gKfHx8UKvVVFdX4+bmxokTJ4BvG+XY+p0bTD0+vb297NmzB2dnZ9zc3Cw2\nGLI1U43P3bt3SUlJwcXFhR07duDv72/tNK0qNzeXY8eOUVJSgq+vL1u3bsXe3h6lUklSUpJ5x0Jb\nJ7WAhBDCRskjICGEsFEyAQghhI2SCUAIIWyUTABCCGGjZAIQQggbJROA+KtcvnwZrVYLQF5eHk1N\nTT/ElJWVWZQNmAlqtZq2trbf5hYVFUVRUdFPY9LT09Hr9ej1etLT06eMaWtr4+LFi/8q38kuXLhA\ne3v7jLYpZj/5PwDxV+nu7jbXgPlfysnJ+W2MVqvl1KlTREZG/uN+jEYj586do6am5h+3MZXU1FTS\n0tKIiIjAyclpRtsWs5dMAMIq9Ho9FRUVODo6MjAwgEqlwtXVldu3bwNQVVWFQqGgvb2dM2fOYDKZ\n8Pb2pqCgAIVCgUqlYufOnXR2dvLp0yeKioowGAzodDq6urpYuHAhAHfu3KGuro6hoSEyMzNJSkoy\n59DY2EhXVxenT58Gvq0MnJ2d2b9/vzmmrKyMR48e8fr1a9LS0vDz86O0tJTR0VEMBgOHDx8mNjaW\nvLw8QkNDCQ0NJTs7G39/f548eYKHhwdqtZra2loeP35Mfn4+R48e5ePHj5w/f57R0VGMRiMnT54k\nODj4t+PW0tKCUqk0j9vBgwfx9fXl+fPnBAQEsH79epqbmxkZGaG8vJwVK1agUqmIj4/n3r17ODg4\nkJWVRXV1Nf39/eTm5hIXF4eTkxMbNmxAo9GQkJAwk79qMYvJIyBhNT09PeTn53P16lUuXbrEggUL\naGpqYtWqVVy/fp2hoSGOHz9OeXk5Go2G4OBgCgoKzOe7u7tz5coVkpOTqaysJDw8HJVKxaFDh8wl\nNYxGI42NjVRWVlJaWmrRf1xcHPfv3+f9+/fAtzv0Xbt2/ZCn0Wjkxo0bpKamUltbS2FhIc3NzRQW\nFqJWq3+If/r0KRkZGWi1WubNm4dGoyE7O5vAwEAKCwuJioqivr6eiooKWlpa2Lt3L1VVVdMaM51O\nR0hIiPn42bNn7Nu3j2vXrtHd3c3g4CANDQ1s376dhoYGc5xCoaCpqYkVK1ZQVVVFdXU1xcXFFv0q\nlUp0Ot208hB/B1kBCKtZuXIlS5YsAWD+/PnmKqheXl4YDAZ6e3tZt24dy5YtAyApKcnigjVxkff3\n9+fWrVtT9hETE2PeFOTdu3cWn82dO5fo6GhaW1vx9vbG29ubRYsW/dDG95uGFBcX097ezs2bN+np\n6bEoLjbBw8ODgIAAc24jIyMWn8+ZM4fy8nJ0Oh0vXrzgwYMHFsXLfqW/v98iR4VCYe5r8eLFFmM4\nMDBgjpsoH+3l5YWnpycODg7mcZ6wdOlS+vv7p5WH+DvICkBYjaOjo8Xx5Br/JpPJ4nh8fJwvX76Y\nj52dnQHMZX+nMtHmz2ISEhLQarVoNBp27949ZYyLi4v559TUVHp7ewkMDCQzM3PK+Im8JvqdXG3l\nw4cPJCYmMjAwQEhIyE+/7J2KnZ2dxbhNfl7/s30Svj/n+7LJk8/91ViKv49MAGLWCgoKoqenx3wn\n29DQQFhY2C/Psbe35+vXr9PuQ6lU8ubNG/R6PbGxsb+MHR4epq+vj5ycHDZt2kRbW9sf9TWhr68P\nOzs7MjMzCQsLo7W1ddrtLF++3OLOfiYNDg7i4+PzX2lbzE7yCEjMWgqFgoKCArKzsxkbG8PLy8tc\nOfVnwsPDKSkp+aM671u2bGF4ePi3b7+4u7uTmJhIfHw8Dg4ObNy4kdHR0T/efGX16tWsWbOGbdu2\nYWdnR2RkJA8fPpzWuZs3b0av1xMdHf1HfU6HXq8nJiZmxtsVs5dUAxU2a3x8nLGxMTIyMjhy5Ahr\n1661dkq/9fnzZ1JSUmhoaPjhEdq/YTQaSU5Opr6+Xl4DtSHyCEjYrLdv3xIREUFQUND/xcUfvn2/\ncODAAerq6ma03ZqaGrKysuTib2NkBSCEEDZKVgBCCGGjZAIQQggbJROAEELYKJkAhBDCRskEIIQQ\nNuo/JbZb7mXk05YAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x224fa26f5f8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# now to take a statistically relevant number of permutation samples\n",
    "for i in range(50):\n",
    "    # generate sample with above function\n",
    "    perm_sample_1, perm_sample_2 = permutation_sample(mtemp, ftemp)\n",
    "    \n",
    "    #compute ECDFs\n",
    "    x_1, y_1 = ecdf(perm_sample_1)\n",
    "    x_2, y_2 = ecdf(perm_sample_2)\n",
    "    \n",
    "    # Plot ECDFs of permutation sample\n",
    "    _ = plt.plot(x_1, y_1, marker='.', linestyle='none',\n",
    "                 color='blue', alpha=0.02)\n",
    "    _ = plt.plot(x_2, y_2, marker='.', linestyle='none',\n",
    "                 color='magenta', alpha=0.02)\n",
    "    \n",
    "\n",
    "_ = plt.plot(xm, ym, marker='.', linestyle='none', color='blue')\n",
    "_ = plt.plot(xf, yf, marker='.', linestyle='none', color='magenta')\n",
    "\n",
    "# Label axes, set margin, and show plot\n",
    "plt.margins(0.02)\n",
    "_ = plt.xlabel('monthly rainfall (mm)')\n",
    "_ = plt.ylabel('ECDF')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Most of the ECDFs from the permutation samples do not overlap the observed data, suggesting that the hypothesis is not supported.\n",
    "\n",
    "The empirical difference of the means is the difference of means of the male and female observations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# the test statistic is the difference of means, so that function will be used\n",
    "def diff_of_means(data_1, data_2):\n",
    "    \"\"\"Difference in means of two arrays.\"\"\"\n",
    "\n",
    "    # The difference of means of data_1, data_2: diff\n",
    "    diff = np.mean(data_1) - np.mean(data_2)\n",
    "\n",
    "    return diff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-0.289230769230727\n"
     ]
    }
   ],
   "source": [
    "empirical_diff_means = diff_of_means(mtemp, ftemp)\n",
    "print(empirical_diff_means)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A distribution of permutation replicates of this test statistic can be used to identify the probability of the test statistic exhibiting of a mean difference at least as extreme as what was observed. In other words, how likely is it that the observed difference of mean temperatures between males and females would be -0.289 degrees F supposing the null hypothesis that the true population mean difference is 0 degrees."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def draw_perm_reps(data_1, data_2, func, size=1):\n",
    "    \"\"\"Generate multiple permutation replicates.\"\"\"\n",
    "\n",
    "    # Initialize array of replicates: perm_replicates\n",
    "    perm_replicates = np.empty(size)\n",
    "\n",
    "    for i in range(size):\n",
    "        # Generate permutation sample\n",
    "        perm_sample_1, perm_sample_2 = permutation_sample(data_1, data_2)\n",
    "\n",
    "        # Compute the test statistic\n",
    "        perm_replicates[i] = func(perm_sample_1, perm_sample_2)\n",
    "\n",
    "    return perm_replicates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By repeating the drawing of permuted replicates, a distribution of replicates is created"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "permutation_replicates = draw_perm_reps(mtemp, ftemp, diff_of_means, 10000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The visualization of the permutation replicates is shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAETCAYAAAArjI32AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzs3Xd4lFX68PHvTHoyk0ZCIIRQIoHQ\niwgISBFBOgoYQCnKKrqir7oi4A8LLAtiB1fWDrvoQgB3EWRVJFSRFhCQHkMSSCFMeiY9M+f9I2Qk\nQhokmcxwf66Li8xT7/PM5M6Z85znHI1SSiGEEMJmaa0dgBBCiFsjiVwIIWycJHIhhLBxksiFEMLG\nSSIXQggbJ4lcCCFsnCTyG0hISCAsLIyxY8da/o0ZM4aNGzdaJZ4NGzbw1Vdf1Wi7tWvX8sknn9R1\naBVKTk5m1KhRjB07ll9++aXcuqlTp9K2bVsuXbpUbvnBgwdp27Ytn3/+eZ3HN2nSJMaOHcuIESPK\nvdd/+ctf6vzcteHixYv8v//3/+rs+PHx8eWuy+jRo5k0aRLff/+9ZZv33nuPzZs3A6WfvYEDB/LE\nE0+QmJjIyJEjGTt2LCdOnKizGMXvHK0dQEPl6urKN998Y3mdkpLCqFGj6NixI+3atavXWI4cOUKb\nNm1qtN3kyZPrOqxKHTx4ED8/P1avXn3D9YGBgXzzzTfMnj3bsmzTpk34+fnVS3zr1q0DSv9ojx49\nutx7bQsSEhKIi4ur03N4eHiUuy6XLl1ixowZODo6MmTIEJ5//nnLuk2bNjFnzhxGjhzJxo0badKk\nSb38QRalJJFXU0BAAC1atCAuLo527dqxYcMG1q5di9lsxtvbm1deeYWQkBDmzZtHZmYmly5dYuDA\ngaSlpeHq6sr58+dJS0tj8ODBeHt7s3PnTgwGA4sXL6ZPnz7MmzePNm3aMHPmTADL6+DgYHbs2MG+\nfftwdXVl2LBhvPrqq6SlpWEwGGjWrBnvv/8+R48eLbddeno6GRkZvPrqq0RHR7No0SIyMzPRaDQ8\n9thjjBs3joMHD/Lee+/RvHlzoqOjKSkpYeHChfTo0YOoqCjeeOMNzGYzALNmzWLYsGHXXZeIiAjW\nrFmDVqvFz8+PV155hZSUFN5//31ycnKYOnUqa9asuW6/MWPGsGXLFksiz8/P5+jRo/Tp08eyTUpK\nCosWLSI5OZni4mJGjhzJk08+CcBHH31EZGQkBQUF5OfnM3fuXO677z4++OADEhMTMRgMJCYmEhAQ\nwFtvvUXjxo1r9H5v376djz76iJKSEtzc3Jg3bx5dunThvffe4/LlyyQlJWEwGOjcuTM9evTgm2++\nITExkblz5zJixAjee+89kpKSLNt16NCBxYsX4+HhQXJyMosWLSIlJYXi4mJGjx7NE088QXx8PI8+\n+ijBwcEkJyfz73//m3Xr1rFz504KCwvJz89n/vz59O/fn9dee42UlBQef/xxFixYwPjx44mKigJK\na9Nlrzds2MCmTZvIzc3F29ub1atXExERQUREBGazGV9fX1555RVatWpV5TVp3rw5zzzzDJ9//jlD\nhgzhxRdfpGPHjiQkJHD69Gneeecd0tLS+OKLL8jJyWHGjBmsXr260mt58uRJrly5Qvv27XnjjTdY\nuXIl27dvx2w207x5c1577TX8/f2ZPHkyPXv25MiRIyQnJ9OvXz8WLlyIRqMhMjKS5cuXo5TCw8OD\nRYsWERoaSlRUFO+88w4FBQVotVqeffZZBgwYQEpKCnPnziUrKwuAwYMH88wzz9To89HgKHGdS5cu\nqa5du5ZbdvToUdWzZ0+VlJSkDh48qKZMmaLy8vKUUkrt3btX3X///UoppebOnaumT59u2W/u3Llq\n4sSJqqioSF25ckWFhoaqf/3rX0oppVavXq0effRRy3afffZZuf3KXl/78+rVq9XHH3+slFLKbDar\nP/3pT+rzzz+/brsVK1aohQsXquLiYnXvvfeqH374QSml1OXLl1X//v3V0aNH1YEDB1RYWJg6ffq0\nUkqpzz//XD388MNKKaWmTZumvv32W6WUUmfOnFGvv/76ddfp559/VkOGDFFpaWlKKaW+/vprNXz4\ncGU2m9XXX3+tnnjiiRte30ceeUR99913atSoUerYsWNKKaU2bdqk3njjjXJlmDp1qoqMjFRKKVVQ\nUKCmTp2qtm7dqhISEtTUqVNVfn6+Ukqpb7/9Vo0aNcpS7nvvvVfl5OQopZSaNWuWWr58+Q3jUOrG\n7/Vvv/2mRo8erTIzMy3l79u3ryooKFDvvvuuGjJkiMrJyVF5eXmqe/fu6s0331RKKfX999+r4cOH\nK6WUevfdd9XAgQNVWlqaMplM6tlnn1VvvfWWUkqpKVOmqF27dimllMrPz1cPP/yw+uGHH1RcXJwK\nDQ1VR48eVUopdfHiRTV9+nRVUFBguUZjx45VSim1b98+NWbMGKWUUnFxcapHjx6W+K99vX79etWr\nVy/L9fj555/VI488Yrl2u3btsly7a/3xmGXOnDmjunfvrpRS6i9/+YtatWqVUkqpSZMmqR9//NFy\nzqeeeqpa13LEiBGqpKREKaXUhg0b1AsvvKCKi4uVUkp9+eWXatasWZbjv/DCC8pkMqns7Gx19913\nq8OHD6vLly+rHj16qDNnziillPrf//6nZs2apdLT09XQoUNVYmKiUkqp5ORk1b9/f5WcnKyWL1+u\nFi5cqJRSymg0qmeffdZyfWyV1MgrUFBQwNixYwEwmUz4+Pjw1ltv0bRpU9asWUN8fDyTJk2ybJ+d\nnU1mZiYAPXr0KHesQYMG4eTkhL+/P+7u7vTv3x+A4OBgyz7VNX36dKKioli1ahVxcXFER0fTpUuX\nCrePi4ujsLCQoUOHAqXfLIYOHcrevXvp1asXgYGBhIWFAdC+fXv++9//AjB8+HAWLVrEjh07uPvu\nu3nhhReuO/bevXsZMWIEvr6+ADz44IP87W9/IyEhoVplGTt2LJs3b6ZLly5s2rSJ+fPn88UXXwCQ\nl5fH4cOHycrKYvny5ZZlZ8+eZcSIEbz55pts2bKF+Ph4jh8/Tm5uruW4d911FzqdzlKmsppXde3b\nt4+UlBSmTZtmWabRaLh48SIAd999t+X4/v7+5d7Pa881fPhwy7WZMGEC77zzDk8++SRHjx7l3Xff\n5d1337WU68yZM7Rt2xYnJyfL+9m8eXOWLFnC5s2biY+P55dffiEvL69GZQFo166dJd5du3YRGxtL\neHi4ZX1GRgY5OTno9fpqHc/Nza3a567qWnbt2hUHBwdLbKdPn2b8+PEAmM1mioqKLPsNGjQIrVaL\nXq+nefPmZGZmcuXKFcLCwizNncOHD2f48OFERkZiMBh46qmnLPtrtVrOnz/PPffcw6xZs0hISODu\nu+/mpZdeslwfWyWJvAJ/bCO/ltlsZuzYscyZM8fy+sqVK3h5eQHg7u5ebntnZ+dyrx0dr7/sGo0G\ndc2wN8XFxTc891tvvcWJEycYP348vXr1oqSkpNx+f2QymdBoNOWWKaUoKSmxlPNGMUyaNIlBgwax\nb98+9u7dy9///ne+//57XFxcyl2HP7r22FUZPXo048ePZ8aMGRiNRkJDQ8sdWynFunXrLIkjPT0d\nFxcXTp06xZ///GdmzJhB37596dmzJwsXLrTsW1GZqstsNtOvXz/eeecdy7Lk5GQCAgKA6r2ff1xu\nNpvRarWYTCag9OZg2XHS09NxdXXFYDDg6uqKVlvaB+HXX39l9uzZPProo/Tr148ePXqwZMmS685T\n1Wfn2s+jyWRi/PjxlvZtk8mEwWCodhL/9ddfy71PVanqWv4xtieffJKHHnoIgMLCQrKzsy3r//i+\nwvXX3mw2c/78ecxmM6GhoZZ7IVDaVOfr64uTkxORkZHs37+fAwcOMGHCBL744gtLhcYWSa+Vm9Cv\nXz+2bt3KlStXgNIeItOnT7+lY/r4+HDy5Emg9AN36NAhyzoHBwdLcvzpp5+YPn0648aNo1GjRvz8\n88+W5HDtdmVat26No6Mj27Ztsxz7hx9+4O677640nkmTJnHmzBkefPBB/vrXv5KdnY3BYCi3Tf/+\n/fnf//5Heno6AF9//TXe3t60aNGiWmUOCAigbdu2vPzyy5ZvP2V0Oh1du3Zl1apVQOk3nsmTJxMZ\nGcnhw4fp2LEjjz76KHfddReRkZGWa1Abevfuzd69e4mNjQUgMjKScePGUVhYWKPjbN++HaPRiMlk\nYuPGjQwePBgvLy86dOhguQmclZVFeHg4u3btum7/Q4cO0aVLF2bMmMGdd97J9u3bb/hee3p6UlBQ\nwIULFwD49ttvK4ypf//+bNmyhdTUVAC++uorHnvssWqVJyYmho8++ohHH320upegRteyX79+rF+/\nHqPRCJT2ipk/f36lx+/atSvnz58nJiYGgG3btjF//ny6detGTEwMR44cAeDUqVMMGzaMtLQ0li1b\nxqeffsp9993HggULaNWqFdHR0dUuU0MkNfKb0K9fPx5//HEee+wxNBoNOp2Ov//979fVfGti6tSp\nvPjiiwwbNoygoCB69+5tWXfPPffwxhtvAPD000/z5ptvsnz5cpycnOjevbvla+q125VxcnJi5cqV\nLF68mA8++ACTycTTTz9N7969OXjwYIXxvPjiiyxZsoT3338fjUbD7NmzCQoKKrdN3759mTFjBtOn\nT7fcOPv4448tNcrqGDt2LC+//DIffPDBdevefvtt/vrXvzJ69GiKiooYNWoUY8aMITU1lW3btjF8\n+HDMZjODBg0iKyvLkgBuVbt27Xj99dd57rnnUErh6OjIypUra9SkAODr68uf/vQnMjIy6NWrF48/\n/jhQmqAWLVrEli1bKCoqYty4cYwYMYL4+Phy+48ePZrt27czYsQIzGYzAwcOJCMjg7y8PNq0aYNW\nqyU8PJyIiAheeOEFZs6cia+v7w1vSpcZMGAAM2bMYMaMGWg0Gjw9PW947QFyc3Mtf2C1Wi0uLi68\n9NJLlqak6qjJtZw8eTJXrlyxNPsEBgbe8BvItRo3bsybb77JnDlzMJlM6PV63n77bfz8/FixYgVL\nly6lqKgIpRRvv/02TZo0YcaMGcybN49Ro0bh5ORE+/btuf/++6tdpoZIo2r6vVMIUaX33nuPvLw8\n/u///s/aoYjbgDStCCGEjZMauRBC2DipkQshhI2TRC6EEDauWr1W0tLSePDBB/niiy8ICQmxLN+x\nYwcffvghjo6OjB8/3tL/szIGQw4+Pu5kZNT8wQZbIeWzbVI+22XPZfP3r7ivf5WJvLi4mFdffbVc\nZ/yy5UuXLmXjxo24ubkxefJkBg0ahL+/f5UBOTo6VCNs2yXls21SPttlz2WrTJVNK8uWLWPSpEnX\nDToUExNDcHAwXl5eODs7WwZaEkIIUb8qrZH/5z//wdfXl/79+183trXRaCz3WK+Hh0e1Hsjw8Sl9\nJLeyrwn2wCrla9my9P86Ht4U5P2zdfZcPnsuW0UqTeRff/01Go2G/fv3c+bMGebOncs//vEP/P39\n0el05QYqys3NrdZ4DRkZefj76zEYcm49+gbKWuXzNZf2JE2v43PL+2fb7Ll89l62ilSayK+dlWbq\n1Km8/vrrljbwkJAQ4uPjyczMxN3dnaioKMtY2kIIIepPjcda2bJlC3l5eYSHhzNv3jxmzpyJUorx\n48dbRjQTQghRf6qdyMtmebm2++HgwYMZPHhw7UclhBCi2uSBICGEsHGSyIUQwsZJIhdCCBsnE0sI\nUYVdxxKvWzawazMrRCLEjUmNXAghbJwkciGEsHGSyIUQwsZJIhdCCBsniVwIIWycJHIhhLBxksiF\nEMLGSSIXQggbJ4lcCCFsnCRyIYSwcZLIhRDCxkkiF0IIGyeJXAghbJwkciGEsHGSyIUQwsZVOR65\nyWRiwYIFxMbG4uDgwNKlSwkODrasX7VqFRs3bsTX1xeAhQsX0rp167qLWAghRDlVJvKdO3cCsG7d\nOg4ePMjSpUv5xz/+YVl/6tQpli1bRseOHesuSiGEEBWqMpEPGTKEgQMHApCUlISfn1+59adOneKT\nTz7BYDAwcOBAZs2aVSeBCiGEuLFqTfXm6OjI3Llz+fHHH1mxYkW5dSNHjmTKlCnodDpmz57Nzp07\nGTRoUIXH8vFxB8DfX38LYTd8VimfVlNv576d3j+9zrXS9bbI1uOvjD2XrSIapZSq7sYGg4GHHnqI\nrVu34u7ujlIKo9GIXl964b766isyMzN5+umnKzlGDv7+egyGnFuPvoGyVvl8e5Q2b6UfOVmn57nd\n3j97m7PTnt8/ey9bRarstbJp0yY+/vhjANzc3NBoNDg4OABgNBoZNWoUubm5KKU4ePCgtJULIUQ9\nq7JpZejQocyfP5+HH36YkpISXn75ZbZt20ZeXh7h4eE8//zzTJs2DWdnZ/r06cOAAQPqI24hhBBX\nVZnI3d3dWb58eYXrx40bx7hx42o1KCGEENVXrZudQtgye2vjFuKP5MlOIYSwcZLIhRDCxkkiF0II\nGyeJXAghbJwkciGEsHGSyIUQwsZJIhdCCBsniVwIIWycJHIhhLBxksiFEMLGSSIXQggbJ2OtCPEH\nJ2NSWbX5JBevGHFxdkCr0XBnO3+aNvKwdmhC3JAkciGuys4tYtX/znA8Jg2AQD8PSkxmDJn5bI9K\noHf7ANo097ZylEJcTxK5EEByWi7vrT9OalYBHUMaMbZvS0ICvQCI2BHNrl+S2H8qhcJiEx1bN7Jy\ntEKUJ23k4rYXnZDJkjVHSM0qYEzflix5qq8liQME+LozvHcw7q6OHItOJctYZMVohbieJHJxWzt/\nKZN3I45TUGRi5sgwxvVvjUajuW47Tw9n7gprjFnB4bMp1GCqWyHqnCRycdv6LSGL9zYcp8Rk5qlx\nHenbqWml2zdvrKNpI3eSUvM49ltqPUUpRNUkkYvb0sWUHN7bcIziYjOzxnSge6h/lftoNBp6hjVG\no4F1kdEUl5jrIVIhqlZlIjeZTMyfP59Jkybx8MMPc/HixXLrd+zYwfjx4wkPD2f9+vV1FqgQtSUn\nr4j31h8nv9DEn0aHcWe7xtXe11vnQttgbwyZBUSdvVKHUQpRfVUm8p07dwKwbt06nn32WZYuXWpZ\nV1xczNKlS/niiy9Ys2YNERERGAyGuotWiFtUUFTC9qgEsnKLmHxvG3q3b1LjY4S18EEDbD9yqfYD\nFOImVNn9cMiQIQwcOBCApKQk/Pz8LOtiYmIIDg7Gy6v0Dn+PHj2Iiopi+PDhdROtELfAbFbsOZZM\nTl4xHVr54uSkLTcxc3UnZNa7O9O1jR+/RKcSk5RVroeLENZQrX7kjo6OzJ07lx9//JEVK1ZYlhuN\nRvR6veW1h4cHRqOx0mP5+LgD4O+vr3Q7W2eV8mk19XZuW3r/9DpXAH46lsjl9DxaBXoyoHvQdb1T\nri3TtT+X7X+t8YND+SU6lZ9+TaF3l6A6irzu2NL7V1P2XLaKVPuBoGXLlvHiiy/y0EMPsXXrVtzd\n3dHpdOTm5lq2yc3NLZfYbyQjIw9/fz0GQ87NR93AWat8vubSLnHpdXxuW3v/cowFXEjK5vhvqXh5\nONO7fQDG3MLrttvw41mgNHHnGAsqPWb3EF8C/Tz46Xgi4/q2wEvnUiex1wVbe/9qwt7LVpEq28g3\nbdrExx9/DICbmxsajQYHBwcAQkJCiI+PJzMzk6KiIqKioujWrVsthS1E7TDmF3Pg1GWcHLUM6t4M\nJ8db76yl0WgY0iMIk1mx+1hSLUQpxM2r8hM9dOhQTp8+zcMPP8zMmTN5+eWX2bZtGxERETg5OTFv\n3jxmzpzJpEmTGD9+PAEBAfURtxDVopRi/8nLlJgUd4U1xtPDudaO3btDAC7ODuw9kYTZLA8ICeup\nsmnF3d2d5cuXV7h+8ODBDB48uFaDEqK27DmeRHJaHs38PWgd6Fmrx3Z1dqRXWGP2HE/mdHw6HVvJ\nGCzCOuSBIGG3MnIKidjxG06OWvp0CLjho/e3qn+XQAD2Hk+u9WMLUV2SyIXd2rT3AgVFJnqE+uPu\n6lQn52jd1JNm/h4cPW8gJ08G0xLWIYlc2KUEg5Gffk2mmZ8HdwTVXT9vjUZD/86BmMylbfFCWIMk\ncmGXNu6KQSmYMDAErbb2m1Su1adDAI4OGvaeSJZREYVVSCIXdudsfAYnYtJoF+xN55C6vwGpd3em\nWxt/ElNzuZCUXefnE+KPJJELu7Ppp1gAJg66o05ucN7IPVdveu45Ln3KRf2TRC7sym+JWZy/lEnH\n1r60alq73Q0rE9bSh0aerhw6c4X8wpJ6O68QIIlc2JnvDsQDMLJ3i3o9r1ajoX/nphQWmzgsw9uK\neiaJXNiNpNRcfolOpXWgJ6FWmO2+b6emaIC90rwi6pkkcmE3vjtYWhsf0btFvbWNX6uRlysdWvsS\nk5RNgqHyUUCFqE2SyIVdyDQWsv/kZTw9nMnKLWTXsUTLv/o0oEvpmOY7jiTU63nF7U0SubALu35J\nxKyuzt5jhdp4mW5t/PDzcuXnk5cx5hdbLQ5xe5FELmxeicnMrmNJODlqa31grJrSakuHty0qMbPr\nl/r9NiBuX5LIhc07fPYK2blFtAnyqpWxxm9V/y6BuDo7EHk0gRKT2drhiNuA9T/1Qtyi7VEJaIC2\nwfXfU+VG3Fwc6de5KVnGIg6fka6Iou5Ve6o3IRqimKQsYpOz6XqHH3r32ps0oioV3UQtm8B5yJ3N\niTySwHcHL9KrQwBaK7bbC/snNXJh07ZHlfYOuffOhjUBcmNvN+4KCyDBYOR4dKq1wxF2ThK5sFkZ\nOYVEnb1CM38P2rfwsXY41xl1d0s0wOZ9cTIqoqhTksiFzdr5SwIms2JIjyCrdjmsSDM/D+5s15j4\nlBxOxKRZOxxhxypN5MXFxcyZM4cpU6YwYcIEIiMjy61ftWoVI0eOZOrUqUydOpULFy7UabBClCkq\nNrHrlyQ8XB3p3aGJtcOp0Oi7WwKweV+s1MpFnan0ZufmzZvx9vbmrbfeIiMjgwceeIB7773Xsv7U\nqVMsW7aMjh071nmgQlzrwOkUjPnFjOjdAhcnB2uHU6GgxjqCA3TEJufw1Y/nCWqsA36/KSpEbag0\nkd9///0MGzbM8trBofwvzKlTp/jkk08wGAwMHDiQWbNmVXlCHx93APz99TcTr82wSvmuzoRTH+e2\n5vunlGLXsSS0Wg0ThrTF38cNAL3OtdbOcbPHutF1ubtzIBd/PM+JC2m0a9UIjUZj9c+/tc9fl+y5\nbBWpNJF7eHgAYDQaefbZZ3nuuefKrR85ciRTpkxBp9Mxe/Zsdu7cyaBBgyo9YUZGHv7+egyGnFsM\nveGyVvl8zaVf3dPr+NzWfv/OxGcQl5zNXWGNoaTEEkuOsaBWjq/Xud70sW50XZwdNLRqqic2OYdT\nMam0aKJnw49nb7h/fdTUrf3+1SV7L1tFquxHnpyczNNPP82UKVMYPXq0ZblSiunTp6PXlx58wIAB\nnD59uspELsSt2h51CSjtq93QVNS/vMsdfsRdzuFYdCrNA3TSr1zUqkpvdqampvLYY48xZ84cJkyY\nUG6d0Whk1KhR5ObmopTi4MGD0lYu6tyVzHyORafSqqmeECuPq1ITnh7OhDTzIiu3iFiZ11PUskpr\n5B999BHZ2dmsXLmSlStXAjBx4kTy8/MJDw/n+eefZ9q0aTg7O9OnTx8GDBhQL0GL21dkVAKK0tp4\nQ+xyWJnOIY24kJjN8d/SaNnUEwetbcUvGq5KE/mCBQtYsGBBhevHjRvHuHHjaj0oIW4kv7CEn35N\nwkvnTM92ja0dTo3p3JwIbe7F2YuZxCRkEdpAxoYRtk/GWhENRlXjl+w6lkh+oYn7e7Xgp1+T6zO0\nWtMppBHRCVmciEkjpJknDg7yTJ64dfIpEjahqNjED4cu4erswODuttsH283FkXYtfMgrLOHcpUxr\nhyPshCRyYRN++jWZ7NwiBncPwsPVydrh3JKOrXxxctRy8kK6jFcuaoUkctHglZjMfHfgIk6OWu7r\n2fC6HNaUi7MD7Vr4UFBkIvpSlrXDEXZAErlo8A6eTiEtu4B7ugTi5VF/Y47XpbAW3jg6aDgZm47J\nLLVycWskkYsGzawUW/fH46DVcP9dwdYOp9a4OjsS2tyb/MISYhKkX7m4NZLIRYN2KcXI5fQ8+nRo\nQiOv2htLpSFo39IXrba0Vm42y8iI4uZJIhcNllKKXy+koQFG9Glh7XBqnburI22CvDDmFxN32T7H\nBxH1QxK5aLCSUvNIzy7kznaNaeLrbu1w6kT7lj5ogDNxGTJeubhpkshFg/XrhdJZdUbaYW28jN7d\nmeYBOtKyCzBk5ls7HGGjJJGLBsmQkc+VjHya+XkQHGDf40uHXZ1v9ExchpUjEbZKHtEXDdKpuHQA\nOrTyrfDRfXvR2McNX08XLqYYSc3Mx8/bzdohCRsjNXLR4GTnFnExxUgjT1cCfO0/qWk0GsJa+KCA\nyKMJ1g5H2CBJ5KLBOW2pjfvY3FC1N6tlUz2uzg7s+/UyxSXygJCoGUnkokHJLyzht8RsdG5Odt82\nfi0HrZbWgZ4Y84v5Jdpg7XCEjZFELhqU6IQszGZF+5Y+aG+ziRfaBHkBsOd4kpUjEbZGErloMMxK\nEX0pE0cHDSHNvKwdTr3z0rnQJsiL03EZ0hVR1IgkctFgJBlyyS0ooXWgJ06Ot+dH854ugQDsPWGb\nE2cI65Duh6LBOH91ooU2zW/fKdDyCktwctQSeSQBb50zWq3GMkOSEBWptNpTXFzMnDlzmDJlChMm\nTCAyMrLc+h07djB+/HjCw8NZv359nQYq7FtqVj4Jhlz8vFxp5Glfg2PVhJOjllZN9eQXlnA5Pc/a\n4QgbUWmNfPPmzXh7e/PWW2+RkZHBAw88wL333guUJvmlS5eyceNG3NzcmDx5MoMGDcLf379eAhf2\nZc/x0qaE0Nu4Nl6mVVNPzl/KIjY5m0A/D2uHI2xApYn8/vvvZ9iwYZbXDg4Olp9jYmIIDg7Gy6v0\nplSPHj2Iiopi+PDhlZ7Qx6d08CN/f/vuWmaV8l3t5VEf567Nc5jNigOnLuPsqKXjHf4Non1cr7Pe\ntwKdhwu6k5e5mGJkyF3OdfKkCQwiAAAgAElEQVR+2vPvnz2XrSKVJnIPj9LagNFo5Nlnn+W5556z\nrDMajej1+nLbGo3GKk+YkZGHv78eg8F+h+20Vvl8r45pnV7H567t8p2NzyA1q4A7grwoKCiioNaO\nfHP0OldyjNaNokWAnlOx6ZyNTaNX29r9lmvPv3/2XraKVFn1SU5OZtq0aYwdO5bRo0dblut0OnJz\ncy2vc3NzyyV2Iapr/6nLALRu6mnlSBqO1oGlv0uxyTJ7kKhapYk8NTWVxx57jDlz5jBhwoRy60JC\nQoiPjyczM5OioiKioqLo1q1bnQYr7E9xiYmoc1fw0bvcFuOqVJe3zgVvnTMJV3LJKyi2djiigau0\naeWjjz4iOzublStXsnLlSgAmTpxIfn4+4eHhzJs3j5kzZ6KUYvz48QQEBNRL0MJ+HP8tjfxCEwO7\nNbttxlWpDo1GQ6umnvwSnUrUOYOlf7kQN1JpIl+wYAELFiyocP3gwYMZPHhwrQcl7F/Z0LQ7j5b+\n7+hg/RucDU3Lpnp+iU7l8NkrkshFpeS3R1hNYZGJRIMRH70LPnoXa4fT4OjdnWnk5cqZuAxy8oqs\nHY5owCSRC6u5eCUHs4JWTeUmeUVaNtFjVoqj52VERFExSeTCauKvzhzfookk8oqUXZvDZ69YORLR\nkEkiF1ZRUGQiOS2PRp6u6N2drR1Og6Vzc6JVU0/OxmeSLc0rogKSyIVVXLqSg1LQoonO2qE0eD3b\nNZbmFVEpSeTCKqRZpfrubFf6ZOfhM9K8Im5MErmod8b8YmlWqQE/LzdCAj05ezGD7FxpXhHXk0Qu\n6t3R8wZpVqmhnu0aoxQckeYVcQOSyEW9i7raA0OaVarvznaNATh8JsXKkYiGSBK5qFfG/GJOx2VI\ns0oN+Xq6ckczL85dyiTLWGjtcEQDI4lc1Kuj5w2YlaKFPARUY9K8IioiiVzUq7IHW1oESPt4TZU1\nr0TJw0HiDySRi3pjzC/mTFwGLZvopVnlJvjoXWgT5MW5i9K8IsqTRC7qTVmzSs+wxtYOxWbd2a4x\nCog6J80r4neSyEW9KWtWubOtJPKbdWfbxmiQsVdEeZLIRb3Izi3iTFwGrZrq8feWmYBuVlnzSvSl\nTDJypHlFlJJELurF4bNXMCtFr/ZNrB2KzesZFoACjpyTWrkoJYlc1Iv9py6j0UAvaR+/ZT3a+kvz\niihHErmocykZeVxIyqZ9S1+8dDIT0K3y1rkQ2tyb6IQsaV4RQBVzdpY5fvw4b7/9NmvWrCm3fNWq\nVWzcuBFfX18AFi5cSOvWrWs/SmHTDpwqfay8d3uZnPtmlM1veq2eYY05dymTqHNXuO/O5laISjQk\nVSbyTz/9lM2bN+Pmdv0NqlOnTrFs2TI6duxYJ8EJ26eU4sCpyzg7auke6m/tcOxGj1B/vvrxPIfP\nSiIX1UjkwcHBfPDBB7z00kvXrTt16hSffPIJBoOBgQMHMmvWrCpP6OPjDoC/v30/om2V8mk19Xbu\n6p7j/MUMUjLyuadrM4KDfCzL9TrXugqtVjT0+O5o5UenED9O/JaKxskRvxr2BLLn3z97LltFqkzk\nw4YNIyEh4YbrRo4cyZQpU9DpdMyePZudO3cyaNCgSo+XkZGHv78egyHn5iK2AdYqn69ZAZBex+eu\nSfk27/4NgO5tGpXbJ8dYUCex1Qa9zrVBxwdgMOTQpbUvJ35L5YefYxnas/q1cnv+/bP3slXkpm92\nKqWYPn06vr6+ODs7M2DAAE6fPn2zhxN2qKCohAOnU/D1dKFjq0bWDsfudG/bGI0GDp+VoW1vdzed\nyI1GI6NGjSI3NxelFAcPHpS2clHOoTNXKCwycU/nQLRXm31E7fHycKZdsA8xidmkZTXsbxCibtU4\nkW/ZsoWIiAj0ej3PP/8806ZNY8qUKdxxxx0MGDCgLmIUNmr3sUQ0GujXuam1Q7FbPctGRJSHg25r\n1ep+GBQUxPr16wEYPXq0Zfm4ceMYN25c3UQmbNrFlBxik3PoEtIIX8+GfePQFpV1SSwoLkEDRB5J\nwMXZgYFdm1k3MGEV8kCQqBO7jyUBMEASS51ydXakSSN3UrMKyMmTiZlvV9WqkQtRHWW1xKJiE3tP\nJOHu6kh6jrTd1rVWTT1JTssjNtk+e2uIqkmNXNS68wlZlJgU7Vr4yE3OehDcRIeDVsOFxCyUUtYO\nR1iBJHJRq8xmxdn4DBwdNIQGeVk7nNuCs6MDzRvryM4rllr5bUqaVkStik/JIa+ghLbB3jg7OQA3\nHitE1K7WgZ7EXc5h/8nLtA70tHY4op5JjVzUGqUUZ+IyAAhr4VPF1qI2Bfp54OrswMEzKZSYzNYO\nR9QzqZGLWnMlI5/UrAKaN9bh6SGTK9cnrVZDy6Z6zsZnsi4ymqDGOss66ZJo/6RGLmrNiZg0ADq2\n9rVyJLen1oGl9yRikrKtHImob5LIRa2IScoiOS2Ppo3cZU5OK2nk6YKXhzOXrhgpKjZZOxxRjySR\ni1rx7b44ADqFyOBY1qLRaGgV6InZrIhPkd4rtxNJ5OKWxV/O4XhMGo193Ajwkdq4NbVuWtpj5YI0\nr9xWJJGLW/bt/jgAOrVuhEYjDwBZk87diQAfN1LS8zHmF1s7HFFPpNeKuCllfcMzcwo5cs5AIy9X\nAv3crRyVgNI+5SkZ+cQmZUtT121CauTilvx6obSnSucQqY03FC2a6NFqNVxIypZH9m8TksjFTcvO\nLSIuOQcfvQtB/h7WDkdc5ezkQHBjHVm5RRgy860djqgHksjFTTt5IR0FdGrtK7XxBqZN89I+5dGX\nsqwciagPksjFTTHmFxOTlIWnhzPBTW6/Wcsbuia+7ujdnYi7nENegdz0tHeSyMVNOR2bjlKltXGt\n1MYbHI1Gwx3NvDCZFQdOy+TM9q5aifz48eNMnTr1uuU7duxg/PjxhIeHW6aCE/YvM6eQ6IQsPFwd\nadVURtprqEKaeaHRwJ5jSXLT085V2f3w008/ZfPmzbi5lX/Qo7i4mKVLl7Jx40bc3NyYPHkygwYN\nwt/fv86CFQ3D5r0xmMyKDq19ZeKIBszd1ZEgfx0XrxiJTc6R4W3tWJU18uDgYD744IPrlsfExBAc\nHIyXlxfOzs706NGDqKioOglSNBx5BcVs3ReLq7MDdzSTiSMaurbB3gBsP3LJypGIulRljXzYsGEk\nJCRct9xoNKLX/36Ty8PDA6PRWOUJfXxKHxrx97fvG2RWKd/V2nFdnjti+znyCkro06kpPl72+wCQ\nXudq7RBqRaiHC+cuZXL4zBWeHN+FRl6l36zt+ffPnstWkZt+slOn05Gbm2t5nZubWy6xVyQjIw9/\nfz0Gg/0O6mOt8vmaS9tB0+vo3IVFJjbtisHDzYkWAR7kGO1zYmW9ztWuyjaoWzP+9f05Nm4/z4P3\ntLbr3z97L1tFbrrXSkhICPHx8WRmZlJUVERUVBTdunW72cMJG7DneBLG/GJG9WuFs6ODtcMR1dSn\nQxM8XB3Z9UuiDG9rp2qcyLds2UJERAROTk7MmzePmTNnMmnSJMaPH09AQEBdxCgagBKTme8PXcTZ\nScvofq2tHY6oARcnBwZ2a4Yxv1i6ItqpajWtBAUFWboXjh492rJ88ODBDB48uG4iEw1C2eBY0QmZ\nZOQUEtbCh/2/Jls5KlFTg7sH8f3Bi3x38CLjBrWxdjiilskDQaJKZqU4eSEdrQY6tJJJlW2Rj96F\nfp2bkpKex56rf5yF/ZBELqoUfzmHnLxiQpp54e7qZO1wxE0a2bsFDloNET+ew2yWB4TsiSRyUSl1\ntTauQSZVtnV+3m707dSUREMuB89IW7k9kUQuKpVoyCUjp5CWTfXo3Z2tHY64RaP6lNbKN++Lw2Q2\nWzscUUtkhiBRIaUUJ2JKJ47o2FpmmrFVu/7QJt62hQ+nY9PZcyyJQd2DrBSVqE1SIxcVOncxk9Ss\nAoIa6/DRu1g7HFFLenVogpODlv/ujSVXhri1C5LIRYW27o8DSoeqFfbD3dWJTiG+GPOL2bIvztrh\niFogiVzcUGxyNqfiMmji646/t1vVOwibEtbSB39vVyKPJJCcllv1DqJBk0Qubmjr/ngAOoVIbdwe\nOWi1PDSoDSazYs0P52S8chsniVxcJzE1l6PnDbRq6kkTX/sd4fB21z3Uj653+HH2YiZ7T8jTurZM\nErm4zpZ9sUBpVzWZVNl+aTQaHhkaiquzAxE7fiMjp9DaIYmbJIlclHPpipFDZ64QHKCjaxs/a4cj\n6pivpysTB91BfmEJX26TJhZbJYlclLNp7wUAHryntdTGbxMDugbStrk3v0Snsu/Xy9YOR9wEeSBI\nWMQmZ/NLdCohzTzpJA8A2b1rHxTq0MqXC8nZfLX9PKHB3jSWnko2RWrkwuK/e67WxvtLbfx2o3N3\nolf7xhQWmfhsy2l5fN/GSCIXABz7LZWTsemEtfAhrKV0ObwdtWrqyV1hjfktMYtvfoqzdjiiBiSR\nC4pLTKzdfh6tRsOUITLpwO1Ko9EwbVhb/Lxc+fbnOH69kGbtkEQ1SSIXfH/wIobMAobcGUQzf521\nwxFW5O7qxJ8f6Iijg4ZPt5wmPdt+JqG2Z3Kz8zaXmpnP1v3xeHo44+vlct1IeeL2Uvb+92jrz8HT\nV1j276P87fHeODpIna8hqzKRm81mXn/9dc6dO4ezszOLFy+mRYsWlvWLFy/m6NGjeHh4ALBy5Ur0\nen3dRSxqjdms+Ozb0xSVmJk+/A4KZYZ1cVVoc29S0vOJu5zD17tjCB8sTW4NWZWJfPv27RQVFRER\nEcGxY8d44403+Mc//mFZf+rUKT777DN8feUGma3534F4zidkcWdbf3q3D2D38SRrhyQaCI1GQ5+O\nTUjPKeSHQ5coLDYRHFC+gjawazMrRSf+qMrvS0eOHKF///4AdO3alZMnT1rWmc1m4uPjefXVV5k0\naRIbN26su0hFrYpNzuabn2Lx0bsw7f520t1QXMfJUcuAroE4aDXs+/Uy2blF1g5JVKDKGrnRaESn\n+/0GmIODAyUlJTg6OpKXl8cjjzzCo48+islkYtq0aXTs2JF27dpVeDwfn9JBmPz97bv5xSrl02qq\nde4sYyHvbziOyazo37UZ5xKzAdDrXKt9qppsa4ukfL9vN6hHENsPX2L3sSQmDG6Ds5MD0HB/hxtq\nXHWpykSu0+nIzf19vGKz2YyjY+lubm5uTJs2DTe30qfAevfuzdmzZytN5BkZefj76zEYcm419gbL\nWuXzvTozenol5y4xmXl77S/k5BXT5Y5GeHs4kWOsWc8Evc61xvvYEilfeYGN3Alr4cOZ+Ay++zmW\nQd2bodFoGuTvsD3nlsr+QFXZtNK9e3f27NkDwLFjxwgNDbWsi4uLY8qUKZhMJoqLizl69CgdOnSo\nhZBFXVBK8eW2c5xPyKJFEz2dQ+QxfFE9Pdr606SROwmGXI7/Jv3LG5oqa+T33Xcf+/btY9KkSSil\nWLJkCatWrSI4OJh7772X0aNH89BDD+Hk5MTYsWNp00bubjdUm/fFsed4Mi0C9PTt1ETaxUW1abUa\n7ukSyP/2x3MiJk3mcG1gNKqex600GHLs+usPWLFppUdHANKP/H5Duqxf8Jm4DA6fvYLOzYn7ewXj\n7nrzjxBI04Ntu5XyZeQU8N2BiwC8OqMnQQ3sATJ7zi231LQibN9vCVkcPnsFNxcH7usZdEtJXNze\nfPSu9O3UlBKTYsXGE9KTpYGQRG7nzl/K5OeTl3F20jLkziD07s7WDknYuLL7K6lZBaz4+oQ8SNYA\nSCK3Yz9GXeLAqRRcnR0YdldzfPT23aVO1J8udzSiT4cmXEjK5tMtpzGbZWYha5JEbqe+OxjP2u3R\nuLk4MFSSuKhlGo2GR0e0o12wN0fPG/jXDzJNnDVJIrdDm/fFsmFnDD56F4bdFYy3TnoYiNrn6KBl\n9oOdCQ7Qsed4EhE7fpNkbiWSyO2ICQ2fdXmQTXtj8fNyZd7D3fH0kDZxUXfcXR15IbwrTRu5s+3w\nJf6z54IkcyuQ7gt2orDIxLK7Z3KoWWea+Xnw/ENd8PWU5hRRd64d8rhvpyb8cOgSW/fHk1dYwsP3\nhaKV5xTqjSRyO5CVW8SKjceJbdaZjsln6XjvKE7I7C6iHrm7lj6fsP9kCjuPJpKbX8zMkWE4OTpY\nO7TbgjSt2Lik1Fz+9q8oYpNzGBCzn3m7PrQMaiREfXJzcWTuw924o5kXh85cYemXR2WGoXoiidyG\nnYhJZcmaI6RmFTCuXyuePLAGR7P06RXW4+HqxJzJXenbqQlxl3NYuPowJ+XbYZ2TRG6DzErxzU+x\nLN9wgqISM38aFcaYfq2QFknREDg5OvDYiDAevi+UvIIS3l1/nH/9cI6CohJrh2a3pI3cxmTkFPL5\n1tOcjsvAw9WRgd2aUVRiZtexREZYOzhx27v2BqiDg4b7ewez70Qyu35J5PCZFO5s15jgAB0ajUZm\nGKpFkshthFKKqHMG/vX9WXILSgjy9+DuTk1xdZb2cNFwNfJ0ZeTdLTj+WxqnY9PZfSyJpo3c6R7q\nb+3Q7IokchuQmpnPlz+e50RMGs6OWqYODQUNMgytsAkOWi3dQ/0JCfTi8NkrJKXmsnV/PJfT8xjT\ntxXNGzesERRtkSTyBsyYX8x3B+OJjEqgqMRMWAsfpg5rSxNf93JfYYWwBV46Z+7t0YzktDyORady\n5JyBI+cMdGzly9CezWnfylf6nt8kSeQNUJaxkB1HE9l+JIH8whJ89C5MHxhC7/YBUgsXNk2j0RDo\n50HTRu746l35/uBFTsamczI2HT8vV/p1akqvDgEEXJ3bV1SPJPIGQinF+l2/8VtCFrFJOZiVwsXJ\ngTvb+hMa7E1hsYndx5OsHaYQtUKj0ZBhLKRXhwBCgrw4fymTuORsNv0Uy6afYmneWEfXO/xo39KH\n1oFeODlKB7vKSCK3IqUUcZdzOHrewOGzV7iSkQ+A3t1JPsDituHn5YqfVxN6tmvMxZQc4i7nkGgw\ncumKkS0/x+HooKFtsA/tW/jQprk3LQJ08sToH0gir0fFJWYSU43EJedw7lImZy9mkGUsnWHF2VFL\nq6Z67gjyoomvuzShiNuOk6OWkGZehDTzoqjYxOX0PC6n5ZGcnsep2HROxaYD4KDV0LyxjtaBnlf/\neRHg43Zb/85UmcjNZjOvv/46586dw9nZmcWLF9OiRQvL+vXr17Nu3TocHR156qmnGDRoUJ0G3NCY\nzYq8whLyCkvILyghr6CY6OQcklKySc8uJD27gLTsAtKzC0nLLsB0zQD8nh7O9OkQQPdQfzq2asT+\n05etWBIhGg5nJweCA/QEB5TOU5lfWMLltDwMmfmkZhVYau47jpbe9PdwdSQ4QE9IkDc+Hk6l7fB+\nHujdnG6LBF9lIt++fTtFRUVERERw7Ngx3njjDf7xj38AYDAYWLNmDV9//TWFhYVMmTKFvn374uxc\nN0On5hYUU1JixqxKmyXMSv3+s1mhVOlTj0qVJljz1eUms8JkMmNSCpPp92UlZnPpzyZlWVdcYqbY\nZKa4xEzJ1f9L/5nILzSRV1DMlcx8ikrMFBeXblsdbi4O+Hq64KN3xdfThcY+bnh5OKPRaMjJL5Yk\nLkQl3FwcaRXoSatATwBMJjPpOYWkZhZgyMonNbOAM/EZnInPKLefs6MWb50LXjpnvHQueHk44+bi\ngIvT7/+cnLRoNZrSf9qy/0Gr0aDRaNBqQKu9+vM16zWUvcayXKPVoKX89hoNaCjthllX8+VWedQj\nR47Qv39/ALp27crJk7/P0H7ixAm6deuGs7Mzzs7OBAcHc/bsWTp37lzrge44msCX287X+nFrSkPp\nV0BnJwd07k44O2lxdnQo97/ewwVlNuPh6oSHmyPuro44aKWtW4ja4uCgxd/bDX9vN8LwAUqbLkvM\nkGTIIdNYRFZuEXkFxeTkF2HIzKchjJI+c2QYfTs1rfXjVpnIjUYjOt3vHfYdHBwoKSnB0dERo9GI\nXq+3rPPw8MBoNFZ6PH9/fbn/qyt8WBjhw8JqtM9t53Lp18yJVg5DCFG/qqwm6nQ6cnNzLa/NZjOO\njo43XJebm1susQshhKh7VSby7t27s2fPHgCOHTtGaGioZV3nzp05cuQIhYWF5OTkEBMTU269EEKI\nuqdRVUywV9Zr5fz58yilWLJkCXv27CE4OJh7772X9evXExERgVKKWbNmMWzYsPqKXQghBNVI5EII\nIRo26UohhBA2ThK5EELYOEnkQghh4+olkRcUFPDMM88wZcoUHn/8cdLT06/b5j//+Q8TJ07kwQcf\n5MMPP6yPsGpNdcoHkJ+fz9ixYy29gGxFdcq3bNkywsPDGT9+POvXr7dClDVjNpt59dVXCQ8PZ+rU\nqcTHx5dbv379eh588EEeeughdu7caaUob15V5Vu9ejUTJ05k4sSJ/P3vf7dSlDevqvKVbfOnP/2J\ntWvXWiHC+lUviXzt2rWEhoby73//m3HjxrFy5cpy6y9evMjatWtZs2YNGzdupLi4mOLi4voIrVZU\nVb4yixYtsslxH6oq34EDB7h48SIRERGsXbuWTz/9lKysLCtFWz3XDj3xl7/8hTfeeMOyrmzoiXXr\n1vH555/z7rvvUlRUZMVoa66y8l26dInNmzezbt06IiIi+Omnnzh79qwVo625yspX5v3332/wn8Pa\nUi+J/NrH/O+55x72799fbv3PP/9Mx44dmTt3Lo888gjdu3fHycmpPkKrFVWVD+Dzzz+nW7dutGvX\nrr7Du2VVla9bt24sWbLE8tpkMlkeGmuoqjv0hF6vtww9YUsqK1+TJk347LPPcHBwQKvVUlJSgouL\ni7VCvSmVlQ/g+++/R6PRcM8991gjvHpX679tGzZs4J///Ge5ZY0aNbI88enh4UFOTk659RkZGURF\nRbF27VoKCwuZPHkyGzduxNPTs7bDu2U3U779+/cTHx/PokWLOHr0aL3FejNupnwuLi64uLhQXFzM\nvHnzCA8Px8PDo95ivhm1PfREQ1NZ+ZycnPD19UUpxZtvvkn79u1p1aqVFaOtucrKd/78eb799ltW\nrFhhc820N6vWE3lZu9u1Zs+ebXmUPzc397oE7e3tzV133YVOp0On0xESEkJcXFydDL51q26mfBs3\nbiQxMZGpU6dy4cIFTp06hb+/P2FhDW/smJspH0BWVhbPPvssd911F7NmzaqXWG+FvQ89UVn5AAoL\nC3n55Zfx8PDgtddes0aIt6Sy8m3atImUlBSmT59OYmIiTk5ONGvWzK5r5/XStNK9e3d2794NwJ49\ne+jRo8d16w8dOkRhYSF5eXnExMQQHBxcH6HViqrK984777Bu3TrWrFlD//79mTNnToNM4hWpqnwF\nBQXMmDGD8ePH8/TTT1sjxBqz96EnKiufUoo///nPtG3blkWLFuHgYHuz7VRWvpdeeokNGzawZs0a\nHnjgAWbMmGHXSRzq6cnO/Px85s6di8FgwMnJiXfeeQd/f3/efPNN7r//fjp37szq1avZvHkzSimm\nT5/OuHHj6jqsWlOd8pWZN28eI0aMsKkPVlXlO3r0KH//+9/L/XFasmQJzZs3t2LUlbP3oScqK5/Z\nbOaFF16ga9eulu1feOEFunXrZsWIa6aq96/MBx98gJ+fH5MnT7ZitHVPHtEXQggbJw8ECSGEjZNE\nLoQQNk4SuRBC2DhJ5EIIYeMkkQshhI2TRC4alLZt2wKl47uUDXa0fv16+vfvz7Jly9i9ezf9+/fn\nL3/5izXDrNDYsWNvet+pU6dy8ODBWoxG3C4a9oAY4rZ1bb/fb7/9lqVLl9KvXz/mz5/P7NmzCQ8P\nt2J0Ffvmm2+sHYK4DUkiFxw8eJCPPvoIJycnEhISGDx4MO7u7mzfvh2ATz75BD8/P/bs2cOKFSso\nKSkhKCiIv/71r/j4+PDdd9+xatUqCgoKKCoqYsmSJXTv3p2pU6fSqVMnjhw5Qnp6OgsWLGDAgAHl\nzp2QkMCcOXPIy8ujS5culuUffPABABqNhl9//ZWFCxcydepUIiMj2b9/P1qtlrvuuovXX3+dzMxM\nXF1deeWVV2jfvj3z5s0jMzOT+Ph45syZg5+fH0uXLqWgoAAfHx8WLlxI8+bNK4wvMTGR+fPnk56e\njqurK4sXL6Zdu3Zs2rSJf/7zn5jNZjp06MBrr7123WBTbdu25dy5c3zwwQekpKQQHx9PYmIiEydO\n5Kmnniq3bVFREf/3f//HyZMnadasGRkZGZZ1n3zyCd999x0mk4l+/foxZ84cNBoN//rXv/jyyy/R\n6/W0bt2a4OBgnnnmGXr37k3Hjh0xGAxs3LiRVatW3XD/G5VBq9Xy8ssvEx0dDcCUKVN46KGHau8D\nJuqeEre9AwcOqG7duqmkpCSVl5enunbtqtauXauUUmrevHlq9erVKi0tTY0ZM0ZlZmYqpZRau3at\nevnll5XJZFLTpk1TaWlpSimlNmzYoGbNmqWUUuqRRx5RixcvVkopFRkZqR544IHrzv3EE0+o9evX\nK6WU+u9//6tCQ0OVUkqtWLFCrVixwnKcAwcOKKWUmjt3rvr666+VUkqFh4erU6dOKaWUio6OVkOH\nDrVsM3fuXKWUUoWFhWr06NEqMTFRKaXUnj171PTp0yuN7/HHH1dffvmlUkqpXbt2qWeffVadP39e\nTZ48WRUUFCillHr77bfVhx9+eF15ro1/woQJqrCwUKWmpqquXbuqrKysctt+9tln6sUXX1RKKRUb\nG6s6deqkDhw4oHbv3q2eeeYZVVJSokwmk3rhhRfUpk2b1JkzZ9TQoUNVTk6OKigoUBMnTrRco9DQ\nUMs1qmj/ispw8OBB9fjjjyullLp8+bKaM2fOdeUSDZvUyAUAoaGhNG3aFAAfHx/69OkDQGBgINnZ\n2Rw/fpzk5GSmTZsGlD4i7eXlhVar5cMPP2THjh3ExsZy6NAhtNrfb72UDTXapk0bMjMzrzvvoUOH\neOeddwAYM2YMCxYsqFgM1zoAAAO3SURBVFa8ubm5nDx5kvnz51uW5eXlWWq1ZcMixMXFcenSpXK1\n4WtHMrxRfIcPH+bdd98FYMCAAQwYMIAvv/yS+Ph4S021uLiY9u3bVxpjr169cHZ2plGjRnh7e5OT\nk1NuwLFDhw5ZmohatmxpeUR+//79nDhxggcffBAoHcsmMDCQ9PR0Bg0aZBn1b+TIkWRnZ1uOV/aN\npqL9c3JybliGyZMnExsby8yZM7nnnnt46aWXKi2XaHgkkQuA68Z//+NASiaTie7du/PRRx8BpaPn\n5ebmkpuby4QJExgzZgw9e/akbdu2fPXVV5b9ypoeKptQQ10dJUKj0ZT7I1AZs9mMs7NzuTbpy5cv\n4+3tDYCrq6tlu6CgIMt2JpOJ1NTUSuO7dpRApRQxMTGYTCaGDx9u+UOTm5uLyWSqNMZrm100Go2l\nnBUtKzuvyWRi+vTpPProowBkZ2fj4ODAxo0bMZvNFZ6vrMyV7X+jMnh6erJ161b27dvH7t27eeCB\nB9i6dWuDHEZa3Jj0WhHV0qVLF44dO0ZsbCwAK1eu5M033yQuLg6NRsOTTz5Jr169+PHHH6tMcNe6\n++672bx5MwDbtm2jsLCwWvvp9XpatmxpSdD79u3j4Ycfvm671q1bk5WVRVRUFABff/01L774YqXH\nvvPOO9m6dStQOunJK6+8YilbWloaSilef/3168Ztr6k+ffqwZcsWzGYziYmJlrHqe/fuzTfffENu\nbi4lJSU8/fTT/PDDD/Tp04fdu3djNBopKipi27ZtN/wDWdH+FZUhMjKSOXPmMHDgQBYsWIC7uzvJ\nycm3VDZRv6RGLqrF39+fJUuW8Nxzz2E2mwkICOCtt97C09OTsLAwhg8fjkajoV+/fhw5cqTax331\n1VeZM2cOERERdOzYsUYTUrz11lu8/vrrfPbZZzg5OfHee+9dl9icnZ1Zvnw5f/vb3ygsLESn07Fs\n2bIqY1qwYAH//ve/cXNzY/Hixdxxxx3Mnj2b6dOnYzabCQsL44knnqh2rDcyZcoUoqOjGT58OM2a\nNbMMxTp48GDOnj3LQw89hMlkon///jzwwANoNBqmTZtGeHg47u7u+Pj43HBmn8r2v1EZtFot27Zt\nY+TIkbi4uDBmzBhLN1BhG2T0QyFsRGxsLLt372bGjBkAPPXUU/z/9uzYBoAQhgEge2QaqNMyOQvR\n/Qq8pbsJXFmOsvcea623wXjOIocQVTXOOaO7v+tnzvk6Fj9gkQOE8+wECKfIAcIpcoBwihwgnCIH\nCHcBxDdxh5ckO+gAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x224fa335f60>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot the permutation replicates\n",
    "sns.distplot(permutation_replicates)\n",
    "plt.xlabel(\"mean difference in degrees\")\n",
    "\n",
    "# draw a line showing the empirical difference of means. the area to the left is the probability\n",
    "plt.vlines(ymin=0, ymax=4,\n",
    "           x=empirical_diff_means, \n",
    "           color=\"red\")\n",
    "\n",
    "plt.title(\"Permutations of Mean Temperature Differences\")\n",
    "\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The probability of observing a difference in mean temps at least as extreme as the empirical difference of -0.289 is given by the ratio of permutation replicates that meet that extreme value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0121\n"
     ]
    }
   ],
   "source": [
    "p_probability = np.sum(permutation_replicates <= empirical_diff_means) / len(permutation_replicates)\n",
    "print(p_probability)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As the number of permutation replicates increases, the distribution approaches the normal distribution. This means that the probability or p-value can also be ascertained by looking up a z-score representing the number of standard deviations away the test statistic is from the mean of zero."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0120426374914\n"
     ]
    }
   ],
   "source": [
    "# get standard deviation of permutation replicates\n",
    "perm_std = np.std(permutation_replicates)\n",
    "\n",
    "# convert empirical difference to standard deviations\n",
    "z_score = empirical_diff_means / perm_std\n",
    "\n",
    "# look up p value for z_score\n",
    "p_lookup = 1 - stats.norm.sf(z_score)\n",
    "\n",
    "print(p_lookup)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These P values indicate that were there no difference between male and female temperatures, the observation in the original sample would be 1.26% or 1.32% likely to occur. Given a 95% confidence level, the null hypothesis that there is no difference cannot be retained. The significance level of 0.05 is greater than the p value of 0.126"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Conclusion\n",
    "\n",
    "The original problem asked if the commonly accepted normal body temperature of 98.6 degrees F reported by Carl Wunderlich was a statistically correct value. In tests against this sample of 130 individuals, evidence was not found that corroborates 98.6 degrees as a global mean, although it was still found to be in the range of normality. Furthermore, a significant difference between male and female body temperatures was discovered, further challenging the simplicity of a universal figure of normality. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "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.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}