ACIC 2019.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Atlantic Causal Inference Data Challenge 2019\n",
    "[Official Site for the ACIC Data Challenge](https://sites.google.com/view/ACIC2019DataChallenge/data-challenge)\n",
    "\n",
    "Quoting directly from the challenge website, here is a description of the data provided:\n",
    "\n",
    "> Data are in the form (Y, A, V1, ... ,Vp), where \n",
    "> - Y is the outcome (either binary or continuous)\n",
    "> - A is a binary treatment indicator\n",
    "> - V1 through Vp are a mix of continuous, binary, and categorical pre-treatment covariates (no mediators). The number of covariates, p, varies across datasets.\n",
    "\n",
    "The goal of the challenge is:\n",
    "\n",
    "> Provide an estimate of the population average additive treatment effect (ATE) of a binary treatment on a binary or continuous outcome, and a 95% confidence interval. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import statsmodels.api as sm\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "import os\n",
    "\n",
    "pd.options.display.max_columns = 500\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's start by finding an easy dataset to work with by summarizing the target."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>testdataset2.csv</th>\n",
       "      <th>testdataset3.csv</th>\n",
       "      <th>testdataset1.csv</th>\n",
       "      <th>testdataset4.csv</th>\n",
       "      <th>testdataset5.csv</th>\n",
       "      <th>testdataset7.csv</th>\n",
       "      <th>testdataset6.csv</th>\n",
       "      <th>testdataset8.csv</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>500</td>\n",
       "      <td>525</td>\n",
       "      <td>500</td>\n",
       "      <td>500</td>\n",
       "      <td>625</td>\n",
       "      <td>500</td>\n",
       "      <td>475</td>\n",
       "      <td>500</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>17.9758</td>\n",
       "      <td>18.8849</td>\n",
       "      <td>0.3</td>\n",
       "      <td>0.268</td>\n",
       "      <td>0.3504</td>\n",
       "      <td>18.8151</td>\n",
       "      <td>17.5423</td>\n",
       "      <td>0.274</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>5.47596</td>\n",
       "      <td>12.6132</td>\n",
       "      <td>0.458717</td>\n",
       "      <td>0.443361</td>\n",
       "      <td>0.477477</td>\n",
       "      <td>12.4771</td>\n",
       "      <td>4.77874</td>\n",
       "      <td>0.446456</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>9.58147</td>\n",
       "      <td>-8.17043</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>1.09136</td>\n",
       "      <td>4.35973</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>13.9635</td>\n",
       "      <td>9.14433</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>9.45514</td>\n",
       "      <td>14.2522</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>16.9006</td>\n",
       "      <td>17.1066</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>16.3349</td>\n",
       "      <td>16.9655</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>20.6377</td>\n",
       "      <td>26.6696</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>26.4147</td>\n",
       "      <td>19.797</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>45.1715</td>\n",
       "      <td>89.7877</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>81.8532</td>\n",
       "      <td>42.8032</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>shape</th>\n",
       "      <td>(500, 24)</td>\n",
       "      <td>(525, 24)</td>\n",
       "      <td>(500, 24)</td>\n",
       "      <td>(500, 17)</td>\n",
       "      <td>(625, 24)</td>\n",
       "      <td>(500, 24)</td>\n",
       "      <td>(475, 24)</td>\n",
       "      <td>(500, 17)</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      testdataset2.csv testdataset3.csv testdataset1.csv testdataset4.csv  \\\n",
       "count              500              525              500              500   \n",
       "mean           17.9758          18.8849              0.3            0.268   \n",
       "std            5.47596          12.6132         0.458717         0.443361   \n",
       "min            9.58147         -8.17043                0                0   \n",
       "25%            13.9635          9.14433                0                0   \n",
       "50%            16.9006          17.1066                0                0   \n",
       "75%            20.6377          26.6696                1                1   \n",
       "max            45.1715          89.7877                1                1   \n",
       "shape        (500, 24)        (525, 24)        (500, 24)        (500, 17)   \n",
       "\n",
       "      testdataset5.csv testdataset7.csv testdataset6.csv testdataset8.csv  \n",
       "count              625              500              475              500  \n",
       "mean            0.3504          18.8151          17.5423            0.274  \n",
       "std           0.477477          12.4771          4.77874         0.446456  \n",
       "min                  0          1.09136          4.35973                0  \n",
       "25%                  0          9.45514          14.2522                0  \n",
       "50%                  0          16.3349          16.9655                0  \n",
       "75%                  1          26.4147           19.797                1  \n",
       "max                  1          81.8532          42.8032                1  \n",
       "shape        (625, 24)        (500, 24)        (475, 24)        (500, 17)  "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# summarize the low dimensional test datasets\n",
    "summary = pd.DataFrame()\n",
    "for f in os.listdir('TestDatasets_lowD'):\n",
    "    if f.endswith('csv') and not f.endswith('cf.csv'):\n",
    "        df = pd.read_csv(f'TestDatasets_lowD/{f}')\n",
    "        summary[f] = df['Y'].describe()\\\n",
    "                            .append(pd.Series([df.shape], index = ['shape']))\n",
    "        \n",
    "display(summary)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Of the 8 datasets, 4 are binary outcomes and 4 are continuous. Each dataset has roughly the same size (475-625 examples) and either 17 or 24 columns. Let's start with a regression case, so we don't have to deal with the added complexity of the logit (we'll get to that later). With that in mind, `testdataset3.csv` looks good."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Y</th>\n",
       "      <th>A</th>\n",
       "      <th>V1</th>\n",
       "      <th>V2</th>\n",
       "      <th>V3</th>\n",
       "      <th>V4</th>\n",
       "      <th>V5</th>\n",
       "      <th>V6</th>\n",
       "      <th>V7</th>\n",
       "      <th>V8</th>\n",
       "      <th>V9</th>\n",
       "      <th>V10</th>\n",
       "      <th>V11</th>\n",
       "      <th>V12</th>\n",
       "      <th>V13</th>\n",
       "      <th>V14</th>\n",
       "      <th>V15</th>\n",
       "      <th>V16</th>\n",
       "      <th>V17</th>\n",
       "      <th>V18</th>\n",
       "      <th>V19</th>\n",
       "      <th>V20</th>\n",
       "      <th>V21</th>\n",
       "      <th>V22</th>\n",
       "      <th>C</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.000000</td>\n",
       "      <td>525.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>18.884909</td>\n",
       "      <td>0.201905</td>\n",
       "      <td>0.294648</td>\n",
       "      <td>0.094381</td>\n",
       "      <td>0.109962</td>\n",
       "      <td>0.112857</td>\n",
       "      <td>0.090381</td>\n",
       "      <td>0.276038</td>\n",
       "      <td>0.067238</td>\n",
       "      <td>0.534876</td>\n",
       "      <td>0.112438</td>\n",
       "      <td>0.029410</td>\n",
       "      <td>0.077867</td>\n",
       "      <td>0.254495</td>\n",
       "      <td>0.171886</td>\n",
       "      <td>0.209486</td>\n",
       "      <td>1.712362</td>\n",
       "      <td>0.109067</td>\n",
       "      <td>0.827390</td>\n",
       "      <td>0.136914</td>\n",
       "      <td>0.112590</td>\n",
       "      <td>0.109010</td>\n",
       "      <td>0.548267</td>\n",
       "      <td>0.975067</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>12.613198</td>\n",
       "      <td>0.401804</td>\n",
       "      <td>0.522733</td>\n",
       "      <td>0.242423</td>\n",
       "      <td>0.337628</td>\n",
       "      <td>0.347767</td>\n",
       "      <td>0.298808</td>\n",
       "      <td>0.591365</td>\n",
       "      <td>0.206849</td>\n",
       "      <td>0.813346</td>\n",
       "      <td>0.368377</td>\n",
       "      <td>0.166184</td>\n",
       "      <td>0.340710</td>\n",
       "      <td>1.008381</td>\n",
       "      <td>0.534332</td>\n",
       "      <td>0.539255</td>\n",
       "      <td>1.802777</td>\n",
       "      <td>0.532904</td>\n",
       "      <td>1.248497</td>\n",
       "      <td>1.160294</td>\n",
       "      <td>0.395746</td>\n",
       "      <td>0.600167</td>\n",
       "      <td>1.730427</td>\n",
       "      <td>0.792503</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>-8.170434</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>9.144325</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.457425</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>17.106626</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.270000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>1.380000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.160000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.809151</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>26.669562</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.420000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.340000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.760000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.140000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.050000</td>\n",
       "      <td>2.610000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>1.350000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>1.323887</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>89.787655</td>\n",
       "      <td>1.000000</td>\n",
       "      <td>3.220000</td>\n",
       "      <td>2.300000</td>\n",
       "      <td>2.530000</td>\n",
       "      <td>3.330000</td>\n",
       "      <td>3.330000</td>\n",
       "      <td>5.260000</td>\n",
       "      <td>2.000000</td>\n",
       "      <td>6.250000</td>\n",
       "      <td>5.550000</td>\n",
       "      <td>2.320000</td>\n",
       "      <td>2.380000</td>\n",
       "      <td>20.000000</td>\n",
       "      <td>5.060000</td>\n",
       "      <td>4.870000</td>\n",
       "      <td>12.500000</td>\n",
       "      <td>6.320000</td>\n",
       "      <td>10.710000</td>\n",
       "      <td>15.430000</td>\n",
       "      <td>4.760000</td>\n",
       "      <td>12.500000</td>\n",
       "      <td>16.660000</td>\n",
       "      <td>7.005336</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                Y           A          V1          V2          V3          V4          V5          V6          V7          V8          V9         V10         V11         V12         V13         V14         V15         V16         V17         V18         V19         V20         V21         V22      C\n",
       "count  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.000000  525.0\n",
       "mean    18.884909    0.201905    0.294648    0.094381    0.109962    0.112857    0.090381    0.276038    0.067238    0.534876    0.112438    0.029410    0.077867    0.254495    0.171886    0.209486    1.712362    0.109067    0.827390    0.136914    0.112590    0.109010    0.548267    0.975067    1.0\n",
       "std     12.613198    0.401804    0.522733    0.242423    0.337628    0.347767    0.298808    0.591365    0.206849    0.813346    0.368377    0.166184    0.340710    1.008381    0.534332    0.539255    1.802777    0.532904    1.248497    1.160294    0.395746    0.600167    1.730427    0.792503    0.0\n",
       "min     -8.170434    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    1.0\n",
       "25%      9.144325    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.457425    1.0\n",
       "50%     17.106626    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.270000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    1.380000    0.000000    0.160000    0.000000    0.000000    0.000000    0.000000    0.809151    1.0\n",
       "75%     26.669562    0.000000    0.420000    0.000000    0.000000    0.000000    0.000000    0.340000    0.000000    0.760000    0.000000    0.000000    0.000000    0.140000    0.000000    0.050000    2.610000    0.000000    1.350000    0.000000    0.000000    0.000000    0.000000    1.323887    1.0\n",
       "max     89.787655    1.000000    3.220000    2.300000    2.530000    3.330000    3.330000    5.260000    2.000000    6.250000    5.550000    2.320000    2.380000   20.000000    5.060000    4.870000   12.500000    6.320000   10.710000   15.430000    4.760000   12.500000   16.660000    7.005336    1.0"
      ]
     },
     "execution_count": 48,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pd.read_csv('TestDatasets_lowD/testdataset3.csv')\n",
    "df['C'] = 1 # add a constant for regression intercept\n",
    "df.describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Assuming a randomized, controlled trial, we can estimate the treatment effect $Y_i$ as a linear function of $A_i$:\n",
    "\n",
    "$$Y_i = \\beta_0 + \\beta_1 A_i + u_i$$\n",
    "\n",
    "where $A_i = 1$ if the subject was treated and $A_i = 0$ if not."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                            OLS Regression Results                            \n",
      "==============================================================================\n",
      "Dep. Variable:                      Y   R-squared:                       0.082\n",
      "Model:                            OLS   Adj. R-squared:                  0.081\n",
      "Method:                 Least Squares   F-statistic:                     46.98\n",
      "Date:                Sun, 17 Mar 2019   Prob (F-statistic):           2.02e-11\n",
      "Time:                        22:49:29   Log-Likelihood:                -2052.6\n",
      "No. Observations:                 525   AIC:                             4109.\n",
      "Df Residuals:                     523   BIC:                             4118.\n",
      "Df Model:                           1                                         \n",
      "Covariance Type:            nonrobust                                         \n",
      "==============================================================================\n",
      "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "A              9.0124      1.315      6.854      0.000       6.429      11.596\n",
      "C             17.0653      0.591     28.884      0.000      15.905      18.226\n",
      "==============================================================================\n",
      "Omnibus:                      109.617   Durbin-Watson:                   2.001\n",
      "Prob(Omnibus):                  0.000   Jarque-Bera (JB):              280.253\n",
      "Skew:                           1.043   Prob(JB):                     1.39e-61\n",
      "Kurtosis:                       5.909   Cond. No.                         2.61\n",
      "==============================================================================\n",
      "\n",
      "Warnings:\n",
      "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n"
     ]
    }
   ],
   "source": [
    "X = df[['A', 'C']]\n",
    "y = df['Y']\n",
    "\n",
    "model = sm.OLS(y, X)\n",
    "fitted_model = model.fit()\n",
    "\n",
    "print(fitted_model.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Unsurprisingly, the $R^2$ is terrible. Our model estimated an average treatment effect (ATE) of 9.0124. Let's use our model to make predictions and test them against the counterfactuals in `testdataset3_cf.csv`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>ATE</th>\n",
       "      <th>EY1_i</th>\n",
       "      <th>EY0_i</th>\n",
       "      <th>EY1_i_pred</th>\n",
       "      <th>EY0_i_pred</th>\n",
       "      <th>ATE_pred</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>2.090099</td>\n",
       "      <td>8.550162</td>\n",
       "      <td>7.550162</td>\n",
       "      <td>26.077697</td>\n",
       "      <td>17.065253</td>\n",
       "      <td>9.012444</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>2.090099</td>\n",
       "      <td>28.928069</td>\n",
       "      <td>27.928069</td>\n",
       "      <td>26.077697</td>\n",
       "      <td>17.065253</td>\n",
       "      <td>9.012444</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>2.090099</td>\n",
       "      <td>29.485588</td>\n",
       "      <td>30.693588</td>\n",
       "      <td>26.077697</td>\n",
       "      <td>17.065253</td>\n",
       "      <td>9.012444</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>2.090099</td>\n",
       "      <td>14.936482</td>\n",
       "      <td>13.142482</td>\n",
       "      <td>26.077697</td>\n",
       "      <td>17.065253</td>\n",
       "      <td>9.012444</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>2.090099</td>\n",
       "      <td>18.028003</td>\n",
       "      <td>17.028003</td>\n",
       "      <td>26.077697</td>\n",
       "      <td>17.065253</td>\n",
       "      <td>9.012444</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        ATE      EY1_i      EY0_i  EY1_i_pred  EY0_i_pred  ATE_pred\n",
       "0  2.090099   8.550162   7.550162   26.077697   17.065253  9.012444\n",
       "1  2.090099  28.928069  27.928069   26.077697   17.065253  9.012444\n",
       "2  2.090099  29.485588  30.693588   26.077697   17.065253  9.012444\n",
       "3  2.090099  14.936482  13.142482   26.077697   17.065253  9.012444\n",
       "4  2.090099  18.028003  17.028003   26.077697   17.065253  9.012444"
      ]
     },
     "execution_count": 73,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "test_results = pd.read_csv('TestDatasets_lowD/testdataset3_cf.csv')\n",
    "\n",
    "test_results['EY1_i_pred'] = fitted_model.predict(pd.DataFrame({'A':np.ones(df.shape[0]),\n",
    "                                                                'C':np.ones(df.shape[0])}))\n",
    "test_results['EY0_i_pred'] = fitted_model.predict(pd.DataFrame({'A':np.zeros(df.shape[0]),\n",
    "                                                                'C':np.ones(df.shape[0])}))\n",
    "test_results['ATE_pred'] = test_results['EY1_i_pred'] - test_results['EY0_i_pred']\n",
    "\n",
    "test_results.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5,1,'Counterfactual vs. Predicted with Treatment')"
      ]
     },
     "execution_count": 99,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3gAAAEWCAYAAAA0DzVNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzsvX2cnVV57/1bM7OT7ABmAsQXRgOo\nNfRgCpGonOJphXqatijmCQpV2+ppq8/p09MW5KQNrZ8Clpr0pAhtT9XSx0d9KmJQ6BSkNbYFbKWC\nJ2ES0xyJtSCBATEYBkkyJPOyzh/3vSZr1l6v98vst9/38wkke9/7vte93q51Xeu6riWklCCEEEII\nIYQQ0v0MtLsAhBBCCCGEEEKqgQoeIYQQQgghhPQIVPAIIYQQQgghpEeggkcIIYQQQgghPQIVPEII\nIYQQQgjpEajgEUIIIYQQQkiPQAWPeBFCrBJCjAkhnhdC/Ga7y6MQQtwnhPjVdpejKoQQnxZCXJ//\n/T8JIfYt0HOlEOLVJe+xVwjxZs/3PdVWhBCSCmXpwkBZSkgGFbwFQAjxbiHEDiHEISHEU0KIvxNC\nvGkBnlt6wgHw2wDuk1KeJKX80xJlmZt0uxEhxPuEEDN5G/5QCLFLCPHWOp4lpfxnKeWqyDJ9rY4y\npCClPFtKeR8ACCGuFUJ8diGeK4R4sxDiCc/3f5e31yEhxJQQ4pj270/UUJ5fFULcV/V9E57/hG9x\nQEi3Q1lKWZpCv8pSIcR7NFk3KYSY1f59qNJCZ88bysfIGVXfO/L51wshPt2OZ3cyVPBqRgjxQQA3\nAfgIgJcAWAngYwDe3s5yhRBCDOV/PR3A3naWpYP4upTyRADDAD4J4DYhxMnmRVrdkTYipfxZKeWJ\neZvdAuB/qH9LKf+reT3bjZDOhbK0p6AsrREp5S2a7PtZAE9qsu9E83rWc48ipeSfmv4AWAbgEIB3\neq5ZjExoPZn/uQnA4vy79wH4mnG9BPDq/O+fBvDnAO4G8DyABwG8Kv/un/JrD+dluDz//K0AdgGY\nAPAvAH5Mu/d3AfwOgG8COArgHgAzAF7I7/EaABcDGAPwQwCPA7jWKN+b8vtO5N+/D8AHAEwBOJbf\n5y7zXbT3uT7/+3IAXwJwAMCz+d9frl17H4BftdTnaQAmAZysfbYGwDMAGgBeDeCrAJ7LP9sW2Zbz\n2gLACXn51wJ4M4An8rr7HoC/iqjrNQAeytttG4DPa+/+ZgBPaNe+AsAdeV38AMD/BPCjebvM5HU6\nofWnPwawH8DTAD4BoKndayOAp5D1tV8220C77kIAe7R//wOAb2j//hqA9Vq/eQuAn8nbeCov026t\nrf4AwP35+34FwKnavS5BtvCZyK/9UVt/1/tIXv+TAGbzZx0CcJqn/T6t6lf77C152X83b7dPaeXZ\nnZfnawBeq/3mQwAeyd9jL4BL8s9XG+3xTP75Z/P22p5//k/IFqd/lt//WwDO0e7/cgB/nbf1owB+\nXfvuegC35vd8HsC/Anhd/t2teV1M5s/5YLvnP/7hn6r+gLKUspSyVLWVU5Y66nteHWifP5G/wx4A\nx/LPfPLnPwJ4IG+DpwD8KYBG/t2/YP4YuRTH5evV+f2eBPC2vC3/DcBBAL+t3X8AmSz+97w/fR7A\n8vy7V+f3/6W83AcAbNL6hl5XO9s9X3XKn7YXoJf/5IN0GsCQ55oP54PmxQBW5APlD/Lv3oewUDoI\n4A0AhpDtUnzedm3+79cB+D6ANwIYBPDefAAqIfhdZJPoK5BPZDAm/3yyWJ0Pxh9DNvGpyWllPum8\nC5kAOAXAuVpZzQW2Tyidkk8SSwGcBOALAEa1a+eVy7jvPQDer/17K4BP5H+/FcDv5eVfAuBNkW05\n1xZ5Xf9W/q7L8jqZBvBHyIRC01fXABYBeAzAlXk9vSOfnFqEUv7b3QBuRCYI58rs6B83AbgTwMl5\nvd0FYLPWH58G8Nr8Xp8z20C7zxJkwv3U/H2/h2yCPil/v0kAp2j95i35368F8FnjXvchm7Rfk//2\nPgBb8u9eg0wo/Oe8Ln4bwHcALIroI3P1FNF+c7/TPntL3m4fydukCeD1eR29Pq/7X87LrspzGYCX\nIes/70YmUF6Sf/eryFyw9Gd8Nu8Ha/I6/Soywfnu/P5bAPy91ta7kAm5RciE2ncB/FT+/fV5va/L\nr92K+QulJwC8ud3zHv/wT9V/QFlKWUpZqtrKKks99T1XB8bnTwDYiUypayIsf16ft8EQgFcC+DaA\n/6a1owRwhnZ/JV9/L2+bX8vb8bMATkTW518AsDK//r8jU1xH8jr7JI4r+ErB+0T+3euQGU5+JP/+\negCfbvc81Wl/6KJZL6cgs+RPe655D4APSym/L6U8AOA6AL+Y8Iw7pJTfyJ9xC4BzPde+H8BfSCkf\nlFLOSCk/g2yQnK9d86dSysellJO2G0gp75NS7pFSzkopv4lskv9J7V3+QUp5q5RySkr5AynlroR3\n0Z/zAynl7VLKI1LK5wH8ofacEJ9DJhghhBAAfj7/DMgm/9OR7fa8IKVM8bs/XwgxgWyCfheA/0tK\n+Vz+3SyAa6SUR/O689X1+cgmvJvyevoigP/leOYbkFlSN0opD/vKnL/r+wFcKaU8mNfbR/L3BzLl\n5FNSyn+VUh5GJkCsSClfALADwE8gs6x+E5ml8YK8/P8mpfyBv7rm8Skp5bfzurkNx/vp5QDullL+\nvZRyCpnFtAngxxPuXYZpZJbzY3nZPgDgY1LK/5W32/+XX/d6AJBS3ialfCrv/59DJgDXBp5xu5Ry\nLK/TUQCHpJSfk1LOILM4r8mvOx/Ai6SUH8nL8x1kQu7ntXt9VUq5Pf/tX8E/3gnpFShLKUspSzNc\nsrQIfyKlfCK/l1f+5DLxQSnltJTyEQA3I9yPXkCmgE4h25FbAeBGKeWhvM/vQ6boAcD/DeB3pZTj\neZ1dC+AyIYSup1ybt9tDyDxozinx7j0P/W7r5QcAThVCDHkE02nILFCKx/LPYvme9vcjyCwjLk4H\n8F4hxG9ony0ynve472FCiDci23V4bf7bxcgsgkBmrfz3uGL7EUIsRWZp+xlkLiYAcJIQYjBf3Pr4\nIoA/E0KcBuBHkFl+/jn/7reRuTh8QwjxLIAbtEV8iAeklK6A/gP5pKTw1bUEMC5lZnrK0fuAzisA\nPBZY2ChWILPS7szkEwBAILPMIX/2zohnKr6K4y4zX0Xm3vOTyITrVyPKo+Pqp/P6v5RyVgjxODIr\n3kLwtJTymPbv0wG8RwhxpfbZIlUeIcT7kFmLT8+/OxGZZdb7DO3vk5Z/q7o4HcDKfOGjGERmpVWY\n9XhC4NmE9AKUpQWhLJ2j12VpEfQ+6pU/QoizANwA4DxkdTOEzJXZxzNaH1OGDpf8WwngLiHErPa9\nRLYjn/1DyirfvefhDl69fB2ZBWO955oncXyxCGSd/Mn874eRDSQAgBDipSXL8ziAP5RSDmt/lkop\nb9Wuka4f53wOmdvCK6SUy5BtmasZ8HEAr3L8znbfI9DeD4D+flcBWAXgjVLKFyGzfkF7lhMp5QQy\n3/TLkLnC3aoEgJTye1LK90spT0NmMfpYBdnRgNb389X1UwBGhCY5kLW7jceRTbo2Y4z5zGeQTZhn\na89cJo8HVT+FTMiFnqlQQukn8r9/FZlQ+km4hVKo/5jM6/95nbwCwHj+ka+PpD7Lhq3drrO0221C\niFcC+DgyV5NTpJTDAB7G8T5ZtjyPI7Pm6s8+SUr5toLvQkivQFnqvy9l6XH6VZYWQX9GSP78BbK4\n71fn/ej3UZ3sAzLl9z8bz19iKHUx70FyqODVSO5y8PsA/lwIsV4IsVQI0RBC/KwQ4n/kl90K4ENC\niBVCiFPz61Vq3N0AzhZCnCuEWAKPG4CDp5H5Siv+EsB/FUK8UWScIIS4WAhxUsI9TwJwUEr5ghDi\nDcgmfcUtAN4ihLgsT5t7ihBCuQ+YZQEyf+93CyEGhRA/g/nb/Schm2An8uxa1ySUEciE5y8hiz1Q\nLiUQQrxTCPHy/J/PIpsYQlbMIvjq+uvIXAN/M6+nDcjcR2x8A5kw2ZLfY4kQ4oL8u6cBvFwIsQjI\ndr/y594ohHhx/r4jQoh1+fW3AXifEOI/5FbdUJ3+C7KFwRuQBYXvRbaAeiOyxAM2ngZwhuFW4eM2\nABcLIX5KCNFAthg5mj8b8PeRpwGcIoRYFvmsGG4G8OtCiNfn7XaiEOJtQogTkFkLJbIAbyGy84jO\nMsrz8vw9ivB1AMeEEFfl7TwohFgthDgv8ve2MUZI10NZSllKWVo7IflzErKEOoeFED+KTKkHAOS7\ndD9AOfnzCQAfEUKsBAAhxIuFEJdE/lbVVdBo0U9QwasZKeVHAXwQWfa9A8isJP8NWSwOkAWH7kDm\nl70HWTao6/PffhtZ4Pg/IMs6lHpOy7UAPiOEmBBCXCal3IHMr/x/IpuQv4MsuDiF/wfAh4UQzyMT\noLepL6SU+wH8HLJF+kFkQkf5SH8SwH/Iy6Le/beQZVWaQBZzoD4HsgDnJjJL2gMAvpxYzjuRuZQ8\nLaXcrX3+egAPiuwsmDsB/JaU8lFg7pDR9yQ+x4qvrnOXwA35v59FFod2h+M+M8jq6NXIsnk9kV8P\nZAHwewF8TwjxTP7Z7+TPekAI8UNkfWdVfq+/Q1av9+TX3BN4h8PI+uNezY3x68jcXL7v+JlyMfqB\nEOIh3/3zZ+wD8AvIsko+k7/r27TnOfuIlPJhZIu6R/J+leKO5SrPg8h26D6OrG2+nZcPeczAn+L4\nQuEszHdR+Xtk4/RpIUSM1dF89jSy8fMGZLF9zyCzmr4o8hYfAXBdXhdXpD6fkE6GspSyFJSltREh\nf65CluDm+fzzbcYtrgHwubxfbihQhI8i65v/mI+Jf0Ee+x7BNmRuuweFEN8o8OyeRMx3XSaEEEII\nIYQQ0q1wB48QQgghhBBCegQqeIQQQgghhBDSI1DBI4QQQgghhJAegQoeIYQQQgghhPQIXXHQ+amn\nnirPOOOMdheDEEJIzezcufMZKeWKdpejW6B8JISQ/iFWRnaFgnfGGWdgx44d7S4GIYSQmhFCPNbu\nMnQTlI+EENI/xMpIumgSQgghhBBCSI9ABY8QQgghhBBCegQqeIQQQgghhBDSI1DBI4QQQgghhJAe\ngQoeIYQQQgghhPQIVPAIIYQQQgghpEeggkcIIYQQQgghPQIVPEIIIYQQQgjpEbrioHNCCCGEEFIt\no2Pj2Lp9H56cmMRpw01sXLcK69eMtLtYhJCSUMEjhBBCCOkzRsfGcfUdezA5NQMAGJ+YxNV37AEA\nKnmEdDlU8EhPQ+skIYQQ0srW7fvmlDvF5NQMtm7fRzlJSJdDBY/0LLROEkIIIXaenJhM+pwQ0j0w\nyQrpWXzWSUIIIaSfOW24mfQ5IaR7oIJHehZaJwkhhBA7G9etQrMxOO+zZmMQG9etalOJCCFVQQWP\n9Cy0ThJCCCF21q8ZweYNqzEy3IQAMDLcxOYNqxnCQEgPwBg80rNsXLdqXgweQOskIYQQoli/ZoQK\nHSE9CBU80rMoocUsmoQQQgghpF+ggkd6GlonCSGEEEJIP8EYPEIIIYQQQgjpEajgEUIIIYQQQkiP\nQAWPEEIIIYQQQnoEKniEEEIIIYQQ0iMwyQohhBBCCCGkFKNj48xc3iFQwSOEEEIIIYQUZnRsfN7Z\nw+MTk7j6jj0AQCWvDVDBI4SQSGidJIQQQlrZun3fnHKnmJyawdbt+ygn2wAVPEIIiYDWSUIIIcTO\nkxOTSZ+TemGSFUIIicBnnSSEEEL6mdOGm0mfk3qhgkcIIRHQOkkIIYTY2bhuFZqNwXmfNRuD2Lhu\nVZtK1N9QwSOEkAhonSSEEELsrF8zgs0bVmNkuAkBYGS4ic0bVjOEoU0wBo8QQiLYuG7VvBg8gNZJ\nQgghRLF+zQgVug6BCh4hhESghBazaBJCCCGkk6GCRwghkdA6SQghhJBOp9YYPCHElUKIvUKIfxVC\n3CqEWCKEOFMI8aAQ4t+EENuEEIvqLAMhhBBCCCGE9Au1KXhCiBEAvwlgrZTytQAGAfw8gD8CcKOU\n8kcAPAvgV+oqAyGEEEIIIYT0E3Vn0RwC0BRCDAFYCuApABcB+GL+/WcArK+5DIQQQgghhBDSF9Sm\n4EkpxwH8MYD9yBS75wDsBDAhpZzOL3sCgDWgRQjxASHEDiHEjgMHDtRVTEIIIaSroHwkhBDio04X\nzeUA3g7gTACnATgBwM9aLpW230spb5ZSrpVSrl2xYkVdxSSEEEK6CspHQgghPup00XwLgEellAek\nlFMA7gDw4wCGc5dNAHg5gCdrLAMhhBBCCCGE9A11Knj7AZwvhFgqhBAAfgrA/wZwL4B35Ne8F8Df\n1FgGQgghhBBCCOkb6ozBexBZMpWHAOzJn3UzgN8B8EEhxHcAnALgk3WVgRBCCCGEEEL6iVoPOpdS\nXgPgGuPjRwC8oc7nEkIIIYQQQkg/UvcxCYQQQgghhBBCFggqeIQQQgghhBDSI1DBI4QQQgghhJAe\ngQoeIYQQQgghhPQIVPAIIYQQQgghpEeggkcIIYQQQgghPQIVPEIIIYQQQgjpEajgEUIIIYQQQkiP\nQAWPEEIIIYQQQnoEKniEEEIIIYQQ0iNQwSOEEEIIIYSQHoEKHiGEEEIIIYT0CFTwCCGEEEIIIaRH\noIJHCCGEEEIIIT0CFTxCCCGEEEII6RGo4BFCCCGEEEJIj0AFjxBCCCGEEEJ6BCp4hBBCCCGEENIj\nUMEjhBBCCCGEkB6BCh4hhBBCCCGE9AhU8AghhBBCCCGkR6CCRwghhBBCCCE9wlC7C0AIIYQQQki7\nGR0bx9bt+/DkxCROG25i47pVWL9mpN3FIiQZKniEEEIIIaSvGR0bx9V37MHk1AwAYHxiElffsQcA\nqOSRroMKHiF9Bi2UhBBCyHy2bt83p9wpJqdmsHX7PspI0nVQwSOkj6CFkhBCCGnlyYnJpM8J6WSY\nZIWQPsJnoSSEEEL6ldOGm0mfE9LJUMEjpI+ghZIQQghpZeO6VWg2Bud91mwMYuO6VW0qESHFoYJH\nSB9BCyUhhBDSyvo1I9i8YTVGhpsQAEaGm9i8YTXDF0hXwhg8QvqIjetWzYvBA2ihJIQQQoBMyaNC\nR3qBWnfwhBDDQogvCiEeFkJ8SwjxH4UQJwsh/l4I8W/5/5fXWQZCyHFooSSEEEII6W3q3sH7EwBf\nllK+QwixCMBSAL8L4B+llFuEEJsAbALwOzWXgxCSQwslIYQQQkjvUtsOnhDiRQB+AsAnAUBKeUxK\nOQHg7QA+k1/2GQDr6yoDIYQQQgghhPQTdbpovhLAAQCfEkKMCSH+XyHECQBeIqV8CgDy/7/Y9mMh\nxAeEEDuEEDsOHDhQYzEJIYSQ7oHykRBCiI86FbwhAK8D8HEp5RoAh5G5Y0YhpbxZSrlWSrl2xYoV\ndZWREEII6SooHwkhhPioU8F7AsATUsoH839/EZnC97QQ4mUAkP//+zWWgRBCCCGEEEL6htoUPCnl\n9wA8LoRQ+dd/CsD/BnAngPfmn70XwN/UVQZCCCGEEEII6SfqzqL5GwBuyTNoPgLgvyBTKm8TQvwK\ngP0A3llzGQghhBBCCCGkL6hVwZNS7gKw1vLVT9X5XEIIIYTUx+jYOLZu34cnJyZx2nATG9et4vEr\npC2wL5JOoxP6ZN07eIQQQgjpIUbHxnH1HXswOTUDABifmMTVd+wBAC6syYLCvkg6jU7pk3UmWSGE\nEEJIj7F1+765xYticmoGW7fva1OJSL/Cvkg6jU7pk9zBI4QQQkg0T05MJn1OSF20sy92ghse6Tw6\nZX6kgkcIIYSQaE4bbmLcslg5bbjZhtLUAxfv1VJXfbarL3aKGx7pPDplfqSLJiGEEEKi2bhuFZqN\nwXmfNRuD2LhuleMX3YVavI9PTELi+OJ9dGy83UXrSuqsz3b1xU5xwyOdR6fMj1Tw+pjRsXFcsOUe\nnLnpblyw5R4KL0IIIUHWrxnB5g2rMTLchAAwMtzE5g2re2bngov3aqmzPtvVFzvFDY90Hp0yP9JF\n00Gvu2fQvYAQUpRenx9JmPVrRnq2zbl4r5a667MdfbFT3PBIZ9IJ8yN38Cz0g3sGLZSEkCL0w/xI\n+hvXIp2L92L0Yn12ihseIS6o4FnoB+WHFkpCSBH6YX4k/cvo2DgOH51u+ZyL9+L0ojLUKW54hLiI\nctEUQrwGwMcBvERK+VohxI8BuERKeX2tpWsT/aD80L2gu6BLHOkU+mF+JMXo9nnKDF1QLF/awDVv\nO7ur3qWTUPXWzX3DRie44RHiIjYG7y8BbATwFwAgpfymEOJzAHpSwet15YcWyu6C8ZKkk+j1+ZEU\noxfmKdvuNAAsXTTUNe/QqVAZIp1OtxuoTGJdNJdKKb9hfNaqIfQIPneCbs88qYTwxOTUvM+XL23Q\nvaBDoUsc6SR60d2KlKcX5inuThPSnZRdm/dibHnsDt4zQohXAZAAIIR4B4CnaitVm3G5EwCghZIs\nOFx0kE6iV92tSDl6YZ7i7jQhGd20m1WF94DPQNWp7x0iVsH7dQA3AzhLCDEO4FEAv1BbqToAmzvB\nBVvu6foO0AtCuN/gooN0mrCluxUx6YV5auO6VS0xeNydrpZOm8tIK93mbl2FctaLa+MoF00p5SNS\nyrcAWAHgLCnlm6SU3621ZB1IL3SAXkxX3OvQJa6/6UXXEdJ79MI8xcyI9cK5rDvoNnfrKtbmvbg2\njs2iOQzglwCcAWBICAEAkFL+Zm0l60BooSTtgC5x/U0vuo6Q3qNX5inuTtcH57LuoNs2M6pYm/fi\n2jjWRfNvATwAYA+A2fqK09n0QgfoFSHcb3DR0b90m7Al/QvnKeIjZS6jK2f76LbNjCrW5r24No5V\n8JZIKT9Ya0m6gF7pABTChHQP3SZsCSHERuxc1m0xYL1Gt21mVLU277W1cayC91dCiPcD+BKAo+pD\nKeXBWkrVwfRaByDVQ8sjqZJuE7aEEGIjdi6jK2d76cbNDK7NW4lV8I4B2Arg95AflZD//5V1FIqQ\nboWWR1I13ShsCSHEJHYuo1t6+6HC1P3EKngfBPBqKeUzdRaGkG5mdGwcV922GzNSzvuclkdSFgpb\nQshCU4c3SsxcRrd0QsoTq+DtBXCkzoIQ0s2onTtTuVOYlke6cZJ2wH5HuhH224Wnnd4odEsnpDyx\nCt4MgF1CiHsxPwavr45JIMSFLWZAR7c80o2TtAP2O9JpxChu7LftoZ1xcHRLJ71GO4xUsQreaP6H\n1ACtk92PLzbAtDwygJy0A/Y70knEKm7st+2h3XFwdEsnvUK7jFRRCp6U8jNCiEUAXpN/tE9KOVVb\nqfoIWid7A1fMwKAQ2Lxh9by2bLfgJBn9ZlhhvyOdRKzixn7bHhgHR0g1tMtINRBzkRDizQD+DcCf\nA/gYgG8LIX6itlL1Eb6GJ93DxnWr0GwMzvus2RjEDZed0zKAXQKSgnPhUIaV8YlJSBw3rIyOjbe7\naLXBfkc6iVjFrWy/HR0bxwVb7sGZm+7GBVvu6ekxXiUumcY4OPYpYsfVL9plpIp10bwBwE9LKfcB\ngBDiNQBuBXBeXQXrF2id7A1SYgYYQN5++tHti/2OdBKxO0Sx/da2Iw+AHjIFYRycHXpddT7t8M7x\n9Yt27YbHKngNpdwBgJTy20KIRk1l6ivoBtE7xMYMdLrg7AfXxX40rHR6vyP9RaziFtNvXYurJY2B\nvjPkVAnj4FrpR+NgN9EuBdzXL9plXI1V8HYIIT4J4K/yf78HwM56itRf0KrePtqpyHSq4OwX62S/\nGlY6td+R/iPF4BDqt67FlSuzcacYcvrBmNZr1G0crKJP9HO/apcC7usX7TKuxip4vwbg1wH8JgAB\n4J+QxeKRktCq3h76RZFJpV+skzSsELJwuBac6o/6/sptu+Ys3inzTeriuhMMOTYZdOW2Xbhi2y6M\naC6mXBt0FiHjoE+5CileVaxL+n1t0y7vnFC/aIdxNVbBGwLwJ1LKjwKAEGIQwOKYH+bX7gAwLqV8\nqxDiTACfB3AygIcA/KKU8lhyyXsIWtUXnn5RZFJxTYLjE5MYHRsvXTedYlmkYYWQ8lRxjl0VC1LX\n4mq42cDR6dm2G3Js9WSTQTL///jEJDZ+YTcggKkZOfdZPy3UOxWfcdDXl4FwPGgV6xLXPa69c29f\nyLuF8s4xx/SFZ63A7TvH2z7X6MQqeP8I4C0ADuX/bgL4CoAfj/jtbwH4FoAX5f/+IwA3Sik/L4T4\nBIBfAfDx6BKTnqJdC/5+jMGKwTU5Aii9uKjSslhFv6FhhZDiVHWOXRWLWtei+9pLzp4rQ7sWtq56\ncrmPKqZmZctnvW6ErHM9kHpv364zYO9TF2y5x5sVPdTPq1iXuK6dmJzCxGR2ulkvGwsWwjvHNqZv\n3zmOS88bwb0PH+gYJTpWwVsipVTKHaSUh4QQS0M/EkK8HMDFAP4QwAeFEALARQDenV/yGQDXggpe\nX9JOV4J+jcEKYZscFWUXF1Xtmva7C8pC4MpI2A8WYBJHVefYVbGoDe3It7OfuuppUAjMyFYlLoSq\nl6LKUKd4UdjKVde8nnrv0PUu42CRvqx/51qXLGvG5zT0GWl1etVYsBDeOa4xfe/DB3D/posqe05Z\nYhW8w0KI10kpHwIAIcR5AGJm35sA/DaAk/J/nwJgQko5nf/7CQDWWhdCfADABwBg5cqVkcWsliqD\nXccnJucm9JEC9+rUSTmE7/3b6SbJGCw7qt6v2LbL+n2ZHc6qdk3pXlsvtsVNyF2sW+enbqXd8nF0\nbNy5iLSdY+czplVlbOvUHXnX/DYjJZqNweBOnslpw83CytBCGceKzAd1zuup9y5allBfDvXzjetW\nYeMXdrfs3h4+Nh0dIuEz0pr0qsdS3XNBt3iARR10DuAKAF8QQvyzEOKfAWwD8Bu+Hwgh3grg+1JK\nPdumsFxqNWFJKW+WUq6VUq5dsWJFZDGro4qDkPV7AJiz1qXeqxsOZbYd8Bh6/9gFQh2sXzOCzRtW\nY2S4CQFgZLiJzRtWd+QCYaFZv2YEIzUcil3VQdu+OMF+OHy27kN2bYubqVk5p9wp1IKnG+anXqOd\n8lG1twvbOXa+A7M76UDtqsbmqEUsAAAgAElEQVSWfp8BYVv2HJc5aq41r2oMCDQG53+q6sWngPgo\n+rsUis4HdS6aU+9dtCy+vhzTz9evGcGJS1r3XaZmZHQb2dY2y5fadwC7wWOpEw+Vr2otUzexO3jf\nBHAWgFXI5qGHEVYOLwBwiRDi5wAsQRaDdxOAYSHEUL6L93IATxYpeN3UFexa5F6hsrTbep5yBpFe\nfpeLykINkk61+HYCdexwVnVPl5VU4LiFtFfdNhfCAp8a78Ed1f7CJ9eKnGPXKQmPqhpb5n1sMk7V\nky6DUtyiryzoYbEQOw917X6VIfXeRcsS05dD/XziyJT13qkuy77snEB3eCx1ajhGt3iAxSp4X5dS\nvg7Av6oPhBAPAXid6wdSyqsBXJ1f+2YA/11K+R4hxBcAvANZJs33AvibYkWvl9SJ0DY5F51sU8rS\nCQMg9Qwihc1FpRMHSVW0WxFPwSWoAOCCLfcUegd13bV37p0L9l7SiHUiOI5tchVodQXoRSVjIZSp\n2BgOdW23uKuQavC1q8sLImRM6wRjW1Vjy6UADwqBWSmd86arDmyfFVVAFiL2vMzuV12L5tR7lymL\nry/H9PM62qhTjCipdKrxsFvq06vgCSFeiixGrimEWIPjXgQvAhBMsuLgdwB8XghxPYAxAJ8seJ9a\nSRlkLiVreGkDzzqsMa57pZalEwZA0YWcHotX1SDpVCWqExRxW5l8dRWyAtreIab+j07Pzv392SNT\nyfVgm1zb6e67kCyEMmVb3DQGxLwYPGC+uxgTFvUOoTHsGm8jw80Fm8vqmOerGluu62elxKNbLk4u\nl42iCshC7DxUuft14Vkr5s5HLNPOqQvydi7g62qjTjCipNLJxsNuqM/QDt46AO9D5kr5Ue3z5wH8\nbuxDpJT3Abgv//sjAN6QUMa2kDLIXErW4qEBZxB1yoD1laWoq0aVpJxBpLC5qJRlIZWo1AVGJyji\nOkXqKsZVeCHO+VH306+/YMs9faFkLIQF3rd76+rz3eCuQsLEjOF2uyfVNc9XNbbaOUZD778QiktV\nu19Vt3PqWqNdC/hu2R1aCJjtvBxeBU9K+RkAnxFCXCqlvH2BytQRpAwylzL13OQUbrz83NJZNH1l\n6QTrecwZRGWziMawUEpUEcHTaZaoInUVeoeYe9ZVD+1edC4UC/WeKe5iXJD0DjFjuN3tXdc8X9XY\navcYret3KfcHyvePMu3cKZ48Rcuh2kj9/sptu7B1+76+m1fLjqVO6QftIjYG77VCiLPND6WUH664\nPB2B2SluvPzcwsHBVU2mrvu0a2Fr1pHvgMd2u4pWrUQVETydZokqUleud5Bw76CZ96yrHtq96CxK\nqgDq1PfsBncVEiZ2Xmhne9c1z1c1tjp1jC4kRfqHORcWdbvvlHCIsuWoMumPqtdlzQaEyBK5lOmX\nC6U4lRlLndIP2kmsgndI+/sSAG8F8K3qi9N+inSKdu4epAyAqgalrY5u3zne9mMGFkqJKrLAcPWR\nC89aUThpSRmK1JXvfJ3xiUlrohPznq57HEk458dFtykZRQVQt70n6R46zRBlw1XGASFw5qa7S8dq\n1WmQ7RYWeufDNhfGyBMbnRIOUbYcVbyHWa8quRlQXZbYuhWnomOpU/pBO4lS8KSUN+j/FkL8MYA7\naylRmynSKaq22BWx6se4LFQ1KBdi4BQRMLGKdlnhVWQRZOsjF561ArfvHK+kTVLfqYhRQn8H106e\nKZRt5/wA8zNpAsWSrQDd7YJBAUQ6jW5wdXYZifRzVq/ctgtXbNtVWzjAQrOQ81w7dj5sc2GMPLFR\n1w5vTBvo11gPeE4oRxXv4TvSBKguS2wnHdul8J2X6zMmd0r5qyB2B89kKYBXVlmQTqHooKrKYlfX\n5FrlYrJuV8gyOxuAX9Guon6LLoLMPnLBlnsqaZMi71TEKKFPfC4ksmx6vnuq2FFdwQPKWye7zQWj\n0+IyCekG90KzjAOW81TVv7ptTnCdhRc7z6nfx8S8uxay7TA8uea8GHliUscudGwGaZeHS5FyVPEe\nMbJEHbcVO+Y7/dguhc/N1zd+OqX8VRCl4Akh9uD4nDkA4MUA/qCuQrWTdruo1GUdqXIxWXcdlREw\nIUW7CuFV1SKoqjYp+k4pRolY4TUy3MT9my4K3q8u62Q37YCVHUe9ZGkknUM73Qtj+7RexjM33e29\nZ7fMCa7F5ZLGQNQ85zpgXd/RVErfcLOBw8em54490Rey7TA8+Y7eiJEnwHzltsjOn+t+LiOC2Qah\n3bLUclSxmx5zpumyZiNJqenEY7ts84YvpMRVppjyu+aoTpTHsTt4bwWwHMB/AjAM4G+llDtrK1Ub\nabeLSl3WkSqVsrrrqE4BU9W9q1gEVdUmCyGQqxZedVonF8IVpwrKjKNeszQSUrRPxyxiu2FX3LW4\ndM275jv55milmiglxfSeUM/aun1fW4zcVWRL1H+vu3cWcdN1Kcsmehv4+pgAkmVJkfwK+s7tcLOB\nY9NhmS0EkpSyTju2yzVvbN6wGps3rMYVCWUKrSlcz9rx2MHKwm2qJFbBezuA9wO4A1lf/ZQQ4i+l\nlH9WW8naRN0uKkUPkY21jrjuX6VSVncduepgWbMxl5CkaDaodu3QxlqYYtrEvNfw0gaePdIqsKt8\np6qFV53WyQEhCidsWUjFqcw4cs0FV922u5ayElI3Ra3/Pku9Qs2FnWhlV6Qugs35vSoD6I2Xn7vg\nRu6yawpXDF/KDmDofjb0NqhiF9KkSH4FnxKvMygENm9YnayU+drKFZ9f1VrENn5988b9my5KKlNo\nfeh61q0PPh7c4W0HsQrerwI4X0p5GACEEH8E4OsAek7BA+pzUSl6iGxjQODIsWnrIh4IWxf0+1cl\n3Op043HVweFj03OTlpkNyhZUX6VSVYaQhSk1Ds68V2NAoDEo5txtYt4pdaFTtfCqoj/6ki1csW0X\nrrtrL65529lJ7hML7WJSdBy5hO+MlF6FtJMXuKS3CfW9MvHvwPHkTy73vA+N7sEtD+zv2Pg81xw7\n3Gzg6PRsUGa5DH2pZSgbn13meImi7VC1N0fM78w2KJPkDahWuY1hVsrCSpmrrVLXVyn9xrWOCu1w\np5QpdK1P7vrK0C5iFTwBQK/FmfyzrqBTFjUxi0dzcl2W+8r7Ju6QdUHdv52xFSnYBIxPwQVag+pd\nW+ZFlKqyhCxMKc+23WtqNnPHOGHxUKmJEljYo0DK9kf126tu222dYFVmzhT3iU5MfGKbv3xuaS6F\nNDZRQLfEF5CFpWgfGB0bb8maa+t7Zbwr9LnEtYDWlTvF5NQMruiQQ6Rdc+y1l2RHEOsueGqMA5gb\nn4demC71fH0+LxOf3Q7FuWrPHNf9BoXArJRzWbC3bs8OIVf9LLS2sNXVxi/sBgSs8ZBl8iuEUOfX\nmtm8gYXx7krtN6511KAlPhI43vYpZQpd6+sXvjK0CyEdmue8i4T4IID3Avjr/KP1AD4tpbypxrLN\nsXbtWrljx45Cv7Ulhyjjm12GMzfdbU2dKwA8uuVi6298B0gD2UDcvGE1ADh9jX33X2iKLhLOCATS\nmwgB2Lp2GXeJohRpdxPdx95Gyr1cfSpUN526yHfVr8I1+dret2jd1IVt/mo2BnHpeSMtQlnH1h9C\n75byLDXv1NH+QoidUsq1ld+4RykjH2Nx9Y1QHwglZ9LH1ejYODZ+cfc8T4TGoMDWd5wDoNxuf0iO\nAtmYec/5K3H9+tXR960a3xzra4OQbLDNj40BgROXDJU+9NpVtylGx7IU7Z/qt6HMpeb9bM+Lqc+Y\nfqhQ69NQv0+5pw01x9/78IEFle2pstYn55uNwdrlk81QpZ7VqTIy9hy8jwoh7gPwJmTzxX+RUo6V\nK+LC4PLNBqqzNMUufItYmXzWmRFjMnLRbiuCwhegeu/DB+YFCJuuli4h5cJltyhi7WrHuXnm82Nj\nTGKI3aWyvXe7lBxbALmKwbRlONNJcZ+I2amsS9FNiS+49+ED2LxhtXP30tYfQu3ebfEFpDpCfbqo\n63LIfaylT5pdWaJwAoOYM8nMR9/ywH6sPf3k2sZz6L6+nTNfG7jGtgBw4+XntsyfVRq3Xc+emJya\nWwwvxGHYQPiIpJAiF+vp4/KkUV5GrvdNWX+YLojmekmVy7YDZ8O1hlLyZKFle6q3jC9MJEYRjiVW\n4QeA5Usbc2Ega08/ueMM4NHn4EkpHwLwUI1lqYXQgNKFVJEJOWWbuYibW0zsk+08NYXIy3TBlnva\n3uFcAkp3ndFTO6t63Lp9X5Jy5yNV2W3nuXmK0CKpquQket2Uee8qFaCYAHKfcgekuU+EFgp1uSMV\niS9Qz4vtW6F277b4AlINMX26qOty6Ht9DG7dvg9Ts/P72tSsLGRgiD3WxUTm5QgpBlWuC2LxtYFv\nbFcVmuGqh5gspkA1RiFfW/jes8gRFL7wiZi5z/a+sXUFYM4V17ynGUd6+87xuR04lxFU1dWV23ZZ\n11LtmMtjkqTp7b2s2XDmGqiyj8f2EwBYumgoqv+1i6IHnXcNsSmUi07IKdnsigQvu5SDC89aMZdR\n0re8rTOgPFXwuSYRV/lD1slUiviVV5F0o0i768Ts4qYEwMf43Bd9b1eMwXV37S3kClQ0gFzn/Fcu\nx0P7n4tWsIta0W0xb7FtvhDxBT5Dw+jYuHMntFPjC0g1XHfX3mCfLuqF4JO/AsCFZ62Y+3eVBoYy\n84Z+36rXBWUUHF8b1J08zFcPMVlMFWVkecgDyLdzZ/N0mJyKP4LCJFZRM+9jq6sBALPG70yXQx1b\nHOmXdj+FExYPQQB46bIl3iRinXL2qi9JmmrXbd94fM7oMzE5hQFku2Zl3YpduMZt0X7SbnpewYtN\noVx0Qk7NZpeq5dsWcLFb8iZVulUVEXwp1iuFzzqpdjHXfPgr1gQsqTEAtsmrqqQbZaw7sRksY11Q\ndItf1RntUl1XimbWS+G7P5isLLFOintryvjwzSO2+AJ94Rbbt1zKIJD1EdtC2hdfUGfmWbIwjI6N\nB7MzA8W9EHzyVwK4fef4nEtklQkMyswb+n1jlN+U55cpl68NyhgRYxbtoURh5rNdSdHKGIViPIBs\n8sU1t/kIlTNWqTXvo7eTyvhqKncCmLcjF0OsO2wZQ0AVu9JmX7v0vBHnDr0tKdIssvCbunJKpI7P\nTjdy9ryCt37NCHY8dtDaiYDyBzSmZLMrav0wF3A+l8wQZlnNwFHdp9hHEYXYNrmEYutirJPXvO1s\nZ/axlInHNnktxBlzIWLjwlJcUEI+965+rTJvufpuiusK0Kp8mrt9VaT+Vu6MVbiWxu5kxI4P9UzX\nGKg6vsBWD675RJ2T1KnxBaQ8ahza0Pt0UQXCXNCa6GPCNc8VMTC4xqntyAEd/b6xym/K808bbgYT\nqfhcEAF3GxQxIsYu2kMKq/lsV9KTMkahWA8gvU+FdnJjj6AwMdtCZTuPOarI5Vqvv8+9Dx8otF5S\nuNZiRcaxL8FbyqaBra/dvnPcqXy73jN0vp9Z7pT5KmXe6AYjZ88reKNj4y2dyJZFs+jWdciSoyal\nKn3yfcJFHTr91HOTmLWMkEFx/HSL0bFxbPzC7nlxD88emcLGL4YPSy5ioUzdjYy1TpZ1gVS/tS3I\nFw8NBHdQ6ibm/apwLQj5uyt8fTd2l3Z8YjJqt68xIDAgYO3LOs3GIJY0Bgor47HjM9YCGjM+QnFC\nVccX6M/V+5KrvdQ5SUBnxheQ8vjma7NPF+0D6neuLHi6kgDY57kYA0NMvI5+5IC6TsUpqb9fmR+b\ncOSY+9iBIusCFVrhmmcAe8IPvW6qHoexhqgYQ6dtd6bKzIwpHkCqT/n6t60/lDG8l3XJN8tf1nvL\n9e4pfSgmljV21ys1DKEMRdfbrnFbpp+0k55X8GydSil3+g5G0a1r1cChbHahiTRmcoix9qt3ch0r\noJfxurv2tgS1A9l5LCGrTJGYDNc7KuHty/IVmpTKCL7RsXGn4HhucqolC5m+A7VQA9z1fqHjE1yY\n7WROiBOTU2gMCCx3CHaX5S7WdWVQiCjBMDUrIRwnborcUuPLdBV74HusdTKkbIfGqJlUwlVPodjK\nMmeSmYLPZRXudPcTUh7XPN5sDGDr9vlnfBWZ6/R+6orxNHcKbc8Jze+++csWr2OTrea48JG6w6Oe\n71sHqL/bvqtLzsQaomxn7A0OiHm7nWb9ffaB/Vi+tIEbLz8XAAr1p1ijo4nqUz63Xz2FfRX1q/dR\nVW7X+8bsAJvvfuTYNG55YD+WNRtY0hiY69d1uMMqYmJZY5+TEobQGBCYnpVWubR8aSP4rKIhV7Gb\nCd1Czyt4sTtNZXaBXFvuehC5rxyxBxDHWPsVIx5rlzrc0uf6FpqAUhXi0DvWOXBizhZyobKQAWF3\nwrIWHXNCFyLbUXUpvTHWtVjXAtdu2tJF2bk+sZm3VNls58XoqPeJUUydxj2HL37sGC5qnfQp2ylj\n1DXGBBA8j7CoN4DL4GUqed3gfkLKY5vH1eJKjc2i3iauDLg6VfUz3/w19vs/Xej3LpqNgUI7PAAK\nhYLUmcghJpPh1u2t2U0BYEb7zFV/zx6Zsh7kfcW2Xbjurr3ecJBUo6PO+MQkzth0N5Y2Bqw7uXWd\nT2Yrt238+HYjbTu9uiydmJxCszGIGy8/d25zICaEo8jaNtT3UsavL5eAHoag3F2tZzcOClzztrML\nlztmLPWSt8pAuwtQNy7rgu3z9WtGcP+mi/Dolou9KXJtrF8zgkvPG4G+2aCCyEfHxr3lCFn1gLC1\n35ywNq5bhWZj0Hr9+MQkbnlgv/d9QlaZ9WtGsHnDaowMNyEcZdCJecc6UJPf+ESWbVRNtqNj485y\nKfTJy+dOaLtv2XJOTE7NCTHz6IiYsqvyX3vJ2VHt5JsQXX1BxeOZ77x+zQhOWOy3HalJ3dVHdQYd\nW3hlx3CV1snQ/fR6Hx0bxwVb7ona5Yt9TuxY8sWxxI5l0jvY5vETlwy17JIUmatd42FQiMr7Wdmk\nJilK1PSsLDzP+9YBKWuVqnDNwSpJ3OjYeDBbKeCvv6lZad11e/bIlFdmumTuDyeno3ZxAODI1Cxm\nZiWWL22U7nNq3j5z091Wuecrtzl+XPW+fGkDmzesxr0PH/DKJtOLaPOG1RhuHq+TJY3jS/vQGsiH\nr+/56tJWV7Z31sMQlNw+YXHr/ANk88bWd5wT1XYLPZZi+8ZC0/M7eHWmDzatIkcsVgc1EH3liLHq\npVr79R1JV7IMF41BEVU/KZaOOrKKxRDaqvc9X5+8UpKHFBEesRbk2LLbXFx9pKbgVhQ50NWMrfTt\n9tWZxbFK66TvfvoYTd3lS3lOTB+NzchK+of1a0aw/jNbgV2ZHHrgkR+4L/7yKdH3vcFzn1e/+ETs\nPziJYx+bwUNDg1h5chOnnrg4+t4mt++fwLHp1jG1aGgQeGBr4d+7WHTrILByOKmMAPA3h47ikQOH\nMavtZg4IgVeuOAEArN+tOGkxHrp5CsemZ7CogrrSWQ/gTYeO4jvfPwzbqmDRLYO4HfDXzZdPSa4/\n8xm2uvT1H5Eb/WSB2K2hvxjAM6csTarDZw4dxWkHDuMGvW0+LvDMihNa7uMrtxo/qt73H5xsbdd7\ngJf67mG535sOHcXNZt/5ZFa+lQcn8Snb2HDUu46vv5564mLgy62/cdXVipMW468OHcP0bJY3dGhw\nAGecshSn/nN8/S3aNogHIsaBd5x9ef5vnnG1QwTPHDqK7z5zBC+dncUN+hcfAx4A7Pc791zgppui\n7l+WnlfwqkjAYSPFZ98VMJuS4KVIzFsouN2GEIi2kqRQpPxVEFoMu8pl7hnFBngXVVhTfhcqe5HF\nemwK7thYtVDsA4C5cxwHHDt0sVkci7qf+No05nzB2PstazbmvasroDz2mWXGUp0GL9IbLBoadCtL\nkTxz6Kjzu6GBgXmLr2PTM3jkwGEAKKy4rDy5aV3QrTw5Tr64fj/rGKt6/aQsENXnvuv175YvbeDA\n80crrStbmb7z/UPO93z1i090fq/6hK3+YnEphq5+CGSK3dDAAAYGB+bqKlbBnJ6ZxXe+fwjf+f4h\na/3b2nP/wcmWd5uVEvsPTra0Q+z4OfXExc42jHkf/X6+8rnuc2x6Bg/tnyjdX01cZXn6hy/M+2x6\nZhbPvzAdXX+qzOr/vnEQW+5nDEUwZXyZv3WVt+rxmkLPK3iAPduRWnAVVfhSfPZ1f3bbc2IWXSkL\nM3PB68qAZYu7qcs1y5dVrGxb+AgthtUOqjlEJTBPadm4bhWucOy06uhtXUU5XdeqMtl2g44cm04u\nQ0xwsc9YYCqorvZWyl0oNgeIy+JYJh7NV8ZQQgezngDg8NHWRASNAYHDx6bndihD2cJUBj/fOCij\npKUYvKo81JZ0OJpFeb8jpmfzhtV43Zq4hGBv33KP03DmkkdldpFPBfA1S7neENlfXb93GbVUWV3x\nT7455NT8T8x3FzjqseyOu837yNcmN43uaTmXTO8Tqv5Csdc2XO9i64c6AvNjsF11FUJvL1d7xpbB\nVW69rmIIvbt5v0sdclllVffVS5n+asNVFhdmIiTAns8idJSQWf6YcrvmqZjx5fqtq4zt8JARRba4\nF5q1a9fKHTt2VHKvIhOyDVeWShehZ6Rk0QxdYwua14OcVXmqTmMcwiy/LeVv1UpmTHu72tKcvF0H\nqpukvoN5FqGPgbxgszLb4Tr/lcvx0P4JTE6Zx6XGn2mYQsqCw9VfYwVxzKToutdgbn2PSbSSeiZQ\nzPgCgBMWDWLRkP3oBpNUY0vdypfLlbSOPmUihNgppVxb2wN6jCrlo8LVv1z9/8QlQ/MWaT5jmGux\nZlss1/EOqffwyQ/X/DPcbGDXNeHkLiFcBrUydZWyRtDnoNj6jMmc6nqG7V6uDOWmfLAd+xSLuper\nPYWwJ/tyyaiq+p7tOA/b/XxyOSazdRUKSNGM3jq6ETjmKJ/Q70L17vNuuylPYlPktyZVz22xMrLv\nFLzQAjV2YL7q6r+1TjoCcE5qC6HF+wTOCYuHOsoaX5d1EqhucjQFSMwRAPrvQ8qFS7Fb2hjA4sbg\nvCyazcaAVZHz0Q6FOUTMxBh7zyrvFUOKlTj2UFrXdVXNF0UWG773tJ0jWiVU8NKoQ8FzEdP/GwPC\nucgeyb0P6o4DTZmnfMps6Pge3/zjOqIhhTpkZN1rBFP2mgeAp84fHxrdg89aksL9wvkrcf361S3P\n1mVqzBmqqkyPbrnY254LnY0zhVB/DylfZRWQ1LWRD1vfjpl3XBnDfW3ku2+Z35pUvfaPlZF94aKp\nU/a4AoXLKiUBp09u3QlFfM94bnKqEotildSVeMWWXrnZGMR7zl+Jex8+0OIGF+vyZnNv8w1wX/8J\nTYjLT1jcMiG86uq/jayB41R9jlIVMa2++LyYXbeYe+lUWQcpfdO3rtDf1VX+KuYL25wWk6Lc92z1\nXkXT55PuJfbsShsCx8+PSz1iJzX+1pVc6wpj7nfJ/B2PHZznXaLO6jLnJd/4VTv3asx98LZdePcb\nW5USH3WENtS5Rkg9j9D2e1OpdmVRvvfhAy2f2cJxYnb1QufmAcAJi4YKKcBVeWeFrt+8YXUwxMKl\nlJTNheALWYo1dCps/dMVSqNjM5SHZP/Gdauc/SPmtzFKbTtj3PtOwfPFZLmEwlW37W45sNJ1hpfP\nQrkQBwf7zrU5c9Pd8/yc2x1fU1fiFVc76jEEtsVpTH2Y1w4GXFBck0QohtM2yYXit1LuZZIiYFzx\ncLEUiX1zlS92kq3KuJISK+lCALjhsuOJjFxC1xyzuiU2dtz6zqfyKWex71n3QcyksyjT/yXm97VY\ndz+X0RVoPZtUfecb7/p1Lllx64OPt8y3tr4eG5sNZDtJaifKp+SZY9wMpTBDG1INLVXJXdtcVOY8\nQte5iS65F3umGRDO1KzWRL72LKIAhzYNbF48ofZ03XPzhtXBXaK6kmzFGARjOW24ae1bRX0NfWVb\nv2YE19211xlGEfot0DqP2T5rl3zsOxdN31Z2yEKgXwvYrZCh71IaOrSYcx2MHbKYxPjb+8px4Vkr\nKond87UFUHyQpPhGp26dF3FF0N0fYv3UbeVyuQWH0O/lShBSRX9NIUVRiXU/8cV8VOnuaItfSY37\n+K7mDhPTp3zHRfjaKTQWfDEksf28XfEFJGMhXTTLuGIVGYM+F0XA7erp+s687smJyaSFo62vx8Zm\nKwaFwL9v/jnrdzGupWXdNqtws/+QI+lKSjISk9QkKUXn9JDscbVn1f03ZJx0Pa+K9i+jgLgU+1TD\nz9LGQMsZiT45t6QRF89uEqoXn4zs1OOD6KLpwLdbE9NJlRVPNbpvoJQdRCHLj+kKoZA4vi1u22FK\n2Y62lUP3h0+1HoaskzaFoyoLpQ2fhSbWQhlTHnW/mMWRAHDhWStaPn/XG19hjUUYHBCYcSgYunXO\n1aeWNAa8ZwUWxSdIUnYBQ2cZ6vdyLV6qcJFQ7zM5NdMSk7PjsYMtCx5fbJ2OOSfZlNSUXQWd0Fhw\n9f8Yy7f+DNIfmH3VFmPlMiCmjMGQISy0Q/cL569sWSTa7uFzF7cZimx9/Zq3nZ2k9PqMdK657to7\n987VfdnQhrJu9q6YOH1eNPHNEUWTcxw5Nt3i4RBDSPbY2rOoDPG1VREvHt/nsfVXxgPHtYawKWUh\nlp+weG5NFVpjTU7NYPHQgNeIYCOm3XxHZXX78UF9p+AB7g6e6u5V1lUNSI8hUIu50OSggpiLnK+m\nE6PQ2GIbbNgmh9t3jrdYDi/Yck8phcPWjq7Ftk3w+NwmUpU7fYKJVQ4lgG3feBx3f/OpeXELyq1H\nLfQHhcC73viKeWfE+RLKuPqUq0xl3BrLHF8QWw7b52UXL7HZA82YnPVrWs/qc2WJtQkNfS4505HV\ntYirUmhO8y28VJn0BZgt22e3C0GShi3OqYibUkqmTpNlzQaef2HaOSZu3zk+Zzx0LXzVM22Ledcu\ngmvs6u9rU3p1XDFlgKxCe+AAACAASURBVHssT0xOzR19U4WLpW3tEhsrdotFuVPMSGlNRmLWm29O\niUWPcawyFriKOHOFr61C8tXVnj6FpMgRTSm41hD3PnxgXgxgKHsqcPx8aLO8V3pcZG+8/NzoZ8Rm\nfHatF99z/squDz2oTcETQrwCwP8P4KUAZgHcLKX8EyHEyQC2ATgDwHcBXCalfLaucqQQY0kHwhNp\n7OLWd53PShPryhCbYlZhe6+URX5oog0praFnlrFQxi62fYsLn4VSZSCzZVwD0t1PpmalVYBdv361\nNX4jZiJKVdjK7MzEtnVsOVIWNEUNL77x6Hsf9b3qazdq6ZWLHNDuFOCOVN0hJQ0Arr7jmy0ZWFPO\nzguVmfQ+rrZ3jbdQkgjXWLvurr1e5U6dLRmKfb734QPOs+pE/syt2/c5jwtKHbutGZdbxxyQeWK4\n8MlrNXfWEUcVu2bZun1fWBmT/gyi5rOqCBKyyZYyc1UVxnvAH/Pm27X0tacr4Yh5dm8d+NZmMUZK\nnVQF9rThZtBbR2fpoqGouqhSoe806tzBmwZwlZTyISHESQB2CiH+HsD7APyjlHKLEGITgE0AfqfG\nciRRhbtXzOLWdb6Lus430cdavVwWyhQXGtehtC58i/hY14I6LJSjY+P40u6n5urBZd0J7bK5LJTX\nXmK3FFWVPrgKl0lXvbrSC5dZMJRV0s340hircFl849bXd3Vhay6MXAsF34LKldnLtp6NrwfR8q9L\nz0tfxFS18CHdRZU78oDfFdHnEjwy7D6U20T3tFHPNGWny4tE/c6UIcpQZ97DrAv12w+N7mnxuPAl\nWPEl+Ri3vE9VC9LYBHOxWVR9SVViPVmUoTR2h298YnLOZbNsIpoi+BRK1+dFzhpdv2bE2UfqztS+\nrNmwjs9lzca8f4c2FpSB5YIt97T03dSs5ql14TNU9Rq1KXhSyqcAPJX//XkhxLcAjAB4O4A355d9\nBsB96AAFr8jgdOFb3PrOPtOvu/Hyc72KgR5n5+LIsWkAsKbPjX2vIjl4XO8f61pQtYXSpmS94DhP\nLmqXLWCh1CkSs+ei7OTtckV46zkv81qri1BGSS+barsoPiXOlTUXaB2DMcq4T5m8f9NF3sxeilgX\nFNuzJOxpxgmxUXRH3iVXfa6ILlTCg5jdAWD+XKMWcDZPitj38O06ue7h8rhwsX6NO6ufLierXpC6\n2kMZoJWCFGvw9cmqGDnmSwzmU/BlXlYzHhpw7/Lp6zFzTk053N1nAPHtcBeRuy55VHc8tMu7+Nj0\nzLxjO2weU2q9GmMcAeKzmrt2Ql3hN6mGqtA5mZ2867cgMXhCiDMArAHwIICX5MofpJRPCSFe7PjN\nBwB8AABWrlxZa/mKDE5f47oWt8uajajdHNUxFw+1Jr/Q0ePsljUbODY9gyOa4qJSobvS59qsk+b7\nPBdIsGBjWbNhvZfPtUAPIq/aQpmyOAkdewDEpX0uGjjuQwJWi1cs69eMtCQDkcjiVdaefnLhbFG2\nsVBGSbe5aMWm2i6DzwBhE1g+QouYkAt2zCLK5YJitked5+yR9rCQ8hEotiPvk6tFjltQc0fMb11z\nTVHPghhDnRq7ZWN/zd0QRZ0ueDF1mpLowqdohJ5ltp1tJzW0jnJJcL2dbWfkPXtkChu/uHvu32b/\n3fiF3bjurr0thsaiBpAYRb1q+VqGCYdcOjI1iyN53apdcZvrs21NZKunFANGSl3E7lS7vPhc52Ta\nFEXbmY4jC6wI1q7gCSFOBHA7gCuklD8UngBjHSnlzQBuBrI00PWVMN06GTo82NXhhEBwYlSHmMYo\ngmYK1wu23DM3yGLew/c+IUHsSvOr4iOURczs+K7tdD2IXF1b1SBIEeqxxxCkLGxCmBlPhz1B+mWz\nlh45Nm21bF512+7oe5r3d53J4zt41Xe/ImfSVIHPAKEHkccsTEOWVJ8yWSYDrK09fAmGfIaqbrBQ\n9it1yEffcThF4tF9ctUlI12p0IebDa+HR2NA4MQlQ8Ed/qKeBTFzjz52fQu+UFIZ3y5mSIkseqRR\nbII5PdGFa24xMzfblBPXkVSDQgSPazANwCmdX1eer7trrzWj+NSMxBXbdjkzkNti430GO1/ik1B7\nuVxNY+Vr1XN4rGFmcmoGX9r9FE5YPF/FKBu64ULfDPF5tsTuVAPHdwdtc9gtD+5v8W7T19quMx0X\nwlVYp1YFTwjRQKbc3SKlvCP/+GkhxMvy3buXAfh+nWXwUTQVs8uap++Y2QafKzuQQk1uMdZCm4Wi\nSutkSBC7XFht7hN6x/e5ui20ddIm1H3lC/1WEWq/mMWIr2+6knuEgtl97zUjZaGJJ+RqmGqdVG7F\nNup2PwnFNvjcvHQaAyKYwjsl06sLW3243DFNbIYkXfgA5Y4qId1F6Dgcm3IX2jEIJWQA7IcE2+TN\ntZecPffvMh4eLkXmyLHpeQtxc14KuSbaxq6KKbQt+EIJnHzPcSmRZY40Uv8OhY+oRBdAOH7MZ/xz\nzXOzUka1o24ATklgdjhvZwBBL4kYY+/k1IxTGVS46jymvXyupiH5WkXcrE0BjfVkmZicajH0u8ZR\nUdluO4/RFX6jnhOzU63Woq45zNU11PW+MV1FPoVY6syiKQB8EsC3pJQf1b66E8B7AWzJ//83dZXB\nxHSD8KUyBtydzqcw+QZfKGuSslz5FEEBOAVarCKT4sLlE6Y265Cr7KrOfEHkKkYxRXDHXB8r1H3X\n6oQslL7+MSgEtr7znODgVgLMdQineWyDSvax47GDczEfqfF/RSaeWKOCK5V6rAIKLMyZNDGxDb4+\nonZfQym8beMqxWUt1QVNRy3CfMq5+rvtO+7w9R4piS9mpYxq75A8MseAMijG7EwU9fBwKTLKOKsw\n56XGgN/zyLWKUJ4pRRI4mbiUSOV9EdOGSiG5YtuuFncxtWPhUvBijvzR3cZ971w0hizFXXFAAIeP\nGa7+M3JufquSUEZXm0H2yLHpYHvFuJq6CHmm+eIPXUdF2VwvY5Me+Vx8bWuxEOrIDteYsLlcxu5U\nq53XVLms+m/oNwsVHlHnDt4FAH4RwB4hhFrR/y4yxe42IcSvANgP4J01lmGOFDcIwG+dDDW6q/Fc\nncvcUva5RfripGIUmVQXLsB9Zo7NPz2Uyn39GncQuRmjGLI4mdab0GLaJ9R1Iec7tHpk2O1aExOI\nHmuhVLj6wqAQ1p2aWx7Yj7Wnn4z1a0YKxf+pvhu7eA8t4nxnCtoOWHehu2jpVK1kxPjzu/qIcsM2\njUYuxdkcV7GWaJ8LSoxAUouwkDHG913VWRVJe4ldcMxKiUe3XBx1bWgs+XZ4YuOBi4x/lyLjM25M\nzUqnbBsOnMnnU+J8xh119E7IAKS8L1KMeYB9zPr6ge46GWPY811jSyIX2hFODQcoMr/Vhc0gW4aY\nePxQsj9X/KEZX6YzOXX8+BFFSkiKcvGNWYuF2LrdfWSHyyXSNCr5ztK7+o70A9w3rluF0bHxoDdO\n3d5ICiGLpElcYNauXSt37NhR/AZXXIGH7rwPx6bjGmnR0CBWntzEqScutn7/zKGjeOTAYcw66m7R\n0CBet3IYzxw6iv0HJ3FsembungBaPjOfY7v/gBB45YoTnGXSf/vdZ45genb+NrX6vXp2iNDzHto/\nEV2f5r1c7zcwIDA907q9rupT55lDR/Gd7x/yPtesX1eZzfuXvW5oYKCl/n3vomP2meVLGzjw/NGW\nunL1Pf0ZDzxyEKmnDKk6i+1/vr4KwDtOYinybN84sY1L/frQ90Ba/1e8+sUnzhsDtrkhpr6GBgew\n9vTlzneLucf5rzzF288BeMeA67f7X/EaXHb/7d5nhxBC7JRSri11kz6itHxEvHFhuNnArmtaEx0V\nyTTnembIkKk/06YkhGK4ADg9I9Q+nWv0xCQYsd0zNZbd9h6hNopJEGZDb9PYNom5LnRNqnKe2l9c\n1xetpzKUeaZPWTD7iV6nLuVlJLDL5DJk6OV5dMvFLR5xQmAu7MS1qzcoBG647BynN1vs2AfcY9iG\n674h5VSNz5gsrmochcapAOadlVuEWBm5IFk0O4EU5c63AAcwt0hzKVIrT262LLSOTc/gkQOH8coV\nJ0Tf/7gyli3o9x+cnPe967f7D05mR8trqN/76mHR0KB3UasTv7htXXCb76ee51LYrAvJg+HFiKpz\n9UxXmc3Py143PTuLl7xoCZ7+4QvzPld9w4Wtzxx4fhYrTlqMZ49MzasrX1se/zxdoCxf2sD+g5Mt\nCoLqP2af8PXVWSkLKXdDAwMYGBDBvugrp1kmQFoVZrOPqP+HDCmpyh2Aueeov9vmhhgjjM0Iomht\nj1aUAudS5F3Kpv6d696HjrpjKEnnEuu6dNjiSlUkCzUQtwsUOr4oNTlaaPHrc7HSF3s+zxfbPV0x\nt+MT/sPWdUJtNCNlIQVUT3AW48EwOjaOw5ZxruJ6VRKY0PmlqW62vkQmJq4yAvGJ1BTmO/hwnSmb\n2ib6b1XfcMXjKxddYL5rsS9u1hcGFKoelZzL9IhrNgbnFBeX4hTabR6fmPTGrpvlKJOQDAifpafG\npy+zJjA/Tti3QywAvOf8lQvm5dIfCt5NN+Gql4YtlMoa8roEVzDb52/wJGKItVCcCuBrBS2Ul3qs\nkzHun6Nj4/hlz9kwAHBVpMVXWXts73eq8dnbE+rM9Y42lGXFVWbz/lVdt9vRN1ykvP/XxsadmcjU\n9bFtpOMTRr62tPXV8LMGAIhCFnjA3wdc7+FalPnGpW2MuyyQNgFvPgdwLyDv33TR3Lg4w3Pely8W\nV40t3w7H69aMzLVbSwwpjscLm+mdVf/9tQ9/xWrFXNoYwC87S006FeV2rB/KPTQocHR6vjFBxTDp\nfa5oivgY926f4piSVMyV1U5HVz58icViEi3ZfgekH7auo7676rbdzp0Z39zkQ7VVKIGNa+G+fGkD\nF//Yy+a5s1V9fqkv87CZIMdWxgEBWJJmzt3DnoxqAJs3/FiUS5++yE+RF7obri/rqWvHakZK51oA\nOB43q3bZrty2y/sePlR/Do15X1+dnJrx7miqswxDLpuxRinA7xKpXLZdYybk5mm2ky+s5obLwvkX\nqqQ/FDyUS61cxEJZRTrYohZK17B1WRLNuAjf2TCpwaopvsYp55mkWG+UhTL2/lVdV5WF0vV5FsM2\nfwGmnu+zYPrwtae+8LIJsFQr5fSsxOWvf3lUGm9XeWLjExWuseFKZ+0a+y7ffCGAS88bmZcNzXyO\nC7Odh5sNZ6xwjBA0F5WqXlSckeqfvt0YtStgtotrbXBkajY5WJ60n9Gxcdy+c3xu0TUjJWam7Y1s\n9tOisi40f7rk37V37g3KORPX/ORLGuNzH4yR47YkJi7lMDaBkfo+pICee91XgnkGdPT38cktX3KV\nex8+YI1drOr80gvPWmGdV83zAV1l9Cl3rrj5JY3BefXhU3B1I7hu/L9y2y7nbua1l9hjqW341jw+\nVW1GStxkxDwWUe70d4yJb/RdF7PbHDIS2RQtW5bPUFI81861ryy+MRLrcr0Q9I2CF9K6fRSxUBY9\nb0enjIXSxHe0gb4b6bIOmpZb8z62rKSpB2+mtFGK9UbdU+3QxGRoq/K6WJY5FvTmwbfH23q+cqes\nqKFU10XQFUebwlPEBWVqRs4L2NYFop5p05zAlULoEppF3WFsypJr7Ktz8WzB4rfvHMfyQFp1G+bc\ncO0lZ7cYW0xihWBsQpTYue45T/9aqBTQpDpSDDRmP3XJulAiiND86ZJ/evp1k9Tssq6kMSHjnG+x\n7VvMjY6NO38Xm8AoRu74xqfrfdSziyi2MYmZyqAMEDHPSH3eacNNd18z5vBYmW9zYSy7m5m65lEI\nZGf+hTxali/N1hmu+DJdSY9d34ZiT1U9Fs0Wahuna08/2bm+DcnB0E6e7wxKndhz+eqmbxS8Mtn2\nikxqsb7sttTxsXECOj4B7bIk6qjO77PsmLscNuu/bbdgx2MHW3Zp9Pc0LZS+dtHrbEljwHmdiWqr\n2F21qq+LQTgycZufu9paSkRnfPIdKmyiJigg3eVCuaAUWdRs/MJuQMtKaZ4TNDE5hYG8fLrQLOKi\npN5DTyHtu486QsSVjc+VDtqFbWFqLiZihGDs7qpNaYtZfCpCx6uQ7iK2zWI9GRQuY4LZT21JB1JT\nlJtyLuZeRbLZ+bwjfIs5Nce5UGWJGa822asvPkPn9pkcOTaND43usR6qrZ6nyuhS5l0yoIqMgSED\nhP4MX1ZSW3ycT2bYyh6qe9ecW3Y3M+Si60LCf+bfdzUDhy++7EOje+ZcuAWAwQGBmVm/Qd+1Dr7w\nrBVRiUuK9B3Xmsw1rq7Ytmsuzu7+TRd53a/V5755zXxf37l8dRO/Qu5iVKWP54sk1TjqwMsQrk6m\nLJS2+6xfM4LNG1ZjZLgJgUz4mNmOzDJt/MJubPzi7rnPUg6XdQloAUQdOB1rwfXV2/o1I7jwrBUQ\nmJ+m9rMP7Pe+Z2x7mHVm7mD5WKi0tEUZHRt3TsKmFdFn2Q61od4Xr3nb2Wg2BoNle2FqFjseO+g1\nACiXCx0lGO7fdNFc7JmJb1EzNSuDwe2zyBTbGy8/F0AWX3D46DQag/5zq1yoFNKqn7lQ5Xa1xXOT\nU/PGvw9zbtBZv2YE92+6CI9uuThYh655LkYRi118KjauW+V8r04fa6QVV5sNNxtOGabQZZ0N3S0Y\niJfHG9etipqfgLCcs90r1cNElX3jF3e3GHWWL23gpsvPxdjv/3Sye6Mq/4VnrQCQblC21eehF+xz\n4PKlDfzC+SsxbHiFPHtkCrc8sN+pWCp8bZKyXknFZ4Awn+Fq62svOdu5JivaP2x1f+W2XbUZv9av\nGcENl50TPS5CmOPetW7d8dhBfPaB/XNtLAHMzEqcsGgwem5Q16nwhlB/rarvKHx1r89BsfOOGhtK\nwT9z09246rbdwTG0kPTFDl7RIHBFqoUyxjrpWtDaiDlctqyFMnbiUVmbbIdIug6eNLG9Z6g9fO6j\nIXT3wpRd3IU6yDl1cZ1q2Va4EomEdrwmp2bmLHcuBoXwZoIL7WiXEXwTk1NOd5hnj0zNSxQSOpT1\ntOFm0Nihl9vVFgNCzI2RGy8/t5K00EVjllzoiqFvbLl2F11nAVYplMnC4OpbKTFCPvTxbXMX88mV\nKiz9VbnTX3fXXqfRyRc/B/jjbyUy74u1p5+cLMtdawkzgYf+vvc+fKBFSY3xEIhxZYtZr4Qw69C1\nIzkoRItiEdPW6js9Hjn0Gxu2uq/j/DOzPi49byQok3V8u5e2/mrKJD1Tp84LU7PBczHNHbULttyT\n3F/LrMVCuSkUag2qh/KE1llmiIqrPdrl1dIXCl7ZhCehSc107YqJd0lp8JjDZVMSlNhIURpch0jG\nDCIfIQtlymRmThRAfByS/szY68vgUyhS3B58LpeuvqAm31AMZ6juZ6T0ZoILCdCiSqsi1h3G954x\nKaRNNzCX8cccI7akLFXHqKbMKWZMpa99zTb1nX9UlxGE1Euob+kuy2bcjpk90YZuTHDNUaHDiUfH\nxnHtnXtbfhc7jqpwp3eVXX3uczUPodYRqbLc50VgO7PQ9xsbZiyl+uPK6hizXjEx5xQ9nn98YhKN\nAWGNt/bJmyJy3fUb17gInXeWavxyGQfMMt++cxzveuMrrOPO5jrpyu5pu7dtnePz3PHFpdlkhWsM\nufprmbVYyiHsQGsoT+isPV9CN50BIaKPf6iSvlDwyu5u6Z3UhfrOZUU3rZMpC9qYcpa1UG5ctyqY\n0MGGntWszAIdSLNQurBNZlu378ORY9NJu7hld31D6H3KV+MutwdVRt9krYgJ8i3q368Tk/DD50Jl\ny3I7C8wTVikod0uznjZvWD0vVlQ/BkAZKmJ328y2sMXN6klZXMHfsePWF/8Rm/p6ZLg5FwMRGrMj\nw82W5+nzhNotLXtwK2k/sYtivYc9e2TKmTFWYe4yx2DLKll0bltIT4wUzxwbKrYXiJflRdY3rt+4\njgtIicdL3aWyzSkmrh0eAPMUDN8xA0CxzOQ+5SKU7n9kuBm9IwXYla0sW7Y70Zcvj4P5XPMdbbtp\ntvrwvacrLs2WaMZHyvpPX3f6xkjK2tFWhlAypdh7u4xXddMXCl7IImYmB1HnAM3IbFIxs0PaUB3D\n1RliLPo2Uqz8ZSyU69eM4Lq79iZn/gP8Wc1sNAbEvOQZgD+VrW8RajvqAmidKF2ozEhVZAuLJdaq\nZC6udXxtnbqY8VlPUylaP7ZFzYVnrcC2bzwOXy01BgVOXDxk7bfLmg2rwNy8YbXXLTLVgq7vgroO\nTFULN7MtqrROxih3wvF+NmwuPEDrAnBqVuLaO/dSwetRXBn4YjAVsJT5Qb/WtVBbumgoqNxV6Ynh\nO7qkaDZhHZUx2WfIMef1C89akewq7ZrjLj1vBF/a/ZT1Hc1Fdegg8xhGPee5mpg7PKZiaCbhsrV1\nqlwPKYS+Odflfm8r98Yv7MaiIbsi5+pTLpkCxPdt13ubCfXe9cZXBA05qryqblIN86k71Pq6swpP\nudQEUksaA06vKeWm7DL4VrVREKIvFDyfRcy1SFL/j1Fc9I7hs3QoTOuL6faiWOj0qmYyjyoxd0qA\n1sX81u37cMW2XS2HwLqsisBxVzy9nmxWKR91WicV+iI5ZqdFZZqKTcurMPu6GWNgK5cvfs2F69DY\n04abhS3mNl99m/XbjPEA7OdCCdHquumaXG0xDjHn85muay5SrZMxAiDVOqnKEfM7FdsCxC1cqz6W\ng7QP0+BTpm3NDHJFPVeKJB6JCanw/dY29n1Hl0xO+Q9wjuHwsemWcyR9SiqQxe7pTxTIzuKM8diw\n7f74jiPQF9VlU/+r94qtLXMOvfbOvcHdUbOtU+V6qN+NeHZCbQqLS6GdmpWYOpY+l5fFNx71NdH1\n6zNZoOL+Qjt6IddGHTPsITb+0sQMlUoNGbJlZTfDs/S16LNHpoKuw2duutv6rIWKyesLBQ9o9eXf\nuj07byvWrcmGAFomtdh76dYXV1rWkHWyanyD3bWgj8V2YLLL/cZ8TOixpoJWZPAo66TvDLki1knb\nIsPXR1SfMg/sHJ+YxBXbduH3/noPjhybcQrT0GLAFOgul6Kli4awdJH7eIMXLbEHbV941oqW51+5\nbRd2PHbQeT6NC1c7umI8zHu7YunMXVtbnfniCRUfGt0zz3LuatUi1smYM3dS+3lMjKHiRc1MNBRR\nIkn3kupWpXAZ4cwFtssV2+fRAaQtzGM8JELx3qEYLddOve0AZ9v7uTDPmwX8RiD1dx2JLIlKCNvu\nT6pxtEzq/7I7PLF9U2/rVO+MUL+z3U8AeM/5KwG0uo+ayngMAsBQyZ1SF74dKnPsXr9+9ZyiB8B7\nnECKcqfvctrGX0qqfz0LdqrcevbIlHUX0LdODyWHqXqjIJW+OCZBx0xrW1S5GxRZ3ImZmtmVKtqk\njHWyDHpKV/OIh43r3KnPyyh3Cl0o6aQuIgcshdTv7Uv5baaH1vG5mvpS2duISbVve8ajWy7G/Zsu\nwr0PH7DWyeFjM97U4j5/9dQU+j4BYh4DoOrHVm4J4LMP7G85GuPKbbvwodH47KG+z9evOX6cgBqT\nvklUr7/QAkqhj51zr/sKPhuRMRawx1H63kUvo6+tU4WEslDG/E4Ju9j+qw7JJd1NEYW+MSjw4686\n2fm9mYXRnDe2vvMcbH3HOd7jGFLS2Me8Q5Eddf0dXHJelV1/l8vf8AqcsCjelm7Kfd/6oOq1Q5Hf\n1fWsXC+29ofYI66A+W1t63+++TnU72z3u/Hyc7H29JNb5K3tGIoYZP6f5UsbUWVOQZXfha+NYo8T\ncGEbv7bxl3KKXKyHigvXGhXwJzMy1x6KlHmrDvpmB09RlUV6Rkqrth8T31LGOqko4gYXY510WSar\nwjZIUgWES9lU9/Gl/N66vfVg6hApqewVqf3M7BMxdWJzNfL5q9t+7zucdv2aEfzeX+/BYYvrSLMx\nYO1/vt0h04KtFL+7v/mU1RW5bGbY0FhU9RezSCq6s2GLozRd4GJ2x21tHRtLp1BKmy3+17YDE+ty\n1hgUuOZtZ0eVgXQesQmfbKgwAl/yFFOOxcQPm67lKYlHQnNnkR1183Pf3GR6C9l2eCTiDwYPrQ+q\n3CEoks24jmc1BgW2vuMcb+KTGFxHvMQqRzH9LnYntIyNvOwh6T58+ReGPYY7s2587yeAeRmX1d+v\n1A4ZL+p9pUjxUPHhKkORdXrKvFUHfafgFelAjYFsIjYXYbZFlytZhC+eJ3UhWzRwPCbex+VT7sNV\nPzaWWXbQyqbI1+8D+AdVkcFfl1XTdV7Q6Nh4tOuw+ZzUurS5FOl9rzE4AFjSnByZmsWRidbsWUXa\n0ucaARSfHE3/eRu+RD7Lmo3kDJU6tjFcVFFUZdXvo8ZzStyPK/uaLzmUzeVMT2yk4mdtZ5iRzqaI\nO1NjQGDrO+cvvn3zaoxBJkamxS7Mff3ZjPeJ/a1NSQXCc5NN5qrsijFyf3RsHIePTreUR7+ujBHM\nJNVoVMezQrkHfMbTAZHN21Ue25KiECpS1wwDIjvewOfGm3LP1A0Al/hQn7vuFyNjY9wwy6wfFCru\nNCY7tA+XwlbU4Fyk/1RFzyt4RQM2zaQgLgFmG3SpDZq6kHUpauowStfvYqyToQnezDKqJ02JEQzC\n4l7p8mN3TXW+QzsVrjZYKAtl6DkCwA2XtVooU8/8M8vmmoRc2Z5U+7n63nORCogyFKixkmqtNA0N\n5rgtmoY/FOd6mmOh1RgQOHxsek4BS1Xuhh0WyjIeBKqtbYmhfOPFxJZ9zVU/of6xkOdFkuop1B8t\nc7hrvhtuNqL6QWyyIVd6edOgajtzMsatLWURFyPnfTI3JPd9yrdyt656h8C8ny2rspprQsqyeodQ\nyv7UsvsUnY9e1hnHtaSsM5RCC/gVpdh1SOycHLNz/9zkVDCu37fui3XD1NcPMTG6NlTcaaqRwsSX\nLVuVvx27cUXo694OaQAAG09JREFUaQXP1jEbtgAuA5tLnmvgVRUsmaIUuiY4l9uoIsY6qX5nSzYS\nIyRDE4YtU6dr4AB266Tr0M6Y+ksd/I1BUchC6XuOCsK2lbds4HlqXZqWOJMUQTU+MYmrbtttbfvG\ngAhmPFP9OiSgirgnh1KJL9ZSVKt4stQjQ1S7qvgLW/mLuqCYZ4mVcf+xzVmxLmcmZbKAkvZTpD/a\nEoH43OKB8I5CETdp/RBx/VDs23eOB7Pg2g5gVrs+sRl0YwjJ3NSxpTC9HsqONV/7FM2KbGuvK7bt\nwnV37Z3boStSdled+o4VKktqHcSsMwaFaDHyKvkWa2SwlStmTo7duXfFtPmS/ChcBoAiRg/9M5es\nU/d1edG5jv/QCRmkzHvrruQLdd5mCj2t4LkyBArh3pIG7B3wwrNWWM8BufCsFVFliW38mOt8i27f\n4irWOqkm3dQOq0/Wvh2TWEIWk6I7OuqeIcWlzDEV5nNsB2rb8C24ho2FiOs+oUVDygSUuruq73SZ\nlt4djx1sUbJ0VN8ICZTUHaPRsXFnKnHzfkCW2j1Gybadwah2C13ld43d0A6cblipIk7BpGpr+kKl\ngSbl8MkS35lvKlud3l90xUiPsbn2zr3zdoGKHkvjkucmyg3ZFTftc5OOzaAbS5k44tAYqsqQEhOb\nX+QZLgXV5ZIfS9nYbCBNYSvipRCzzpiV0im/1W9NJceWldMsl0t26f0pxpAscNwwGrqf7beu8VfU\n6BGztvR5/ly/fnXLsUamwVcZpFy4+sKOxw5a20IvdzsQssSZLQvF2rVr5Y4dO5J/5zuLw3cKvW0H\nr6gLE2C3lth2w2zXxRzkbSIA3Hj5uVblIhQPWASX20zonUODrSoh68LVPwRgTcPvokrLja+fmX7s\nC2UtMp9lc4FyYSt3aHfY1y4+661LoPjqFLAnKQjFtdmsrwpf+W+8/FznAcPqjCGT4WZj3gG/vvTU\nPnxl1nH1LZe1OLU9Qgghdkop1xb6cR9SVD4CYUu+y/DgcpNPOTtR7yMuGakrjSkrFd8cHjN+yvRf\nk6Jzdew4HxlulpIDsTInhPmeobK7YtCLPMv2e988FrMeU5Spn//T3v3HWlLWdxx/f/fuZbkLhgVc\nG7r8NJJFzcpu3eIaTCu0sFQrXSmIBFKtJsSkRqFkmyUl/DAmi9lYbJOm1bamTeMPqNIrinGhIGlL\nqgK9C8uqq2ixcjFdf7Aq2Qtedp/+MTN358yZ55lnzpnz4875vJKb3XvuOXOeeebHM8/3+TU7N88N\ndz1eel+vk7+hyXqKfGVX/vti16mbmZ5ihVE60Vqo/AztW+xzbmyFO0tn2eRhoeNat5IfqqzH5HmT\nYsvIVrfghW8y5ad3vkteTB/l+YMLbP/c48HoZGw3Jl+EMusqlm1752Ub2HnZBu+N44SZ6a4xOtnn\nm4xOgj+ikaUxdnxB2Qx+g+7m1cusSEVNj0GKHXg/zHFPZRG1bE27qkK8GOWLaR0OTXrSS4tRL58p\nm1wkzxd9hfB5FWop23zGSaWLKGeLH2efKwuIxAilOVM3QukrVHudeEGGq6qlwdH9IDkzPYVZd9es\nmG5befnrr+y6OPPkmWCLf0joHh7TutxrC3TVZBR1xHTzM44+YJd1gYzRRCt82X2jSv65pG75VZWn\noTKybrfyXvMnNJ6+7j2yTrf8qonTIH7ohe/cm56yNND7TOnfQz3bitd6Ns6z+Jybf2/o83W6pha3\nE1MehtZGzviCwaPuydLqFrxeZgjLouW9LpaYiYmW5Fva6kQos237Ihm+CTXK0lZHkxH8mOhkMX+a\nbqmqG8kr+3wT0bmy7Yb2uamIaxM2fei+Rs614piYX7yw2DUr6/SUcfyqld7JYnzf40tjVQRy+9b1\nXH/XntLu3HUjlLHnlS+tZS0mdSt5vjTn8943W2goQlnVg6EuteDV008LXl4oql9sJfJNpJSNcI85\nL0MtyjfN7i0dEhGjahKQqLLH4I6aE3bcNLu3dJxv2XVfd8hGWVAndP3XmQilqjyJSWuvvQrKvq8J\noX0KPWs9XdLqG9q3UP6GPnfNltM7Fg2vEtvilk9T3V5ldcxMrwCsVi84H18+xfY4yVT1yKrbkl52\nPftUteA13eNKLXjUG2+VyWYN7Ge2O+isuYeiJb2sOxcaTBqa8bMsbbF8EbGY/t69pqHYEpktjv3o\nD35W6+bo0+u4o5iozvzBBc7acW9PF3NVZKnpcU/9DKJ//oXuKbwzsVHK2KUDFg87nOvuXh0agB46\nTod+9RKvOeVlpddmFoFcacZi4cZdNfFOr+cVlE9CBOV5kj3AQXXE3JdHxULMF4kMRSh7baGQ8VKn\n+3PVpGOxLTi+2f18435iZGeqryUgpmXMOdj+ufCs1HlZmmN6ovgmi7nti/u6uqjlr618Za+qC3lV\nHuSFeo3E9hZpoqVi/uBCR0+Ffh6GQ2VkaBz07Nx813eFzpdQ/oby5POPzbP5jJOi9yu2xS1mYqx8\nWntd93hhMbz8eHYsY/YvtETP9Xfu4bo790QFKkI9Z+r2evJdz2VC3UPrXEOD0OoWvLzYCEhWmNWJ\nmJTJRx/6jZaUbTvUdz12jEGdG6evZcFX0KyZmea4VStrt0JlYlsiBzWGsE7/7yr5KG4TkZyYFrxi\n5cY3YUyob39V/oaOYZ0JaupEf2NbdZtogYfywqc4Ji7WICLhWX4Uu3auAE5Y7V8PKrbrSUYteOOr\nqRa8Oi3PvuvrxNXTvPV1p3Q97ITWSi2WFYdy3bXKVLXCFFl6Q/ONJyUw6VpVBL5qbA50jweMucbr\n5nuVqtZ730RgsWWNrydLkUFwTdGyKfF7GY9fNWdCaCmfsnKv6jjXmbch9Bmf2bl5b5p7GcsYc972\nK+a4hfar7vZCY/vqPkfWGf+azWvRzzVUl1rwCmIjIIfS8S6+98cuKlwWnWzqgqrqux4ToawTRZid\nm/deIL6FkPNriJVFKMsm6sgqFlOWNP1XFWKx+xAqnIsPuaFt1inI8hYWD3Nd5GxyMarG6c3OzXc9\n8D93aLE0Ih3q21+VvlCE8oWKCF/sdoqycWxV+dVEC7x3qY90baDYY+Y7z8qik75j6wt2LI01Kqz+\nMjVl3gp23YfEUITygnPWah28lqjT8py9Vjyvnzu06F2mwNez5ODCYkdZEWIQHBpQJrtdl80OCUnw\n1yebLTRmXKpPcTxgzP3ON3ao1/taWV4V9ysry/PHvKq3SGicWVFoaEnGNytq3fH4Vcu+hFquyu5h\n2fniC/yX5VPVc1idci+U5iPO1Z4Urr/umdXBd4g7brt2748O1FRtr2psn0/Z8Y4NwAB9X0ODtGLg\n3zAmtm9dz8z0VMdr0yuM1dOdWZBN4XvBOWu73g+waqUxPdX5NDW9wihbXm9h8TC33rOP829/cKlw\ny9bYasrC4mGuv2sPZ+24l/Nvf3DpwXPnZRuWWiFCn921ez+zc/Ocf/uDHdsAll4P3QzXrZlZ+i5L\nfz/+2JVdi1Jmk8U4Otcqyn/u6i2nMzM9VasCleWxz02ze7n+zj3Mpw/s2cV80+xebrx7b2kLRn6y\ngCwPztxxL9fduad25S7v4MJiV77kvytW/vhmeZePbO3avb+0oFw87Ljuzj0dxzh2Ku4yoYkM6uxX\n7KQ2dQamDzIyCclNPcvDkKww9XWxhKPnZPHazR/bW9722q77UVbBuuGux7uvt8OOG+56vOuahriH\nxCkzjKR15djpFXzqa//LqpUrOHH1dEe6vvrtHwcn25DlZdumdTy840L+5/a38vCOCyu7eR23qjtG\nnF+mIL+dJtaMzdYPLSvPq1e4LT83Q+kKrQP2ma//MGqq+eI9KzYfyu7NvT4UTll37lQtRwP+tIaW\ntCmTv3dn97g66u53VRkZ+2xUVJUfZWkoy/vQtnx8aQ5tp+zZLnTM1sxMdz3f5k2ZecujMk0M1anz\n/vz967hV3c+hPsXj7TtmkDzDZ+dSE9fQIE1MBa/sgt91xbmceNyqrvdmBdTOyzawZqazQnZo8Qg4\nOh50dl1xrreLx8GFxY7KRVU0YWZ6qnYl0Dk6Ki/Zg+LDOy6sLPSyz/gqQFUPylmkIl+Y+8YR5ZU9\nBJQ9LMbIWlSKQuMiqgrnfOQ2trLQy7HLvquu0INY1fby50nMTca3vbKHrJjPFV1wztrK87RYQIfM\nzs1HPeyFVBULsZWY2Aeg/PbKjm3Z/StrVQuNmyveF6D6uBjw0Xecyx1XbuTFl44sBWYOLizywuIR\n7rhy41K6tA7eZKtz/KvuFyFG58QUZdfD1VtO7ynN27euT7oGFmRjbX37WBXsM45WSIvfF5MPZffm\nXh8Ky9Iac+zK0pqvrIWuc1/lCpLjV1XJyutlv0NlZMwxiD2HQ4HHbZvW8dF3nFvrMz51vzv//JIv\nB3zPMwbsueVidl1+rjcN2UzMxeuv10ps3eNa5/39VB5D1/bczRfXap2re9yaNDFdNKF88gpft5Fs\n8oBdu/d3Rd8XjzhWH7OSuZuPjsNpovtlNm4JyteQ66VZvKpratYdsrgN35pceWtmpksfuGO7wxYv\njn4eCn3dWUJTCYf4Irc+oWNXpelITkz+Z+dJTHdeX/qy/PZ1W43Zr9m57kXIi2Jmf8x3J6vT7aMf\nTU+3HhOdLI4TiD3P8veFqvMjeygNLdhedX8ZRnRSRq/O8S/rAuobb1c1hjvbXvG1e5/4UWUZWUxb\nWXfT/BhiX9keGq4RGr8cM2SjrOUP4oZf+NJTFHPsqrrt1l2XtDhb8vSUdbSy+MbgNf0wHHMMYs/h\nqnFvvXymie34Wpd8521+ofGqSZSy92Xf7Ru/W3Xc6pzPdc+D2OfQ/Psz6wLndcx31LmGBmmiKnhl\nqg5QbISy1xsvHI30FWeGLFs8PGbWo2L0wJeu0DpfVRWgmekpbr30taV/i82L4g207gWZV7c7S6hw\nzm4kVbORZtspm8o3tsI/iMJr+9b1peupFWVBDPCvr1aVvuzzvdzcs+8NnSehGTJ9Y79Cx316hVXm\nS6yYSkydc7pupajXCOX2reuD95HsPhQbnez12MvyV/f4FytlvgfDWy+NX8ct75a3vTZY9vjSFhrX\n69vHugsrl32fb5Krspa/7HPgv1/XqSDFHrte8ibmnn1wYZHpFcaJhcmgsv0b9MNw6BjUOYfrfFe/\n6mwn1PJcNRN1L9c11D9uVedzneU+iupet73sfxPX0CBNfAWv6gDFRiiHEZ0EuO2L+2pFKIsXUHGW\nn6ajk/nvDE1IEnuxFJlnxjNfdxbfdMhXveG00gHyMZHb/Hbylbtii9I1W04vnU3u+GNXemc3bIJv\nAoSifMTOtw91b9J1C+WqLj6+7YT6voeiytu3ru9popyi2EpM1YxtdbeX12uEctumdd77SD5COe7R\nSRm9fo9/0+dP2UQLZvR1vw2lcfMZJ/WV9l5bhEL369jtNZH3dbZRds8u6w2V3+4wtPUeVlUOhva3\n3/Oyjn6fP0Lb9e1D1XUbu//jfu5MzDIJIaGTqp9Fi/tdSDt2m3l1t+9LYz/Ryaq0hqbQDy3sWjdd\noejoh7dtqLyZhPK62OoaysfibHLDvvgHcR42pdcphEOLmt5x5cbg/vqmUy5GvvOV8X4eFKsWTK2z\npEReP9duzDkxqvNGyyTUM+jyUaQJVQtRS7PGudyX/oz1Mglmdgnwl8AU8PfOudtHkY5MKPLQTw19\nELX7piOUyy06WSddVd9fFXGqav3Mf9bXopRNJDNK4xxl6rWLX6h1Kea4l/099Jl+fHjbhqVztuo8\nqmPQEcpxPm9Ehq3J1oVJpPG6w6X7twy9Bc/MpoDvABcBzwCPAFc5577p+4wilDLuFJ3sXS8PTopO\ntpda8OpR+Th4ut/0T3ko0oxxbsE7D3jKOfd9ADP7LPAHgLeCJ8OhCGXvFJ3sXS999xWdFJFhCY35\n1T0nju7ZIsM1igreOuCHud+fAd5QfJOZXQtcC3D66XHr20jvQrMS6gZcTbMJDt+oZqYSGTWVj8Ol\n9R6boXu2yPCMYqHzshURu3q3Oec+4Zzb7JzbvHbt2iEka7KFIpRSbdum7oV31fVERAZB5eNw+Xpi\nqIeGiIyrUbTgPQOclvv9VODZEaRDchSh7J+ikyIi7aMeGiKy3IyiBe8R4GwzO8vMjgHeCdwzgnRI\njiKUIiIi3dRDQ0SWm6G34DnnXjKz9wO7SZZJ+KRzbt+w0yGdFKEUEREppx4aIrKcjGQdPOfcl4Ev\nj+K7pZxmuBIRERERWf5GUsGT8aQIpYiIiIjI8jaKMXgiIiIiIiIyAKrgiYiIiIiItIS6aIqIiMjY\nm52b1zhxEZEIquCJiIjIWJudm++Y6Xn+4AI33r0XQJU8EZECVfDGiKKTIiIi3Xbt3t+xjA/AwuJh\ndu3er3JSRKRAFbwxoeikiIhIuWcPLtR6XURkkmmSlTERik6KiIhMsl9fM1PrdRGRSaYK3phQdFJE\nRKTc9q3rmZme6nhtZnqK7VvXjyhFIiLjSxW8MaHopIiISLltm9ax87INrFszgwHr1syw87INGsIg\nIlJCY/DGxPat6zvG4IGikyIiIpltm9apQiciEkEVvDGRFVqaRVNERERERHqlCt4YUXRSRERERET6\noTF4IiIiIiIiLaEKnoiIiIiISEuogiciIiIiItISquCJiIiIiIi0hCp4IiIiIiIiLaEKnoiIiIiI\nSEuYc27UaahkZj8GftDnZl4O/KSB5Cx3yoeE8iGhfEgoHxLjkA9nOOfWjjgNy0ZD5SOMx7EfB8qH\nhPIhoXxIKB8S45APUWXksqjgNcHMHnXObR51OkZN+ZBQPiSUDwnlQ0L5MLl07BPKh4TyIaF8SCgf\nEsspH9RFU0REREREpCVUwRMREREREWmJSargfWLUCRgTyoeE8iGhfEgoHxLKh8mlY59QPiSUDwnl\nQ0L5kFg2+TAxY/BERERERETabpJa8ERERERERFpNFTwREREREZGWmIgKnpldYmb7zewpM9sx6vQM\nkpl90swOmNmTuddOMrP7zey76b8npq+bmf1Vmi9PmNlvjC7lzTGz08zsq2b2LTPbZ2YfTF+ftHw4\n1sy+YWaPp/lwW/r6WWb29TQf7jSzY9LXV6W/P5X+/cxRpr9pZjZlZnNm9qX090nNh6fNbK+Z7TGz\nR9PXJurakKNUPk7eua8yMqEyspPKyHaVj62v4JnZFPDXwO8BrwGuMrPXjDZVA/WPwCWF13YADzjn\nzgYeSH+HJE/OTn+uBf5mSGkctJeAG5xzrwa2AH+SHvNJy4cXgQudc+cCG4FLzGwL8BHgjjQfngPe\nm77/vcBzzrlXAXek72uTDwLfyv0+qfkAcIFzbmNuPZ9JuzYElY+pSTz3VUYmVEZ2UhmZaEf56Jxr\n9Q/wRmB37vcbgRtHna4B7/OZwJO53/cDp6T/PwXYn/7/48BVZe9r0w/wBeCiSc4HYDXw38AbgJ8A\nK9PXl64PYDfwxvT/K9P32ajT3tD+n0pyY74Q+BJgk5gP6T49Dby88NrEXhuT/KPyUed+um8qI1VG\nqox07SofW9+CB6wDfpj7/Zn0tUnya865HwGk/74ifb31eZN2HdgEfJ0JzIe0y8Ue4ABwP/A94KBz\n7qX0Lfl9XcqH9O8/B04ebooH5mPAnwFH0t9PZjLzAcAB95nZY2Z2bfraxF0bAuj4woSf+yojVUam\nVEYmWlM+rhx1AobASl7T2hCJVueNmR0PfB64zjn3C7Oy3U3eWvJaK/LBOXcY2Ghma4B/BV5d9rb0\n31bmg5n9PnDAOfeYmb05e7nkra3Oh5zznXPPmtkrgPvN7NuB97Y9Lyadjq9f6/NGZaTKSFAZWdCa\n8nESWvCeAU7L/X4q8OyI0jIq/2dmpwCk/x5IX29t3pjZNEnB9Snn3N3pyxOXDxnn3EHgIZLxFmvM\nLAvu5Pd1KR/Sv58A/Gy4KR2I84FLzexp4LMkXVA+xuTlAwDOuWfTfw+QPNCcxwRfGxNOx3dCz32V\nkZ1URqqMhHaVj5NQwXsEODudDegY4J3APSNO07DdA7wr/f+7SPrbZ6//UToT0Bbg51kz9HJmSRjy\nH4BvOef+IvenScuHtWlUEjObAX6XZAD1V4HL07cV8yHLn8uBB13asXw5c87d6Jw71Tl3Jsn1/6Bz\n7momLB8AzOw4M3tZ9n/gYuBJJuzakCUqHyfw3FcZmVAZmVAZmWhd+TjqQYDD+AHeAnyHpG/1n486\nPQPe188APwIWSaIL7yXpG/0A8N3035PS9xrJDGrfA/YCm0ed/oby4E0kzeRPAHvSn7dMYD68DphL\n8+FJ4Ob09VcC3wCeAv4FWJW+fmz6+1Pp31856n0YQJ68GfjSpOZDus+Ppz/7svvhpF0b+uk4J1Q+\nTti5rzJyKR9URnbnycSWkW0rHy1NpIiIiIiIiCxzk9BFU0REREREZCKogiciIiIiItISquCJiIiI\niIi0hCp4IiIiIiIiLaEKnoiIiIiISEuogifSImb2/KjTICIiMo5URsqkUAVPRERERESkJVTBE2mI\nmV1jZt8wsz1m9nEzO8PMvmtmLzezFWb2H2Z2cfreWTN7zMz2mdm1uW08b2YfSf/2b2Z2npk9ZGbf\nN7NL0/e828y+YGZfMbP9ZnaLJz3bzewRM3vCzG4bTi6IiIh0UxkpMjyq4Ik0wMxeDVwJnO+c2wgc\nBn4b+Ajwt8ANwDedc/elH3mPc+71wGbgA2Z2cvr6ccBD6d9+CXwYuAh4O/Ch3FeeB1wNbASuMLPN\nhfRcDJydvm8j8Hoz+61m91pERKSaykiR4Vo56gSItMTvAK8HHjEzgBnggHPuVjO7AngfSSGS+YCZ\nvT39/2kkBc1PgV8BX0lf3wu86JxbNLO9wJm5z9/vnPspgJndDbwJeDT394vTn7n09+PT7/j3/ndV\nRESkFpWRIkOkCp5IMwz4J+fcjR0vmq0GTk1/PR74pZm9Gfhd4I3OuUNm9hBwbPqeReecS/9/BHgR\nwDl3xMzy16ujU/F3A3Y65z7e+y6JiIg0QmWkyBCpi6ZIMx4ALjezVwCY2UlmdgZJ95NPATcDf5e+\n9wTgubTgOgfY0sP3XZR+xwywDXi48PfdwHvM7Pg0PeuytImIiAyZykiRIVILnkgDnHPfNLObgPvM\nbAWwCPwp8JskYw4Om9kfmtkfA58G3mdmTwD7ga/18JX/Cfwz8Crg0865fNcTnHP3pWMe/ivtDvM8\ncA1woLc9FBER6Y3KSJHhsqMt3SKyHJjZu4HNzrn3jzotIiIi40RlpIi6aIqIiIiIiLSGWvBERERE\nRERaQi14IiIiIiIiLaEKnoiIiIiISEuogiciIiIiItISquCJiIiIiIi0hCp4IiIiIiIiLfH/++8U\nuizFuWEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10e9f61d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.rcParams[\"figure.figsize\"] = (15,4)\n",
    "\n",
    "f, (ax1, ax2) = plt.subplots(1, 2, sharey=True)\n",
    "\n",
    "ax1.scatter(test_results.index, test_results['EY0_i'])\n",
    "ax1.plot(test_results.index, test_results['EY0_i_pred'], color='red')\n",
    "ax1.set_xlabel('example')\n",
    "ax1.set_ylabel('outcome')\n",
    "ax1.set_title('Counterfactual vs. Predicted without Treatment')\n",
    "\n",
    "ax2.scatter(test_results.index, test_results['EY1_i'])\n",
    "ax2.plot(test_results.index, test_results['EY1_i_pred'], color='red')\n",
    "ax2.set_xlabel('example')\n",
    "ax2.set_title('Counterfactual vs. Predicted with Treatment')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The population ATE (from `testdataset3_cf.csv`) is 2.09 while the predicted ATE is 9.01. This huge overstatement of the effect and the poor fit of the simple regression model indicate that the assumption of a randomized, controlled trial is likely false. Of course, this is expected given that we were provided with 22 features to use in addition to the treatment indicator.\n",
    "\n",
    "This experiment was a helpful detour in understanding the causal inference process because we will use this basic model in assessing more complicated regression models. Next let's take a look at a multivariate regression."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The challenge description provides a few helpful details to aid our work:\n",
    "\n",
    "> The causal parameter of interest is the population average additive treatment effect (ATE), $E[Y(1) - Y(0)]$. Causal assumptions of consistency and strong ignorability are guaranteed by the contest organizers, therefore the target statistical estimand equals $E[E(Y \\mid A = 1, V) - E(Y \\mid A = 0, V)]$, where $Y$ is the outcome, $A$ is a binary treatment indicator, and $V$ is a pre-treatment covariate vector containing all true confounders, and possibly additional variables related to either the outcome, treatment assignment, or neither. "
   ]
  },
  {
   "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": 2
}