GB213-review-in-Python.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# GB213 Review in Python\n",
    "\n",
    "## We're not covering everything\n",
    "\n",
    "We're omitting basic probability issues like experiments, sample spaces, discrete probabilities, combinations, and permutations.  At the end we'll provide links to example topics.  But everything else we'll cover at least briefly.\n",
    " \n",
    "We begin by importing the necessary modules."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import scipy.stats as stats\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Discrete Random Variables\n",
    "\n",
    "(For continuous random variables, [see further below](#continuous-random-variables).\n",
    "\n",
    "Discrete random variables taken on a finite number of different values.  For example, a Bernoulli trial is either 0 or 1 (usually meaning failure and success, respectively).  You can create random variables using `scipy.stats` as follows.\n",
    "\n",
    "### Creating them"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "b1 = stats.bernoulli( 0.25 )  # probability of success\n",
    "b2 = stats.binom( 10, 0.5 )   # number of trials, prob. of success on each"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Computing probabilities from a Discrete Random Variable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.75, 0.25)"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b1.pmf( 0 ), b1.pmf( 1 )  # stands for \"probability mass function\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The same code works for any random variable, not just `b1`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generating values from a Discrete Random Variable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0, 1, 1, 0, 0, 1, 0, 0, 0, 0])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b1.rvs( 10 )  # asks for 10 random values (rvs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The same code works for any random variable, not just `b1`.\n",
    "\n",
    "### Computing statistics about a Discrete Random Variable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.25, 0.1875, 0.4330127018922193)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b1.mean(), b1.var(), b1.std()  # mean, variance, standard deviation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The same code works for any random variable, not just `b1`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plotting the Distribution of a Discrete Random Variable\n",
    "\n",
    "Here's a function you can use to plot (almost) any discrete probability distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_discrete_distribution ( rv ):\n",
    "    xmin, xmax = rv.ppf( 0.0001 ), rv.ppf( 0.9999 )\n",
    "    x = np.arange( xmin, xmax+1 )\n",
    "    y = rv.pmf( x )\n",
    "    plt.plot( x, y, 'o' )\n",
    "    plt.vlines( x, 0, y )\n",
    "    plt.ylim( bottom=0 )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example use:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAASjUlEQVR4nO3df2xd533f8fdnlJywDVo5sTZMlB0pqKvGm7ZwUJ1uwZJhcSIHK2xBSFAFyOAOBowV9dY5mwZrBWzA/UNJNWzeH95moXEXbG1d1xE0YVDLGnG0/VE4kxylUWWXiKxmNslsVuHI+xHOlpTv/uCVd81S4aF5ySs+fL8Aguc85znnfh/cez88eu65R6kqJEnt+nPDLkCStLIMeklqnEEvSY0z6CWpcQa9JDVuw7ALmO+mm26qbdu2DbsMSVpTnn/++T+tqs0Lbbvugn7btm2cOnVq2GVI0pqS5L9da5tTN5LUOINekhrXKeiT3JlkMsm5JA8usP3zSV5I8q0kX03y/r5tV5J8s/dzbJDFS5IWt+gcfZIR4DHgE8AUcDLJsap6oa/baWBXVX0/yS8Avwr8XG/bbFV9aMB1S5I66nJGfztwrqrOV9WbwJPA3f0dquprVfX93upzwNbBlilJeqe6BP0Y8Erf+lSv7VruBX63b/3dSU4leS7JnoV2SHJfr8+pCxcudChJktRVl8srs0Dbgre8TPI5YBfwsb7mW6pqJskHgGeTnKmql952sKrDwGGAXbt2eTtNXXeOnp7m0MQkMxdn2bJplP27d7Bn/Ied70jXjy5BPwXc3Le+FZiZ3ynJHcAvAx+rqjeutlfVTO/3+SQngHHgpfn7S9ero6enOXDkDLOXrgAwfXGWA0fOABj2WhO6TN2cBG5Nsj3JDcA+4G1XzyQZBx4H7qqqV/vab0zyrt7yTcBHgP4PcaXr3qGJybdC/qrZS1c4NDE5pIqkpVn0jL6qLie5H5gARoAnqupskkeAU1V1DDgEvAf4nSQAL1fVXcAHgceT/IC5PypfmHe1jnTdm7k4u6R26XrT6RYIVXUcOD6v7aG+5Tuusd8fADuXU6A0bFs2jTK9QKhv2TQ6hGqkpfObsdIi9u/ewejGkbe1jW4cYf/uHUOqSFoag15axJ7xMQ7u3cnIG69DFWObRjm4d6cfxGrNuO7uXildj/aMj/HoA4cBOHHixHCLkZbIM3pJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mN2zDsAqSujp6e5tDEJDMXZ9myaZT9u3ewZ3xs2GWtqPU4Zg2eQa814ejpaQ4cOcPspSsATF+c5cCRMwDNBt96HLNWhlM3WhMOTUy+FXhXzV66wqGJySFVtPLW45i1Mgx6rQkzF2eX1N6C9ThmrYxOQZ/kziSTSc4leXCB7Z9P8kKSbyX5apL39227J8m3ez/3DLJ4rR9bNo0uqb0F63HMWhmLBn2SEeAx4FPAbcBnk9w2r9tpYFdV/RXgaeBXe/u+F3gY+DBwO/BwkhsHV77Wi/27dzC6ceRtbaMbR9i/e8eQKlp563HMWhldzuhvB85V1fmqehN4Eri7v0NVfa2qvt9bfQ7Y2lveDTxTVa9V1feAZ4A7B1O61pM942Mc3LuTkTdehyrGNo1ycO/Opj+UXI9j1sroctXNGPBK3/oUc2fo13Iv8Ls/ZF9fpXpH9oyP8egDhwE4ceLEcItZJetxzBq8LkGfBdpqwY7J54BdwMeWsm+S+4D7AG655ZYOJUmSuuoydTMF3Ny3vhWYmd8pyR3ALwN3VdUbS9m3qg5X1a6q2rV58+autUuSOugS9CeBW5NsT3IDsA841t8hyTjwOHMh/2rfpgngk0lu7H0I+8lemyRplSw6dVNVl5Pcz1xAjwBPVNXZJI8Ap6rqGHAIeA/wO0kAXq6qu6rqtSS/wtwfC4BHquq1FRmJJGlBnW6BUFXHgePz2h7qW77jh+z7BPDEOy1QkrQ8fjNWkhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIa1ynok9yZZDLJuSQPLrD9o0m+keRykk/P23YlyTd7P8cGVbgkqZsNi3VIMgI8BnwCmAJOJjlWVS/0dXsZ+HngnyxwiNmq+tAAapUkvQOLBj1wO3Cuqs4DJHkSuBt4K+ir6ju9bT9YgRolScvQZepmDHilb32q19bVu5OcSvJckj1Lqk6StGxdzuizQFst4TFuqaqZJB8Ank1ypqpeetsDJPcB9wHccsstSzi0JGkxXc7op4Cb+9a3AjNdH6CqZnq/zwMngPEF+hyuql1VtWvz5s1dDy1J6qBL0J8Ebk2yPckNwD6g09UzSW5M8q7e8k3AR+ib25ckrbxFg76qLgP3AxPAi8BTVXU2ySNJ7gJI8tNJpoDPAI8nOdvb/YPAqSR/CHwN+MK8q3UkSSusyxw9VXUcOD6v7aG+5ZPMTenM3+8PgJ3LrFGStAx+M1aSGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqXKf/M1bqd/T0NIcmJpm5OMuWTaPs372DPeNjwy5LA+bz3A6DXkty9PQ0B46cYfbSFQCmL85y4MgZAEOgIT7PbXHqRktyaGLyrTf/VbOXrnBoYnJIFWkl+Dy3xaDXksxcnF1Su9Ymn+e2GPRaki2bRpfUrrXJ57ktBr2WZP/uHYxuHHlb2+jGEfbv3jGkirQSfJ7bYtBrSfaMj3Fw705G3ngdqhjbNMrBvTv9gK4xPs9t8aobLdme8TEefeAwACdOnBhuMVoxPs/t8Ixekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNa5T0Ce5M8lkknNJHlxg+0eTfCPJ5SSfnrftniTf7v3cM6jCJUndLBr0SUaAx4BPAbcBn01y27xuLwM/D/zmvH3fCzwMfBi4HXg4yY3LL1uS1FWXM/rbgXNVdb6q3gSeBO7u71BV36mqbwE/mLfvbuCZqnqtqr4HPAPcOYC6JUkddQn6MeCVvvWpXlsXy9lXkjQAXYI+C7RVx+N32jfJfUlOJTl14cKFjoeWJHXRJeingJv71rcCMx2P32nfqjpcVbuqatfmzZs7HlqS1EWXoD8J3Jpke5IbgH3AsY7HnwA+meTG3oewn+y1SZJWyaJBX1WXgfuZC+gXgaeq6mySR5LcBZDkp5NMAZ8BHk9ytrfva8CvMPfH4iTwSK9NkrRKOt2muKqOA8fntT3Ut3ySuWmZhfZ9AnhiGTVKkpbBb8ZKUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJalynoE9yZ5LJJOeSPLjA9ncl+e3e9q8n2dZr35ZkNsk3ez//drDlS5IWs2GxDklGgMeATwBTwMkkx6rqhb5u9wLfq6qfSLIP+CLwc71tL1XVhwZctySpoy5n9LcD56rqfFW9CTwJ3D2vz93Al3vLTwMfT5LBlSlJeqe6BP0Y8Erf+lSvbcE+VXUZeB14X2/b9iSnk/znJH9zoQdIcl+SU0lOXbhwYUkDkCT9cF2CfqEz8+rY57vALVU1Dnwe+M0kP/ZnOlYdrqpdVbVr8+bNHUqSJHW16Bw9c2fwN/etbwVmrtFnKskG4MeB16qqgDcAqur5JC8BPwmcWm7hgqOnpzk0McnMxVm2bBpl/+4d7Bmf/48taW3xdT14Xc7oTwK3Jtme5AZgH3BsXp9jwD295U8Dz1ZVJdnc+zCXJB8AbgXOD6b09e3o6WkOHDnD9MVZCpi+OMuBI2c4enp62KVJ75iv65WxaND35tzvByaAF4GnqupskkeS3NXr9iXgfUnOMTdFc/USzI8C30ryh8x9SPv3q+q1QQ9iPTo0McnspStva5u9dIVDE5NDqkhaPl/XK6PL1A1VdRw4Pq/tob7l/wt8ZoH9vgJ8ZZk1agEzF2eX1C6tBb6uV4bfjF2jtmwaXVK7tBb4ul4ZBv0atX/3DkY3jrytbXTjCPt37xhSRdLy+bpeGQb9GrVnfIyDe3cy8sbrUMXYplEO7t3p1Qla03xdr4xOc/S6Pu0ZH+PRBw4DcOLEieEWIw2Ir+vB84xekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEbhl3AWnf09DSHJiaZuTjLlk2j7N+9gz3jY8MuS9I70Or72aBfhqOnpzlw5Ayzl64AMH1xlgNHzgA08eKQ1pOW389O3SzDoYnJt14UV81eusKhickhVSTpnWr5/WzQL8PMxdkltUu6frX8fjbol2HLptEltUu6frX8fjbol2H/7h2Mbhx5W9voxhH2794xpIokvVMtv58N+mXYMz7Gwb07GXnjdahibNMoB/fuXPMf3EjrUcvv52auuhnWZVF7xsd49IHDAJw4cWLFH0/SyhnW+3ml86vTGX2SO5NMJjmX5MEFtr8ryW/3tn89yba+bQd67ZNJdg+s8j5XL4uavjhL8f8vizp6enolHk6SBmY18mvRoE8yAjwGfAq4DfhsktvmdbsX+F5V/QTwL4Ev9va9DdgH/CXgTuBf9443UC1fFiWpbauRX13O6G8HzlXV+ap6E3gSuHten7uBL/eWnwY+niS99ier6o2q+hPgXO94A9XyZVGS2rYa+dVljn4MeKVvfQr48LX6VNXlJK8D7+u1Pzdv3z8z8ZTkPuC+3ur/TrKkP2UbN2/bmZENN8xvryuX38wXf/bMUo61HHN/24bipiR/OowHHtKYhzZecMyraF2MeYD59f5rbegS9AuNtjr26bIvVXUYONyhlkUlOVVVuwZxrLVivY15vY0XHPN6sVJj7jJ1MwXc3Le+FZi5Vp8kG4AfB17ruK8kaQV1CfqTwK1Jtie5gbkPV4/N63MMuKe3/Gng2aqqXvu+3lU524Fbgf86mNIlSV0sOnXTm3O/H5gARoAnqupskkeAU1V1DPgS8O+TnGPuTH5fb9+zSZ4CXgAuA79YVVcWfKDBGcgU0Bqz3sa83sYLjnm9WJExZ+7EW5LUKm+BIEmNM+glqXHNBP1it2loTZKbk3wtyYtJzib5pWHXtFqSjCQ5neQ/DbuW1ZBkU5Knk/xx7/n+68OuaaUleaD3uv6jJL+V5N3DrmnQkjyR5NUkf9TX9t4kzyT5du/3jYN4rCaCvuNtGlpzGfjHVfVB4GeAX1wHY77ql4AXh13EKvpXwO9V1U8Bf5XGx55kDPiHwK6q+svMXQSyb7hVrYh/x9ytYfo9CHy1qm4FvtpbX7Ymgp5ut2loSlV9t6q+0Vv+X8y9+df+/VQXkWQr8HeAXxt2LashyY8BH2Xuyjaq6s2qujjcqlbFBmC0972cH6HB799U1X9h7irFfv23k/kysGcQj9VK0C90m4bmQ++q3t1Cx4GvD7eSVfEo8E+BHwy7kFXyAeAC8Ou96apfS/Kjwy5qJVXVNPDPgZeB7wKvV9XvD7eqVfMXquq7MHcyB/z5QRy0laDvdKuFFiV5D/AV4B9V1f8cdj0rKcnPAq9W1fPDrmUVbQD+GvBvqmoc+D8M6J/z16vevPTdwHZgC/CjST433KrWtlaCfl3eaiHJRuZC/jeq6siw61kFHwHuSvId5qbn/naS/zDcklbcFDBVVVf/tfY0c8HfsjuAP6mqC1V1CTgC/I0h17Ra/keSvwjQ+/3qIA7aStB3uU1DU3q3gf4S8GJV/Yth17MaqupAVW2tqm3MPcfPVlXTZ3pV9d+BV5Jc/Y9LP87cN81b9jLwM0l+pPc6/ziNfwDdp/92MvcA/3EQB23ivxK81m0ahlzWSvsI8HeBM0m+2Wv7Z1V1fIg1aWX8A+A3eicx54G/N+R6VlRVfT3J08A3mLu67DQN3g4hyW8Bf4u52zFPAQ8DXwCeSnIvc3/wPjOQx/IWCJLUtlambiRJ12DQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMb9P4/uFQ+8nHS8AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_discrete_distribution( b2 )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Continuous Random Variables\n",
    "\n",
    "(For discrete random variables, [see further above](#discrete-random-variables).\n",
    "\n",
    "Continuous random variables take on an infinite number of different values, sometimes in a certain range (like the uniform distribution on $[0,1]$, for example) and sometimes over the whole real number line (like the normal distribution, for example).\n",
    "\n",
    "### Creating them"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for uniform on the interval [a,b]: use loc=a, scale=b-a\n",
    "u = stats.uniform( loc=10, scale=2 )\n",
    "# for normal use loc=mean, scale=standard deviation\n",
    "n = stats.norm( loc=100, scale=5 )\n",
    "# for t, same as normal, plus df=degrees of freedom\n",
    "t = stats.t( df=15, loc=100, scale=5 )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Computing probabilities from a Continuous Random Variable\n",
    "\n",
    "For a continuous random variable, you cannot compute the probability that it will equal a precise number, because such a probability is always zero.  But you can compute the probability that the value falls within a certain interval on the number line.\n",
    "\n",
    "To do so for an interval $[a,b]$, compute the total probability accumulated up to $a$ and subtract it from that up to $b$, as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.3413447460685429"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a, b = 95, 100  # or any values\n",
    "n.cdf( b ) - n.cdf( a )  # probability of being in that interval"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The same code works for any continuous random variable, not just `n`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generating values from a Continuous Random Variable"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([100.86564343,  99.08962212, 101.15723854,  97.48893396,\n",
       "       107.10855389,  99.65035029, 104.57998766,  98.63143601,\n",
       "       101.24020389, 106.12612087])"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "n.rvs( 10 )  # same as for discrete random variables"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plotting the Distribution of a Continuous Random Variable\n",
    "\n",
    "Here's a function you can use to plot the center 99.98% of any continuous probability distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_continuous_distribution ( rv ):\n",
    "    xmin, xmax = rv.ppf( 0.0001 ), rv.ppf( 0.9999 )\n",
    "    x = np.linspace( xmin, xmax, 100 )\n",
    "    y = rv.pdf( x )\n",
    "    plt.plot( x, y )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Example use:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAD4CAYAAADlwTGnAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXxV9Z3/8dcnNxsJJJAVCAECYQs7REBE3Bcqiygqtlbb0moXp7+p019Hp2Nrnen8xum0ztQ609Haat1RUaOguFB3WcIWdglbEhJIQiCQfbmf3x/36KQxgUu4ybk39/N8PPLg3HO/N/edQ3I/53zP+Z6vqCrGGGPCT4TbAYwxxrjDCoAxxoQpKwDGGBOmrAAYY0yYsgJgjDFhKtLtAO2lpKTo8OHD3Y5hjDEhZePGjZWqmno2rwm6AjB8+HDy8/PdjmGMMSFFRA6d7WusC8gYY8KUFQBjjAlTVgCMMSZMWQEwxpgw5VcBEJGrRWSPiBSKyN0dPB8jIs87z68TkeHO+igReUJEtonILhG5J7DxjTHGdNUZC4CIeICHgXlADnCziOS0a7YMOK6q2cCDwAPO+huAGFWdCEwH7vi8OBhjjHGXP0cAM4BCVd2vqk3Ac8Cidm0WAU84yy8Cl4mIAArEi0gk0AdoAk4GJLkxxphz4s84gAyguM3jEmBmZ21UtUVEqoFkfMVgEVAGxAE/UtWqcw1tTE9pavHy/mcVbCs58cW6hD5RzJs4iIz+fVxMZsy586cASAfr2k8i0FmbGUArMBgYAHwoIu+o6v6/erHI7cDtAEOHDvUjkjHda19FDY9/fJDXC0o5XtcMgDi/5arwy1W7mJmVxI25mVw7JYOIiI7+BIwJbv4UgBIgs83jIUBpJ21KnO6eRKAK+Crwpqo2A+Ui8jGQC/xVAVDVR4BHAHJzc22GGuMaVeWpdUX88+s7Abhy/ECum5rBnFEpRHl8PaZFx+p4ZcthXt58mLuWb+XlzYf59Q2TSUuIdTO6MWfNn3MAG4BRIpIlItHAUiCvXZs84DZneQmwRn1TjRUBl4pPPDAL2B2Y6MYEVlVtE9/580bufWU7s0Yk8+HfX8JDN0/lkrFpX3z4AwxNjuOHl41izd9dxL8snsiGg1Vc/Z8f8u6uoy6mN+bsnbEAqGoLcCewGtgFLFfVHSJyv4gsdJo9BiSLSCFwF/D5paIPA32B7fgKyZ9UtSDAP4Mx56ziVCNL/vsTPvisgnvn5/Cnb5xHWr/T79GLCF+dOZTX/2YO6QmxLHsin+c3FPVQYmPOnQTbnMC5ublqN4MzPam6rpmlj67lQGUNf/7WTGZkJZ3192hobuX2Jzfy0d4KHrp5GtdMGtQNSY3pnIhsVNXcs3mNjQQ2Ya2uqYVvPr6ewvJT/M/Xc7v04Q8QG+Xh97dMY9rQAfzt85v5y57yACc1JvCsAJiw5fUqdz6zmS3FJ/jt0qlcNPqsbqX+JXHRkfzxm+cxOr0f331yI9sPVwcoqTHdwwqACVt//PgAa3aX8/MF45k3MTBdNgmxUfz5WzMYEBfND5/dTF1TS0C+rzHdwQqACUvbD1fzwJu7uSInnVvPHxbQ753cN4bf3DSZA8dq+UXezoB+b2MCyQqACTt1TS388LnNJMVH88D1kxAJ/CCu2SNT+N5FI3k+v5iVBWUB//7GBIIVABN2/un1XRyorOXBG6eQFB/dbe/zoytGMzmzP/esKKD0RH23vY8xXWUFwISVjYeqeHZ9Ed+ek8Xs7JRufa8oTwS/XTqFplYvv1y5q1vfy5iusAJgwkarV/l53g4GJsTyt5eP7pH3HJYcz/cuymbltjI+Lqzskfc0xl9WAEzYeH5DMdsPn+QfrhlHfIw/t8EKjDsuGkFmUh/uy9tBc6u3x97XmDOxAmDCwom6Jn61ejczs5JY0MOjdGOjPNx7TQ57y2v486eHevS9jTkdKwAmLPz6rc+orm/mvoXju+WqnzO5IieduaNT+Y+3P6PiVGOPv78xHbECYHq9A5W1PL3uELfMGsa4QQmuZBARfr4gh/rmVn63Zq8rGYxpzwqA6fUeencv0ZER/M2lo1zNMTK1LzfkDuHZ9cV2WagJClYATK9WWF7DK1sOc+v5w0ntF+N2HH5wSTaK8vBfCt2OYowVANO7/fbdvcRGebhj7gi3owAwZEAcN+Zmsjy/mJLjdW7HMWHOCoDptT47eorXCkq5bfZwkvu6v/f/uR9cko0gdhRgXOdXARCRq0Vkj4gUisjdHTwfIyLPO8+vE5HhzvqviciWNl9eEZkS2B/BmI795zt7iYvycPuFwbH3/7nB/ftw84xMXsgvoeiYHQUY95yxAIiIB9/UjvOAHOBmEclp12wZcFxVs4EHgQcAVPVpVZ2iqlOArwMHVXVLIH8AYzqyr6KGldvK+OYFWQzoxvv9dNX3L8kmIkL47/f3uR3FhDF/jgBmAIWqul9Vm4DngEXt2iwCnnCWXwQuky9fbH0z8Oy5hDXGX3/48ADRkRF844LhbkfpUHpCLNdNzWDFphKO1di4AOMOfwpABlDc5nGJs67DNs4k8tVAcrs2N9FJARCR20UkX0TyKyoq/MltTKcqaxpZsamE66cNISWI+v7b+/aFWTS2eHlyrY0ONu7wpwB0NGyy/Uzyp20jIjOBOlXd3tEbqOojqpqrqrmpqec2LZ8xT356iMYWL8vmZLkd5bSy0/px6dg0/vzpIRqaW92OY8KQPwWgBMhs83gIUNpZGxGJBBKBqjbPL8W6f0wPqG9q5cm1h7h8XBrZaX3djnNG37lwBFW1TazYdNjtKCYM+VMANgCjRCRLRKLxfZjntWuTB9zmLC8B1qiqAohIBHADvnMHxnSrlzaVUFXbxHeC7MqfzswakcTEjET+8OF+vN72B9bGdK8zFgCnT/9OYDWwC1iuqjtE5H4RWeg0ewxIFpFC4C6g7aWic4ESVd0f2OjG/DWvV3nsowNMGpLIjKwkt+P4RUT49oVZ7K+s5d3d5W7HMWHGr3EAqrpKVUer6khV/aWz7meqmucsN6jqDaqaraoz2n7Yq+p7qjqre+Ib87/e31vBgcpals3JcuWOn131lYmDGJwYyxOfHHQ7igkzNhLY9BpPrz1ESt8Y5k3o2fv9n6soTwQ3zxjKR4WV7K+ocTuOCSNWAEyvUHK8jnd3l7P0vEyiI0Pv1/qmGZlERghPrytyO4oJI6H3l2JMB55dX4QAN88c6naULknrF8tVEwby4sYS6pvsklDTM6wAmJDX1OLl+Q3FXDo2nYz+fdyO02VfnzWM6vpmXitof5W1Md3DCoAJeW/uOEJlTRO3zArNvf/PzcxKYlRaX562kcGmh1gBMCHvqU8PMTQpjrmjQnsUuYhwy6xhbC2ppqDkhNtxTBiwAmBC2t6jp1h/sIqvzRxKREToXPrZmcXTMoiL9vD0WjsZbLqfFQAT0p7fUExkhLBk+hC3owREQmwU8ycN4vWCUmobW9yOY3o5KwAmZDW1eHl582EuH5ceVDN+nasbczOpbWpl5bYyt6OYXs4KgAlZa3aXc6y2iRvP6x17/5+bPmwAI1LjeSG/+MyNjTkHVgBMyHohv5j0hJiQP/nbnohwY24mGw4et5HBpltZATAh6ejJBv6yp5zrpw0h0tP7fo2vm5aBJ0JYnl/idhTTi/W+vxwTFl7aVIJX4YbczDM3DkFp/WK5ZEwaL20qoaXV63Yc00tZATAhR1V5Ib+EGcOTyEqJdztOt7kxdwgVpxp5/zObJtV0DysAJuTkHzrOgcpabsjtXSd/27tkbBopfWNYbieDTTexAmBCzopNJfSJ8vCViaF12+ezFeWJ4Nopg1mzu5zjtU1uxzG9kF8FQESuFpE9IlIoInd38HyMiDzvPL9ORIa3eW6SiHwqIjtEZJuIxAYuvgk3Dc2tvF5QxrwJA4mPiXQ7TrdbPC2D5lbldRsTYLrBGQuAiHiAh4F5QA5ws4jktGu2DDiuqtnAg8ADzmsjgaeA76rqeOBioDlg6U3YWbO7nFMNLSyeluF2lB6RMyiBsQP78fImuxrIBJ4/RwAzgEJV3a+qTfgmd1/Urs0i4Aln+UXgMvHNyXclUKCqWwFU9Ziq2s3OTZet2FRCekIMs0emuB2lR4gIi6dmsKnoBAcqa92OY3oZfwpABtD2LFSJs67DNs4k8tVAMjAaUBFZLSKbROQnHb2BiNwuIvkikl9RYVc8mI4dq2nkvT0VXDvFd418uLh2agYRAi9vPux2FNPL+FMAOvpLUz/bRAJzgK85/y4Wkcu+1FD1EVXNVdXc1NTeNarTBM5rW0tp8SrXTevdV/+0l54QywXZKby8uQTV9n96xnSdPwWgBGg72mYI0H7Koi/aOP3+iUCVs/59Va1U1TpgFTDtXEOb8PTy5sPkDEpgzMB+bkfpcddNy6C4qp78Q8fdjmJ6EX8KwAZglIhkiUg0sBTIa9cmD7jNWV4CrFHfrspqYJKIxDmF4SJgZ2Cim3BSWF7D1pJqrguTk7/tXTV+IHHRHlbYyWATQGcsAE6f/p34Psx3ActVdYeI3C8iC51mjwHJIlII3AXc7bz2OPAbfEVkC7BJVVcG/scwvd2rWw4TIbBwymC3o7giLjqSq8cPZGVBGY0tdh2FCQy/LqRW1VX4um/arvtZm+UG4IZOXvsUvktBjekSVeXVLaVckJ1CWr/wHUaycMpgVmw+zPt7Krhy/EC345hewEYCm6C3pfgERVV1LJwcnnv/n7sgO4Xk+Ghe3dr+FJwxXWMFwAS9V7eUEh0ZwVUTwnuvN8oTwTWTBvHOzqPU2HSRJgCsAJig1tLq5fWCMi4bm0ZCbJTbcVy3aMpgGlu8vLXjiNtRTC9gBcAEtU/2HaOyppFFYXryt71pQwcwZEAfXtli3UDm3FkBMEHt1S2l9IuN5OIxaW5HCQoiwsLJg/m4sJKKU41uxzEhzgqACVoNza2s3nGEeRMGEhvlcTtO0Lh2agatXmWV3SHUnCMrACZordldTk1jC4umhOfgr86MTu/H2IH9eGWL3RvInBsrACZovba1lJS+Mcwakex2lKCzcMpgNhedoLiqzu0oJoRZATBB6VRDM2t2lzN/0qCwuvOnvxZM8p0UX2ndQOYcWAEwQemdXUdpbPGyYHLvnvaxqzKT4piS2Z/XbFCYOQdWAExQem1rGRn9+zA1c4DbUYLWgsmD2VF6kn0VNW5HMSHKCoAJOsdrm/jgswrmTxpEhHX/dGr+pEGIYEcBpsusAJig8+aOI7R4lQVhfu+fM0lPiGVmVhKvbS21iWJMl1gBMEHnta2ljEiJZ/zgBLejBL0Fkwezr6KWXWWn3I5iQpAVABNUyk828On+Y8yfPBgR6/45k3kTfFdJvVZg3UDm7FkBMEFl1bYyVGHBJLv6xx9J8dHMyU6xbiDTJX4VABG5WkT2iEihiNzdwfMxIvK88/w6ERnurB8uIvUissX5+n1g45ve5rWCMsYO7Meo9PCb97erFkweTMnxerYUn3A7igkxZywAIuIBHgbmATnAzSKS067ZMuC4qmYDDwIPtHlun6pOcb6+G6DcphcqPVHPxkPH7eTvWbpyfDrRngheL7BBYebs+HMEMAMoVNX9qtoEPAcsatdmEfCEs/wicJlYB645SyudD7D51v1zVhJio5g7OpWVBWV4vdYNZPznTwHIAIrbPC5x1nXYxplEvhr4/AYuWSKyWUTeF5ELO3oDEbldRPJFJL+iouKsfgDTe7xeUMrEjESGJce7HSXkLJg8iCMnG9hYdNztKCaE+FMAOtqTb7+b0VmbMmCoqk4F7gKeEZEvXdunqo+oaq6q5qampvoRyfQ2h47VsrWk2m790EWXj0snNirCBoWZs+JPASgBMts8HgK0/y37oo2IRAKJQJWqNqrqMQBV3QjsA0afa2jT+3zef33NJOv/74r4mEguHZvGqm1HaLVuIOMnfwrABmCUiGSJSDSwFMhr1yYPuM1ZXgKsUVUVkVTnJDIiMgIYBewPTHTTm7xeUMa0of3J6N/H7Sgha/6kwVTWNLJu/zG3o5gQccYC4PTp3wmsBnYBy1V1h4jcLyILnWaPAckiUoivq+fzS0XnAgUishXfyeHvqmpVoH8IE9oKy2vYVXbSrv45R5eMSSM+2mODwozfIv1ppKqrgFXt1v2szXIDcEMHr3sJeOkcM5pe7vWCUkTgKxOt//9c9In2cHlOOm9sP8L9iyYQ5bFxnub07DfEuEpVeW1rKTOGJ5GeEOt2nJA3f9JgTtQ183FhpdtRTAiwAmBctfvIKfZV1Fr3T4DMHZ1Cv9hIXttqg8LMmVkBMK56vaAUT4Qwb8JAt6P0CjGRHq4aP5C3dhyhsaXV7TgmyFkBMK7xdf+UMXtkMsl9Y9yO02ssmDyYU40tvL/HBlWa07MCYFyz7XA1RVV1X0xwbgJj9shkBsRF2b2BzBlZATCueW1rKVEe4arx1v0TSFGeCOZNHMQ7u45S32TdQKZzVgCMK7xeZWVBGXNHpZIYF+V2nF5n/qRB1DW1smZ3udtRTBCzAmBcsanoOKXVDXb1TzeZmZVMar8YuzeQOS0rAMYVr20tJSYygstz0t2O0it5IoRrJg5izZ5yTjU0ux3HBCkrAKbHtbR6WbmtjEvHptE3xq/B6KYLFkweRFOLl7d3HnU7iglSVgBMj1u7v4rKmiYWWvdPt5o2dAAZ/fuQZ91AphNWAEyPy9t6mL4xkVwyNs3tKL2aiLBg8mA+2ltJVW2T23FMELICYHpUY0srb24/wpU56cRGedyO0+stmDyIFq/yxnYbE2C+zAqA6VEffFbJyYYWFkyx7p+ekDMogZGp8eRtsW4g82VWAEyPyttayoC4KOZkp7gdJSyICAsnZ7D+YBVHqhvcjmOCjBUA02Pqmlp4Z+dR5k0cZPeq70ELJg9C1XfjPWPa8uuvUESuFpE9IlIoInd38HyMiDzvPL9ORIa3e36oiNSIyI8DE9uEond2lVPf3Gr3/ulhI1L7MiEjwQaFmS85YwFw5vR9GJgH5AA3i0hOu2bLgOOqmg08CDzQ7vkHgTfOPa4JZXlbSklPiGFGVpLbUcLOwsmD2VpSzcHKWrejmCDizxHADKBQVferahPwHLCoXZtFwBPO8ovAZSIiACJyLb6J4HcEJrIJRSfqmnj/s3IWTBqMJ0LcjhN2FkwejAg2JsD8FX8KQAZQ3OZxibOuwzbOJPLV+CaJjwf+HvjF6d5ARG4XkXwRya+osHuY90arth2huVW5dmr7Xx3TEwYl9mFmVhKvbDmMqrodxwQJfwpAR7tr7X+DOmvzC+BBVa053Ruo6iOqmququampqX5EMqHmlS2HGZEaz/jBCW5HCVuLpmSwv6KW7YdPuh3FBAl/CkAJkNnm8RCg/XHkF21EJBJIBKqAmcC/ichB4G+BfxCRO88xswkxpSfqWX+gimunZOD0DBoXfGXCIKI8wqtbDrsdxQQJfwrABmCUiGSJSDSwFMhr1yYPuM1ZXgKsUZ8LVXW4qg4H/gP4F1X9XYCymxDxeb/zIhv85arEuCguHpNG3tZSWr3WDWT8KABOn/6dwGpgF7BcVXeIyP0istBp9hi+Pv9C4C7gS5eKmvD16pZSpmT2Z1hyvNtRwt61UzIoP9XIuv3H3I5igoBf9+JV1VXAqnbrftZmuQG44Qzf474u5DMh7rOjp9hVdpL7FrS/cti44bJxvltwv7LlMLNtNHbYs+GYplu9svkwnghhvt36OSjERnm4avxA3th2hIZmmy843FkBMN3G61Ve3VLKnOwUUvrGuB3HOK6dOphTjS02X7CxAmC6z/qDVRw+Uc910+za/2Aye2QK6QkxrNhkVwOFOysAptus2FRCfLSHK3MGuh3FtOGJEK6dksF7e8o5VtPodhzjIisApls0NLeyatsRvjJxEH2ibeKXYHPdtCG0eNVuEBfmrACYbvHWzqPUNLaw2Lp/gtKYgf3IGZTAis3WDRTOrACYbrFiUwmDE2OZlZXsdhTTieumZVBQUk1h+Sm3oxiXWAEwAVd+qoEP91ayeFoGEXbnz6C1cIrvzqx2Mjh8WQEwAZe3xXergcVTh7gdxZxGWr9Y5o5K4eXNh/HarSHCkhUAE3AvbTrM5CGJZKf1dTuKOYPF04ZQVt3Ap3ZriLBkBcAE1PbD1ewqO8mS6bb3HwquzEmnX2wkL+QXn7mx6XWsAJiAenFjCdGRESycbFf/hILYKA+Lpgzmje1HONnQ7HYc08OsAJiAaWhu5eXNh7lq/EAS46LcjmP8dMP0TBpbvDYmIAxZATAB886uo1TXN3NjrnX/hJJJQxIZk96P5fklbkcxPcwKgAmYF/JLyOjfh9kj7TbDoUREuCF3CFuLT/DZURsTEE6sAJiAKD1Rzwd7K7h+WgYeu/Y/5CyemkFkhNjJ4DDjVwEQkatFZI+IFIrIl2b7EpEYEXneeX6diAx31s8QkS3O11YRWRzY+CZYrNhUgiosmZ555sYm6CT3jeGycWm8vPkwza1et+OYHnLGAiAiHuBhYB6QA9wsIu2nd1oGHFfVbOBB4AFn/XYgV1WnAFcD/+NMGm96Ea9XeWFjCeePSGZocpzbcUwX3ZibSWVNE+/usnkCwoU/RwAzgEJV3a+qTcBzwKJ2bRYBTzjLLwKXiYioap0zpzBALGDDDXuhT/cf49CxOm46z/b+Q9lFo1NJT4jhuQ1FbkcxPcSfApABtO0YLHHWddjG+cCvBpIBRGSmiOwAtgHfbVMQviAit4tIvojkV1RUnP1PYVz1zLoiBsRFcfUEu+9/KIv0RHBTbibvf1ZBcVWd23FMD/CnAHR0Rq/9nnynbVR1naqOB84D7hGR2C81VH1EVXNVNTc1NdWPSCZYVJxqZPWOI1w/bQixUXbf/1B304yhACy3k8FhwZ8CUAK0PbYfArQfMfJFG6ePPxGoattAVXcBtcCEroY1weeFjcW0eJWlzgeHCW0Z/ftw8ehUnt9QbCeDw4A/BWADMEpEskQkGlgK5LVrkwfc5iwvAdaoqjqviQQQkWHAGOBgQJIb13m9ynPri5mZlWQ3futFvjpzGOWnGu1kcBg4YwFw+uzvBFYDu4DlqrpDRO4XkYVOs8eAZBEpBO4CPr9UdA6wVUS2AC8D31fVykD/EMYdH++rpKiqjq/OtL3/3uSSMakMTIjl2fV2Mri38+uSTFVdBaxqt+5nbZYbgBs6eN2TwJPnmNEEKTv52ztFeiK48bxMHlqzl+KqOjKT7NLe3spGApsuOXqygbd3HmXJ9CHERNrJ395m6XmZCPCMHQX0alYATJc8vfYQrarcMmuY21FMNxjcvw+Xj0vnufVFNDS3uh3HdBMrAOasNba08sz6Ii4Zk8aw5Hi345hu8o3Zwzle10ye3Sa617ICYM7aqm1lVNY0cdvs4W5HMd3o/JHJjE7vyxOfHETVBvH3RlYAzFl7/JNDjEiJ58Jsu+1zbyYi3Hr+cHaUnmTjoeNuxzHdwAqAOStbik+wtfgEt54/jAi77XOvt3hqBv1iI3n8k4NuRzHdwAqAOStPfHKQ+GgP19uk72EhPiaSG3MzeXP7EY6ebHA7jgkwKwDGb+UnG1hZUMaS6UPoF2tz/oaLW88fRqsqT356yO0oJsCsABi/PfHpQZq9Xr5xQZbbUUwPGpYcz+Xj0nlq3SHqmr50M18TwqwAGL/UNrbw1NoirsoZSFaKXfoZbu6YO4ITdc28YBPH9ypWAIxflucXU13fzHfmjnA7inHB9GEDmDq0P3/4aD+tXrsktLewAmDOqKXVy2MfHSB32ACmDxvgdhzjAhHhjrkjKK6q583tR9yOYwLECoA5ozd3HKHkeL3t/Ye5K3IGMiw5jkc+2GcDw3oJKwDmtFSVRz/YT1aK70SgCV+eCOHbc7LYWlLNhoM2MKw3sAJgTuuTfcfYWlLNty/MwmMDv8LekumZJMVH81/vFbodxQSAXwVARK4WkT0iUigid3fwfIyIPO88v05EhjvrrxCRjSKyzfn30sDGN93tt+/uJT0hhuun2cAvA32iPSybk8V7eyooKDnhdhxzjs5YAETEAzwMzANygJtFJKdds2XAcVXNBh4EHnDWVwILVHUivikjbXKYELJu/zHWHajiuxeNtAnfzRduPX8YCbGRPLTGjgJCnT9HADOAQlXdr6pNwHPAonZtFgFPOMsvApeJiKjqZlX9/F6yO4BYEYkJRHDT/R5aU0hK3xhutgnfTRv9YqP41pws3t55lJ2lJ92OY86BPwUgAyhu87jEWddhG2cO4WoguV2b64HNqtrY/g1E5HYRyReR/IqKCn+zm2608dBxPiqs5I65I2zv33zJN2dn0S8mkt/9Za/bUcw58KcAdHTmr/01YKdtIyLj8XUL3dHRG6jqI6qaq6q5qampfkQy3e2hNXtJio/ma7Ns7998WWJcFLfNHs6qbUf47Ogpt+OYLvKnAJQAmW0eDwHaTxH0RRsRiQQSgSrn8RDgZeBWVd13roFN99tcdJz39lTw7QuziIuOdDuOCVLL5mQRF+3hP9+1o4BQ5U8B2ACMEpEsEYkGlgJ57drk4TvJC7AEWKOqKiL9gZXAPar6caBCm+6jqjzw5m6S46O59fzhbscxQWxAfDTL5mSxsqCM7Yer3Y5juuCMBcDp078TWA3sApar6g4RuV9EFjrNHgOSRaQQuAv4/FLRO4Fs4F4R2eJ8pQX8pzAB8+HeStbur+JvLs2mb4zt/ZvT+87cEfSPi+LfVu9xO4rpAgm2Id25ubman5/vdoyw5PUqC373EdX1zbz7dxcRE2knf82ZPfrBfn65ahfPfGcms0faNKFuEZGNqpp7Nq+xkcDmCyu3lbGj9CR3XTHaPvyN375+/jAGJcbywJt77B5BIcYKgAGgudXLr9/aw9iB/Vg0pf1VvsZ0LjbKw48uH83W4hOs3nHU7TjmLFgBMAA8vfYQB4/V8eMrx9g9f8xZu25aBiNT43ngzd00tXjdjmP8ZAXAUFXbxG/e/ow52SlcNs7O0ZuzF+mJ4B/n53CgspbHPzngdhzjJysAhl+/tYfaplZ+viAHEdv7N11zyZg0Lh2bxm/fLaT8VIPbcYwfrACEuZ2lJ3l2fdP4+8kAABDHSURBVBFfnzWMUen93I5jQty983NobGnlV2/aZaGhwApAGFNV7nttB4l9ovjR5aPdjmN6gayUeL51QRYvbCxha7HdLjrYWQEIY3lbS1l/oIofXzWGxLgot+OYXuLOS7NJ6RvDz17dbhPIBzkrAGGqqraJX7y2k8mZ/Vl6nt3wzQROv9go7p0/jq0l1TzxyUG345jTsAIQpv555U5O1jfzwPUT7bJPE3ALJw/mkjGp/PtbeyiuqnM7jumEFYAw9P5nFazYdJjvXTySsQMT3I5jeiER4Z8XTwTgp69stxHCQcoKQJipbWzhH1ZsY0RqPD+4JNvtOKYXy+jfh59cNYYPPqvglS2H3Y5jOmAFIMz8vzd2cfhEPf963SSb6ct0u6+fP5xpQ/tzX95Oyqrr3Y5j2rECEEbe2XmUp9YW8Z0Ls5iRleR2HBMGPBHCv98wmaYWL3+3fCteuyooqFgBCBPlpxr4yUsF5AxK4MdXjXE7jgkjI1L78vMFOXyy7xiPfrjf7TimDSsAYcDrVX78QgF1TS389uYpdqtn0+NuOi+Tq8cP5N/f2mOzhwURvwqAiFwtIntEpFBE7u7g+RgRed55fp2IDHfWJ4vIX0SkRkR+F9joxl+PfXSADz6r4B+vySE7zW73YHqeiPCv108kOT6GHz67mVMNzW5HMvhRAETEAzwMzANygJtFJKdds2XAcVXNBh4EHnDWNwD3Aj8OWGJzVj4prORf39zNvAkD+dpMG/Bl3NM/Lpr/WDqFQ1V1dj4gSPhzBDADKFTV/araBDwHLGrXZhHwhLP8InCZiIiq1qrqR/gKgelhJcfr+MEzmxiREs+vbphsd/o0rps1IpmffmUcb+08ysN/KXQ7TtjzpwBkAMVtHpc46zps40wiXw0k+xtCRG4XkXwRya+oqPD3ZeY0GppbuePJjbR4lUduzbUJ3k3Q+OYFw1k8NYPfvPMZa3bbDGJu8qcAdLTb2P7YzZ82nVLVR1Q1V1VzU1NT/X2Z6YTXq/zdC1vZUXqS/7hpClkp8W5HMuYLIsK/LJ5IzqAE/s+zW9h95KTbkcKWPwWgBMhs83gIUNpZGxGJBBKBqkAENGdHVfmnlTtZWVDGPfPGctm4dLcjGfMlfaI9PHJrLvExkdz2x/UcPmGDxNzgTwHYAIwSkSwRiQaWAnnt2uQBtznLS4A1ajf/cMUjH+znTx8f5FsXZHH73BFuxzGmUxn9+/D4t86jrqmV2/64nhN1TW5HCjtnLABOn/6dwGpgF7BcVXeIyP0istBp9hiQLCKFwF3AF5eKishB4DfAN0SkpIMriEyAvLixhP/3xm7mTxrEP14zzk76mqA3dmACj96aS9GxOpY9kU9tY4vbkcKKBNuOem5urubn57sdI+S8kF/MT14q4IKRKTz2jVwb7GVCyhvbyvjBM5uYPmwAf/rmDLtooQtEZKOq5p7Na2wkcC/w3PoifvJSAXOyU3j0VvvwN6Fn3sRBPHTzNDYVneDWx9Zx0gaK9QgrACHuz58e5O4V25g7KpVHb82lT7R9+JvQdM2kQTz81akUlFTz9cfWU1Vr5wS6mxWAEOX1Kr9cuZOfvbqDy8el88it0+32zibkXT1hEP99y3R2lZ3kuv/6mAOVtW5H6tWsAISg+qZWvv/0Jh798AC3nj+M398yzbp9TK9xRU46z35nJicbWlj8Xx+z/oBdUd5drACEmKJjddzwP5+weucR7p2fwy8WjifSY/+NpneZPiyJl78/m6S4aG75wzqeXHvIppXsBvbJEUJWbSvjmt9+SNGxOh79ei7L5mTZpZ6m1xqWHM+K789m1shk7n1lO3faXUQDzgpACKhpbOGnL2/j+09vYkRaX1b+8EIuz7ERvqb36x8XzePfOI+fXD2GN7cfYf5DH7HxkHUJBYoVgCD33p5yrnrwA55ZX8S352Txwh3nk5kU53YsY3pMRITw/Yuzee72WbS0Kkt+/yn35e2wQWMBYKMtglRZdT0PvLGbV7aUMjI1nhe/ez7Th9k8viZ8nTc8idU/msuv3tzNE58e5O2dR7l3fg5XjU+3rtAuspHAQaamsYXfv7ePP3y0H68X7rhoBD+4JNsu8TSmjfyDVdyzYht7y2s4b/gAfnpNDlMy+7sdy1VdGQlsBSBIVNc389TaQ/zxowMcq21i4eTB/N+rxlh3jzGdaGn1sjy/hN+8vYfKmiauyEnn+xePZOrQAW5Hc0VXCoB1AbmsuKqOp9cV8fTaQ5xqbOGi0an86IrRYb83Y8yZRHoi+OrMoSycMpg/fOi7C+7bO48ye2Qyy+ZkcfGYNDwR1jV0OnYE4ILGllb+sruCZ9cX8cHeCgTfvVC+d9FIJmQkuh3PmJBU09jCs+uK+MNH+zl6spHBibHcdN5QluQOIaN/H7fjdTvrAgpiDc2tfLKvkpUFR3hr5xFONbQwMCGWG8/L5KbzMsPiF9SYntDc6uWdnUd5Zn0RH+6tBGDq0P7MnzSYq8anM2RA7+xWtQIQRFq9yp4jp1i7/xgf7K1g7f5jNDR76RcbyVXjB3LNpEFcmJ1io3iN6UbFVXW8VlDKyoIydpT6pp4cmRrPxWPSmJOdwrRhA0jsE+VyysCwAuASVaW0uoHth6vZUXqSrcUn2FR0nFMNvuuUR6TEM3d0KhePSWX2yBSiI+1D35iedqCyljW7y3lvTznrDlTR1OJFBMak92PasAFMGJzIhIwERqf3C8mr7rqtAIjI1cB/Ah7gD6r6r+2ejwH+DEwHjgE3qepB57l7gGVAK/BDVV19uvcK1gLQ0uqloqaR0hP1FFfVU1RVx8Fjtewrr2FfRS01zqCUCIHstL5MH5bEjKwB5A5Lsit5jAky9U2tbC4+Tv7B42w4WMWW4hNf7LBFCGQmxZGd2pcRqfEMTYpjSFIcmQPiGJgYG7ST1XTLVUAi4gEeBq7AN/n7BhHJU9WdbZotA46raraILAUeAG5ypn9cCowHBgPviMhoVW09m5CB4PUqjS1eGppbqW9upa6plbqmFmobW6lpbKGmsZmT9S1U1zdzoq6Z43VNVNY0UlXbRPmpRiprGmlfK9MTYshO68v10zLITuvL+IxExg1MsHvyGxPk+kR7mD0yhdkjUwDfUXxxVT07SqvZdeSUs2NXw4eFlTS1eP/qtfHRHtITY0mJjyEpPpqkvtEMiIsisY/vq19sFH1jIomPiSQ+xkNcVCR9oj3ERkUQE+khyiNBM3DNn1I2AyhU1f0AIvIcsAhoWwAWAfc5yy8CvxPfT7gIeE5VG4EDzpzBM4BPAxP/f+0qO8nfPLuZllYvza1Ki9dLU4vz5azzV9+YSPrHRZHcN4aBCbFMGJxIemIsAxNiGZQYS2ZSHEMG9AnJw0RjzJeJCEOT4xiaHMe8iYO+WO/1KhU1jRRX1XH4RD1Hqhs4crKB8pONHKttZF9FDesPNlFd30yr17/PmAiB6MgIoj0RREd6iPYIkZ4IIj3CZWPT+Ok1PTdtuj8FIAMobvO4BJjZWRtVbRGRaiDZWb+23Wsz2r+BiNwO3A4wdOhQf7P/lbhoD6PT+xIZ4duQURERvo0cGUGUJ4LYqAhiozzEREYQF+0hNspDXLSvQveLiSI+xkNinygS+kQRZSdmjTH47kOUnhBLekIsp+tbUVVqm1o5Udfk61FoaOFUYwv1Ta3UN7VS19xKY3MrDc2tNDT7dkqbWrw0tnhpafXS4lWaW70MTOzZqwH9KQAdHau0L3WdtfHntajqI8Aj4DsH4EemLxmWHM9/fW16V15qjDHnREToGxMZtOcHOuPPrm4JkNnm8RCgtLM2IhIJJAJVfr7WGGOMC/wpABuAUSKSJSLR+E7q5rVrkwfc5iwvAdao7/KiPGCpiMSISBYwClgfmOjGGGPOxRmPV5w+/TuB1fguA/2jqu4QkfuBfFXNAx4DnnRO8lbhKxI47ZbjO2HcAvzAjSuAjDHGfJkNBDPGmF6gK+MA7HIXY4wJU1YAjDEmTFkBMMaYMGUFwBhjwlTQnQQWkQrgUJtVKUClS3H8ZRkDwzIGhmU8d8GeD76ccZiqpp7NNwi6AtCeiOSf7ZntnmYZA8MyBoZlPHfBng8Ck9G6gIwxJkxZATDGmDAVCgXgEbcD+MEyBoZlDAzLeO6CPR8EIGPQnwMwxhjTPULhCMAYY0w3sAJgjDFhKqgKgIj8SER2iMh2EXlWRGKd21CvE5G9IvK8c0vqYMv4uIgcEJEtztcUlzP+HyffDhH5W2ddkoi87WzHt0VkQBBmvE9EDrfZjl/p4Ux/FJFyEdneZl2H2018fisihSJSICLTgjDjxSJS3WZ7/szFjDc4/9deEclt1/4eZzvuEZGrgi2jiAwXkfo22/H3Lmb8lYjsdn7nXhaR/m2eO/vtqKpB8YVvqsgDQB/n8XLgG86/S511vwe+F4QZHweWuL0NnUwTgO1AHL7bfb+Dbx6GfwPudtrcDTwQhBnvA37sYq65wDRge5t1HW434CvAG/hmvZsFrAvCjBcDrwfJdhwHjAHeA3LbrM8BtgIxQBawD/AEWcbhbdu5vB2vBCKd5Qfa/F93aTsG1REAvg+DPs6sYnFAGXApvonmAZ4ArnUp2+faZwy2Gc7GAWtVtU5VW4D3gcXAInzbD9zfjp1ldJWqfoBvPou2Ottui4A/q89aoL+IDKKbnWVGV3SUUVV3qeqeDpovAp5T1UZVPQAUAjOCLKMrOsn4lvM3A7751oc4y13ajkFTAFT1MPDvQBG+D/5qYCNwos0P3OGk8j2lo4yq+pbz9C+dw7IHRSTGrYz49qznikiyiMTh21PNBNJVtQzA+TctCDMC3Olsxz+63U3l6Gy7ZQDFbdq5+bt5uv/b80Vkq4i8ISLj3Yl3WsG0HU8nS0Q2i8j7InKh22Ec38J3FApd3I5BUwCcP/ZF+A5fBgPxwLwOmrp23WpHGUXkFuAeYCxwHpAE/L1bGVV1F75Dw7eBN/EdFrac9kU97DQZ/xsYCUzBV2B/7VZGP0gH64LtmupN+O4PMxl4CHjF5TwdCYXtWAYMVdWpwF3AMyKS4GYgEfkpvr+Zpz9f1UGzM27HoCkAwOXAAVWtUNVmYAUwG9+h9edTV7o9qXyHGVW1zOkKaAT+RA8cwp6Oqj6mqtNUdS6+Q8i9wNHPuyicf8uDLaOqHlXVVlX1Ao/i8nZ0dLbdSvjfoxZw93ezw4yqelJVa5zlVUCUiKS4lLEzwbQdO+R0qxxzljfi618f7VYeEbkNmA98TZ0TAHRxOwZTASgCZolInIgIcBm+uYT/gm+iefBNPP+qS/mg44y72vzxCb7+1+2n+R7dTkTSnH+HAtcBzwJ5+LYfuL8dO8zYrg99MS5vR0dn2y0PuNW5GmgWvu7AMjcC0klGERno/E4iIjPw/b0fcyVh5/KApSISIyJZ+C4GWO9ypr8iIqki4nGWR+DLuN+lLFfj62FYqKp1bZ7q2nbs6TPbZzjr/QtgN74//CfxndEe4fwghcALQEwQZlwDbHPWPQX0dTnjh/iK51bgMmddMvAuvqOBd4GkIMz4pLMdC5xf6EE9nOlZfIf7zfj2qJZ1tt3wHXI/jG9vcBttrhoJoox3AjucbbwW39GqWxkXO8uNwFFgdZv2P3W24x5gXrBlBK5vsx03AQtczFiIr69/i/P1+3PZjnYrCGOMCVPB1AVkjDGmB1kBMMaYMGUFwBhjwpQVAGOMCVNWAIwxJkxZATDGmDBlBcAYY8LU/wepwKhWIynfQAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_continuous_distribution( n )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Confidence Intervals\n",
    "\n",
    "Recall from GB213 that certain assumptions about normality must hold in order for you to do statistical inference.  We do not cover those here; refer to your GB213 text or notes.\n",
    "\n",
    "Here we cover a confidence interval for the sample mean using confidence level $\\alpha$, which must be between 0 and 1 (typically 0.95)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(70.29847811072423, 350.0544630657464)"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "α = 0.95\n",
    "# normally you'd have data; for this example, I make some up:\n",
    "data = [ 435,542,435,4,54,43,5,43,543,5,432,43,36,7,876,65,5 ]\n",
    "est_mean = np.mean( data )  # estimate for the population mean\n",
    "sem = stats.sem( data )     # standard error for the sample mean\n",
    "# margin of error:\n",
    "moe = sem * stats.t.ppf( ( 1 + α ) / 2, len( data ) - 1 )\n",
    "( est_mean - moe, est_mean + moe )  # confidence interval"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Hypothesis Testing\n",
    "\n",
    "Again, in GB213 you learned what assumptions must hold in order to do a hypothesis test, which I do not review here.\n",
    "\n",
    "Let $H_0$ be the null hypothesis, the currently held belief.  Let $H_a$ be the alternative, which would result in some change in our beliefs or actions.\n",
    "\n",
    "We assume some chosen value $0\\leq\\alpha\\leq1$, which is the probability of a Type I error (false positive, finding we should reject $H_0$ when it's actually true).\n",
    "\n",
    "### Two-sided test for $H_0:\\mu=\\bar{x}$\n",
    "\n",
    "Say we have a population whose mean $\\mu$ is known to be 10.  We take a sample $x_1,\\ldots,x_n$ and compute its mean, $\\bar{x}$.  We then ask whether this sample is significantly different from the population at large, that is, is $\\mu=\\bar{x}$?  We can do a two-sided test of $H_0:\\mu=\\bar{x}$ as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.05, 0.35845634462296455, False)"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "α = 0.05\n",
    "μ = 10\n",
    "sample = [ 9, 12, 14, 8, 13 ]\n",
    "t_statistic, p_value = stats.ttest_1samp( sample, μ )\n",
    "reject_H0 = p_value < α\n",
    "α, p_value, reject_H0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The output above says that the data does NOT give us enough information to reject the null hypothesis.  So we should continue to assume that the sample is like the population, and $\\mu=\\bar{x}$.\n",
    "\n",
    "### Two-sided test for $H_0: \\bar{x}_1=\\bar{x}_2$\n",
    "\n",
    "What if we had wanted to do a test for whether two independent samples had the same mean?  We can ask that question as follows.  (Here we assume they have equal variances, but you can turn that assumption off with a third parameter to `ttest_ind`.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.05, 0.02815503832602318, True)"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "α = 0.05\n",
    "sample1 = [ 6, 9, 7, 10, 10, 9 ]\n",
    "sample2 = [ 12, 14, 10, 17, 9 ]\n",
    "t_statistics, p_value = stats.ttest_ind( sample1, sample2 )\n",
    "reject_H0 = p_value < α\n",
    "α, p_value, reject_H0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The output above says that the two samples DO give us enough information to reject the null hypothesis.  So the data suggest that the two samples have different means."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Linear Regression\n",
    "\n",
    "### Creating a linear model of data\n",
    "\n",
    "Normally you would have data that you wanted to model.  But in this example notebook, I have to make up some data first."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.DataFrame( {\n",
    "    \"height\" : [ 393, 453, 553, 679, 729, 748, 817 ],  # completely made up\n",
    "    \"width\"  : [  24,  25,  27,  36,  55,  68,  84 ]   # also totally pretend\n",
    "} )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As with all the content of this document, the assumptions required to make the technique applicable are not covered in detail, but in this case we at least review them briefly.  To ensure that linear regression is applicable, one should verify:\n",
    "\n",
    " 1. We have two columns of numerical data of the same length.\n",
    " 1. We have made a scatter plot and observed a seeming linear relationship.\n",
    " 1. We know that there is no autocorrelation.\n",
    " 1. We will check later that the residuals are normally distributed.\n",
    " 1. We will check later that the residuals are homoscedastic.\n",
    "\n",
    "To create a linear model, use `scipy` as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "LinregressResult(slope=0.1327195637885226, intercept=-37.32141898334582, rvalue=0.8949574425541466, pvalue=0.006486043236692156, stderr=0.029588975845594334)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model = stats.linregress( df.height, df.width )\n",
    "model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A linear model is usually written like so:\n",
    "\n",
    "$$ y = \\beta_0 + \\beta_1 x $$\n",
    "\n",
    "The slope is $\\beta_1$ and the intercept is $\\beta_0$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-37.32141898334582, 0.1327195637885226)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "β0 = model.intercept\n",
    "β1 = model.slope\n",
    "β0, β1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From the output above, our model would therefore be the following (with some rounding for simplicity):\n",
    "\n",
    "$$ y = -37.32 + 0.132x $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To know how good it is, we often ask about the $R^2$ value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.8949574425541466, 0.8009488239830586)"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "R = model.rvalue\n",
    "R, R**2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this case, $R^2$ would be approximately $0.895^2$, or about $0.801$.  Thus our model explains about 80.1% of the variability in the data."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualizing the model\n",
    "\n",
    "The Seaborn visualization package provides a handy tool for making scatterplots with linear models overlaid.  The light blue shading is a confidence band we will not cover."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWAAAAFgCAYAAACFYaNMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deXRc5Z3u++9bs0ZrsOVRsgYMDkOYbDMYbCekE9Kh4WQmIQkQgjHp7sPpvn1Oku51c/pk3V4rffqsdLi3G4NJIGQgJpCJTncS6BAPjB6YDQZjybZky7aseSrV9N4/9pZsjAfZVu1dw/NZS0tVu0qqd1uqx1vvfvfvZ6y1iIiI9wJ+D0BEpFgpgEVEfKIAFhHxiQJYRMQnCmAREZ+E/B7AVLj22mvt7373O7+HISJyPOZYGwviCPjQoUN+D0FE5JQVRACLiOQjBbCIiE8UwCIiPlEAi4j4RAEsIuITBbCIiE8UwCIiPlEAi4j4RAEsIuKTgrgUWURkqq3bfpD7NrTS3jtCfXUpdyxrZsXCuil9DR0Bi4gcZd32g3zz8W0cHIxTVRLm4GCcbz6+jXXbD07p6yiARUSOct+GVsJBQ2kkhDHO53DQcN+G1il9HQWwiMhR2ntHKAkH37WtJByko3dkSl9HASwicpT66lJGk+l3bRtNpplXXTqlr6MAFhE5yh3LmkmmLSOJFNY6n5Npyx3Lmqf0dRTAIiJHWbGwjm9dfx51FTH6R5PUVcT41vXnTfkqCC1DExE5hhUL67i4oZpIKEBJJHjyLzgNOgIWETlKKp2hs3+U7uExLDZrr6MjYBGRIwyNpTg0OEbGZi94xymARUQAay2HhhIMxpOevaavUxDGmL8yxmwzxrxujPmpMSZmjGkyxrxgjNlhjHnEGBPxc4wiUvgyGUtnf9zT8AUfA9gYMxf4r8Aia+35QBC4EfhH4J+ttQuAXuA2v8YoIoUvnbF0DsSJH7Xu1wt+n4QLASXGmBBQCnQCHwQecx9/CPgvPo1NRApcOmPp7B9lzIfwBR8D2Fq7F/g/wB6c4O0HtgJ91tqU+7QOYO6xvt4Ys9IYs8UYs6Wrq8uLIYtIAUmlM+zrGyWRyvg2Bj+nIKqBG4AmYA5QBnz0GE895qlIa+0aa+0ia+2iGTNmZG+gIlJwnGVmcZJp/8IX/F0F8SGgzVrbBWCM+QVwJVBljAm5R8HzgH0+jlFECkw8mebgwBipjL/hC/7OAe8BLjfGlBpjDHAN8AbwR+BT7nNuBn7t0/hEpMAMjaXo7I/nRPiCv3PAL+CcbHsReM0dyxrga8BfG2PeAWqB7/s1RhEpHN1DYxwciGM9uMBisny9EMNa+z+B/3nU5lZgiQ/DEZEClM5YDg7GGU34s9LhRHQlnIgUrFya7z0WBbCIFKSBeJLuoUROTTkcTQEsIgWne2iM/lFvLys+HQpgESkY1lq6BscYGkud/Mk5QAEsIgUhk7EcyNGTbcejABaRvJdKZ9g/EPf1suLToQAWkbyW6ysdTkQBLCJ5q38kSc9Ibq90OBEFsIjknWQ6Q9fgmC81fKeSAlhE8spAPEnPUMKTnm3ZpgAWkbyQzjhLzEYS+bHEbDIUwCKS80YTaboG8/NE24kogEUkp/UOJ+gdSfg9jKxQAItITkqlMxwsgBNtJ6IAFpGcMzyW4tDQGOlM/p9oOxG/uyKLiExIZywHB+IcGIjnTPh2DY5l7XsrgEUkJwzGk3T0juRMIZ2e4QTfefJtrvnOerbvH8jKa2gKQkR8lUxn6B5K5MzysrFkmsde7ODhF9oZdeefH3x6F//4qfdP+WspgEXEN/2jSXqHc+Oiioy1/HH7Qe7f2MZBd9qhMhbiLz54FrcubcrKayqARcRziVSGQ0O5s8Lh9b393LNuJ9v3DwIQChg+fvFcvnB5AwtmVhAOZme2VgEsIp6x1tI3kqRvNJkTBXT29o1y/4ZWNuw4NLFt2dnTuf3qZuZWlWT99RXAIuKJeDLNoaGxnKjZOxhP8uPn9/DLl/aScldbLJxVwZ3LW7hg3jTPxqEAFpGsstbSM5ygfzTJptYe1m5up3NglNmVJdy4uJ4lzTWejSWVzvD4K/v44XO7GYg7J/3qKqLcfnUzH1g4g4Axno0FFMAikkWjCeeoN5nOsKm1h7uf2kEoYKiMhegeHuPup3ZwFwuyHsLWWp7d2c19G1rp6B0FoDQS5PNLGvjkJXOJhoNZff3jUQCLyJTLZCzdwwkG44c7E6/d3E4oYChxw64kHGQ0mWbt5vasBvDbBwZZvW4nr3T0AxAw8LELZnPzlY3UlEWy9rqToQAWkSk1PJaieyjxnsplnQOjVMbeHTmxcID9A6NZGUfX4BgPPNPGE9sOMH66b0ljNXcsb6FpellWXvNUKYBFZEqkM5buoeO3hJ9dWUL38NjEETBAPJlhVuXUrjYYTaR5ZHM7j2xpZ8w94dc0vYxVy5tZ3OjdfPNkKIBF5IwNxpP0DCdOWL/hxsX13P3UDkaTaWLhAPFkhlTGcuPi+ikZQzpjeWLbfh54Zhfdw075yurSMLcubeSj588mGPD2BNtkKIBF5LSdymXES5pruIsFrN3czv6BUWZN4SqIF3f3snr9TnZ2DQMQCQX49KXz+NySekojuRtzuTsyEclZ1lrnMuKRU7ugYklzzZSecNvTPcK9G3byfGvPxLYPva+O265qYmZlbMpeJ1sUwCJySkYSzkm2ZNq/Cyr6RhI89Nxu/u2VfYzPelwwt5I7V7SwcFalb+M6VQpgEZmURCpD9/AYown/6jckUhl+8dJefvL8bobdccypirFyWTNXnzUd4/GFFGdKASwiJ2StpXckSb+P9Rustax/u4v7N7bR2R8HoDwa4ouXN3DDRXOJhPKztLkCWESOa2gsRe+wv9MNb+wb4J51O3mj0ymKHgwYbrhoDl+8fD7TSsK+jWsqKIBF5D3GUmm6hxK+lovc3x/n/o2t/PGtroltS8+qZeXVzdTXlPo2rqmkABaRCZmMpXfEKZzjl6GxFA+/sIefv9hBMu1MeSyoK+erK1q4sL7Kt3FlgwJYRIDjX0LslXTG8ptXO/nBs7sm/gOYXh7hK1c18aFzZ3peqcwLCmCRIpdKZ+geTjDsUzNMay0vtPVw3/pWdveMAE6NiM8tbuDTi+YR86lSmRcUwCJFbCCepGfIv55sO7uGuHd9K1t39wJggI+eP4tblzZSWx71ZUxeUgCLFCG/1/R2D43x4DO7+O3r+ycqlV3aUMWq5S201JX7MiY/KIBFiojfPdniyTSPbungp5v3EE86c80NNaWsWt7MZU01eXchxZlSAIsUiSO7U3gtYy3/+eZBvrexlUNDTqWyaSVhbrmykeven5uVyrygABYpcKl0hp7hxHHr9GbbKx19rF63k7cPDAEQDho+eck8Pn9ZA+XR4o4gX/feGFMFfA84H7DAl4G3gEeARmAX8Blrba9PQxTJW9ZaBkZT9I74c5Kto3eENRvaePqdwy3fV5w9g9uXNTF7WvZbvucDv//7uRv4nbX2U8aYCFAK/C3wB2vtt40xXwe+DnzNz0GK5Jt4Mk3XoD/TDQOjSX70/G5+9fK+iQLt586u4M4VLZw3x7uW7/nAtwA2xlQCy4BbAKy1CSBhjLkBWOE+7SFgHQpgkUlJpTP0jCQYins/3ZBMZ/jVy/v48fO7GXRff1ZljNuvbmLFOTOK7gTbZPh5BNwMdAEPGmMuBLYCdwEzrbWdANbaTmNM3bG+2BizElgJ0NDQ4M2IRXJY/0jSl+kGay1Pv9PNmg2t7O1zGmyWRYLcdFkDn7hkXt5WKvOCnwEcAi4B/tJa+4Ix5m6c6YZJsdauAdYALFq0yJ9V5CI5IJ50VjckUt5PN7y1f5B71u3ktb2HW77/2YVzuPmK+VSV+tvyPR/4GcAdQIe19gX3/mM4AXzAGDPbPfqdDRz0bYQiOSydsfQMJxiMe1845+BAnO893cZ/vnn47Xl5cw13LGtmfm1utHzPB74FsLV2vzGm3RhzjrX2LeAa4A3342bg2+7nX/s1RpFcNZkuxNkwkkjx003tPLq1Y+KIu2VGGXcub+GS+dWejqUQ+L0K4i+Bn7grIFqBW4EA8DNjzG3AHuDTPo5PJKckUhkODY15Xqc3nbH89vVOHnxmF70jzhF3TVmELy9t5CPnzSraCynOlK8BbK19GVh0jIeu8XosIrlsvE7vQDzl+SXEm3f1cO/6VtoOOS3fo6EAn11Uz2cX11MSKdxKZV7w+whYRE5iIJ6k14fphrZDw9y3oZVNbYdbvn/43JncdlUTMyoKv1KZFxTAIjlqLJXm0FCCMY+nG3pHEvzg2V38+6udEy3fL6qfxqrlLZw9s8LTsRQ6BbBIjslkLD0jCQY8bguUSGV4bGsHD2/aw4hbpnJedQkrr25m6Vm1upAiCxTAIjnEj9UN1lqe2t7F955u5cDAGACVsRBfvGI+1184h3BQF1JkiwJYJAf41YX49b393LNuJ9v3DwIQChg+fvFcvnB5AxWx/G75ng8UwCI+8qsL8b6+Ue7f2Mb6tw+3fF+2YDq3L2tmbtWpVSrb1NrD2s3tdA6MMruyhBsX17OkuWaqh1yQFMAiPhmMJ+kdTnrahXgonnIrle2daPl+zswKvrqihQvmnXqlsk2tPdz91A5CAUNlLET38Bh3P7WDu1igEJ4EBbCIx/y4mCKVzvBvr3by0LO7GHArldVVRPnK1U18cGHdabd8X7u5nVDAUOJ2Li4JBxlNplm7uV0BPAkKYBGP+HExhbWWZ3c6lcrae51KZSXhIJ+/rJ5PXTKP6Bm2fO8cGKUy9u4YiYUD7B8YPaPvWywUwCIeGBpL0TOU8HS6YceBQVavb+Xl9j7AqVT2pxfM5pYrG6kpm5pKZbMrS+geHps4AgaIJzPMqlTHi8lQAItkkR/t37sGx3jgmTae2HZgouX74sZqVi1voWn61FYqu3FxPXc/tYPRZJpYOEA8mSGVsdy4uH5KX6dQKYBFssCP9u+jyTSPbG7nZ5vbibuVyhprS1m1vIUlTdmZj13SXMNdLGDt5nb2D4wyS6sgTokCWGSKed3+PZ2xPPHGAR54uo3uYafle3VpmFuXNvLR87Pf8n1Jc40C9zQpgEWmiB/t31/a08vqda2803W45funL53H55Y0UFbkLd/zgX5CImfIWkv/aJK+kaRn/dj2dI9w34ZWnmvtnth2zcI6bru6iVmVMU/GIGdOASxyBkYSKbqHEp5NN/SPJHnouV08/sq+iUplF8ytZNXyFt43u9KTMcjUUQCLnIakO90w7NF0QyKV4Zcv7eXHL+xmeMxZUTGnKsbKZc1cfdZ0VSrLUwpgkVPg9eoGay3r3z7E/Rtb6eyPA1AeDfHFyxu44aK5avme5xTAIpPk9XTDm50D3LNuJ9v2DQAQDBhuuHAOX7xiPtNKVKmsECiARU4ilc7Q7eF0w/6BON/b2MZT2w+3fF/aUsvKZc3U15R6MgbxhgJY5DjGVzf0jngz3TA0luLhF/bw8xc7JiqVnVVXzp3Lm7m4QS3f/RAOBrK6jloBLHIMXl5Mkc5Y/v21Tn7wzC763LrAteURvnJVE39y7szTrlQmpy5gDCWRoPMRDma9G4gCWOQIXq5usNayyW35vrt7BIBYKMCNS+r59KL6dxW4keyJhAKUhIOURkLEwgFPV5QogEVwjkL7PCwV2do1xOr1rWzd3QuAAT5y3iy+fFUj08vV8j2bjjzKLQ0HCfnY804BLEXN66vYeoYTPPjMLn77+uGW7xc3VHHn8hbOqivP+usXq0goQGkkRGkkSDTk7VHuiSiApWh5uawsnkzz2NYOfrqpnVG3E0Z9dQmrlrdweXNNzgRCoQi6XTrG53L9PMo9EQWwFB0vl5VlrOUPbx7kexvb6Bo63PL9lisbue79s3M2GPJRNOxMKZREgsTyZP5cASxFw1rLwGiK3pGEJ9MNr3T0ce+6Vt464LR8DwcNn7h4LjddNp/ymN56ZyoYMBNHuKWRUNbLbmaDfgukKAyPpegZ9ma6oaN3hPs3trFxx6GJbSvOnsHty5qYPU2tes5EPh7lnogCWAraWCpN91DCkw7EA6NJfvT8bn798j5S7hm2982u4M7lLZw/99RbvouzYqE06hzhloSDeXmUeyIKYClI6YylZzjBYDyZ9ddKpjP8+uV9/Oj53Qy6Ld9nVka5/epmPnDODJ1gO0XBgKE0EqIs6kwvFPK/nwJYCoqX87zWWp5+x2n5vrfPacNeFgly02UNfOKSeapUdgpCgQCl0SBlkRAlkfyfWpgsBbAUDC+Xlb21f5DV63fyakc/4LR8/7P3z+HmK+dTVTo1Ld8LXSgQoCwapCwaKoj53NOhAJa8F0+m6Rn2Zp734ECc7z+ziyffODCx7bKmGu5Y3kxj7dS2fC9E4WCA0khxh+6RFMCSt8ZSaXqHk4wksr+edySRYu3mdn62pYOE2/K9eXoZq5Y3s6hRHYFPJBwMUBZ15nSjIYXukRTAkneS6Qy9IwmG4tkP3nTG8rvX9/PAM230jjgn9GrKItx6ZSPXnj+r4M7KT5VIKEBZJERZNKS58BNQAEveSKUz9I0mGfSoYM6WXT3cu6GV1q5hAKKhAJ9ZNI8bFzcU1YmiyYqGg5S50wvZLuNYKBTAkvMyGUvfaJKBUW8K5uzqHua+9a280NYzse1Pzp3JbUsbqVPL93eJhYPukW7u1lvIZQpgyVmZjFOpbCCeJJ3JfvD2jiR46Nnd/ObVwy3fL5w3jTtXtHD2zIqsv36+KIk4F0aURRS6Z0oBLDknk7EMxL0rEZlIZXhsawcPb9rDSMJZSTG3qoSVy5q56qzagr4QYDKMcSqLja/T1bz31FEAS86w1jIQT9E3kvDkiNdayx/f6uL+ja0cGHAqlVXEQnzpivlcf+Gcop7HHC90UxoJURoOElDoZoUCWHKCl8VyALbt62f1up280elUKgsFDP/l4jl84bL5VBZhy3djDLFwgNJwiFgkoOViHlEAi6/iyTS9IwlGE9m/iAKgs3+U+ze0se7troltV501nZXLmphXXVwt38enFsrcqQUd5XpPASy+SKSctbxeFEUHGIqn+MkLu/nFS3snWr6fM7OCVSuauXBelSdjyBXRcJDyaIjyqOZz/eZ7ABtjgsAWYK+19jpjTBOwFqgBXgS+aK1N+DlGmTqpdIYejy6iGH+937zayQ+e3cWA+5ozyqN85eomrnlfXVG1fI+EAtRVxHRhRA7xPYCBu4A3gUr3/j8C/2ytXWuMuRe4DVjt1+BkamQyll4Puw5ba3m+tYd71++kvdepVBYLB/j8kgY+dem8oqtDEA4GmD2tREe8OcbXADbGzAM+BvwD8NfGWe/zQeDz7lMeAv4eBXBec7oOe7OyAeCdg0Pcu34nL+7pA5xKZdeeP4tbr2yktghbvjvhG1P45iC/j4C/C/wPYHyVey3QZ60d//u0A5h7rC80xqwEVgI0NDRkeZhyOuLJNIeGxiaK12TboaExHnh6F7/ftp/xqF80v5pVy5tpnlGcLd8DxjBrWkwXTOQo3wLYGHMdcNBau9UYs2J88zGeeszDJmvtGmANwKJFi7w5tJJJSaUz9AwnGPLoBNtoMs2jW9pZu6mduBv282tLWbW8mcuaaj0ZQ66qLosU9XrmXOfnEfBS4HpjzJ8CMZw54O8CVcaYkHsUPA/Y5+MY5RSM12zoH016Ms+bsZYn3zjA955uo3vIOU9bVRLm1qWN/OkFs4v+T+6SSJBpRbimOZ/4FsDW2m8A3wBwj4D/xlp7kzHmUeBTOCshbgZ+7dcYZfIG4kn6hpOkMt5MN7zc3sc963byzsEhwGn5/qlL5/H5JQ2URf2eWfNfwBhmFOF8d77Jxd/UrwFrjTH/D/AS8H2fxyMnMJJwrmDzap53T88Iaza08uzO7olt1yys47arm5ilSmUTassjmvfNAzkRwNbadcA693YrsMTP8cjJedkGCJyVFD98bjePv7JvYjXFeXMq+eqKFt43u/IkX11cyqIhKmKaesgHORHAkj+8voItkcrwq5f38uPn90yc1Js9LcbKZc0sWzC96CuVHS0YMEzX1EPeUADLpKTSGXpHkgyNeXchxYYdh1izoZXO/jgAZdEgX7hsPh+/eK6u5jqO6eXRoj/5mE8UwHJCXq9sAHizc4DV63by+r4BwLmQ4voL53DzFY1MK9Wf1kfb1NrD2s3tHBiM01hbxh3LmlmxsM7vYckkKIDlmKx1ulH0j3rTjQJg/0Cc729s4w/bD05su6K5ljuWNdNQW1yVyiZrU2sPdz+1g3DQUFMa5uBgnG8+vo1vgUI4DyiA5T28XlI2PJbi4U17eGxrx0SlsrNmlLNqRTOXNFR7MoZ8tXZzO6GAoSIWJmAMpZEAI4kU921oVQDnAQWwTPB6SVk6Y/mP1zp58Jld9I06Ld9ryyJ8+aomPnzuTM1lTkLnwCjVJeF3VXUrCQfp6B3xcVQyWQpg8XxJGcCmNqdS2a5uJyhioQCfXVzPZxbXU1JklcpOlzGGhppSeoYTHDkzPppMF11x+XylAC5iXi8pA2jtGuLe9a1s2d0LOMU/PnLeLL58VaOWT52C8dq+f77iLL75+DZGEilKwkFGk2mSacsdy5r9HqJMwqQC2BhzNvDfgflHfo219oNZGpdk0fiSssF40rPX7BlO8OAzu/jt650TLd8vbqjizuUtnFVXnJXKTkc4GKC6LEK5e7n1ioV1fAu4b0MrHb0jzKsu1SqIPDLZI+BHgXuB+wHv/k6VKeXHkrKxZJrHXuzg4RfaGXWnOOqrS1i1vIXLm2t0IcUkGGMoiwQpi4YojQTf82+2YmGdAjdPTTaAU9ZaFUXPU9ZaBkZT9I16VxQ9Yy1/ePMg33+6jYODTsv3yliIW65s5Lr3z1adguMIBQJEQgHCQUMk5NyOBAP6j6pAnTCAjTE17s1/M8Z8FfglMDb+uLW2J4tjkykwGE/SN5L0rN07wKsdfaxe18pbB5yW7+Gg4RMXz+Wmy+ZTHtNph2MJBwPUlkcojejfp5ic7Ke9Facg+vh/v//9iMcsoJn+HDWaSNM97F03CoC9vaOs2djKxh2HJrYtP3sGt1/dxJyqEs/GkU+CAUNVaYTKWEhHuUXohAFsrW0CMMbErLXxIx8zxqj2Xw4aSzlLykYT3k3VD8aT/Oj53fzqpX2k3CmO982u4M7lLZw/d5pn48gn4WCAyliYyhIFbzGb7N87zwKXTGKb+CSZztDrYRug8dd8/JV9/Oi53RMt3+sqoqxc1swHzplR9MFijCEUMO6crjOvGw46c7oBXWQinHwOeBZOU8wSY8zFHJ6KqAS00jsHpDOWPg/bvYNzUu/Znd3u0ien5XtpJDjR8r1YK5WFgwEqYiE3bJ3ALfb/hOTETnYE/BHgFpzebN85Yvsg8LdZGpNMwnixnL6RJBmPghfg7QODrF63k1c6+gGnUtnHLpjNLUsbqS6NeDaOXBIMGKpKIppOkFN2sjngh4CHjDGftNb+3KMxyQlYaxkcS3laLAega3CM7z/dxpNvHJhoU72kqYY7ljXTNL3Ms3HkkvF53IpYSFMKclpONgXx18e6Pc5a+52jt0n2DI85xXK8XFI2mkjzyOZ2HtnSzpi7oqJpehmrljezuLHmJF9duGrKIlQV6RG/TJ2TTUFUuJ/PARYDj7v3/wzYkK1Bybv5USwnnbH8ftt+HnhmFz3DTsv36tIwty5t4qPnzyrqSmW1ZVEVhpcpcbIpiP8FYIx5ArjEWjvo3v97nMuTJYv8KJYDsHV3L6vX76S1axhwCr98ZtE8blxcX/QXCtSWR5lWovCVqTHZd1MDkDjifgJonPLRCOBPsRyA3d3D3LehledbD1/g+CfnzuS2pY3UqeU70yuiVKrbsEyhyQbwj4BNxphf4lwB93Hgh1kbVZHyo1gOQN9Igoee3c2/vbpvolLZBXOnceeKZhbOUst3gBkVUbV6lyk3qQC21v6DMea3wNXupluttS9lb1jFJZOxThsgj5eUJVIZfvFiBz95YQ/D7pVzc6pi3LGshavOqtWSKpyLKWZURCfKP4pMpZOtgqi01g64RXl2uR/jj9WoGM+Z87r/GjhL2da91cX9G9vYP+BcYV4RC/HFy+dzw0VzCKtSGeAsM6urjBINqUOHZMfJ/lt/GLiOw0V54PDVcCrGcwa87r82btu+flav28kbnU6lsmDAcMNFc/jS5fOp9Pnk0nh79c6BUWZXlnDj4nqWNPuz1K00EqKuIqr1vZJVJ1sFcZ1782mcZWcbrbXbsz6qAhZPpukd8bZYDkBn/yjf29jGH9/qmth21VnTWbmsKSf6h423Vw8FDJWxEN3DY9z91A7uYoHnIVwRCzOjQu2RJPsmO7H1IHAV8P8ZY5qBl3DC+O6sjazAjKXS9I0kPV9SNjSW4uEX9vDzFw+3fF9QV85XV7RwYX2Vp2M5kfH26uMNOcf7m63d3O5pAOsCC/HSZE/CPWWMWY9zMcYHgFXA+YAC+CQSqQx9I95WKQPnQorfvLqPHzy7m3635fv08ghfubqZD72v7l1tzHNB58AolUcVa4+FA+wfGPVsDFrjK16bbFPOPwBlwHPARmCxtfZgNgeW78ZSafpHkwzFvQ1eay0vtPVw7/pW9vS4Ld/DAT63pIFPXzqPWI62fJ9dWUL38Ni7WtLHkxlmVXpTyL22TOEr3pvsFMSrwKU4R739QJ8x5jlrrXeHJ3nAWsvQWIqBeIqxM7hs+HRPRu08OMTq9Tt5cU8f4FQqu/b8Wdx6ZSO1Od7y/cbF9dz91A5Gk2li4QDxZIZUxnLj4vqsv3ZNWUSXFosvzKks+DfGlAO3An8DzLLW5sS7etGiRXbLli2+vf5YKs1gPMVQPHXG63iPPBl1ZBDd9cHjn4zqHhrjgWd28bvX908sVbm0oYpVK1pomZE/Ld/H/+PZPzDKLI9WQVSXRqgu05yvZN0x5/wmOwXxFzgXYVwK7AYewJmKKGrxpHNibSQxddMMp3IyajSZ5tEt7azd3E486Sxnm19Tyh3Lm7msKf9avi9prvH0hFuVwld8NtkpiBKcguxbrbXeTmpm2brtB1O74DQAABI2SURBVLlvQyvtvSPUV5dyx7JmViysO+HXjKXS9A5PbfCOm8zJqIy1PPnGAb7/dBuHhpwSHVUlYW52W74Xc6WyyaoqjVCj8BWfTXYVxD9leyB+WLf9IN98fBvhoKGqJMzBwTjffHwb34JjhnA8mWZgNJnVFQ0nOxn1cnsfq9ftZMfBIcBp+f7JS+bx+csadLnsJE0rCSt8JScU9Tv2vg2thINmosRiaSTESCLFfRtaJwI4k7EMxlMMxJOeFEI/3smoaxbW8X//6nWe2dk98dwPnDOD269uZtY0VSqbrMqScM6fkJTiUdQB3N47QtVRS49KwkE6ekcmCuT0jyZJZ7wrkLOkuYa7WDBxMmp6WZTK0jDffWrHxDjOnV3JV1e0cO4cVSo7FRWxMNMVvpJDijqA66tLOTgYf1eR8ZFEipmVMdp7RzwN3iMtaa7h4vlV/Oplp+X7UKcz5TGrMsbKZU0sP1st30+VLi+WXFTUAXzHsma++fg2hseSRENBRhIpEmnLpy6Z51v4WmvZuOMQaza2sq/PqVRWFgly02UNfOKS4m35fibKYyGFr+Skog7gFQvr+BbwL398h47eEc/Wnh7P9v0DrF63k9f2DgDOhRR/duEcbr5ivuoTnKbyaIi6Cs2RS24q6gAGJ4QvrK+idyRx8idnyYGBON9/uo3/fPPw1d1XNNdyx7JmGmr9r1SWr8qjIbVSkpxW9AHsp5FEip9uaufRrR0TdYFbZpRx54oWLmmo9nl0+a0sqmkHyX0KYB+kM5bfvt7Jg8/sonfEqVRWWxbhy1c18eFzZ+pCijNU7oavTlRKrvMtgI0x9TiNPWcBGWCNtfZut/3RIzhdl3cBn7HW9vo1zqm2eVcPq9ftZFe3U6ksGgrw2cX1fHZRPSWR3KxUlk/KY5rzlfzh5xFwCvi/rLUvGmMqgK3GmCeBW4A/WGu/bYz5OvB14Gs+jnNKtB0a5t71O9m8y/m/xAAfPm8mX17apD+Vp4iWmkm+8S2ArbWdQKd7e9AY8yYwF7gBWOE+7SFgHXkcwD3DCR56dhf//lrnRMv3i+qncefyFhbMrPB3cAVkmq5wkzyUE3PAxphG4GLgBWCmG85YazuNMSeujJOjxpJpHnuxg59uamfE7f82r7qEO5Y1c2WLWr5PJbURknzlewC7NYZ/Dvw3a+3AZIPJGLMSWAnQ0NCQvQGeooy1/HH7Qe7f2MbBwTEAKmMhvnTFfK6/cA4htXyfUmojJPnM1wA2xoRxwvcn1tpfuJsPGGNmu0e/s4Fjtj6y1q4B1oBTkN2TAZ/Eax393LN+J2/td1q+hwKGj188ly9c3kBFTCExlYwxTC+P6N9V8pqfqyAM8H3gTWvtd4546HHgZuDb7udf+zC8U7K3b5T7N7ay4e1DE9uWLZjO7cuamVvlTU+zYmKMYWZl9F01PETykZ+/wUuBLwKvGWNedrf9LU7w/swYcxuwB/i0T+M7qcF4kh8/v4dfvrSXlHuGbeGsCu5c3sIF86b5PLrCFDCGmZUxLdmTguDnKoinOU6fJOAaL8dyqlLpDI+/0skPn9vFgNv1uK4iyu1XN/OBhTNyruV7oQgGnPDN1c7OIqdKf8OdAmstz+7s5r4NrXT0Oi2CSiNBPr+kgU9eMpeogiFrQoEAs6bFVA1OCooCeJLePjDIvet38nJ7P+BUKvvYBbO5+cpGtbfJsnDQCd+wVpBIgVEAn0TX4BgPPNPGE9sOTLR8X9JYzR3LW2iaXubr2IpBJBRgVmVMy/ekICmAj2M0keaRze08sqWdMbdSWWNtKXeuaGFxoz/1gotNNBxkVmVMxYmkYCmAj5LOWJ7Ytp8HntlF97BTI7i6NMytSxv56Plq+e6VkkiQmRUxAvr3lgKmAD7Ci7t7uWf9Tlq7hgHnz99PXzqPzy2p15pTD5VFQ9SpnKQUAaUKTqWy//377Tzf2jOx7UPvq+O2q5qYqY4KnlI5SSkmRR3A1lq+9Zs3+OGzu0lb5xTbBXMruXNFCwtnqeW71ypL1DZeiktRB7AxhrFUhrS1zKmKsXJZM1efNV1/+vqgujRCtZbzSZEp6gAG+KsPnc2caTH+5NyZWmfqk9qyKNNKVVRHik/RJ86Miig3XTZf4euTGRUKXyleRX8ELP4wxlBXEaUsql9BKV767RfPqaKZiEMBLJ5SRTORwxTA4plQIMDMaVGiIYWvCCiAxSOqaCbyXgpgybpwMMDsaapoJnI0BbBkVSQUYPa0EhUxEjkGBbBkTTQcZHalKpqJHI8CWLIi5tbyVfiKHJ8m5WTKRRW+IpOiAJYpNd5CSOErcnIKYJkyzmoHnXATmSwFsEyJ8aVmCl+RyVMAyxkLBZyLLLTOV+TU6B0jZyQUCDC7Sle4iZwOvWvktAUDRpcXi5wBvXPktASME76RkH6FRE6X3j1yysbDV1XNRM6MAlhOiXHDV/V8Rc6cAlgmzRjDLBVTF5kyCmCZFGMMMyujaiMkMoUUwDIpMyqilEZUu0lkKimA5aRmVEQpV/dikSmnAJYTml4RpSIW9nsYIgVJASzHVVsWpVLhK5I1CmA5ppqyCNNKFb4i2aQAlveoKo1QVRrxexgiBU8BLO8yrSRMTZnCV8QLCmCZUFkSprY86vcwRIqGAlgAKI+FmK7wFfGUAlgoj4aoq4j5PQyRoqMALnJl0RAzKnTkK+IHBXARK42EqKuIYoz6uIn4IWcD2BhzrTHmLWPMO8aYr/s9nkITCweZWanwFfFTTgawMSYI/CvwUeBc4HPGmHP9HVXhiIaDzKqMKXxFfJaTAQwsAd6x1rZaaxPAWuAGn8dUECKhALMrYwTUPl7Ed7kawHOB9iPud7jbJhhjVhpjthhjtnR1dXk6uHwVCQWYPa1E4SuSI3I1gI+VEPZdd6xdY61dZK1dNGPGDI+Glb/Gwzeo8BXJGbkawB1A/RH35wH7fBpL3lP4iuSmXA3gzcACY0yTMSYC3Ag87vOY8lI0HGSOwlckJ+VkmwNrbcoY8xfA74Eg8IC1dpvPw8o70XBQJ9xEclhOBjCAtfY/gP/wexz5KuYuNVP4iuSuXJ2CkDOg8BXJDzl7BCynpySiiyxE8oWOgAuIwlckv+gIuECURkKq7SCSZ3QEXAAUviL5SUfAea4sqpKSIvlKR8B5TOErkt90BJynyt1OFgpfkfylI+A8pPAVKQw6As4z5TE10BQpFDoCziMKX5HCoiPgPFERC6t7sUiB0RFwHlD4ihQmHQHnuMqSMNPLFb4ihUgBnMOmlYSpVfiKFCxNQeQoha9I4dMRcA6qKo1QUxbxexgikmU6As4xCl+R4qEj4BxSXRqhWuErUjQUwDmipixCVanCV6SYaAoiByh8RYqTjoB9VlsWZVpp2O9hiIgPFMA+qi2PMq1E4StSrDQF4ROFr4joCNgH0yuiVMYUviLFTkfAHlP4isg4BbCHFL4iciQFsEdqyxW+IvJuCmAP1JbphJuIvJcCOMu0zldEjkcBnEU1ZRGFr4gclwI4S3R5sYicjAI4C6pLFb4icnIK4ClWpZKSIjJJCuAppGLqInIqFMBTZFpJWOErIqdEATwF1EBTRE6HAvgMVSp8ReQ0KYDPQGVJmOkKXxE5TQrg01QRU/iKyJlRAJ+G8liIGRUKXxE5MwrgU1QeC1FXEfN7GCJSABTAp6A8qvAVkamjAJ6k8miIukqFr4hMHV8C2BjzT8aY7caYV40xvzTGVB3x2DeMMe8YY94yxnzEj/EdrSyqOV8RmXp+HQE/CZxvrX0/8DbwDQBjzLnAjcB5wLXAPcaYoE9jBJzwrauIYozxcxgiUoB8CWBr7RPW2pR793lgnnv7BmCttXbMWtsGvAMs8WOMAKURha+IZE8uzAF/Gfite3su0H7EYx3utvcwxqw0xmwxxmzp6uqa8kGVRkLMrFT4ikj2hLL1jY0x/wnMOsZDf2et/bX7nL8DUsBPxr/sGM+3x/r+1to1wBqARYsWHfM5p6skElT4ikjWZS2ArbUfOtHjxpibgeuAa6y14wHaAdQf8bR5wL7sjPDYSiJBZlXGFL4iknV+rYK4FvgacL21duSIhx4HbjTGRI0xTcACYJNX41L4ioiXsnYEfBL/AkSBJ92we95au8pau80Y8zPgDZypiT+31qa9GFAsHGRmhcJXRLzjSwBba886wWP/APyDh8MhGg4wrSRMIKDwFRHv+HUEnFNKI/pnEBHv5cIyNBGRoqQAFhHxiQJYRMQnCmAREZ8ogEVEfKIAFhHxiQJYRMQnCmAREZ8ogEVEfKIAFhHxiQJYRMQnCmAREZ8ogEVEfGION6PIX8aYLmD3FH/b6cChKf6e+UD7XXyKdd+93O9D1tprj95YEAGcDcaYLdbaRX6Pw2va7+JTrPueC/utKQgREZ8ogEVEfKIAPr41fg/AJ9rv4lOs++77fmsOWETEJzoCFhHxiQJYRMQnRR3AxpigMeYlY8xv3PtNxpgXjDE7jDGPGGMi7vaoe/8d9/FGP8d9Jowxu4wxrxljXjbGbHG31RhjnnT3+0ljTLW73Rhj/l93v181xlzi7+jPjDGmyhjzmDFmuzHmTWPMFYW+78aYc9yf9fjHgDHmvxX6fgMYY/7KGLPNGPO6MeanxphYrr3HizqAgbuAN4+4/4/AP1trFwC9wG3u9tuAXmvtWcA/u8/LZx+w1l50xBrIrwN/cPf7D+59gI8CC9yPlcBqz0c6te4GfmetXQhciPOzL+h9t9a+5f6sLwIuBUaAX1Lg+22MmQv8V2CRtfZ8IAjcSK69x621RfkBzMP5xfsg8BvA4FwVE3IfvwL4vXv798AV7u2Q+zzj9z6c5n7vAqYfte0tYLZ7ezbwlnv7PuBzx3pevn0AlUDb0T+3Ytj3I/bhw8AzxbDfwFygHahx37O/AT6Sa+/xYj4C/i7wP4CMe78W6LPWptz7HTg/RDj8w8R9vN99fj6ywBPGmK3GmJXutpnW2k4A93Odu31iv11H/pvkm2agC3jQnXb6njGmjOLY93E3Aj91bxf0fltr9wL/B9gDdOK8Z7eSY+/xogxgY8x1wEFr7dYjNx/jqXYSj+WbpdbaS3D+1PxzY8yyEzy3kPY7BFwCrLbWXgwMc/jP7mMppH3Hneu8Hnj0ZE89xra82293TvsGoAmYA5Th/M4fzdf3eFEGMLAUuN4YswtYizMN8V2gyhgTcp8zD9jn3u4A6gHcx6cBPV4OeKpYa/e5nw/izAUuAQ4YY2YDuJ8Puk+f2G/Xkf8m+aYD6LDWvuDefwwnkIth38EJnxettQfc+4W+3x8C2qy1XdbaJPAL4Epy7D1elAFsrf2GtXaetbYR58+yp6y1NwF/BD7lPu1m4Nfu7cfd+7iPP2XdyaJ8YowpM8ZUjN/GmRN8nXfv39H7/SX3zPjlQP/4n635xlq7H2g3xpzjbroGeIMi2HfX5zg8/QCFv997gMuNMaXGGMPhn3duvcf9niz3+wNYAfzGvd0MbALewflTLepuj7n333Efb/Z73Ke5r83AK+7HNuDv3O21OCckd7ifa9ztBvhXYCfwGs4ZZd/34wz2/yJgC/Aq8Cuguhj2HSgFuoFpR2wrhv3+X8B2nIOMHwHRXHuP61JkERGfFOUUhIhILlAAi4j4RAEsIuITBbCIiE8UwCIiPlEAS0EyxjQaY14/heevMsZ86STPucUY8y/HeexvT3WMIgpgEcBae6+19odn8C0UwHLKFMBSyILGmPvdmrBPGGNKjDEtxpjfucWINhpjFgIYY/7eGPM37u3Fbi3c54wx/3TUkfQc9+t3GGP+t/v8bwMlbr3dn3i/m5KvFMBSyBYA/2qtPQ/oAz6J04jxL621lwJ/A9xzjK97EFhlrb0CSB/12EXAZ4ELgM8aY+qttV8HRq1Td/emLO2LFKDQyZ8ikrfarLUvu7e3Ao04BVkedcoDAM7lqROMMVVAhbX2WXfTw8B1RzzlD9bafve5bwDzeXf5RpFJUwBLIRs74nYamIlTD/aiE3zNscoSnuh76j0kp01TEFJMBoA2Y8ynYaL/2YVHPsFa2wsMupXAwKmWNxlJY0x46oYqxUABLMXmJuA2Y8x4RbgbjvGc24A1xpjncI6I+yfxfdcAr+oknJwKVUMTOYoxptxaO+Te/jpOT7S7fB6WFCDNX4m818eMMd/AeX/sBm7xdzhSqHQELCLiE80Bi4j4RAEsIuITBbCIiE8UwCIiPlEAi4j45P8HENkXu+Z2GfQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 360x360 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.lmplot( x='height', y='width', data=df )\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Other Topics\n",
    "\n",
    "### ANOVA\n",
    "\n",
    "Analysis of variance is an optional topic your GB213 class may or may not have covered, depending on scheduling and instructor choices.  If you covered it in GB213 and would like to see how to do it in Python, check out [the Scipy documentation for `f_oneway`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.f_oneway.html).\n",
    "\n",
    "### $\\chi^2$ Tests\n",
    "\n",
    "Chi-squared ($\\chi^2$) tests are another optional GB213 topic that your class may or may not have covered.  If you are familiar with it and would like to see how to do it in Python, check out [the Scipy documentation for `chisquare`](https://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.chisquare.html)."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}