aula15-regressao-logistica.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Aula 15 &mdash; Regressão Logística\n",
    "\n",
    "Renato Vimieiro\n",
    "\n",
    "rv2 {em} cin.ufpe.br\n",
    "\n",
    "maio 2017"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introdução\n",
    "\n",
    "Nas aulas anteriores nos dedicamos ao estudo do problema de aprendizado supervisionado de predição de valores para uma variável alvo (dependente) a partir dos dados (variáveis independentes). Vimos que o modelo linear ajustado pelo método dos mínimos quadrados, apesar de ser uma simplificação bem drástica da realidade (que obviamente não é linear), fornecia bastante informações acerca dos dados como um todo, mas sobretudo da relação entre as variáveis. Embora o modelo fosse bastante útil, ele é restrito à predição de valores para variáveis quantitativas. Surge, então, a dúvida sobre a aplicação, ou melhor, o ajuste de modelos lineares para predição de variáveis categóricas. Esse problema é conhecido como **classificação**.\n",
    "\n",
    "Para entendermos os conceitos relacionados ao problema, vamos usar um conjunto de dados real sobre câncer de próstata publicado em 1989. A referência do estudo é:\n",
    "\n",
    "*Stamey, T.A., Kabalin, J.N., McNeal, J.E., Johnstone, I.M., Freiha, F., Redwine, E.A. and Yang, N. (1989)\n",
    "Prostate specific antigen in the diagnosis and treatment of adenocarcinoma of the prostate: II. radical prostatectomy treated patients, Journal of Urology 141(5), 1076–1083*\n",
    "\n",
    "Esse conjunto de dados é amplamente usado para o ensino de diversas técnicas. Frequentemente ele está nos pacotes de R e Python como exemplo. A ideia do estudo era avaliar a correlação entre o valor do PSA e outros dados clínicos de pacientes prestes a terem a próstata retirada (prostatectomia). O conjunto possui as seguintes variáveis:\n",
    "\n",
    "- **lcavol** log do volume do câncer \n",
    "- **lweight** log do peso da próstata\n",
    "- **age** idade do paciente\n",
    "- **lbph** log da quantidade de tumor benigno\n",
    "- **svi** variável indicadora se houve invasão da vesícula seminal \n",
    "- **lcp** log da penetração do tumor na cápsula prostática (tecido que contém a próstata)\n",
    "- **gleason** escore de Gleason\n",
    "- **pgg45** percentual de escores de Gleason 4 ou 5 (pgg45)\n",
    "- **lpsa** log do valor do PSA\n",
    "\n",
    "O escore de Gleason é uma métrica usada pelos patologistas para avaliar a gravidade e o prognóstico das amostras estudadas. Em geral, um escore de Gleason entre 5-6 indica que há grandes chances do câncer ser tratado e curado, ou seja está associado a um bom prognóstico, enquanto escores entre 8-10 indicam que o tumor já é ou vai se tornar bastante agressivo e há poucas chances do paciente ser curado (mau prognóstico). No caso do escore 7, o comportamento do tumor é mais heterogêneo, podendo o paciente ter bom ou mau prognóstico. No entanto, ele é mais associado a um mau prognóstico que a um bom.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import seaborn as sns\n",
    "from seaborn import color_palette\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib\n",
    "matplotlib.style.use('ggplot')\n",
    "\n",
    "import statsmodels.formula.api as smf\n",
    "import statsmodels.api as sm\n",
    "\n",
    "from sklearn.linear_model import LogisticRegression, LogisticRegressionCV\n",
    "from sklearn.metrics import accuracy_score, recall_score\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "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>lcavol</th>\n",
       "      <th>lweight</th>\n",
       "      <th>age</th>\n",
       "      <th>lbph</th>\n",
       "      <th>svi</th>\n",
       "      <th>lcp</th>\n",
       "      <th>gleason</th>\n",
       "      <th>pgg45</th>\n",
       "      <th>lpsa</th>\n",
       "      <th>train</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>-0.579818</td>\n",
       "      <td>2.769459</td>\n",
       "      <td>50</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>0</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "      <td>-0.430783</td>\n",
       "      <td>T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>-0.994252</td>\n",
       "      <td>3.319626</td>\n",
       "      <td>58</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>0</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "      <td>-0.162519</td>\n",
       "      <td>T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>-0.510826</td>\n",
       "      <td>2.691243</td>\n",
       "      <td>74</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>0</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>7</td>\n",
       "      <td>20</td>\n",
       "      <td>-0.162519</td>\n",
       "      <td>T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>-1.203973</td>\n",
       "      <td>3.282789</td>\n",
       "      <td>58</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>0</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "      <td>-0.162519</td>\n",
       "      <td>T</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>0.751416</td>\n",
       "      <td>3.432373</td>\n",
       "      <td>62</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>0</td>\n",
       "      <td>-1.386294</td>\n",
       "      <td>6</td>\n",
       "      <td>0</td>\n",
       "      <td>0.371564</td>\n",
       "      <td>T</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "     lcavol   lweight  age      lbph  svi       lcp  gleason  pgg45      lpsa  \\\n",
       "1 -0.579818  2.769459   50 -1.386294    0 -1.386294        6      0 -0.430783   \n",
       "2 -0.994252  3.319626   58 -1.386294    0 -1.386294        6      0 -0.162519   \n",
       "3 -0.510826  2.691243   74 -1.386294    0 -1.386294        7     20 -0.162519   \n",
       "4 -1.203973  3.282789   58 -1.386294    0 -1.386294        6      0 -0.162519   \n",
       "5  0.751416  3.432373   62 -1.386294    0 -1.386294        6      0  0.371564   \n",
       "\n",
       "  train  \n",
       "1     T  \n",
       "2     T  \n",
       "3     T  \n",
       "4     T  \n",
       "5     T  "
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "prostate = pd.read_csv(\"https://web.stanford.edu/~hastie/ElemStatLearn/datasets/prostate.data\",\n",
    "                       header=0,index_col=0,sep='\\t')\n",
    "prostate.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "No estudo original, o objetivo era modelar o nível de PSA dos pacientes. Contudo, aqui vamos modificar o objetivo para predizer (classificar) se um paciente tem um bom ou mau prognóstico baseado em seus dados clínicos. Para isso vamos criar uma variável categórica *prognóstico* que terá os valores *bom* e *mau* baseados nos escores de Gleason dos pacientes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "mau    62\n",
       "bom    35\n",
       "Name: prognostico, dtype: int64"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "prostate['prognostico'] = prostate.gleason.map(lambda x: x >=7 and 'mau' or 'bom')\n",
    "prostate.prognostico.value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAw4AAAFpCAYAAAA1LH+GAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xt4XHd97/u3LMtxAlbSSBSiJDS04dI0u7VrH5I03RRQL2kbg6H0WwilJe0+9ubAKbceuqH7AKUPBVIKhEJ37ZJwKSHsLxcXnKawqYAG2uIgxd7k2iYFb5LI2a5lgsxJ7EiTOX/MCCTF1ozimVlrRu/X88zjmTVLP39Ga5ak7/wuq69arSJJkiRJS1lVdABJkiRJ5WfhIEmSJKkhCwdJkiRJDVk4SJIkSWrIwkGSJElSQxYOkiRJkhqycJAkSZLUkIWDJEmSpIZWFx1gvojoB8aBezPz0qLzSJIkSaopW4/DK4Hbiw4hSZIkaaHS9DhExFnArwJvBV7TxJdU25tIWpH6ig7QJM9/qfU8/6WVq6nzvzSFA/Ae4HXAuma/YHJysn1ppBVmZGSk6AjL4vkvtY7nv7RyLef8L0XhEBGXAgcycyIinrnEfluBrQCZyfDwcIcSSpIE09PTXHvttVx22WWsW9f051yS1BP6qtXie/wi4m3AS4BZYC0wCHw6M39ziS+r+omD1Dr1Txy6ZqiC57+KsHPnTm688UYuuOACtmzZUnSclvH8l1au5Zz/pZgcnZmvz8yzMvMc4IXAFxsUDZIkddT09DQTExNUq1XGx8c5fPhw0ZEkqaNKMVRJKqtdu3Y9Yizt1NQUAENDQwu2j4yMsHnz5o5lk9RZY2NjzPXSV6tVxsbGeqrXQe23nN8p4O+VIvn7/9hKVzhk5peBLxccQzquhx56qOgIkgqwd+9eKpUKAJVKhT179lg46IT5O6V7eKxKWDhIZXKsTxC2b98OwLZt2zodR1KB1q9fz+7du7//eMOGDQWmUTfyd0r38FgdWynmOEiSVHYXXHDBko8lqddZOEiS1ITdu3fT11dbeKSvr29B74MkrQQWDpIkNWHv3r0LJkfv2bOn4ESS1FkWDpIkNWH9+vX09/cD0N/f7xwHSSuOhYMkSU0YHR1dMFRpdHS04ESS1FkWDpIkNWFwcJCNGzfS19fHpk2bWLduXdGRJKmjXI5VkqQmjY6OcuDAAXsbWiQizgY+AjweqAI7MvPKYlNJOh57HCRJatLg4CDbtm2zt6F1ZoHXZuZ5wIXAyyPivIIzSToOCwdJklSIzNyfmTfV7x8GbgfOLDaVpOOxcJCkgk1PT7N9+3YOHz5cdBSpMBFxDrAB8AIZUkk5x0FS20VEPzAO3JuZly567qXAnwL31je9LzM/0NmExRobG2Pfvn2MjY2xZcuWouNIHRcRjwU+BbwqM6eP8fxWYCtAZjI8PNzhhO0zMDAA0FOvqVd5rCwcJHXGK6kNQRg8zvP/PTNf0cE8pTE9Pc3ExATVapXx8XFGR0cdP68VJSIGqBUN12Tmp4+1T2buAHbUH1YPHjzYqXhtNzMzA0AvvaZe1avHamRkpOl9LRwKsGvXLiYnJxdsm5qaAmBoaOgR+4+MjLB58+aOZJNaLSLOAn4VeCvwmoLjlM7Y2NiCqxHb66CVJCL6gKuA2zPzXUXnkbQ0C4eSeOihh4qOILXLe4DXAUt9jP5rEfEM4F+BV2fm3R1JVgJ79+6lUqkAUKlU2LNnj4WDVpKLgZcAN0fE3vq2N2Tm9QVmknQcFg4FOFbvwfbt2wHYtm1bp+NIbRMRlwIHMnMiIp55nN12Addm5tGI2AZ8GHj2MdrqyTHOF110EV/5yleoVCr09/dz0UUX9cxrkxrJzK8CfUXnkNQcCwdJ7XQx8JyI+BVgLTAYER/NzN+c2yEzp+bt/wHgimM11KtjnH/mZ36Gr371qwD09fVx8cUX99z4WZXfcsY4S1q5SlE4RMRa4AbgJGqZPpmZbyo2laQTlZmvB14PUO9x+P35RUN9+xmZub/+8DnUJlGvGIODg2zcuJEbb7yRTZs2OTFaklRapSgcgKPAszPze/XVFb4aEX+XmV8rOpik1ouItwDjmflZ4Pci4jnUriB7CHhpkdmKMDo6yoEDBxgdHS06iuZZzkIWLmIhaSUoReGQmVXge/WHA/VbtbhEklotM78MfLl+/43ztn+/V2IlWOqP0Y997GOP2N8/SMvFhSwkrWSlKBzg+xeImgDOBd6fmV45UtKK4B+j5eRCFpK0UGkKh8ysAOsj4jRgZ0Scn5m3zN+nV1dVAa9G2E08VjoR/jEqSepWpSkc5mTm/RHxJeAS4JZFz/XkqirQu1cj7EW9eqxcVUWSJC1lVdEBACLicfWeBiLiZOAXgDuKTSVJkiRpTll6HM4APlyf57AKyMy8ruBMkiRJkupKUThk5jeADUXnkCRJ0spxrJXujmf//tolh+bmpTXSi6vilaJwkCRJkjptcnKSyW/v4wmnNP6TeKBSAeDhg/c03Pe+B2ZPOFsZWThIkiRpxXrCKau5/Gmnt7TND95xqKXtlUUpJkdLkiRJKjd7HCQc4yhJaq1mf6/4O0XdxMJBwjGOkqTWmpycZN/d97J63dCS+1WqtcEf99x/pGGbs4enWpJNerQsHKQ6xzhKklpp9bohTr/g0pa1d2i3K9WrWM5xkCRJktSQhYMkSZKkhiwcJEmSJDVk4SBJkiSpIQsHSZIkSQ1ZOEiSJElqyMJBkiRJUkNex0FS20VEPzAO3JuZly567iTgI8BGYAr4jczc1/GQkiRpSfY4SOqEVwK3H+e53wW+k5nnAu8G3tGxVJIkqWkWDpLaKiLOAn4V+MBxdnku8OH6/U8CoxHR14lskiSpeRYOktrtPcDrgIeP8/yZwN0AmTkLfBcY6kw0SZLUrFLMcYiIs6mNcX48UAV2ZOaVxaaSdKIi4lLgQGZORMQzT7CtrcBWgMxkeHi4BQnLYWBgAKCnXlOv8lhJWslKUTgAs8BrM/OmiFgHTETEFzLztqKDSTohFwPPiYhfAdYCgxHx0cz8zXn73AucDdwTEauBU6lNkl4gM3cAO+oPqwcPHmxv8g6amZkBoJdeU6/q1WM1MjJS2P8dEVcDcx8ynF9YEEkNlWKoUmbuz8yb6vcPU5tEeWaxqSSdqMx8fWaelZnnAC8EvrioaAD4LPDb9fsvqO9T7WBMScX6EHBJ0SEkNVaKwmG+iDgH2ADsLjiKpDaJiLdExHPqD68ChiLiLuA1wH8pLpmkTsvMG4BDReeQ1FhZhioBEBGPBT4FvCozp4/xvGOc1RYDAwMcbWPbHlfIzC8DX67ff+O87UeAXy8mlSRJalZpCoeIGKBWNFyTmZ8+1j6OcVa7zH3/29V2NxzX441xjojPZOZzj7H905n5/LYHk9ps165dTE5ONrXv/v37Adi+fXtT+4+MjLB58+ZHnU013fjBYe0DwSNtabcbXn+38IPD5SlF4VBfs/0q4PbMfFfReSQt8KzjbH9mJ0NI7TI5Ocnkt/fxhFMa/0ocqFQAePjgPQ33ve+B2RPOpppu/OCwXR9IdcuHUd3CDw6XtzhCKQoHaiuvvAS4OSL21re9ITOvLzCTtKJFxFvqd9fMuz/nR4H/1eFIUts84ZTVXP6001va5gfvcNi+pN5SisIhM78KeKVYqVzOrv+7at59qF1r5W7gzZ0OJKn3RMS11HowhyPiHuBNmXlVsakkHcuShUNENLXqUmYe74qwUleYmpri6AOzLf+E8L4HZjlp6hGXJOgKmXk5QET8U2b+VdF5JPWmzHxR0RkkNadRj8MstU8Xj6ev/nx/yxJJKpXM/KuIOBV4KvDYRc99sZhUkiSp0xoVDk/qSAqpYENDQzxcfbAtY5xXDQ21tM1Oi4iXAu8Hvgc8MO+pKrW5DpIkaQVYsnDIzEdMfqwPX3o88L8doiStCG8FXpCZf1d0EEmSVJymJ0dHxCDwPuCF9a+biYiPA7+Xmd9tUz5JxVsN/I+iQ3SjZq8P4LUBJEndYDmrKr0XeAxwPrVlGH+E2ieR7wV+u/XRJJXEO4D/GhF/bC/j8kxOTrLv7ntZvW7p4WqVam0dinvub3yxqNnD3TnZXt3NeU6SYHmFwyXAj2bm3Bjnf42Iy4F/a30sSUWKiLv5wcIIfcATgNdFxIK/WjPziZ3O1m1Wrxvi9AsubVl7h3Zf17K2pGY4z0nSnOUUDkeAx7Hwok/D0LYrdUsqzm8WHUDqFJdjbsh5TpKA5RUOHwC+EBHv4gdDlV7NDy4BL6lHZOY/FJ1BUmk4z0kSsLzC4a3AJHAZMFK/fwVwdRtySSqJiHjLcZ46CtwDfC4z/3cHI0kt5XLMDTnPSRKwvMJhVWZejYWCtNI8BXgecCNwN3A28HRgF7AZ+IuI+LXM/FxxESW1kvOcJB3LcgqH+yLiE8A1mfmP7QrUa1yOUT1gFfDCzNw5tyEingtclpkXRsRvA28HLByk3uE8J0mPsJzC4ReBFwHXRkQF+Djwscy8uS3JeoTLMaoH/BK1c3++64C/rt//KPDnHU0kqa2c5yTpWJouHDJzD7CHWlflz1H7Q+KLEbE/M3+yXQF7gcsxqsv9G/AyaheAnPOf+cFSzMMsXKJRUg+JiDXAf6X2e39ujuPHgbdmZuNPuyT1jOX0OMx3B3A78G3gya2LI6mE/hPw6Yj4A+Be4EygAjy//vxTgf+3oGyS2u+/UTvPf48frKr4Bmo/C36nwFySOqzpwiEiTgN+jdqqShdSW5rtHcBn2xNN6qz7mlzH/dCRCgCnr+1vqs2RE05WrMy8KSKeDFwEnAHsB/45M2fqz98A3HCsr42ItfXnTqL28+aTmfmmRfu8FPhTakUJwPsy8wNteCmSHp0twI9l5v31x7dFxG7gLiwc1OW8jsvyLKfHYRL4J+BjwK/N+wEidb2Rkeb/vJ+pT2RfNXxG43aX2XZZ1YuEYxYHDRwFnp2Z34uIAeCrEfF3mfm1Rfv998x8xQkHLZmpqSlmHzjS0uGFs4enmKqsbVl7UhPuA04B5v/eP5nahwiSVpDlFA4/lplt+yEREVcDlwIHMvP8dv0/0rEsZ4WquZWvtm3b1q44hYuI2zPzx+v35y/LuECjpRgzswp8r/5woH47ZluSSuuvgc9FxJ9Tu3bL2cDLgY9ExLPndsrMLxaUT3rUvI7L8ixncvT++gSpp1KbDNk377lW/LD4ELXJlx9pQVuSTsz/Oe/+CS3LGBH9wARwLvD+zNx9jN1+LSKeAfwr8OrMvPtE/s+yGBoa4sH+Iy1fHGHoNHsc1FFzn5K8YdH2/1y/Qe0DgR/tWCJJhVjOHIefBT5BbazyIDANrKN2QagT/mGRmTdExDkn2o6kE5eZX513/4SWZczMCrC+Pk9qZ0Scn5m3zNtlF3BtZh6NiG3Ah4FnL24nIrYCW+ttMjw8fCKxOmJgYABo/aIzAwMDXfH6u8XAwABH29h2tx+rzHxS0RkklcNyhiq9G7giM98dEd/JzNMj4o24DKPU0yLiJOCN1JZiHMrMUyPiF4GnZOb7lv7qH8jM+yPiS8AlwC3zts+fPfYB4IrjfP0OYEf9YfXgwYPLeyEFmJmZaVu73fD6u0W7jtNc291wrBrNxar3HF5IberWvcDu+ocCklaQ5RQOTwGuXLTt7cC3gHe2LNES/MRxYbvd8Pp7Ue2YspK+/++mtuzii4G/q2+7tb59ycIhIh4HzNSLhpOBX6C2Gtv8fc6YN3/qOdSWepZUEhHxk8DfAGupzXE4CzgSEc/PzL2FhpPUUcspHL5LbYjS/cD+iDgPmAIe245gx+Injgvb7YbX34vmjmmvff+X+MTxecC5mfn/RcTDAJl5b0Sc2USzZwAfrn9auar2pXldRLwFGM/MzwK/FxHPAWaBQ8BLT/ClSGqtq4H3A+/KzGpE9AGvBq4CNhaaTFJHLadw+DTwK9SWY70a+BIwA3yyDbkklcdDLPpZUe9JaLhAdWZ+A9hwjO1vnHf/9cDrTzym9Oh5HZclPQV4T32VNOrFw5XAmwtNJanjlrOq0qvm3X9nRHyN2uToz7ciSERcCzwTGI6Ie4A3ZeZVrWhb0gn5BLVeg1dDbWgR8B7g44WmklrE67g0dD21YYQ7523bDPxtMXEkFWU5qyqdCTyQmd+B2qorEXE68ARqF4c7IZn5ohNtQ1JbvIHavISbqV0E6k7gr4A/KjKU1Cpex+WRIuKv+cE1V/qBj0fEBLWVFM+mNkTpMwXFk1SQ5QxV+htql5b/zrxtZ1JbBeWCVoaSVB6Z+RC18cyvrg9ROjg3ZEFSz7pr0eP5SyjfRotGG0jqLstaVSkzb56/ITNvjointTiTpBKJiHcBXwZuyMx/LziOpM74StEBJJXPcgqHf4+IczPz+59CRMS5NDFBUlJX+x7wWmpDFe4E/qF+s5CQetdV1IYq9S2xT0uuFh0Rl1Bb7r0f+EBmvv1E25TUHsspHK4GPhURfwh8E/gx4I+pDVWS1KPmVkCqXwjuQmqrq11NbSnmxkvLSOo6nbpadH2p5vdTu8bLPcDXI+KzmXlbJ/5/ScuznMLh7dSWX30ntYlR36b2icSftSGXpJKIiMcCFwM/R23lsydSG9/8DwXGktQbng7clZnfBIiIjwPPpTaPQlLJLFk4RMSzF22aqN/6+MFqC88Cvtj6aJJK4jvAPuC9wO9k5h3Fxukus4enOLT7uiX3qTzwXQD6Tzm1qfY4rZlr70ld4UxqKzXNuQcXXJFKq1GPw/GuozBXNMwVECc8xlFSab0ZeAa1ZVkvjYi5OQ43ZmZ7Lo3eI5pdw3//g7XF6s44bW3jnU87s1euDSA1LSK2AlsBMpPh4eGCEzU2MDAAHGlLu93w+rvFwMAAR9vYdq8dqyULh06NcZRUXpn5VuCt9bHIPw08j9oFofqpzXPQcTR7fYCVcm0A6RjupTb8ec5Z9W0LZOYOYEf9YfXgwYMdiHZiZmba87nKzMwM3fD6u0W7jtNc291wrJbzYdRy5jjoUZiammL2gSMNhyosx+zhKaYqTXwyKbVA/UKPP1e/PQt4KrUhi85xkHSivg48OSKeRK1geCFwWbGRJB2PhYOkRu4BbgRuAF4D/FNmPlhsJEm9IDNnI+IV1BZc6AeuzsxbC44l6TgsHNpsaGiIB/uPcPoFl7aszUO7r2OombHQ0qO0aGGEXwfmFwoXRQQAmenCCJJOSGZeT234o6SSs3CQdCzHWxhhPhdGkCRpBbFwkPQILowgSZIWW1V0AEmSJEnlZ+EgSZIkqSELB0mSJEkNWThIkiRJaqg0k6Mj4hLgSmrrOH8gM99ecCRJJygi1lK7/sNJ1H7efDIz37Ron5OAjwAbgSngNzJzX4ejSpKkBkrR4xAR/cD7gV8GzgNeFBHnFZtKUgscBZ6dmT8FrAcuiYgLF+3zu8B3MvNc4N3AOzqcUZIkNaEUhQPwdOCuzPxmZj4EfBx4bsGZJJ2gzKxm5vfqDwfqt+qi3Z4LfLh+/5PAaET0dSiiJElqUlkKhzOBu+c9vqe+TVKXi4j+iNgLHAC+kJm7F+3y/fM/M2eB7wJDnU0pSZIaKc0ch2ZExFZgK0BmMjw8XHCixgYGBoAjbWm3G15/L6odU/z+NykzK8D6iDgN2BkR52fmLcttpxvP/2b5nuoeHitJK1lZCod7gbPnPT6rvm2BzNwB7Kg/rB48eLAD0U7MzMxM29rthtffi+aOaa99/0dGRtrafmbeHxFfAi4B5hcOc+f/PRGxGjiV2iTpxV/fded/s3r1PdWLevVYtfv8l9QbylI4fB14ckQ8idofES8ELis2kqQTFRGPA2bqRcPJwC/wyMnPnwV+G/hn4AXAFzNz8TwISZJUsFLMcaiPa34F8Hng9tqmvLXYVJJa4AzgSxHxDWofEHwhM6+LiLdExHPq+1wFDEXEXcBrgP9SUFZJkrSEsvQ4kJnXA9cXnUNS62TmN4ANx9j+xnn3jwC/3slckiRp+UpTOPSy2cNTHNp93ZL7VB74LgD9p5zaVHuc5qJTkiRJ6hwLhzZrdsLZ/ge/A8AZp61tvPNpZzqRTZIkSR1l4dBmmzdvbmq/7du3A7Bt27Z2xpEkNWnXrl1MTk4u2LZ//37gBz+z54yMjDT9816SupWFgyRJTVqzZk3RESSpMBYOkiQdgz0IkrRQKZZjlSRJklRuFg6SJEmSGnKokrQEJ0eq1ZbzngLfV5Kk8rDHQVqmNWvWOEFSLeV7qntMT0+zfft2Dh8+XHQUSeo4exykJfhJr1rN91R3GxsbY9++fYyNjbFly5ai40hSR9njIElSE6anp5mYmKBarTI+Pm6vg6QVx8JBkqQmjI2NUa1WAahWq4yNjRWcSJI6y6FKkiQ1Ye/evVQqFQAqlQp79uxxuNIJiIhfB94M/Djw9MwcLzaRpEbscZAkqQnr16+nv78fgP7+fjZs2FBwoq53C/B84Iaig0hqjoWDJElNGB0dpa+vD4C+vj5GR0cLTtTdMvP2zPyXonNIap6FgyRJTRgcHGTjxo309fWxadMm1q1bV3QkSeoo5zhIktSk0dFRDhw4YG9DkyLi74EnHOOpP8zMzyyjna3AVoDMZHh4uEUJ22dgYAA40pZ2u+H1d4uBgQGOtrHtXjtWhRcOTo6SJHWLwcFBtm3bVnSMrpGZP9+idnYAO+oPqwcPHmxFs201MzPTtna74fV3i3Ydp7m2u+FYjYyMNL1vGYYqOTlKkiRJKrnCexwy83aAiCg6iqQWi4izgY8AjweqwI7MvHLRPs8EPgN8q77p05n5lk7mlNR5EfE84M+BxwF/GxF7M/OXCo4laQmFFw6Setos8NrMvCki1gETEfGFzLxt0X5fycxLC8gnqSCZuRPYWXQOSc3rSOGwkidHNas2iYqeek1SZu4H9tfvH46I24EzgcWFgyRJKrmOFA4reXJUs+Ym5/TSa1J3Wc7kqEcjIs4BNgC7j/H0RRHxP4FJ4Pcz89a2hpEkScvmUCVJbRcRjwU+BbwqM6cXPX0T8COZ+b2I+BXgb4AnH6ONnu1xlCSpGxReODg5SuptETFArWi4JjM/vfj5+YVEZl4fEX8REcOZeXDRfj3b4ygVrd09jpJ6Q+GFg5OjpN4VEX3AVcDtmfmu4+zzBOB/Z2Y1Ip5ObZnoqQ7GlCRJTSi8cJDU0y4GXgLcHBF769veADwRIDP/EngB8LKImAUeBF6YmdUiwkqNTE9Pc+2113LZZZexbt26ouNIUkdZOEhqm8z8KtDXYJ/3Ae/rTCLpxIyNjbFv3z7GxsbYsmVL0XEkqaPKcOVoSZJKb3p6momJCarVKuPj4xw+fLjoSJLUUfY4FGDXrl1MTk4u2LZ//34Atm/f/oj9R0ZG2Lx5c0eySZKObWxsjGq1NoquWq3a66AlTU1NMTN9mAN//+Hvb6tWZqG6jJGYfX309f/gT7Xq7AxTFYfItdp9D8zywTsOLdh26EiFhx5u7litWdXH6Wv7H9FmLy45YOFQEmvWrCk6giRpCXv37qVSqQBQqVTYs2ePhYOO6zGPeQwPPfTQgm0z1T4ebvKPUYBVq/oYWD1vcMjqk3jMYx7Tqoji+CuK9U1N0bfo+B1P35o1rBoaWtjuEm13MwuHAth7IEndZ/369YyPj1OpVOjv72fDhg1FR1KJvfKVryw6gprg32TL4xwHSZKaMDo6Sl9fba5/X18fo6OjBSeSpM6ycJAkqQmDg4Ns3LiRvr4+Nm3a5HKsklYchypJktSk0dFRDhw4YG+DpBXJwkGSpCYNDg6ybdu2omNIUiEcqiRJkiSpob7qctYTLpeuDS6V2JJXeS4Rz3+p9Tz/pZWrqfO/m3sc+nrtFhETRWfwtuKPVbco+vvke2oF33r4WHWLor9PvqdW8K2Hj1VTurlwkCRJktQhFg6SJEmSGrJwKJcdRQdQ0zxWajXfU93DY6VW8z3VPVb0sermydGSJEmSOsQeB0mSJEkNWTh0QEScExG3FJ1DUud5/ksrl+e/eo2FgyRJkqSGVhcdYAVZHRHXAD8N3Ar8FnAR8E5qx+HrwMsy82hE7AOuBX4ZmAW2Am8DzgX+NDP/svPxV4aIOAf4HPA14GeoHZcPAn8E/DDw4vquVwJrgQeByzPzXyLipcCmzHxFva3rgHdm5pc7+BJUTp7/XcDzX23i+d8FPP+bY49D5zwV+IvM/HFgGngN8CHgNzLzP1D74fGyeft/OzPXA1+p7/cC4EJqb2C117nAnwFPq98uA34W+H3gDcAdwH/MzA3AG4E/KSinuofnf/fw/Feref53D8//BiwcOufuzPzH+v2PAqPAtzLzX+vbPgw8Y97+n63/ezOwOzMPZ+a/A0cj4rSOJF65vpWZN2fmw9Q+HRrLzCq1Y3EOcCrwifq41XcDP1FYUnULz//u4fmvVvP87x6e/w1YOHTO4nVv72+w/9H6vw/Puz/32CFm7bX4+z3/WKwG/hj4UmaeD2ym1mUJtW7l+efUWqQaz//u4fmvVvP87x6e/w1YOHTOEyPiovr9y4Bx4JyIOLe+7SXAPxSSTMt1KnBv/f5L523fB6yPiFURcTbw9A7nUnl5/vcOz38tl+d/71jx57+FQ+f8C/DyiLgd+CFqXVyXU+vyuplaNeukp+5wBfC2iNjDwk9//hH4FnAb8F7gpgKyqZw8/3uH57+Wy/O/d6z4898rR0uSJElqyB4HSZIkSQ1ZOEiSJElqyMJBkiRJUkMWDpIkSZIasnCQJEmS1JCFgyRJkqSGLBwkSZIkNWThIEmSJKkhCwdJkiRJDVk4SJIkSWrIwkGSJElSQxYOkiRJkhqycJAkSZLU0OqiAwBExFrgBuAkapk+mZlvKjaVJEmSpDll6XE4Cjw7M38KWA9cEhEXNviaqjdv3lp+6xZFf5+8eevFmyQtqRQ9DplZBb5XfzhQvzX8ITY5OdnOWNKKMjIyUnSEZfH8l1qn285/ScUoS48DEdEfEXuBA8AXMnN30ZkkSZpvenqa7du3c/jw4aKjSFLH9VWr5eqdjIjTgJ3A/52Ztyx6biuwFSAzNz700EMFJJR605o1awD6is7RpKo9DirCzp07ufHGG7ngggurHwpOAAAV40lEQVTYsmVL0XFapt7j0C3nv6SClK5wAIiINwIPZOY7l9jNPxykFuqyPxw8/9Vx09PTXHHFFczOzrJ69Wr+4A/+gHXr1hUdqyW67PyXVJBSzHGIiMcBM5l5f0ScDPwC8I6CY0ns2rXrEWPpp6amABgaGlqwfWRkhM2bN3csm6TOGhsbY+7Dtmq1ytjYWE/1OkhSI2WZ43AG8KWI+AbwdWpzHK4rOJN0TA899BAOk5NWnr1791KpVACoVCrs2bOn4ESS1Fml6HHIzG8AG4rOIS12rB6E7du3A7Bt27ZOx5FUoJ/4iZ/gpptu+v7j888/v8A0ktR5ZelxkCSp1BbPCSzjHEFJaicLB0mSmnDbbbcteHzrrbcWlESSimHhIElSE9avX8+qVbVfm6tWrWLDBkfYSlpZLBwkSWrC6OjogsJhdHS04ESS1FkWDpIkNWFwcJCNGzfS19fHpk2beuYaDpLUrFKsqiRJUjcYHR3lwIED9jZIWpEsHCRJatLg4KBLMUtasSwcJBUuIvYBh4EKMJuZm4pNJHnleElazMJBUlk8KzMPFh1CWopXjZe0klk4SJJ0DF45XpIWclUlSWVQBf5HRExExNaiw0iSpEeyx0FSGfxsZt4bET8MfCEi7sjMG+bvUC8otgJkJsPDw0Xk1Ao3MDAA4PtP0opk4SCpcJl5b/3fAxGxE3g6cMOifXYAO+oPqwcPOh1CnTczMwNAr73/RkZGio4gqQs4VElSoSLiMRGxbu4+8IvALcWmkiRJi9njIKlojwd2RgTUfiZ9LDM/V2wkSZK0mIWDpEJl5jeBnyo6hyRJWppDlSRJkiQ1ZOEgSZIkqSELB0mSJEkNWThIkiRJasjCQZIkSVJDFg6SJEmSGrJwkCRJktRQKa7jEBFnAx+hdiGoKrAjM68sNpUkSZKkOWXpcZgFXpuZ5wEXAi+PiPMKziRJkiSprhSFQ2buz8yb6vcPA7cDZxabSpIkSdKcUhQO80XEOcAGYHfBUSRJkiTVlWKOw5yIeCzwKeBVmTl9jOe3AlsBMpPh4eEOJ5RgYGAAwPefJElaUUpTOETEALWi4ZrM/PSx9snMHcCO+sPqwYMHOxVP+r6ZmRkAeu39NzIyUnQESZJUYqUYqhQRfcBVwO2Z+a6i80iSJElaqCw9DhcDLwFujoi99W1vyMzrC8wkqYMioh8YB+7NzEuLziNJkhYqReGQmV8F+orOIalQr6S2otpg0UEkSdIjlaJwkLSyRcRZwK8CbwVeU3Ccttq1axeTk5MLtk1NTQEwNDT0iP1HRkbYvHlzR7JJkrQUCwdJZfAe4HXAuqKDFOGhhx4qOoIkSQ1ZOEgqVERcChzIzImIeOYS+/XEcsyXX375I7ZdccUVALzuda/rdBwtk8sxS1rJLBxKYnp6mmuvvZbLLruMdetW5IeuWrkuBp4TEb8CrAUGI+Kjmfmb83fq5eWYe3WJ317Uq8fK5ZglNaMUy7EKxsbG2LdvH2NjY0VHkToqM1+fmWdl5jnAC4EvLi4aJElS8SwcSmB6epqJiQmq1Srj4+McPny46EiSJEnSAg5VKoGxsTGq1SoA1WqVsbExtmzZUnAqqfMy88vAlwuOIUmSjsEehxLYu3cvlUoFgEqlwp49ewpOJEmSJC1kj0MJrF+/nvHxcSqVCv39/WzYsKHoSCvOsdbWP579+/cDsH379qb2dx1+SZLUCywcSmB0dJSJiQkA+vr6GB0dLTjRyjM5Ocnkt/fxhFManxID9d6hhw/e03Df+x6YPeFskiRJZWDhUAKDg4Ns3LiRG2+8kU2bNrkca0GecMpqLn/a6S1t84N3HGppe5IkSUWxcCiJ0dFRDhw4YG+DJEmSSsnCoSQGBwfZtm1b0TEkacVxjpMkNcfCQZK0ojnHSZKaY+EgSVrxnOMkSY15HQdJkiRJDVk4SJIkSWrIwkGSJElSQxYOkiRJkhqycJAkSZLUkKsqSSpURKwFbgBOovYz6ZOZ+aZiU0mSpMXscZBUtKPAszPzp4D1wCURcWHBmSRJ0iKl6XGIiKuBS4EDmXl+0XkkdUZmVoHv1R8O1G/V4hJJkqRjKU3hAHwIeB/wkYJzSOqwiOgHJoBzgfdn5u6CI0mSpEVKUzhk5g0RcU7ROSR1XmZWgPURcRqwMyLOz8xb5u8TEVuBrfX9GR4eLiBpewwMDAD01GvqJgMDAxxtY9seV0m9ojSFg6TuFhGrgf8L+DlgGOibey4zn9FMG5l5f0R8CbgEuGXRczuAHfWH1YMHD7YidinMzMwA0EuvqZvMff/b1XY3HNeRkZGiI0jqAl1VOPTyJ44qlp84tsS7gWdT++P+rcAfAi8DPr7UF0XE44CZetFwMvALwDvanFX6vqmpKY4+MMsH7zjU0nbve2CWk6amWtqmJBWpqwqHXv7EUcXyE8eWfOL4fOCizPx2RPxRZl4ZEZ8HtgNvXuLrzgA+XJ/nsArIzLzuRMNIkqTWWlbhEBG/A7wIGAEmqX2SeHV9VRRJK9spwN31+w9GxCmZeUdEbFjqizLzG8CS+0jtNDQ0xMPVB7n8aae3tN0P3nGIVUNDLW1TkorUdOEQEVcAzwXeA/wv4InA7wNPBV53okEi4lrgmcBwRNwDvCkzrzrRdiV1zO3A/wHcCIwDb46IaeDeQlNJkqSWWE6Pw0uBn87Me+Y2RMTfAjfRgsIhM190om1IKtQrgUr9/muA/wasoz4vSZIkdbflFA6H67fF26ZbF0dSt8rMr8+7fyfw8wXGkSRJLbacwuE9wKcj4u3APcDZwP8DvDsifnRup8z8ZmsjSuoGEfEsYF9mfisizgDeTq0H4g2ZeV+x6SRJ0olaTuFwZf3fZy3aPgq8t36/CvSfaChJXekvgF+q3/+z+r+z1FZCe04hiSRJUss0XThk5qp2BpHU9c6sL8W6mloB8SPAQ9RWYJMkSV1uOasqnQk8kJnfmbfth4CTM9M/DCRNR8TjgfOB2zLzexGxBhgoOFdhdu3axeRk4x+P+/fvB2D79u1NtTsyMsLmzZtPKJskScu1nKFKfwP8DvCdedvOAj4AXNDKUJK60p8DXwfWAK+qb7sYuKOwRAWbnJxk3933snrd0mv5V6q1Dt177j/SsM3Zw16JWJJUjOUUDk/JzJvnb8jMmyPiaS3OJKkLZeY7ImInUMnMf6tvvhf4TwXGKtzqdUOcfsGlLWvv0G4vqi1JKsZyCod/j4hzM/OuuQ0RcS7gx1+S5vwbcGFE/DS1uQ1fy8xKg6+RJEldYDmFw9XApyLiD4FvAj8G/DG1oUqSVriI+ElqQxrXUluy+SzgSEQ8LzP/Z6HhpAbue2CWD95xqOF+h47U6uDT1zZeQPC+B2YZOeFkklQeyykc3g7MAO+kdg2HbwNX8YNlFyWtbFcD7wfelZnViOgDXl3fvrHQZNISRkaa//N+pj6RfdXwGY3bXWbbklR2SxYOEfHsRZsm6rc+atdsgNp1Hb7Y+miSusxTgPdkZhWgXjxcCby50FRSA8tZoWpu5att27a1K44klVajHoerjrN9rmiYKyB+9Dj7qUnT09Nce+21XHbZZaxbt67oONKjcT21C73tnLdtM/C3xcSRJEmttGThkJlP6lSQlW5sbIx9+/YxNjbGli1bio4jPRr9wMcjYgK4m9qQxo3AZyLiI3M7ZeZvzf+iiDgb+AjweGofROzIzCuRJEml4tWgS2B6epqJiQmq1Srj4+McPny46EjSo3EL8CfA54Hb6v/+CXArtdWW5m6LzQKvzczzgAuBl0fEeR1JLEmSmracydFqk7GxMarV2uivarVqr4O6xqJ5UF853n6Zedx5UJm5H9hfv384Im4HzqRWfEiSpJKwcCiBvXv3UqnUlvirVCrs2bPHwqHDpqamONrkcozLcd8Ds5w01dOXOjnePKj5mp4HFRHnABuA3SeQSZIktYGFQwmsX7+e8fFxKpUK/f39bNiwoehIUlNaOQ8qIh4LfAp4VWZOH+P5rcDW+v/L8PBwq/7rthkYGACOtKXdbnj9vah2TPH7L2lFsnAogdHRUSYmJgDo6+tjdHS04EQrz9DQEA9XH+Typ53e0nY/eMchVg0NtbTNXhQRA9SKhmsy89PH2iczdwA76g+rBw8e7FS8R21mZqZt7XbD6+9Fc8e0177/Xm9CUjOcHF0Cg4ODbNy4kb6+PjZt2uRyrFpR6heKuwq4PTPfVXQeSZJ0bPY4lMTo6CgHDhywt0Er0cXAS4CbI2JvfdsbMvP6AjNJkqRFLBxKYnBw0CuRakXKzK9Su5ikJEkqMQsHSWqTqakpZh84wqHd17WszdnDU0xV1rasPUmSmuUcB0mSJEkNlabHISIuAa4E+oEPZObbC44kSSdkaGiIB/uPcPoFl7aszUO7r2PoNHscJEmdV4oeh4joB94P/DJwHvCiiDiv2FSSJEmS5pSicACeDtyVmd/MzIeAjwPPLTiTJEmSpLqyFA5nAnfPe3xPfZskSZKkEijNHIdmRMRWYCtAZjI8PFxwIvWKgYEBjraxbd+rkiSp25WlcLgXOHve47Pq2xbIzB3AjvrD6sGDBzsQTSvBzMxMW9vuhvfqyMhI0REkSVKJlaVw+Drw5Ih4ErWC4YXAZcVGkiRJkjSnFHMcMnMWeAXweeD22qa8tdhUkiRJkuaUpceBzLweuL7oHJIkAezatYvJyckF2/bv3w/A9u3bF2wfGRlh8+bNHcsmSUUoTeEgSVLZrVmzpugIklQYCwdJaqPZw1Mc2n3dkvtUHvguAP2nnNpUe5zmatWdYA+CJC1k4SCpcBFxNXApcCAzzy86T6s0u1LV/ge/A8AZp61tvPNpZ7oCliSpEBYOksrgQ8D7gI8UnKOlmv3Eem68/LZt29oZR5KkE1KKVZUE09PTbN++ncOHDxcdReq4zLwBOFR0DkmSdHwWDiUxNjbGvn37GBsbKzqKJEmS9AgOVSqB6elpJiYmqFarjI+PMzo6yrp164qOJZVKRGwFtgJkJsPDwwUnap2BgQGAnnpNkqTeY+FQAmNjY1SrVQCq1SpjY2Ns2bKl4FRSuWTmDmBH/WH14MGDRcZpqZmZGQB66TWpuzjhXlIzHKpUAnv37qVSqQBQqVTYs2dPwYkkSZKkhSwcSmD9+vX09/cD0N/fz4YNGwpOJHVWRFwL/DPw1Ii4JyJ+t+hMkiRpIYcqlcDo6CgTExMA9PX1MTo6WnAiqbMy80VFZ5AkSUuzx6EEBgcH2bhxI319fWzatMmJ0ZIkSSodexxKYnR0lAMHDtjbIEmSpFKycCiJwcFBrxorSZKk0nKokiRJkqSGLBwkSZIkNWThUBJ33nknr3/967nrrruKjiJJkiQ9goVDSVxzzTVUq1WuueaaoqNIkiRJj2DhUAJ33nknR44cAeDBBx+010GSJEmlY+FQAot7Gex1kCRJUtlYOJTAXG/DnAcffLCgJJIkSdKxWTiUwNq1axc8PvnkkwtKIkmSJB1b4ReAi4hfB94M/Djw9MwcLzZR5734xS/mqquuWvBYkiRJKpMy9DjcAjwfuKHoIEV58pOf/P1eh5NPPplzzz234ESSJEnSQoX3OGTm7QARUXSUQr34xS/m6quvtrdBK1JEXAJcCfQDH8jMtxccSZIkLVKGHgdR63V429veZm+DVpyI6AfeD/wycB7woog4r9hUkiRpsY70OETE3wNPOMZTf5iZn1lGO1uBrQCZyfDwcIsSaqUbGBjgaBvb9r26pKcDd2XmNwEi4uPAc4HbCk0lSZIW6EjhkJk/36J2dgA76g+rBw8ebEWzEjMzM21tuxveqyMjI0X912cCd897fA9wQUFZJEnScRQ+x0GSmtHLPY4DAwMAPfWaetX999/Pjh072LZtG6eeemrRcSSpowovHCLiecCfA48D/jYi9mbmLxUcS1Ln3AucPe/xWfVtC/Ryj+Ncj1cvvaZetXPnTu68804+8YlPsGXLlqLjtEyBPY6SukjhhUNm7gR2Fp1DUmG+Djw5Ip5ErWB4IXBZsZGkR5qenmZiYoJqtcr4+Dijo6OsW7eu6FiS1DGuqiSpUJk5C7wC+Dxwe21T3lpsKumRxsbGqFarAFSrVcbGxgpOJEmdVXiPg1QW9z0wywfvONRwv0NHKgCcvra/qTYdANBYZl4PXF90Dmkpe/fupVKpnf+VSoU9e/b01HAlSWrEwkFieeN7Z/bvB2DV8BmN211m25LKa/369YyPj1OpVOjv72fDhg1FR5KkjrJwkIDNmzc3ve/27dsB2LZtW7viSCqh0dFRJiYmAOjr62N0dLTgRJLUWc5xkCSpCYODg2zcuJG+vj42bdrkxGhJK449DpIkNWl0dJQDBw7Y2yBpReqbWyGiC1UnJyeLztDQrl27vt+1Pefo0aMs5/ve19fHSSedtGDbxo0blzW8Ro/Orl27WPw+21+f43DGGQvnOIyMjHT1ManPxegrOkeTuuL8P5blvKeg+99X6g5ddv5LKog9DtIyrVmzpugI6jG+pyRJ3cAeB0lA133i6PkvtVCXnf+SCuLkaEmSJEkNWThIkiRJasjCQZIkSVJDFg6SJEmSGrJwkCRJktRQV6+qVHQAqQd1y6oqnv9S63XL+S+pIN3c49DXa7eImCg6g7cVf6y6RdHfJ99TK/jWw8dKkpbUzYWDJEmSpA6xcJAkSZLUkIVDuewoOoCa5rFSq/me6h4eK0krUjdPjpYkSZLUIfY4SJIkSWrIwqEDIuKciLil6BySOs/zX5LUKywcJEmSJDW0uugAK8jqiLgG+GngVuC3gIuAd1I7Dl8HXpaZRyNiH3At8MvALLAVeBtwLvCnmfmXnY+/MkTEOcDngK8BP0PtuHwQ+CPgh4EX13e9ElgLPAhcnpn/EhEvBTZl5ivqbV0HvDMzv9zBl6By8vzvAp7/krQ0exw656nAX2TmjwPTwGuADwG/kZn/gdofDy+bt/+3M3M98JX6fi8ALqT2C0ztdS7wZ8DT6rfLgJ8Ffh94A3AH8B8zcwPwRuBPCsqp7uH53z08/yXpOCwcOufuzPzH+v2PAqPAtzLzX+vbPgw8Y97+n63/ezOwOzMPZ+a/A0cj4rSOJF65vpWZN2fmw9Q+HR7LzCq1Y3EOcCrwifq49XcDP1FYUnULz//u4fkvScdh4dA5i9e9vb/B/kfr/z487/7cY4eYtdfi7/f8Y7Ea+GPgS5l5PrCZ2pAFqA0rmX9OrUWq8fzvHp7/knQcFg6d88SIuKh+/zJgHDgnIs6tb3sJ8A+FJNNynQrcW7//0nnb9wHrI2JVRJwNPL3DuVRenv+9w/Nf0opl4dA5/wK8PCJuB36IWhf35dS6vG+m9mmWkx67wxXA2yJiDws//f1H4FvAbcB7gZsKyKZy8vzvHZ7/klYsrxwtSZIkqSF7HCRJkiQ1ZOEgSZIkqSELB0mSJEkNWThIkiRJasjCQZIkSVJDFg6SJEmSGrJwkCRJktSQhYMkSZKkhv5/vlpMUgsc+A0AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10a825128>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axes = plt.subplots(2,3,figsize=(13,6))\n",
    "sns.boxplot(x='prognostico',y='lcavol',data=prostate,ax=axes[0][0],width=.2,palette=color_palette(['#67a9cf','#ef8a62']))\n",
    "sns.boxplot(x='prognostico',y='lweight',data=prostate,ax=axes[0][1],width=.2,palette=color_palette(['#67a9cf','#ef8a62']))\n",
    "sns.boxplot(x='prognostico',y='lbph',data=prostate,ax=axes[0][2],width=.2,palette=color_palette(['#67a9cf','#ef8a62']))\n",
    "sns.boxplot(x='prognostico',y='lcp',data=prostate,ax=axes[1][0],width=.2,palette=color_palette(['#67a9cf','#ef8a62']))\n",
    "sns.boxplot(x='prognostico',y='lpsa',data=prostate,ax=axes[1][1],width=.2,palette=color_palette(['#67a9cf','#ef8a62']))\n",
    "#sns.boxplot(x='prognostico',y=prostate.prognostico.value_counts(),data=prostate,ax=axes[1][2])\n",
    "axes[0][0].xaxis.label.set_visible(False)\n",
    "axes[0][1].xaxis.label.set_visible(False)\n",
    "axes[0][2].xaxis.label.set_visible(False)\n",
    "axes[1][0].xaxis.label.set_visible(False)\n",
    "axes[1][1].xaxis.label.set_visible(False)\n",
    "fig.delaxes(axes[1][2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Conjuntos de treinamento e teste\n",
    "\n",
    "O objetivo de um algoritmo de classificação é ajustar um modelo aos dados de forma a predizer a classe de um objeto. A intenção é obter um modelo que seja geral o suficiente para classificar de forma acurada objetos que não tenham sido usados durante a fase de ajuste (treinamento). Dessa forma, frequentemente dividimos o conjunto de dados em duas partes: treinamento e teste. O primeiro é usado durante a fase de ajuste (treinamento) do modelo, equanto o segundo é usado para avaliar sua qualidade.\n",
    "\n",
    "Se notarmos o conjunto de dados de câncer de próstata usado como exemplo, existe uma variável adicional *train*. Essa variável indica se o objeto será usado para o treinamento ou teste do modelo. Assim podemos dividir os dados em treinamento e teste a partir dessa variável."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(67, 10)\n",
      "mau    42\n",
      "bom    25\n",
      "Name: prognostico, dtype: int64\n",
      "(30, 10)\n",
      "mau    20\n",
      "bom    10\n",
      "Name: prognostico, dtype: int64\n"
     ]
    }
   ],
   "source": [
    "treino = prostate[prostate.train=='T'].drop('train',axis=1)\n",
    "print(treino.shape)\n",
    "print(treino.prognostico.value_counts())\n",
    "teste = prostate[prostate.train=='F'].drop('train',axis=1)\n",
    "print(teste.shape)\n",
    "print(teste.prognostico.value_counts())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "O conjunto de dados original foi dividido em 1/3 teste e 2/3 treinamento (em geral, divide-se o conjunto em 70-30% entre treino e teste). Fortuitamente mantivemos a mesma proporção de objetos das classes bom e mau tanto nos dois conjuntos. Isso nem sempre ocorre. É bom o analista se atentar e refazer a amostragem caso seja de interesse preservar as características."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Validação cruzada\n",
    "\n",
    "Embora seja usado em muitas ocasiões, esse esquema de divisão entre treino e teste é bastante criticado por nem sempre fornecer uma boa avaliação do modelo. Isso é ainda mais crítico quando trabalhamos com poucos dados. Nesses casos, uma divisão enviesada dos conjuntos de treino e teste pode acarretar sobre-estimação ou subestimação da qualidade do modelo. Por exemplo, uma mudança drástica na distribuição das classes nos conjuntos de treino e teste pode levar a grandes erros de predição. Se houver preponderância de uma classe no conjunto de treinamento, o modelo tenderá a valorizar esta classe, ou seja se ajustará melhor a ela negligenciando as outras. Isso é natural uma vez que, ao se ajustar à classe dominante, o erro do modelo tende a cair. Contudo, quando a situação se inverter no conjunto de teste, as taxas de erro serão elevadas, refletindo o fato de que os dados não seguiram a mesma distribuição.\n",
    "\n",
    "Uma alternativa mais recomendada nesses casos é o uso de validação cruzada (cross-validation). Na validação cruzada, o conjunto de dados é dividido em $k$ partições (folds). A cada rodada $k-1$ partições são usadas no treinamento e a partição restante usada no teste. A qualidade do modelo é avaliada como a média dos resultados obtidos em cada rodada. \n",
    "\n",
    "O valor mais usado para $k$ é 10 (10-fold cross-validation). Porém esse valor deve ser ajustado conforme a disponibilidade de dados. A partição dos dados pode ser estratificada (mantém-se a proporção das classes) ou não."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Regressão logística\n",
    "\n",
    "Voltando ao problema da classificação dos pacientes em bom ou mau prognóstico, nossa primeira ideia para resolvê-lo nesse momento seria transformá-lo num problema de regressão linear. Ou seja, poderíamos transformar a variável prognóstico numa variável fictícia como fizemos com as variáveis categóricas introduzidas no modelo. Logo, se o modelo desse um valor próximo de 1, atribuiríamos uma classe, caso contrário atribuiríamos a outra."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>OLS Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>           <td>Ip</td>        <th>  R-squared:         </th> <td>   0.319</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>                   <td>OLS</td>       <th>  Adj. R-squared:    </th> <td>   0.287</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>             <td>Least Squares</td>  <th>  F-statistic:       </th> <td>   9.858</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>             <td>Thu, 03 May 2018</td> <th>  Prob (F-statistic):</th> <td>2.03e-05</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>                 <td>10:32:04</td>     <th>  Log-Likelihood:    </th> <td> -33.505</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>No. Observations:</th>      <td>    67</td>      <th>  AIC:               </th> <td>   75.01</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Residuals:</th>          <td>    63</td>      <th>  BIC:               </th> <td>   83.83</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Df Model:</th>              <td>     3</td>      <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Covariance Type:</th>      <td>nonrobust</td>    <th>                     </th>     <td> </td>   \n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>t</th>      <th>P>|t|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>    0.8802</td> <td>    0.400</td> <td>    2.199</td> <td> 0.032</td> <td>    0.080</td> <td>    1.680</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lpsa</th>      <td>    0.1479</td> <td>    0.068</td> <td>    2.190</td> <td> 0.032</td> <td>    0.013</td> <td>    0.283</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lcavol</th>    <td>    0.1184</td> <td>    0.060</td> <td>    1.967</td> <td> 0.054</td> <td>   -0.002</td> <td>    0.239</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lweight</th>   <td>   -0.2128</td> <td>    0.122</td> <td>   -1.743</td> <td> 0.086</td> <td>   -0.457</td> <td>    0.031</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "  <th>Omnibus:</th>       <td> 3.213</td> <th>  Durbin-Watson:     </th> <td>   2.004</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Prob(Omnibus):</th> <td> 0.201</td> <th>  Jarque-Bera (JB):  </th> <td>   2.736</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Skew:</th>          <td>-0.387</td> <th>  Prob(JB):          </th> <td>   0.255</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Kurtosis:</th>      <td> 2.384</td> <th>  Cond. No.          </th> <td>    40.0</td>\n",
       "</tr>\n",
       "</table><br/><br/>Warnings:<br/>[1] Standard Errors assume that the covariance matrix of the errors is correctly specified."
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                            OLS Regression Results                            \n",
       "==============================================================================\n",
       "Dep. Variable:                     Ip   R-squared:                       0.319\n",
       "Model:                            OLS   Adj. R-squared:                  0.287\n",
       "Method:                 Least Squares   F-statistic:                     9.858\n",
       "Date:                Thu, 03 May 2018   Prob (F-statistic):           2.03e-05\n",
       "Time:                        10:32:04   Log-Likelihood:                -33.505\n",
       "No. Observations:                  67   AIC:                             75.01\n",
       "Df Residuals:                      63   BIC:                             83.83\n",
       "Df Model:                           3                                         \n",
       "Covariance Type:            nonrobust                                         \n",
       "==============================================================================\n",
       "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept      0.8802      0.400      2.199      0.032       0.080       1.680\n",
       "lpsa           0.1479      0.068      2.190      0.032       0.013       0.283\n",
       "lcavol         0.1184      0.060      1.967      0.054      -0.002       0.239\n",
       "lweight       -0.2128      0.122     -1.743      0.086      -0.457       0.031\n",
       "==============================================================================\n",
       "Omnibus:                        3.213   Durbin-Watson:                   2.004\n",
       "Prob(Omnibus):                  0.201   Jarque-Bera (JB):                2.736\n",
       "Skew:                          -0.387   Prob(JB):                        0.255\n",
       "Kurtosis:                       2.384   Cond. No.                         40.0\n",
       "==============================================================================\n",
       "\n",
       "Warnings:\n",
       "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
       "\"\"\""
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "treinoCP = treino.copy()\n",
    "treinoCP['Ip'] = treinoCP.prognostico.map(lambda x: x=='mau' and 1 or 0)\n",
    "model = smf.ols(\"Ip ~ lpsa + lcavol + lweight\",data=treinoCP).fit()\n",
    "model.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x10b7d1978>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAuYAAAENCAYAAABD8qCPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzsvXmUJFd95/u5EZFb7dVdXa3u1i61VhotDQKBLCQLbSAQkiAQqyQY4/Eznjl+MzzPzPHBHOwzh3nvzNh4YGwwaAODfNECMhYIARYCIWRJjYT21gbqRerq2pdcI+K+P25mVVZ1dVVmVmZWZtXvc06dqrwZce8vMisyvvHL36KMMQiCIAiCIAiCsLo4q22AIAiCIAiCIAgizAVBEARBEAShJRBhLgiCIAiCIAgtgAhzQRAEQRAEQWgBRJgLgiAIgiAIQgsgwlwQBEEQBEEQWgAR5oIgCIIgCILQAogwFwRBEARBEIQWQIS5IAiCIAiCILQA3mobsIpIy1NBOBS12gYcBjlfBeFQ5HwVhPahovN1PQtz9u/fX/c5BwYGGB4ervu8rcx6O+a1erxbt25dbROWZC2fr2LHobSKLa1iB8y3ZT2er42ild7jerCWjmetHEs156uEsgiCIAiCIAhCCyDCXBAEQRAEQRBaABHmgiAIgiAIgtACrOsYc0EQ6ovv+zcCVwBDWus3LLHdm4GHgGu11rc3yz5BEARBaGXEYy4IQj25GbhsqQ1833eB/wH8qBkGCYIgCEK7IMJcEIS6obV+ABhdZrM/Ae4AhhpvkSAIgiC0DyLMBUFoGr7vbwOuAv5utW0RBEEQhFZDYswFQWgmfwP8mdY68n1/yQ193/8U8CkArTUDAwN1N8bzvIbMK3asnFaxpVXsgNayRRCExiDCXBAETC7brKXeBNxWFOUDwLt83w+01t9duKHW+qvAV0smNqLJRKs0rxA7DqVVbGkVO2DOllwQcdzRR662OUsSBgbXa9XGpILQuogwF4R1jpmahEKuKWtprY8r/e37/s3A9xcT5YIgLE6mEDGWCThu+U1XlUwmwnUV8bgSgS4IVdAywny5Mmu+718AfA94pTh0p9b688XnLgO+CLjA17TWX2iK0YLQ5pipCZiZhni8LvP5vv9t4AJgwPf9vcBfADEArfXf12URQVinzOQCxjIBZrUNqZAwNGQyRgS6IFRBywhzbJm1LwG3LrHNz7XWV5QPFEuvfRm4GNgLPOL7/t1a62caZaggrAVmRXkd0Vp/qIptr6/r4oKwhpnOhaTdQtuI8nJEoAvrmTCs7qxtmaosFZZZW4xzgBe11i9rrfPAbcCVdTVOENYYZnK87qJcEITGMJkNmMiFq23GirECPSKTjgiDdrzFEITqCAJDNh1VtU8recwr4Vzf958A9gP/WWv9NLAN2FO2zV7gLYvtvJ6qPDST9XbM7X684fgoJhGHxFz4ikokVtEiQRAOx3gmYKZQ3YW91REPurAeCAqGbLb6c7edhPku4Bit9bTv++8Cvgtsr2aC9VTloZmst2Nu5+M1E2OQSR/6RDzOtq1HNd8gQRAWxRjDaCYgu4Y9yyLQhbVKPh+Rz9V27rZMKMtyaK0ntdbTxb/vAWK+7w8A+4ByRXFkcUwQhDIOK8oFQWgp1oMoL0dCXIS1RC5buyiHNvKY+75/BHBAa2183z8He1MxAowD233fPw4ryK8FPrx6lgpCa2GMgYkxyGZW2xRBEJYhMoaRdEC+yoSxtUC5Bz2RUDiueNCF9sFEhmzWVJ3suZCWEeYVlFl7P/BHvu8HQAa4VmttgMD3/U8D92LLJd5YjD0XhHWPFeWjkG1aAyFBEGokjKynfD2K8nLC0JBOG2IxG+KiHBHoQmsTBjae3NTh1G0ZYb5cmTWt9Zew5RQXe+4e4J5G2CUI7YqIckFoH8LIMJwuEKytPM8VUSgYgoLBiyviMRHoQmsSBIZspn4nbssIc0EQ6ocxBsZHISeiXBBanSAyjKySKK+guZ/CNvB7F5AGrtda72qWfQYo5A1BvijQ4wqlRKALrUGtlVeWQoS5IKwxjDEwNgL53GqbIgjCMuSCiNFMQLR60Ss3s3Rzv8uxFdC2Y0sR/x2HKUncSMoFeiyuiLWgQD/wWp6XnsuTngnp6HQ54ZQ4m7dU3lV5JfsfeC3PM09kmZmMMGYcANcF11Mkkop8IaJQ9NN0djucekaSzVvih6y5cdBh/56A9JQNyyiFZrguKBeiAKIFOtT1wBTH4kn7XuUzHBLWkeq071dHp0siZXjt1XDeNo4DJ55mj/fl5wuEgcFxxoknFWFxXeUYeno9TjjFbrfw9Vo4tnHQYWQoOuzjSl/j8tfJ8xQGCAoRyaTLUcfHGBiMVfQ+VYIIc0FYQ4goF4T2IV0IGc+Eq9rNU2v9gO/7xy6xyZXArcWcrl/5vt/n+/4WrfVrzbFwPgbI5w1BYIgnHLwWKbF44LU8Tz2WRTkQiymymYinHsvCTioWfrXuf+C1PI8/nD7kYz8Mbbx+eYUQpWB6MuLxf8tw7IkBe18JZtecngoZGQoWXSMMgcP0uApLuyjILlH4KzNjiCdgcjygMHTo81EEu5/Kz9ppjF03M2Nmx1AwPRnw+MMBKEUspmZfr8f/zd4NxOKOPZ7JgNGDEE/aZOLpqZDRgwHxBCSSTsWvcfl7A/b1MwYSScjlIl542r7w9RLnbVMuURCEpRFRLgjtw3Q+ZGyVRXmFLNbEb9sq2TJLFEE2E5HNRJhV/LqhxEvP5VEOeJ715HueQjl2vNH7v/RcnmBxPX0IjmMFblAwvPx8Yd6aQWGFr2MFu+fzUChUZudiYwoIguJPwcx7vYKCIQjKjqf4mgQFM+/4goCqXuPy96Z0k6OUPQ7XtXPsebmCg6oQ8ZgLwhpgNqZcRLkgtDxTuZDJ3GHcj23Kws7a/f39TV0/FnOKHvTq/Y316Oacy0yTSHjzwmtc15DLRBXNvZL9c5npikQxAEqhFEShIYwMiURsdk0TNeF/sgo7D9m4OFYKmzHGvnezU0chSs2NmShEOTbMxvO8Qx5DZa9x6b0BiMLiHKY4j+viupDPRXX7nxdhLghtzmz1FUn0FISWZzIXMtVeoryiJn4LO2uPjY01wbRDcRwbf+55VByDXo9uzokUZDPBvNCaIDAkU05Fc69k/0QK0pX2jjMGY0A59rXK5ebWVA7Q6ATk0uEtJ9AXqztoDBRtN6bo+S/7qqAUalIaUw5EITiuHVv42G67/GucSEEmXZj1jpduDJQDQRgShoZEwmGp//lNg13LHPAcEsoiCO3OxJiURBSENmAyG7SbKAe4G/i47/vK9/23AhOrFV9eCVFkyGUjMjMRhbyxjosmcMIpcUxkhZ4xNgbeRMwmJDZy/xNOieNV6GaNoqKnOaY4/uTYvDW92Arj9SvYPR6HWAWh2AsTTEtjBvC84k9MzXu9vJi9IZs9nuJr4sXUvOPzPKp6jY8/OUYU2nljcWaTYmMxG8NvIjjqeEn+FIR1j5REFIT2YTwbMJNvvSLlFTT3uwdbKvFFbLnEG1bH0uqIjE3MK+QhlrDJgI1k85Y47Dy0SkilVVVWsv/mLXHOfAtlVVnseCVVWfo2zFUb6ep2OeaEFVZlSdW7KguHVGXp6lm8KsvpZ80f6+rx5lVh6ep22XhCdVVZwtDQ3eOx/fQEe14ukMmEdHbZqiwlT7lUZREEwcYCjo1CobLEIkEQVgdjDGPZkEyh9UQ5VNTczwB/3CRz6k5kIJe1Aj0ed1buFV6CzVuqK49Yz/3L960mNGexNU8+vSYTqueth3+qZMNyx7LY63XI2MLjqfD4gsDGnxtsxZV6iu+lEGEuCG2GiUIYHabiNPwmUkGzko8Af4b90nMK+COt9RPNtVIQmkNkDKOZgFyw+lVD1jtRBNlshJO3TYoaKdCF9qdQsCFRq4HEmAtCG2HC1hXlRW4GLlvi+VeAd2itdwB/yVyymCCsKcLIMDwjorzViCLbqTE9ExHIeyMsQrCKohzEYy4IbYMV5QeLnR5ak+WalWitf1n28FfYCg+CsKYohIbRTIGgNaNXBIoCPWNwXUUYikAXLEHB3ritJuIxF4Q2wARBy4vyGvgk8IPVNkIQ6kkhjBhOiyhvF8LQMD1VIJdtjUZFwupRyEerLspBPOaC0PKYIICx4TUlyn3fvxArzM9bYpt5DUtW2gBkMerRWETsaAytYks1duSCiIPTOXprz/9b2hbXpb+/n0AEZN0pFGz5vESisQmiQmuSz0Xk861xXokwF4QWxgQFGBtZa6L8jcDXgMu11iOH225hw5KVNgBZjHo0FhE7GkOr2FKpHfkwYiQd0EjN3N/fzwv7hrjxsSG+8wdvb9xC6xRjbIKoFyjiCYXjiEBfD+TzrSPKQYS5ILQsplCwnvLFOi20Kb7vHw3cCXxMa717te0RhHrQDFEO8NAro3z1ob1kJWmxoQSBIQwN8SbUPxdWl0I+Ip9rrfNJhLkgtCAmn4PxERp+pQdMFMHPfghvOGPFc1XQrOSzwEbg//i+DxBord+04oUFYZVohijPBRHfeXqEX746BUBPwm3cYgJgvee5bERQsA16xHu+tjDGkMualqzMI8JcEFoMk83AxNihbdMasdbEKOaOW+CV3fDHf7bi+SpoVvLvgH+34oUEoQVohijfO5nj648NcWC6AMBpm1Jcd9Zg4xYU5hGGhsxM0XseF3G+FrAVeaKW/TK6ZYT5ShqT+L7/2+JYiHjghDbGZNMwPtactX7zCOafvw3ZTFPWE4S1RKHBotwYwwO/neSOZ0YJIoOr4Nqzj+TcLTEcJQKxmRggl4sIAvGetzthaEV5E/xeNdMywhzbmORLwK2Heb7UmGTM9/3LsUlhbyl7/kKt9epnCglCjZhspimi3GQzmO/fBk/8mx3o6EJ94IaGrysIawVbErFxonwmH/LNJw7yxOtpADZ1eNywc5CzjtvM2FhzbtyFQyl5z2NxRSymUCLQ24ogMOQyES2syYEWEubSmERYz5hc1oavNHqd376AueNmGB+1A9tPR131MdTGTQ1fWxDWAo0W5S+MZLhp1xDjWVuJ6c3burh2xwCpmLQdaQUMkM8b8nmD5ym8mMLzRKC3OoW8IZdr0diVBbSMMK+ShY1JDPAj3/cN8JVimbVDWE91kZvJejvmeh+vyWYI8xno66vbnIesEQRk7rmdzH3fs7HrsRid7/soifMvQSmFSiQatrYgrBUaKcojY/jB7nHu2T2GARKu4oM7BnjLkV0oCV1pSYLAJg86CivQYxLm0mqEoU3yjNqo9n/bCfPDNCY5T2u9z/f9QeA+3/ef01o/sHDf9VQXuZmst2Ou5/E2I9HTHHwdc/tNsP9VO7D1KNQ1N5AZ3EJmfNyOxeNs23pUw2wQhHankYmeY5mAm3YN8eJoFoAje+J8cucgm7sa1KlIqCuRmfOiu64Nc3E95IZqlQkKpiU6eVZLWwnzwzUm0VrvK/4e8n3/LuAc4BBhLgithMmkGxq+YoyBRx7A/PAOKBRAKTjvEtTvX4Hy2urUF4RVJRdEjGYaI8qfeH2Gbz5+kJmCFRAXHNfDVaduJOaKqGtHwtDWQFcK4nGp5LJatFInz2ppm6vz4RqT+L7fCTha66ni35cAn18lMwWhIszMNExNNG7+6UnMXd+A3U/Zgb4NqGuuRx27vWFrCsJaZCYXMJIO6p4wVggj7nxmlJ/9dhKAzpjDx87cxBuP6KzzSsJqYMxcJZd4QuHKjVZTaOX65JXSMsJ8BY1JNgN3Fcc84Fta6x82/QAEoULM1CTMTDVu/ud+g/nuN2Bm2g6ccQ7qimtRyVTD1hSEtchkLiTtFuouyl+fynPjriH2TuYB2L4xyQ1nDdKXaplLslAnwtCQSc+FuHjSSbRhmMiQybRXPPlitMynQK2NSbTWLwMrb1koCE3ATI5DeqYxc+dzNmzlkZ/bgWQK9Z4Pod745oasJwhrFWMM49mQdCGiv6O+8z60Zwr91Aj50CYNvvukfi7d3ie1ydc4pRAXJwdeqdyivOd1o9WbBlVDywhzQVjLGGNsPHmDmvmYvb+1CZ4jQ3bguJNQ11yH6t3QkPUEYa0SGcNoJiBX56/CM4WIb/3mII/ttzfmG1IeN5w9yAkbknVdR2htIgP5nKGQM3gxRSwulVxWSlCwpRBbuWlQNYgwF4QGY4yxdcNz2frPHUXwwL2Yf/0+RBG4Luqd74W3vRPlSN1jQaiGMDKMpAMKdf4q/JWxLDfuGmIkHQBw1pZOPvLGATribl3XEdoHAxQKhkLBhrnYeuhSyaUagsAQFNo7nnwxRJgLQgMxUQTjI5DP13/usWHM7TfDqy/ZgU1bUB+4AbVFyh4KQrU0QpRHxvDjlya4+7lRIgMxR3HN6Rv5vWO6RYAJs5TCXPIgXvQKsGEr7R9LfjhEmAtCgzBhaEV5oVDfeY2Bx3+F+Rc954U/90LUxe9DxaTusSBUSxgZhtMFgjrGp05kA2759UGeG7bha1u6Y3zy7M1s7Wmtc9T3/cuALwIu8DWt9RcWPH80cAvQV9zmv2it72m6oeuAci+6dBVdnCAw5LJrJ2xlMeS7bkFoACYIYPRg/UV5egbzT1/D3HmrFeXdvaiP/wnOu3wR5YJQA/kwYmimvqL8maE0//1n+2ZF+XnHdPNn521rRVHuAl8GLgdOAz7k+/5pCzb7c0Brrc8CrgX+T3OtXJ8EgU1mTE+H5HPRmvUOV0MYGnKZtS3KQTzmglB3TFCAsREIw/rO+9KzmDtumat/fuqZqCs/gursqus6grBeyAYRo3WsUR5EhrufG+XHL9lzNBVz+MgbBzh7a8ueo+cALxarm+H7/m3AlcAzZdsYoKf4dy+wv6kWrnPKu4p6Xqmr6PrzopeqrqxxTQ6IMBeEumIKeSvK61izyRQKmB9/D375EzsQT6De5cPZ50qcqiDUyEw+ZDxbv5vngzMFbtw1xO/GcwAc15/gE2cPsrEjVrc1GsA2YE/Z473AWxZs8zngR77v/wnQCbyzOaYJCwkCm+ioFHieFeiuu/YTRgv5iHzOrAtRDiLMBaFumFwOJkaoZ99u8/o+WwbxwD47cNRxqPffgNqwqW5r1BPf928ErgCGtNZvWOR5hY1nfReQBq7XWu9qrpXCemcyGzCVr9/N8yN7p/n2kwfJBgYFXLq9j3ef1I+7NhL4PgTcrLX+n77vnwt8w/f9N2it572Avu9/CvgUgNaa/v7+VTC1NlzXayt7y1FKkUg4xBPOrED3PI+BgYFVtmzlhKEhmzZ0dPTRUcd+Aq2OCHNBqAMmm4aJceoV/GaiCB76Kea+70EYgOOg3nE5vONylNvSJdZuBr4E3HqY5y8Hthd/3gL8HYd66Fqex/ZNc9ezoxyYLrC5K8aOwRRPDmVmH1916gZ2bmt8+ELJjj0TOTL5kFLVsJirSLiKo/uSs7YstHkxGx/67Sj/+/6X2TuZJzKgANex95qOUvQnXVIxh0xgZucAuOvZUZ4/mKZc63oKwqI9BluRpD/pohSMZqynui/p0Bn3SBeiefPd+vhB9k89jzGwtTvGdWcNVnwMJf76F3t54NXp2eNIeYpsYCiZGFPgulAIoOQ3V0DSUxzVm+CkjUl2j2QZy+6hP+ly8Qm9ANz30gTD6YCBDo93HNvNk0MZfrVnenb/bT0xTuhP1CTKnz4wM2/+i0/o5fTNnVXPUwX7gPIyTkcWx8r5JHAZgNb6Id/3k8AAMFS+kdb6q8BXiw/N2NhYQwxuBP39/bSTvYvhOJBIOLieYmBggOHh4dU2aUXkchGFvFkT7w3ApsHKrwcizAVhhZj0DEyO12++iTHMnbfAy8/bgQ2bUO+/HnXU8XVbo1ForR/wff/YJTa5ErhVa22AX/m+3+f7/hat9WvNsXDlPLZvmq88egDPUXTFHfZN5nhmKE1f0qU36TGaCfjKowf4Q2ioOC/ZEYQRE9lw3te8YWAohIZ9kzm+8ugBLhrJ8JNXJmdtXszGx/ZN86WHX2c8E8yKV8NcVJYxhqF0gAIGO+1x/u2vXsMYQ7YQsdABvbC0cCGy+4OtOqAUHEyHDKdDNhfn++JD+wkiG/vtKLv+3sk8X3xoP+8+qX/ZYyjx17/Yy/2/m559bID0AoMKxorycgyQCQy/G8/y4kiW7oRDf0eC8WzALb8eQilFKubSEVMcnAn4+q6Ds1+QxV3FxpRLPoTbnhrhWqhKVD99YIbbnhrBVYqOmGI8G9Y0T5U8Amz3ff84rCC/Fvjwgm1eBS4CbvZ9/1QgCRxslEFCbUQRZDIRsbiylbvaFGMMuezaq01eDVKVRRBWgElP11eUP/UY5st/NSfKd74d9X/9t7YQ5RWyWEzrtlWypSbuenYUz1EkPfvVcbpgLyDpgkEpO+45irueHW2KHaX1y1FYb3W6YPAcxfeeG5tn82I23vXsKDOFiKi4/0LKdfdELiLpOWQKEZnAkKsyKsQwP+Jrdr7AkC5EOErhOA6uUjjKiuVKjqHEA6/OebBrIVd0oWcC+54mPEUuNGSDiLgL0/mIsWwwewxdcYfBTg/PdUh4Clcp7ismgFbKfS9N4BbXKq1ZyzzVoLUOgE8D9wLP2iH9tO/7n/d9/73Fzf4T8Ae+7z8BfBsbfrZ+VVOLU8gbJicKZDMR+aLXuV0quoShIZNe36IcxGMuCDVjZqbnKqSsdK5sBvMv/wSPP2wHOjpRV34UddqZdZm/HVkYs9qImMlaYjGHM6/QnYjNxnMGUQ7XUQSRwfPsR2qnaxjOBBXPvRI7hqYLh02KCiNDZyLG61N5OstsXszG4cwrhBVcwE1xXs/zCKOsdX3XgCn7XXrtIpPDAI6jUIBR9ndkrBd9uWMoUQ8d4jq2EoRS4LkekcljDIxkQjIFeyfiOYowMgx0JubZ5TqGsWxQVdzyWHYPnXFvyXk81617LHSxJvk9C8Y+W/b3M8Db67qo0HDmxK2BHLiubVzUinXRjbFVZwr59S3IS4gwF4QaMNOTMD1Vn7l+96Lt4Dk+YgdOPA119cdR3b11mb/FqCSmFTg0ZrURMZO1xGIOpFxGMwWSnv3C0XMUhTAi7joEgY2NyAYRAymv4rlXYodXvClY7JLmOoqZnLW19LvEQhsHUi6jaTvXUqjivEEQ4Dqq5koJCivKFfY1DIIAR9nHUWRs50NjMMbgKEXMXf4YSjjKivOVXObDCDxXYQwEoX1fI5gV5R0xh46Yw0w+JJ0PSJQJnlxg6E+6VcXG9iddxrNLzzMv3nbzmvx8EBpEGBrCjMF1S5VcbOKoKoubcJzmVngxkZltqNTG0Td1R0JZBKFKwvHRuohyE4ZEP/4e5uv/y4pyL4Z6t4/6+KfXqigHuBv4uO/7yvf9twIT7RRfDnDVqRsIIhvWYIyhI2YvZB0xG9uZDSKCyMwmMjbajtL65RjAVdamIDJceUr/PJsXs/GqUzfQGXNwWFzQll8sehMO2SAiFXNIeYpElVcShRXPh8znKTpiDpExRFFEaAyRsYmblRxDifOPXllsf6KYX53yFFEUMZIOZhNZFVZEd8VdDHDR8b2ExpAL7E1ELjCExswmi1bKxSfUZx5BWIowtN7pXM6QzUZk0nM/M9MR6ZmoIQ2NjDGEoQ1TyecispmImZmIfF5E+ULEYy4IFWKMgfFRTEdq5XMdfB1zx82w73d24IgjUR/4BGpwy4rnXk183/82cAEw4Pv+XuAvgBiA1vrvsV+Zvwt4EVsu8YbVsbR2dm7r4g9htjrItp4El53Y/Kos5XaEZvGqLNt6ErO2nLgxtWRFk53buvhvl5zE/77/xcNWZRkoq8qyIeXxybMHgdqrsmxIzVVl2VQ2n63Kkgdjq5yUqrIsdwwl/uCcLaQL+3hkf3r2JiPpQj6kqqosTw9l2DOZm/0WYVOHR2fcYTIX0Ves1HL65k6O6UusuJrK6Zs7uRaaXZVFEOYRRYZ83jY1chwbAqOKd+uuW3lzoygyhIG9EYgiU8/WHmse1c7ZuyvE7N9f/wZma6FMUbWsh2M2UQhjo1DIr6h8kzEGHv0F5ge3QyFvv08872LU778H5a3ifXI8zrY3nAG158s1mjV9voodh1KrLZO5kKncyhsH/fq1Gb75xMHZ0JWLju/lylM34K1ibfLyz543n3IstPD5uvv5Pctv1SKslZJ8JRp5PI4CpyjWlaKYqG3HjLFivFCw3vF6sFbem5NOPgoqPF/FYy4Iy2DCEMaGIQiW33ipeaYnMd/9Jjz/pB3o7Uddcz3quJPqYKUgCOPZgJkVNg7KhxF3PD3Cz39nw9V6kh4fe+MAp29ujQ4nroKuREv3MhDWMJGBaF7VFPu3UnVr47HuEWEuCEtggsCK8nBlHjjz/JOYu74BM8XY9B1vQr3nQ6hUa1zsBaHdGcsEpAsrE+X7p/Lc+NgB9k8VADhlIMV/uGA75KaX2bPxeA5s6IiRCmJrvgW70H6IKK8fLSPMV9LK2/f964A/L276V1rrWxppa/TkY5h774ThAzCwGXXp1Tg7djZyyZqZZ2syZW9rM2kY2Awn77De20WOY3a/fa/azpOl+gmxGGw5atFjjp58zMZNv7bHnqXGMBucluyAi6/Efc+1i28LEE/YtPDSOl7M1ggP7EUSL2btL+TnFnUc2Ho06urrlrY9l4coZNG0NuXYuT1v3vGpU8+YFeXR7qfh3jthZIgRpWDjIJx+Fryy227TPwDnXYJz0unWht1Pwy9+BKMHbeeHUlnFZAp1xbWoM84h2v005hc/WnT/ejBrQ4PmF4RWwBjDaCYgu4Lax8YYHnx1iu88NUIhMjgK3nPyBi4+sZf+jhhjuToaXCWeA11xl46YQ1fCIyuiXBDWNC0jzKmxlbfv+xuwCWZvwqqux3zfv1tr3ZCgpOjJxzDf+ooVcZ3dMDGG+dZXiD78hy0nzufZqhx4fa99YsMmOLAfXngGevugu2/ecQB2v6AAmRkraI2xIjifhaH9hxxz9ORjmJu+CNMT82+dTdGDlcvA9/+JEFDHbl9821zW/nYcu+5CSgJ93kFGsO93mJv/luj6/3B425fCRFDIQZCHvGOP7x//HnPFB3FOPNUK3DtvsfMpx853YJ99DXt67f8qeUoKAAAgAElEQVTB1AR8/zaiK+yNB9+/zdo2PVW8OQA2b0V99I9RfRvsnN+/zWafpTrn7V8P8VzV/LE4dPWseE1BaDbGGEYyAbkViPJ0IeRbTwyz6zX7mbMx5fGJnYMc15+sl5k1EXMUXQmHVLGhkiAIrY8xhmzaMD0VMj0ZMT0VMTMVctLJlc/RMsK81lbe2AoQ92mtRwF8378PuAzboazumHvvtEI3UfzQTiSBrB1vMWE+z9axfcwWLJ0sem8VkE5DT//84wC739TE/OYhxoDj2n36N847ZnPvnVZ8H9YYY7NG7vse5ujjl9+2qgM1kEkvbXul8ygHZmagpw8e+CGceKr1Ouez9jlHAapYxsFAJmNFbTwB+Zzd1hh7k5Eu+/q7owtSnai+Ymm3X/zIiuZ4wj4u378eXu1K5vc86OpGJSWcRmg/ImMYSQfkV5Bk9tJolpt2DTGasTfPZ2/t5MNvHKAjtnox3HFX0Z1w59VrFwSh9QgKVoDPTEVFEW7/XmHka+sI8wo4XCvvilt816OT4MGxYVRXzzwPhnE7MWPDDAwM1NTBr1GU2xqEQTFMhDkPruNCGMx2KywdB4Dq6iEcPmD3KXmjbEs+CAPcjvnH7IwNEy3lnTbF9XLZyratlijEWcr2SlEKwgAnmcJMjNHf38/oxBgmiubCbBas67rF1y/pEo0Ow8zkXKKoF8PZOAixOGZyfLZr3+jEGKqja/7/UdKdXXOlLDX/hs2bUV09OBLfLrQpYWTDV2oV5ZEx3PvCOP+ye4zIWO/0B96wkbcf3b0q3mkFpGIOnXGHuCuCXBBaCWMMmbRherIowqcipidDspnDf/44DnR2O3R1u3T1VHdOt5MwXzH16CQY9Q/AxNicxxysd7TflvZqpXJj82x1vTlBXhSSBHnw4rPdCkvHAdj9SvuU0q1tSz5wPcL0zLxjjvoHYHwU1GEEt8KGlSSSlW1b7fXWce28h7O9UqLIzpVNQ68t02R6+2FqvFjMuWTg3LphGNg1piZgZnoufKezC3r6iJQDZfMBmN5+zNTEnEcbrEe7tz6loRadv5CHjYOMKQ9m0vanjK1bt654XUFoNGFkGE4XCGrM8xzPBty8a4jdIzZ0bmt3nE/uHGRLd7yOVlZOR8yhO+GuahlGQRAstXjBE0lFV8+cCO/qdkl1qppv8ttJmB+ulfc+bDhL+fj9jTJCXXq1jWEmOxceEASoS69u1JI1M8/W7l4YO2if6OmFQgEm89DRYUXlguMw3/oKpDrmkhahKHIju8+CY1aXXm3jxssTM+cZo6zYvvhKZmPMl9y2CjGtFKQ6lra9UkwEqS5bheW8S+zYeZfMxZhH5X0RFaRS9gZgbMS+hmATXT3PJtui7Hj5fKU5v3+bfS4Wt6/Fwm1WQvn88YT14BtQl7+/PvMLwipQCA2jmdpF+ZMH0nzj8SGmiyUVzz+2h6tP27AqXuqkp+hJuMTEQy4ITccYQ2amGAs+FTFTFOFLesFd6Oxy6Opx6ep26Opx6Ox2iS3SfXkltJMwvxv4tO/7t2GTPye01q/5vn8v8N993y99/38J8F8bZYSzYyfRh/+wLaqyHGLrEUfOVWXZvBXOv/TwVVlK+wXBoVVZBrcecszOjp1EN/zHw1dlSaTmV2VZuC3UryrLYrYvVZUFBTEPHM8K6o2D8yqYOCedTnT1dbNVWVAKNm+zVVmeedwmgpaO4ZQzUO/7CGbfq0tWRHFOOt0mijaoaopz0ulEV30Mfl6sDNPC/6eCUAn5MGIkHVBLp/BCaPjesyP89JVJwHqpP3rGJs7c0vyumomiIJeQFUFoDoWCYWZBMub0VGQlwWFIphRd3S6dPY4V4Sv0gldDywjzWlt5a61Hfd//S+CR4lSfLyWCNgpnx86WS/Q8HMvaWhTKVe+3wn0a+RpWOrcxxoa9ZJdIRC3NedLps0mT/f39jB4cwvz47rlKN/GE9UbvfDtKKVTZ9pXMWTeSKUgmIZ7APWIbnHthfecXhFUgG0SMpoOa0k+Gpgt8fdcB9kzYG/oTNiS54axBNnQ09/IXd60gT0hSpyA0hIVe8FJM+HJe8K5u6/kuecG7ul28OnvBq6FlhLnW+kPLPG+APz7MczcCNzbCLmFtYqLIxrnnqy9QHOzfg/n631hPOcCRx6LefwNq42CdrawQVQyp6ehCebHVsUEQGsRMPmQ8W1uZg1/tmeKfnhwmFxoUcPlJfVy+vR+3ifHcngM9CY9UTAS5INSLhV7wkgiPlghzS6bUbBhKZ7cNSUl1NMcLXg0tI8wFoVmYIIDxkbnKKZXuF0Xwq39l4r7v2fAapeCCd6HecTnKXYXyaq4LHZ02vt6RFt3C2mMiG8zGg1dDNoi47TfD/Ns+W7K0L+ly/VmDnDSQqreJi1KqstIRc8RDLggrwESGmWkrwG1FFPt3Lnt4L7jrzq+IUvp7Nb3g1SDCXFhXmHzOesqXuq1ebL/Jccydt8JLz9qB/gHrJT/6+AZYuQyxOHR2oZLNERmC0GyMMRycztUkyl8dz/H1XQc4OGNvvHds7uBjZ26iK974m1dHWUHeHXeb6pUXhLVAIX9oRZT09DThEmVRy73gpWTMVvSCV4MIc2HdYDJpm0xaZQMj8/QuzPe+NduNNPHWC8i/80pUosmdAeMJ6OxGJRLLbysIbUqpcVCnV50oj4zhX1+e4LvPjhIaG0Jy1WkbueDYnoZfpD0HuhOudOkUhAqIIkMmHc0lY07amPBlveClkoTFMJTObgfPW3vnmwhzYV1gpiZhZqq6fXJZzL9o+PVDdiDVibryI3Sd9/t1qTVeEUrZOvSdXajY6tRZFoRmEUSGkWKN8mrqpUzlQm59fIinh2wi9+bOGJ/cOciRvY29iVVAR9yhJ+HiiCAXhEMoecFna4JPRsxMLx0LnupQNhmzx2Hzlh6UkyGZam8veDWIMBfWNMYY6yXPpJffuHy/372EueMmW58c4MRTUVd9HNXT1wArF2E2obMb5clpKqx9skHEWKb6cojPHcxw86+HmMzZBNFzj+rGf8PGhsZ2uwq64i6pmCMhK4JA0Qs+U0rEnGvMs6QX3CuLBT+MF7y/v4OxseqLNLQzcsUX1iwmimw5xFy28n3CEHP/PfCzH9iQF89DXXIVvOUClNOEJC7HgVQndHZKQqewbpjKhbPCulLCyPD958f40YvjGGzDng+/cRNv2tbVGCOxHvKuuEOXeMiFdUw+H80lYxbDUNIVeMFLwrtUknA9ecGrQYS5sCYxYWiTPA/XXXSxfUaGMLffBHt/awc2b0N94AbU5m2NMbIcz4OOLlthpc0/qHzfvwz4IuACX9Naf2HB80cDtwB9xW3+i9b6nqYbKqw6xhjGsiGZQnXx5MPpAjftGuKVoift2L4Enzh7kIHOxpULTcVsyIq3xjzky52vxW184HPYDm1PaK0/3FQjhVVh1gteKklYDEnJ55b2gnd1z0/GXKux4I1ChLmw5jC5HExUXnnFGAOPPYi55ztWyCsFb7sI9c73Nr4ueCxWrLDS0dh1moTv+y7wZeBiYC/wiO/7d2utnynb7M8BrbX+O9/3T8M2Dzu26cYKq0oYGUYzAfklKi4sxmP7p/nHJw6SDex+F5/Qy3tO2dAwwew5io0dHsk1WPawkvPV9/3t2G7ab9daj/m+v0oNG4RGks/NhaGUOmPOTEeYCrzg5WEo4gVfOSLMhTWFSU/D1GTFlVfMzBTmu9+E535jB3r6Udd8HHX8KQ20krVcYeUc4EWt9csAvu/fBlwJlAtzA/QU/+4F9jfVQmHVyYe2k2c1mjwfRHzn6REefNUmcfckXD5+5iZOG2zMTa0CuhIuW3oSjBTWnigvUsn5+gfAl7XWYwBa66GmWynUjSgypKfnPOAzRW/4Ul5wz4POkgDvduzfXQ6ueMEbgghzYU1gjIGpcUhXnuRpdj+FuesbMD1pB96wE/XeD6FS1dSDqJJk0grytVthZRuwp+zxXuAtC7b5HPAj3/f/BFt8452LTeT7/qeATwForRkYGKi7sZ7nNWResePwTOcCMukCPcucAp7r0t/fD8CrYxn+9pcvs2/C5ou8cWsPf/T2Y+lL1f8bLUdBZ9yjO+nhOapl3htoyPtTyfl6EoDv+w9iw10+p7X+YT2NEBpDuRe8vCLKUn6rjk5nXhx4V49DIile8GYiwlxoe6pN8jT5POZHd8LDP7MDiSTqimvhjHMa8+GjFCRTNmSl0aEx7cGHgJu11v/T9/1zgW/4vv8GrfW8L0211l8Fvlp8aIaHh+tuyMDAAI2YV+w4FGMME7mQmQqbBvX39zM6OsoDv5vkjqdHCSKDo+DKUzdw0fG9mOw0Y5XndS9LeWJnGCjGi/f4rfLewHxbtm7d2qxlPWA7cAFwJPCA7/s7tNbj5RstvJEu3VS1A67rtZW9C4lCw9RkgYmxAhPjeZ759UHGx3Lksoc/12JxRW9fnN7+GD2l370xvBYL2Wr396YWRJgLbU21SZ5m/6uY79wEw6/bgWNORF1zPap/Y/2NUwpSHdDRtZ5KHu4Djip7fGRxrJxPApcBaK0f8n0/CQwA8hX5GiWMDGPZgFxQeezKdC7gq48e4InXrUIe6PD4xM5Bju2rf2OvtZrYWQGVnK97gYe11gXgFd/3d2OF+iPlGy28kW5ar4c60N/f37zeFCsknytLxiyriLKcF7y8Nf2hXvACUGCqulYfTaGd3pul2DRYebWodaMWhLWHKeStKA+XL7Nmogh+cR/mJ3fbpFDHQV30Xjjv4vqXQXQUpLrWa8nDR4Dtvu8fh73AXwssrODwKnARcLPv+6cCSeBgU60UmkYt8eQvjmS45Sd7GEkXAHjzti6u3TFAKlbfczXuKnqTLnG3tbyETaSS8/W72G+5bvJ9fwAb2vJyU61ch8yLBZ8MZ8V4Ib9ELHjMVkTZuClFLB7Q2W3FuOuuuxvOtkaEudCWmGwaJsYrSvI04yOY22+G371oBwY2oz7wCdTWo+trlONAZ5ftENqMmuctiNY68H3/08C92HjUG7XWT/u+/3ngUa313cB/Av7B9/0/xSaCXq+1rrKtjNAOTOVCpnIhlb65kTH88IVx/uX5MQxWOH9wxwBvPbKrrmFmngM9Ca/uQr/dqPB8vRe4xPf9Z4AQ+IzWemT1rF5bGGPI54ytglImwtMzy3jBu5z5yZjdc17wteJlXq8oU2H1ijWI2b+//sUgWikesVk0+5jNzJStvLLcdsbAbx7B/PO35+LPzzkfdek1qHjtyZeHfOg5DnR2Q0dnWyfIFGNWW/UA1vT5utbsqCV0ZSwTcPOvh3hhxJ6rx/SnuP7MjWzuql+itKOgO+7SGXcqPldb5b2BRWPMW/Z83f38nuW3ahGaJWSj0DBTqgs+aTtjLucFj8UUnT0lEW7DUDq6lvaCryVhvlaO5aSTj4IKz1fxmAttg4kimByD7PIZXyYzg/nn2+DJR+1AZzfqqo+hTt5RP4OUgo5OO/c69ZALwkKyQcRYJiCqwufzm9dn+MbjB5kpNhq64LgebnjbCcxMTtTFJgV0xh264i7u+osjF5pMuRd8VoAv4wVXClKdzmxjnpIIjyekIsp6Q4S50BaYfM5WXqkknvzl5zF33GJFPMDJO1Dv+yiqq2fpHSullNTZ2b2ekjoFYUmqrboCUAgj7nxmlJ/91n4D1hlz+NiZm3jjEZ3EXYeZOtiV8BS9CY+YxNkKDSAMF8SCF0NSKvaCF8NQlvOCC+sHURVCS2OMsXXGZ6aX3zYoYH7yz/Dgj23seSyOuvz98Kbz6uNxUApSKdzBLajx8eW3F4R1QiE0jGUCClW4yV+fynPjriH2TtqKSts3JrnhrEH6UvW5LLkKepMSRy7Uh1kv+ORca/rpqYjMMl7wxSqiiBdcWIqWEea+718GfBGbgPI1rfUXFjz/18CFxYcdwKDWuq/4XAg8WXzuVa31e5tjtdBITCFvveRBsPy2B/Zjbr8JXt9rB7Ydg3r/DaiBzSs3pNxD7rriJReEMtKFkIlsWHHoijGGh/ZMoZ8aIR8aFPDuk/u5bHsfTh3ESqljZ3cVceSCUM6sF7wkwove8KBw+H1icTUbhtJZFOCdnQ6OeMGFKmkJheH7vgt8GbgYWzP1Ed/379Zaz7YF1lr/adn2fwKcVTZFRmt9ZrPsFRqLMQZmpqyXfJnkZBNF8PD9mB/dZQW8UnD+ZagL341y61CqsKMTurrXY9lDQViSyBjGsyGZQuWhK5lCxLd/c5BH99sglQ0pjxvOHuSEDfWpTR53Ff0pbz3WIxdqwBhDLmt4LZ3hwP5csU29rQt+OA7xghdDUcQLLtSLlhDmwDnAi1rrlwF8378NuBJ45jDbfwj4iybZJjQRExRsGcQKGgaZqQnMnbfCi8V/k/6NqGtuQB1zwsoNSSahswcVk06dgrCQbBAxnqmuNvlvx7LcuGuI4bT9BuyMIzr46Bmb6Iyv/KbXUdAVd+lOyA20sDhhaJgplSScDUUpecEXz2aIxVUxEdMK8M5u8YILjadVhPk2oLy20l7gLYtt6Pv+McBxwE/LhpO+7z8KBMAXtNbfPcy+81oGDwwM1MH0+Xie15B5W5l6HLMxBjMzRTSZga5OoHPJ7XNP/Bsz3/oH61kHEm+9gI5rPo6T6qjdCMdBpTpwOrtQscOXaFuP77EgQDHBMxvOVk+phMgYfvzSBHc/N0pkIOYorjl9I793THddPIydMYfuhFRbESwlL3hJeJdXRDkcyil6wcs84DYWXPIThObTKsK8Gq4Fbtdal5fnOEZrvc/3/eOBn/q+/6TW+qWFOy5sGdyI2rStVPO2Waz0mE02YxM8K4klz2Ux93wHdv3SDqQ6Ue/9MIU3nM1ENgfZXPUGuK5tDBRPoYIIJpaukb5W3+NiXWRBWJRCGDGaCQgq1+RM5gJu+fVBnj2YAWBLV4xP7tzM1p6V1yZ3FPSnPJKeiKf1SskLXp6MOTMZLnkpiSdUMRFzToRvO3IjE5OS0N8IlALHUShl/wbbfDssft1Wup1etx11FqFVhPk+4Kiyx0cWxxbjWuCPywe01vuKv1/2ff9+bPz5IcJcaC1MoQBT45BfPmwFwOx52XbwHC12bz/hFNTV16F6+mozQCno6oaO+nYVFIS1xkzeJnhWc/F8ZijNLb8+yFTe+lDOO7qb95++kfgKhbQCOuIO3VKTfN2w0Ate+p2ZOfx/pFLQ2TU/GdPGgh/6/yehKfVBAcpROI7tu+e6Ctdb/LU1xhxy3S01vIwiMEXxvh5LSLaKMH8E2O77/nFYQX4t8OGFG/m+fwrQDzxUNtYPpLXWOd/3B4C3A/9vU6wWasJEofWQp9OVbR+GmJ/9AH72A3vGeh7q4vfBWy+svbFPMgndffVJEBWENUoY2drk1SR4BpHh7udG+fFLtjlQynP4yBkDnL21a8X2xF1FX1Jqkq9lwsAwPW0936VkzEq84OUe8M5iXXBHbtyagqMgnnTwDiPCF2MxZ1hpzHUBF7yYoqs7RibrWC97YAgCiKrpXtaGtIQw11oHvu9/GrgXWy7xRq31077vfx54VGt9d3HTa4HbtNbl78qpwFd8348ABxtjfrikUWEVMcZAesbGhUeVXejNyJAtg7j3t3Zg8zZbBvGIbbUZEYtBVy8qkahtf0FYJ2QKERPZ6hI8h2cKfH3XEL8btyFlx/cnuOHsQTZ2rCyJ2lHQHXfpkuTONYMxhmymGIoyOVeSMJNewgvuFL3g3ct7wYXmEIuphlekUUrhutYDH09YT3ohbwgDsyZDYFpCmANore8B7lkw9tkFjz+3yH6/BOrYZ11oBNXEkUNRxO/6pY0nzxfjxt92EeqdV9ZWKcV1obsXlUxVv68grCOCyDCRDcgG1V3yHtk3zbd/c5BsYGuTX3JiH1ec3L/icBNJ7mx/wsAUO2KWBLgNRQmX84KXWtMXkzFTneIFbwWUst7sWEytyvvhugo3pTDGEAZQKJjZmPW1QMsIc2FtYoIAJsfnxHUl+8xMY773TXj2CTvQ02djyU84pXoDlIKOLluLXOLIBeGwRMYwnikwNF2oyguVDSL0U8P8ao/tztuTcLn+rEFO2bSym+Ck5zDY6RFzxRvajrzyQm7WG76UF9xxoKOrLBmzWB88Hpf3vdVQCuJxBy+2eChK8+1ReDF7kxCGhqBgf9pdooswFxqCDVuZhumpZZsEzdvvhadtbfLpYmWUN+xEvedDqI6lyycuSqoDunokjlwQlmEmHzKZC+n1gqouansmctz42BAHZmxLxNMHU3z8zMEV1ROPOYqepMtgd4LhnIizduW3Lxya1J9Iqnmt6bu6xQveDjiO9Y63iiBfDNdVxVAXQxiW4tFNNfKjZRBhLtQdU4ojD8PlNy7tU8jb7p2/ut8OJJKoK66FM86p/oPAdaGnX+LIBWEZckHERDakUGUylTGG+1+Z5K5nRwgicBW879QNXHh8L06NF25HQXfCpasODYeE1ae7x6Gzp1yEu8TirSnqhMVRChJVJnWuNkopPA88TxEvhrrk86atEkZFmAt1w+RyMD0BhUJ1++3fYxM8D75mB44+AfX+61H9NTTx6eiwyZ21VmsRhHVALdVWSkznQr7xxEGePGCrKm3q9Pjk2Zs5uq/2G2GJI197vOm8Gr7lFFqGZiR1NpryUJegYMjnItpBn4swF1aMCQLM+Ahks9XtF0Xw4I8xP7nbetcdB/X7V8DvXVq9sPY86O1fsmOn0Bx8378M+CK2wtLXtNZfWGQbH/gctq/EE1rrQ8qjCvXHGMNMPmIqH9Z0gdo9nOGmXw8xkbXfhp1zZBfX7hioucmP50Bf0iMhTYIEoSVwHEgknMPWH29XvJjC9RwKeVvRpZX1eVXC3Pf9PuDdwFZgP3CP1nqsEYYJrY8JQ5iZJMxnqhfl46OYO26G375gBzYOoj7wCdS2Y6ozwnaRgE5J7qyGdDrNM888w+TkJD09PZx22ml1mdf3fRf4MnAxsBd4xPf9u8tLmPq+vx34r8DbtdZjvu8P1mVxYUlm8iFTubCq8oclwsjwgxfG+MHucQyQcBUf3DHAW4/qrtmezrhDb8KV87YCSudrEAR4nle381UQyonFFfF4e3vJl0Ip+y2AF7PivFBoTXlesTD3ff/3gTuB54HfAUcDX/Z9/xqt9U8aZJ/Qgpgogplpm9xpDFQZy21+8wjmn78NWdumm3POR116NSpe5VfhHR3QKcmd1bJ7925uvPFGBgcH2bBhA2NjY9x+++185jOfYceOFVcePQd4UWv9MoDv+7cBVwLlvQX+APhy6aZeaz200kWFw1MIDePZgHyN5cRG0wE3/XqIl0btzffRvXE+cfZmBrtqq03uKuhPiZe8UsrP1yOOOIIDBw7U83wVhLaMJV8JjqNIJBWxmCGft0mirUQ1HvMvAZ/SWuvSgO/7H8B6x2qoYye0G7MNgtLTVSV2zu6fSWO+fxv85hE70NmFuurjqJOrvLjEYraEooSt1MQdd9zBBz/4Qc4666zZsccff5yvf/3r/M3f/M1Kp98G7Cl7vBd4y4JtTgLwff9BbLjL57TWP1zpwsJ8ImOYyoXM5KOav7Z9/LUZvvnEQdLFWPSLju/lylM34NUYC56KWS+5xJJXTvn5OjAwwPDwcD3PV2Gd4yhIdqzPyjiOq0gW66EHARRaJAa9GmG+FbhjwdhdwD/UzxyhFamlY+chc7yy24auTBQjn07egXrfR1FdPdVN1NUtYSsrZGJigjPOOGPe2I4dOyi75240HrAduAA4EnjA9/0dWuvx8o183/8U8CkArTUDAzUkAy9niOc1ZN7VtmM6FzCeKRDzoK+KHDzPdenv7ycfRHzj0b38ePdBAHqSHv/+bcdy1pG9NdnjKOhPxehMVH7JWavvTbVMTU1x0UUX4TjOrC0XXnhhM89XYY2i1rEoL0cpRSwGnufYJNH86pZZrEaYfwP4Y+Bvy8b+CLi1rhYJLYOJIsika/aQA5iggPnJ9+HB+2zYSyyGuuz98Obfq05ci5e8brz5zW/m5z//Oe94xztmxx588EHOP//8eky/Dziq7PGRxbFy9gIPa60LwCu+7+/GCvVHyjfSWn8V+GrxoRkeHq6HffMoeSBXm3rZUWv5wxL9/f088+oBvv7YEPunbB3qkwdSXH/WJnqTEWNj1acUJT1FX9IjEygyU5Xvt9bem1rZuXMnd911F+94xztmbXnggQfqdb5WlKxd3O4a4HbgzVrrR+uyuLBqKAWplIjycpRSxOJlMeirlCRajTA/C/j3vu//P9gL7TZgEHjY9/0HShtprevzaSGsGiYKYWYGMtOs5HsdM/Qa5js3wut77cDWo1HvvwG16YjKJ3GUjSPv7KrZDmE+e/fu5cEHH+SnP/0pvb29TExMMD09zfbt26nDufwIsN33/eOwnxPXAgsrrnwX+BBwk+/7A9jQlpdrOxoBVlb+sIQxhp++MMzND++jEBkcBVec3M8lJ/bVVJtcYbuAdq2g2ZAw/3wtCfN6na+VJGsXt+sG/iPw8MqORmgFXFeRTCqUiPJFKSWJlsegN9ODXo0w/wckbGVNY6usTEE2vTJBbgw8/DPMvXdCULC35udfirrwiuoSNZMp6O6V5M46c+6553LuueceMt7X18dzzz339ZXMrbUOfN//NHAv1gN3o9b6ad/3Pw88qrW+u/jcJb7vPwOEwGe01iMrWXe9YoxhOh8xXWP5wxLpQsi3nhhm12szAGxIeXzi7EGO35CsaT7PsQmecVcSPFdK+fna1dXF9PQ0UJ/zlcqStQH+EvgfwGdWuJ6wyqz1yiv1RBWTREuNigoFQ1hjEn01LCvMi9VYYH5C1zy01j+tm0VC07Ee8mkbR77C20IzNYG561Z4ofi53rfRNgs65sTKJ3GUDVtJdqzIFmE+u3fvBuwFfTEGBgbQWt+y0nW01vcA9ywY+2zZ3wb4v7tjxIIAACAASURBVIs/Qo3kgojxbEBQu5McgJdHs9y4a4jRTADA2Vs6+fAZA3TEarsh7og59CbdmjuACpbFztfe3l48z16263S+Lpus7fv+2cBRWut/8X1fhHmbst4qr9ST8kZFUWTIZRsr0CvxmC93R26A4+tgi9Bk6hWyMjvfs09gvvtNG5MOcOZbUO/+ICqZqnySeMI2ChIved257bbblnze8zy+9KUvNckaoVbyYcRkLiS3whJfkTH86MVxvv/8GJGBmKO47pyjOGugttriCuhLuTULemE+i52vrusSFvN9mnG++r7vAP8LuL6Cbecla/f39zfUtnriul5b2bscC49HAZ3dMVy3/UT5aidfH45cLiSXjWyEQJ1ZVphrrY+r+6rCqmIKeesdz2ZW7CEHMLks0ff+ER79hR1IdaDe82HUjp2VT6IUdPegOiSWvFF89rOfXfL5rVu3NskSoRbCyDCZC2dLF66E8WzALb8+yPPDtpfA1u4Yn9i5mdOP3lRTgqerYGOHR0xCV+rGYudreSJqnc7X5ZK1u4E3APfbZr0cAdzt+/57FyaALkzWruX/aLXo7++v6f++VSk/HoWtvJIP2k+Uw+onXy+FMXONipaTUpsGK9c2VXX+FNobk8vZGPJ8rn5z7nmF8btuhYOv24HjT0ZdfR2qtwrvg+dB7wZUrLaGJYKwljHGMJWPmFlhHHmJpw6kufXxIabzVuCff2wPV5+2oeZ48Lir2JDypDZ5e7JksrbWegKYdVf6vn8/8J+lKkt74Di2TrdUXmkMs0micSvOC3UqsyjCfB1gctmiIM/Xb84whAd+iLn/HltW0fVQl7wP3nohyqniAp/qsAme1ewjCOuEdCFkMhtSj3DGQmi4+7lRfvLyBGBjwT96xibO3FJFofMFdMUdehK1hb4Iq0+FydpCG+J5NnFRzs3Go5RNqI3FjK2DnltZmUUR5msUY8xcDfIgqO/cowcxt98Ee14BwN1yFNHV16GO2Fb5JKqY4JmSBE9BWEg+tPXI83VKMBqaLnDjrgO8OmFvzk/YkOSGswbZ0FHbJcBR0JuUePK1wHLJ2gvGL2iGTcLKiMUckilxdjWbUh1011tZgmjLCPPlmhz4vn898P8xF//2Ja3114rPXQf8eXH8r+pRWaJdmUvonKm5S+dh5zYGdj2EuUfPhcO87SJ6P3Ad49MzlU/kurZai4SuCMI8ckHETCFaUT3yhTy8Z4rbnhwmFxoUcPlJfVy+vb/m0JOYo+hPecTaMJFMENY6rqtIdbikM6ttyfrFcRSpDkU+H5HPVS/OW0KYV9rkAPgnrfWnF+y7AfgL4E3YCjGPFfddO5kcFWCCwIrxzExdKqwcMv/MNObuf4RnHrcD3b2oa65DnXBqsRtnhcI8mYKeXpQjnjZBKJEthAxNF2ru2LnonEHEbU8O8297i3Wvky7XnzXISQNVVElaQCrm0J+U0BVBaEVKMeVyfrYG8biD59omRdXQEsKcypscLMalwH1a69HivvcBlwHfbpCtLYUJCjA9Cdls49Z48RnMnbfClI1N5fSzUe/9MKqjithUx7GCXGqTC8IshdAwmQtIu/m6ivJXx3PcuGuIoZkCADs2d/CxMzbV3IVTungKQmvjKESUtyCOa2+WqqFVhPmyTQ6KXOP7/vnAbuBPtdZ7DrPvosHOC+usNqI2ZrNqbpogIJqexOQDSKXsT73XyOdJ3/1tsvf/wA4kknR+4AYSbzl/3snveu6SNWBVKoXTu2HN1CZv1bqqQvsQRoapfEg6H2GAep29xhj+9ZVJ7npmhNDYDpxXnbaRC47tqfmC7SrY0CFdPAWhVZHqK2uLVhHmlfDPwLf1/9/enYfZVdWJ3v+uvfeZaj6VSgIZGGQmMgbQblRAQEAZBHEZEGXQF/t227a329te2/v07dd+u1/b977etm97+9UWBEWEZRBBRBERBdoBCAYJQ5iRhJCkkkpSc51z9nr/WPtUnarUfKZdp36f56mnqnbt4bdOnb3Pb6+9BmOGtdafAG4G3j3DNuNMHGe1GmNjVnvMTZvPuVk6KzQG+ZTH2bYFu/5G2LHNLTjoLagPXMNg51IG9+wZt+60Y8C2tqNCBQ00Rmycx1Uth4xjXn3WWvpGQvoqNPRhqd7hAt/auIOnd7jGpcubE1y3dhmr21Pz3mc6UHSkZShEIeKq2J5ZasobR1wS85kmOcAYs6vk128AXyrZ9swJ2/6i4hHWmc3n3ZCHgwPVPU4Ywq8ewP7sLigUwPNQZ70P3nne3Gq8PQXtS1Cp+ScFQjSK0FoGciF9w5UZ+nCizd2D3PTEDvYOu1kh/2h1C/qtXaSC+dVyK6A15dMqTVeEiC1JyhtTXBLzaSc5ANBaH2iMiapvuRh4Nvr5PuAftdbFthTvAT5X/ZBrwxYKYwl5FWvIAeze3dg7boZXnncLlixDXX4tatUhc9uR70N2CSqQUVfE4pYPLf0jbrbOKvTJphBa7tncw09f3IPF1XBfcfxSTl05/xl03agrvsziKUSM+b4iLeOUN6RYJOaznOTgU1rri4E8sBu4Jtp2t9b673HJPcAXih1BFzKbz8NALwxWt8nK6PF+/xj2h991TWQATnkH6oLLUck51ngnk9DRKaOuiEVtOO+aqwzlq3fu7hrIceMTO3ilxw1denBHiutOXsbS5vnfELemfDLNgXzYCxFjiaSb0EbO08YUi8QcZp7kwBjzOaaoCTfG3AjcWNUAa2S0DXmVm6yMHm9oEPvD2+D3j7oFTS2o91+FOuaEue+sqQlaO+RiIRalQmijMcgL5Cs7hcB+nnijj+882c1gdKBzDmvn4qM7CebZFlwBHRmfbFOS7gE5f4WII6UglfYIAjlHG1lsEvPFzuZy0F/dYQ/3O+arL2DX3wR7owcMR77VJeWt7XPbkVJuKMSm+T8+F2KhyoeWvmHXXKXaz7ZG8iHrn97FI3/oBaA16XP1SUs5dtn8hyH1FGQzAel5tkcXQlRfEChSKYWSjtgNTxLzOrPDQ66GvDiTZi2Omc9jH7wHHv6payaTSKDO/wCc+q6513Z7Hv6SpajevuoEK0RM5QqWvpECgzVIyAG27hvhxg3b2dbnxiY/ZmmGq09aSltq/pfxwIPOTCDtyYWIKQUk0x6JhCTki4Uk5nVgC4Vols4BN/JJLY+98003DOIb0dDvK1ajLr8OtfSAue+s2MkzlQZJzEVEa30+8BVcf5FvGGO+OMV6HwDWA6caYx6vYYhlGcqH9Fe5/Xgpay0Pv9bLHU/vIhdaPAWXHN3J2Ye145XRbEyGQhQi3jzPdfD0fDlHFxNJzGvIjgzDQF9Nm6uMHttaePQh7H13QC7nmp+88zzUWe9DBfN4GySiTp4NMmmQqAyttQ98FTgXN9nXY1rru40xz0xYrxX4C+C3tY9y7mxxuMOR6rcfL9U/UuA7T+5k45uuz0lXU8B1Jy/jkGx63vuUoRCFiD8/mjFS+mwtPpKY14AdGoD+fsiN1Of4ffuwd34bnt/kFnR0oj5wDeqQI+a3w1Qa2rMoTx5/i/2cBrxojHkZQGt9G3AJ8MyE9f4e+Cfgv9Q2vLkpdujsr8KEQDN5cdcQ33xiOz1D7qnaKSubueK4pWQS8z/v/Kg9+XzHNxdCVF9xJk9JyhcnScyrxIYFN9ThQF/Nm6uMi+PZJ7F33eLasQOc8DbUhR9Cpec5CXhTE6otO/N6YrFaCbxe8vsW4G2lK2itTwZWG2N+pLWOZWJeCF378f6R2rQfLxVay30v7OGezT1YIOkr9FuX8EerW8v6oE76is6MNF0RIs5kfHIhiXmF2eEh7J7dMDxUk/HHp4xjZBj74/Xw+CNuQTqDuugK1PGnzn+nre2oZhl5Rcyf1toDvkw0D8EM614PXA9gjKGrq6vi8QRBMG6/obXsG8ozNJwnkYCO5oofcvI4fJ9sNsvugRH+98Ov8Mx2dyN9cDbDn7/rLaxsn3/TFYDWVEBHZubxySe+HvUUl1jiEgfEKxZRecmkIpmSp1mLnSTmFWBzOTcxz9AAheG2sUl66hXPllex678Ju3a4BYceifrA1aj2zvntUClo70Cl5z8km1g0tgKrS35fFS0ragXeCvxCaw1wAHC31vriiR1AjTFfB74e/Wq7u7srHmxXVxfd3d2E1tI/4tqQ17rJCkA2m+WXz27h2xt30p9zjdjPPLSNS4/pJBEO0tMzv2tKcXzyQsFn1yymRii+HnEQl1jiEgeMj2XFihV1jkZUioxPLkpJYj5P1loY6Hejq+Tz9Q4HiEZ7efg+7IM/gjAE30edcwn88dnzbw/uea5N+lxnABWL1WPAEVrrQ3EJ+TrgyuIfjTF7gdEqP631L4DP1GtUFmvdGOS9dUrIAXKFkJsffZ2fPOdupJsTHh85cSnHH1Bedb2nYElTQFKGQhQitortyT1pYiYikpjP0WhC3t/rkt+YsLt3Yu+4Cf7wsluw7EA3DOKBq+a/0yCAjiXzG7VFLErGmLzW+pPAfbjhEm80xjyttf4C8Lgx5u76RugUR1kZ3jfM3uH69QHZ3jfCDRt2sGWf6xh+eGeaa09eRjZT3jmX8BRLmqQ9uRBxJiOviMlIxjVLNgzduONxS8ithd/9Bvuj28cmKfqjs1DnXopKJOa/42QKOrIoT4ZUE3NjjLkXuHfCsr+dYt0zaxFTqYFcgd5hN+xhNlOfanJrLb95vY/bN3UzUrAoBe87Msv5R3SUNTY5uPHJs5mg7P0IIapHknIxFUnMZ2DDgqshH+iPVUIOYAf6sHfdCs/8zi1obUdd9lHU4ceWt+NME7R1yAVDNJTBXEjvcIFcvdqslMTx3ad28vjWfgCyaZ9PnXEYy5PlN4lrSXq0p+WyLkScyXCIYjpyBZ+CzefcEINDg3UdXWUq9qVnsXfcDL173YJjT0Rd/OHyR01pbkW1tpUfoBAxYK1lMB+O1pDX26t7hrhxww66B1wSfsIBTVx1wlJWLW+lp6dn3vtVQHvapzkpT7jE3Mw0U6/W+i+BjwN5YCdwnTHmtZoH2iA8DzKSlItpSK+gCWxuBLtnF3TvcE1XYpaU21yO8MfrsTf9i0vKkynUpR9Brbu+/KS8tU2SctEw+kcKbO/L0TNY/6Q8tJb7X9zD/3jkDboH8gSeYt1xXVx/yvKyk+liJ09JysVclczUewFwLHCF1nriI9ffAacYY44H1gNfqm2UjUMpSGc8lPT9ENOQGnOi9uPRcIeM1Gd2ztmwb27Bfu+bsOMNt2D1oajLr0V1Li1/5zJGuWgQQ/mQfUP1b7JStG84z7d+t5NndrohDw9oSfCxtctZ2ZYse9+BB0uaEgTyQS/mZ8aZeo0xD5as/xvgqppG2CCUgkzGk9FXxIwWdWJuh4fdcId1ngxoJjYM4dc/x95/FxTy4HmoMy6AMy5A+RWoJWtrRzVJUi4WtlwhZN9wgaF8fM7lZ3YM8K2NO9kXjfxy+kGtfHDNEpJB+Q8rpZOnqIAZZ+qd4GPAj6saUQOSIRHFXCzqxJyeeEwaMR27twf7/Zvh5c1uQedSV0u++tDyd+550N6JSskY5WLhKoSW3uHC6MQ8cVAILT/cvJufvuj6gGQCjytP6GLtisrcALcmPdqkk6eoIa31VcApwBlT/H3cTL3ZbLaG0ZXH94OqxZtIeGSa/Jq2KW+kGWIbqSyzJVf2GLNPbcDe/Z2xmUTXno664HJUqrzpuQEZo1wseIXQ0jdSoH8kJD515NDdn+PGJ3bw6h43fOmh2RTXnryMrqYyhi+NFGfybEpIe3JRETPN1AuA1voc4PPAGcaY4cl2NHGm3nI6M9daNpstq/P1VFIpRSLpMVDjycDjNFttuRqlLHOZqTc2WVk5PcO11gXgqWjVPxhjLq5Z4FVghwbduOQbf+sWNDWjLrkKdeyJlTlAMuVm85zvbKBC1FFoLX0jIf11nK1zKo9v7ePW3+9kKG9RwHsO7+DCo7IVmejHU9CZCUhVoBmMEJFpZ+oF0FqfBHwNON8Ys6P2IS48ngeptIfvS9MVMXexSMxLeoafi2vj9pjW+m5jzDMlqxV7hg9orf8Trmf4h6K/DRpjKpS11pd97UXs+ptgzy634IhjUZd+FNXaXpkDpNOu+Yq0SxULTGgt/SMhfTFMyIfzIWZTN79+vQ+AtpTPNSct4+ilmYrsP+m79uTSyVNU0ixn6v1/gBbge1praIDKr2qSiYNEuWKRmCM9w7GFAvbn98DD97mOqEECdf5lcNoZlTvB0xloz8oFQyxI2/tysUvIAV7fO8yNG3awvT8HwJplGT564jJaU5VpbtKc8GhP17aNqlg8Zpqp1xhzTs2DWqCCQJFKS1IuyhOXxLzcnuFprfXjuGYuXzTG/GCyjWrROcUP/Dnvt7D9DXpv+lcKr7/s9rHqEFqu/iTBgasqFpdqbsHv6KzY/kotts4Zi628cRG3pNxayy9e3cedz+wiH4Kv4JJjOnn3W9orNlKKdPIUYmFIJBWplDQzE+VbcFf8KXqGH2yM2aq1fgvwc631U8aYlyZuW4vOKXPpRGKthccexv5kPeRybqDT088lPPsieoMAKhVfUzMqFUKVOlA0SueM2WrU8s6lc8pi1zdS4JaNO/n99gEAljYFXLd2OQd3VGaEIwW0pX1aZNIgIWIvnfYIElJLLiojLol5WT3DjTFbo+8va61/AZwE7JeYx4nt24f9wS2wOeqz2p5FfeAa1KFHVvZAza0ym6cQFfR89yA3/W4He4bc2OSnrWph3XFdpCvUKVM6eQqxcKRSkpSLyopLYj7vnuFa6ywwYIwZ1lp3AacT8ymD7eansHd+G/p73YLjT0VduA6VaarsgVpaUS2SlAtRCYXQ8uMXevjx83uwQMpXrDuui7etbq3YMQIPOjMJEjKagxCxl0gqEkk5V0VlxSIxL7Nn+DHA17TWIeDh2pg/M+mB6syODGN/cgc89rBbkM6gLroCdfyplT9YazuqWWbzFKISdg/muemJHby4ewiA1e1JPnbycpa1lD82eVHCUyxpCioytKIQoroSCWlTLqojFok5zL9nuDHmV8Bx1Y2ufHbra9jv3Qi7osr+Q45wTVeq0SGzrR3VJEm5EJWwcVs/tzy5k4FoZtF3v6WdS47urGitdirwSDcHFes0KoSoHunoKaopNol5o7JhCA/dh33wHghD8H3U2RfD6edUZ4IfScpFnZUzWVic5Aohdzyzm4de3QdAS9LjIycu47jllW1ylkl4LGtJsmtYknIh4i6V9khIm3JRRXLLV0W2pxt7w5exD9ztkvKlB6Cu/yzqne+pUlLeIUm5qKuSycIuAI4FrtBaHzthteJkYccD64lhn5BtvSN86ZE3RpPyI5ek+ZszVlU8KW9OenRmAhn3WIiYUwoyTZKUi+qTGvMqsNbCk49i77kNhl2bVN52Juq8S1GJZHUO2p6tfOdRIeZuQU8WZq3lV3/oxWzaRS60eAouPCrLew7vqHgzk7aUX7FJiIQQ1eP7inRaoaT/h6gBScwrLBzow5obYNMGt6ClDXXZR1FHrKnOAZWC9g5UWpJyEQvlThZWNwO5Arf+vpsn3ugH3JCF1568jMM60xU9jgI6Mj5NCUnKhYi7REKRSkvjAlE7kphXkH3pOfb+4NuwZ7dbcMwJqEuuqt7oKEq5mvJ0pjr7F6KKppgsrPTvVZ+pN/DdTL3P7+zjXx/Zws6+EQBOO6iD6//4YJqTlb1Eegq6mpOkJyTlcZlNNi5xQHxiiUscEK9YGp3CtSeXMcpFrUliXgE2n8Pefxf86gEsQDKFuuByWHt69dqOKgUdnahUZWvzhChTWZOFlarFTL3tHR3c/tgr3LO5h9C6IQs/+NYlnH5QKyP9vYz0V+5YnoKupoC+vEffhL/FZTbZuMQB8YklLnHA+Fhkpt7q8Tw3m6cn8wmIOpDEvEz2za3Y9d+E7S73CA45nML7P4Jasqx6B/U8yC6pXnt1IeZv3pOF1dreoTz/+/4X2PSmm+jrwNYEH1u7nBWtlT+vAg+WNCUIpI2qELEWBIpUWkmHbFE3kpjPkw1D+M2D2Pt/APm8S5bPuIC291/Jnn37qnfgRALaO1GB/OtE/JQ5WVjNPL19gJs37qBvxI1N/o6DW7l8zRKSfuXbkiZ9RWdGJg4SIu6SSUVSxicXdSbZ3TzYfXuw378ZXnrOLch2oS6/FnXQW1B+FTt0pdKuTXk1hloUokLmO1lYLeRDy13P7uaBl/cC0Jz0ufL4Lk46sLkqx0sHimxGJg4SIs6UUqQzHkEg56moP0nM58huegJ7960wGDU+XfvHqAs+WP223k3NqLaO6h5DiAa2oz/HjRu284e9roPnYdkUnz7rSPxcBRuSl2hOeHRk5BIrRJx5ClpaAkZykpSLeJBPjVmyQ4PYe78Hv/u1W9DUjLr4w6g1J1X/4DKbpxBl+e2WXm77fTfDBYsCzj+ig/cemaWrJUlPT+UT89aUT5uMUS5ErCkF6Sbp5CniRRLzWbCvvYS945vQs8stOPwY1KUfrX4Ntoy8IkRZhvIhtz/VzW+3uHFQ2tM+15y0jKO6qjfEaHvapyUpSbkQcVacydOTvh8iZiQxn4YtFLC/uBd++WOwFoIAdd5lcNoZ1W/nLSOvCFGW1/cOc8OGHezozwFw3PImPnLCUlqqVJMtEwcJsTAoIJ2RpFzEkyTmU7Dd290wiFtfcwsOWOU6eC6vwdixvu+S8iBR/WMJ0WCstTz4yj5+8Owu8qEbqvDSY5Zw5qFtVRsCTQGdTQHpQDpmCxF3qYyHL81XRExJYj6BtRYefwT74/WQG3HPu04/B3X2RbVJlIPAjfJSzdFdhGhQvcMFvr1xJ5t2DACwrDnBx9YuY3V7qmrH9BR0ZgJSkpQLEXvptIy+IuJNEvMStm8f9ge3wOan3IK2LOoDV6PeclRtAkgkoGOJJOVCzMPm7kFuemIHe4cLALx9dQv6rV1VrcUuzuaZqML450LUgtb6fOAruHkHvmGM+eKEv6eAbwFrgV3Ah4wxr9Y6zkpIpRRBQpJyEW+SmEfs85uw3/8W9LtZADnuFNRF61CZ6oxvvJ9kynX0lDHKhZiTQmi5Z3MPP31xDxY3dvi647o4bVVrVY/rKzebZ0IeiYsFSmvtA18FzgW2AI9pre82xjxTstrHgB5jzOFa63XAPwEfqn205UkkFYnkwvp83b5thJeeG2Ggv0BTs89hRydZfuD0/c62bxvhsUe2srt7iDAE5Vna2oP9tt389AAvb86Rz9nRZUq5Lz8Bng/5YSgUxvbtedDU6nHsCW5AimJsobXkh3HHU5DKuGtiU7PPkmUeu3aE9O4tjMaTTntYYHgoxIYKz4PW9unLt33bCM88OcRAr5sUrrnV45gT0qPrz+a1mmyd0nLM9jWutkWfmNuREex9d8CjD7kFqTTqoitQJ5xWuyDSGTdxkExCIsSc7BrI8c0ndvByzzAAB7WnuG7tMpY1V7fZWeC5pDyQzmNiYTsNeNEY8zKA1vo24BKgNDG/BPi76Of1wL9qrZUxxrJAJBKK1AKb0XP7thE2bRhCeS7+ocGQTRuGYC3TJq+bNgwBitxI9O9R0Lcvz6YN4ei2m58e4IWnR7AT/oPWuq9wePKYwhD694Vs+NUAfqBIJBT5nCU3Mn4fg/2WZModd/dOCBKQz0V/B3LD4ej6SlkKQF9vYcrybd82wsbfDrjjKLeTvn0hGx8d5MQoVZvptZrs9dz46CBYSyLpzfo1roXYJOblPE7TWn8Od1dfAD5ljLlvNse0W19zHTy7t7sFhxyB+sDV2B1vYm/8n9DTDdkueMd78I5cM7pd+OCP4D8egJEhSKbh9LPxznrfrMsaPv80PPJTt/8ly1Hv03gdnVOv/9QG7H3fd3GmM+6WdHAAupajzrsM77i1U25b+OFtcP9dMDwIqQycewn+RevcPu+4Cba/AWHB3R5bCzYc+1kpOGAl6rKr8Y5bO34bgPYsANv37XHrl6w7aewzxDtpOff0QCHvmvkcuBqOOs41NZrF/mZrrjHu/vkPKWzbUrHji/l54o0+vvNkN4N5d5E/57B2Lj66s+rJctJXdGYCfEnKxcK3Eni95PctwNumWscYk9da7wWWAN01ibBMQaBIpRdWUg6uFld5jLaHDwLI5y0vPTcyZdJY3GZ4KByt/bYW8nmXHBe3fXlzbt5xKeU+kq21ZDIe/X2Trzcy4saxAMjlXLM/pcbXwIOrhbcW8jlLkFGTlu+l50bI58dq9FFj27z0nLsrmOm1muz1HBxwnx2Zptm/xrUQi8S8nMdpWutjgXXAGmAF8DOt9ZHGmAn//v3Zr3/J3QJ6HurdF8E734N98Vm45zb3jso0Q+9euOc2wgvX4R25xiXlD97r7tr8AHLD8OC9hDCr5Dx8/umx/Te3QX8v9tavEV75iUkTvPCpDdhbv+beMcqDN7e4P3Quhb09025b+OFtcM/tY7GODME9t1PY/gY8sxEGeqPb49B9jR40+tnzYNvr2Jv+hcKZF7hyD/S6OAoF2LXDref7gBpdN7zmU2OJfDH25tZp4520nO65l4t/ZAi2vAIvPAPtHdDaMWP5Z2M+MYap1IzriuoZyYesf3oXj/zBNTtrTfp89KSlrFnWVPVjpwKXlHvydEuIcbTW1wPXAxhjyGazdY4IfF/R3BLM+DQ6CAK6urpqFNXsDA/2kUqNj933LcOD4ZSxFrcZ7M+hRisOLDaEVCoY3baQ34PywM6YJU0i2q0N3es25U6ier7icZQ/Fs/4/Sl3AxHacTEWBUHA8GC0mSLKzN23sGDd32DG12qy19OGBZSKyjHFdvUQi8ScMh6nRctvM8YMA69orV+M9vfrGY8ahtB1AOqD16JWHASAfeSnLtFMRqM4JFMwMuxquI9c42rKFWO3gr7vktT/eABmU2te3H+6yX0PAmDI1dhOktzZ+77v1kmloWeruMiviwAAG59JREFUe5cD7NsLB6ycdlvuv2ssKSf6XsjDYw+522flQZifOlZrXe354IDbVyHvtvE89/Po62ghEbjXc3BgNJ5xsUP0ffJ4Jy9n6Grzkym374EBd9s9MABt2Wn3N1vziVGlM64KogLHF3Pzxr4RbnhiO9t6XY3P0V0Zrj5pKe3p6l/K0lFSLk3ORAPZCqwu+X1VtGyydbZorQOgHffUehxjzNeBr0e/2p6enspHOweep8g0KYaGZz5fu7q66O6O1wOAVAaGBvPjRpDJ5y3pjDdlrMVtPN8lrSpq9qE8GB7Oj27rB4pCfp4tkYotZDzI5/OjTUv2o6IUo+B+tqFl0kunte4B/YQYi7q6ukhl3Me+O44tbobyXJlh5tdqstezmE7l8/kpt6uUFStmP9R2XBLzch6nrQR+M2HblZMdZOIdfepd76H5/R9GJceGUtu9twfV1DL+rirtY/f2kM1m2TUyPFZDXOT5MDJMNpvFD/xpawp27+1BtbbjldyhWb8Z29M96R3azp5uVIsbfzlfyLukGKCQJwiCabfdPjzkYlUTYs2NjDVfme7ktERJeAGG8y4p9vyx52Oj60XNXqJ1vSie0thnKuuk5SzGVtx3PgdecrTsM712szG/GKnY8cXsPfTqPu54ehe50OIpuPjoTs45rL0mtddNCY9sJi6XSyEq5jHgCK31obgEfB1w5YR17gauxlV2XQ78PO7tyz0F6Yxa0DfRhx2dZNOGIfJ5O1r/Z0NGOyxOt00y6TE4UHAf08rVm5Vu+5ajErzw9MiU+5mOta6Ozw8U+bwlkWBcG/OiZNJ9dI8UXEvUfI7ReEoT+eID+kRCTVm+w45OsvG3eXIj43JzEkk1uv5Mr9Vkr2eQcLnMXF7jWlhUnzQT7+hz517Knv4B6B8YXce2Z7G9e8dqzMHVmLdn6enpcctzw2M15jBaq9vT45L36WoKbOdSbN8+QtJjC4eHIDv5HXuY7YK9Pa52tljjDeAH7i5vmm1JpV0TkNKLU1hwSa7nj3WhntgDpEgRNfWJavYL+dGmP+O2K/4crRtG8YyLfYayTlrO4n6L+y7efgfJsTvc6co/C/OJUTW3VOz4cTKXO/p6uO0p9xovaQq47uRlHJpNz7BFZbQkvZrUyAtRa1El1yeB+3D9u240xjyttf4C8Lgx5m7gBuDb0dPo3bjkPbaUaoxZPZcfmIS1cxsxpLjNH16y5PNjo7K0tI0fleWoNa7ZXyVGZQlDixeUPypLS+vU5Vt+YJIT38bYqCxq/1FZZnqtJns915wko7JMp5zHabPZdvbe8R7XBnxkGBJJdytYKLjlAKef7dpaF6Ia57Dg7t5OP3vmfbd1oN77QdemmaGxZjL5POq8yybdRJ132dj6re3QszPaV7tLCqfZlnMvcW3MC/nxsZ76rrE25tMm5srdPja1QGkb8xDG3fZ6xR4dFppaRuMZF/sMZZ2ynJ4/djPS1ASDg+67tTO+drMxnxjt0GDUZr/844u5WbuimSuPX0omUZvOXK0pn7aUzCsgGpcx5l7g3gnL/rbk5yHgg7WOaz4UUVLeIEOYLj9w7kni8gOTrDlu5sqio9Y0cdSaaVeZ1bFmpczjFI81003JTPFMtU69E/GJ4tJVefRxmtY6ibsjv3vCOsXHaTD+cdrdwDqtdSp6HHcE8Oh8A/GOXAMXrnPJ4WC/+x51/ISog+dZ74VEKhotJAVnvXfmjp/tWVRTM95xa1FXfsKNaNLf65ZP03lw3Po2hANWudFJwnDGbf2L1sGFH3IjxxTy7vuFH8L/+F+irv0Lty/Pd7fBQcIlm6U/+z4cuBp1zafwL1o3tg1A4MOSZe7L811yHq1bjGcuZZ20nCsPdjcFybTrbLnqUFeeZStm9drNxnxi9LJLKnZ8MXtXndDFdScvq1lS3p6WpFyIhSSV8fAbJCkXi5eyU9WW1pjW+r3APzP2OO0fSh+naa3TwLeBk4gep5V0Fv08cB2QBz5tjPnxLA5ptz7xWMXLMWlTlvYsKlP9ESPqJY4dZ6qpUcsbNWWJ66eafey5Vyu+08nOVwV0ZHyaErVLyuPynopLHBCfWOISB4yPJe7n6/ObX595rQpKp715z+oZp/9xJTRSeRqlLHM5X2OTmNdB9RNzpaC9A5Vu3KQcGufEma1GLW+lPuirNMX3tIn509v7uf+lvXQP5OlqCjj3sHaAKZdt6x0hH1oSvscBLQnOPaydNcubUUBnU0A6mLxWfsPWPv7t0W10DxRGZxl9+8pmdg2FbO/LsbwlwaXHdLJ2ZctML9Po/u58djfdgwW6Mv7otsXl2/tyNCU8+kfy7BlyvaRWtiX56IlLWbuyhdt/v5O7nuthMB+SCTxOXdE0LpbjlmV4ascgr+8dJlewBB4c1JGeMsaJ7+3SOEr3N11Z57PNZK/LPS/2sqVnoKrHmY3pzveJMVTqmLOJRRLzMcmUIlnGrJ6Ndk1vpPI0SlkkMZ+d6ibmSrma8nSm4seIm0Y5cWarUctbiQ/6aE6C5ymZkwC4onROAq31nwLHG2P+JJqT4FJjzExTfE+ZmD+9vZ/bNu3CV4qk70YBGMwVsNbSlAxGlw2M5FFK4SvYF808p5Tr3Bn4Hlcc18W7DmkjNU1S/j8e2cLAJCOMtqU8upoSDBcs+dDyiVOWzyr5/Nrj2wk8RXMqQf9wjnxoOfvQNh54ZR+BpygUQrb357GAr8Z6drSmfE5cnuHhP7jZPTwFhdB1/2hJKJa1JNk7lGfPUIFMoBiMRjhSCtpSPgnfmzTG0vd2aXwpX7FnKM/eoQIdaZ/2dDBpWeezzVSvSzoR4BNW7TizNdX5PjGGSh5zNrFIYu4kkuXP6tlo1/RGKk+jlGUu52tc2pg3FqWgo3NRJOVCTDA6J4ExZgQozklQ6hLg5ujn9cDZ0ZwE83L/S3vxlSIVuOHRUoFiKO8SutJlwwXLUD5kMO+GW/SjwYUG85bAUzz4yt4pk3KAO5/dzWCUlEeTz43qGwlRSpEOPAJPceezu2eM+85ndxN4bpvSbe96rmd0+d6SqatDC76n8JRiMBfyUJSUB9GyYhXLQN6ilGIgGm1hIG/xlCLwFAoYyNlZxTgxvtH95eyUZZ3PNlMeN+FX9Tjlmur/V81jivESifKTciHiJi6jsjQOz4OOJahUauZ1hWg8FZvie7YzCfYMvU5zcvzEP2E4AgoCf+wSF1o3CK5VdnRWPGUtYWhpSyfYPVSYdjz67sFXJp1Hw+17bGz7Zt/SPZifcWz77sFXaE0lUNHMd0EQ0Oxb3uwdoTlang+Hxx9TqWgQpJDQQtIvzlENxVGSirHkw2F8TzFSsCSiDnEKKISW5lRi0hhLZ0AsjQ8Y3V8+tFOWdT7bTP26TP2aVuI4szXVrJATY6jkMecay2Lk+4pUWpJy0XgkMa8kT+F3LkX19tY7EiEWvNnOJJhN++wZypMqmdGtOA9XvmSGWi+q5vaUohDNQmdxtdADI3m6MsG0j0y7Mj5v7pt8kjtPjc0eN5QPZ9xXcX+7B3OupjVw8xIM5UPSgUf/cG60BjYf2rFjWjtacx6GlkLo2o1PFkvgKUYKIZ7CjTWswFqL77n9TxZj6WPj0vjA1cznCiFJ35uyrPPZZqrXpSWdrOpxZmuqR+kTY6jkMWcTS9znHagm31ekM3FtxSNEeeR2s1I8D7JLpaZcLHZzmZOA6ab4nq1zD2unYC3DeYuNvqcDz7X7LVmW8l1zg0ygCK1rk20tNCc98qHl0mM6pz3Opcd0Upz80zI+QW9JeljrmsrMZl/F/eVDt03ptpccnR1d3l7ymN61I7eE1pJJeLzrINeOOR8tK6YpTYFy7eujESqaAkVo7WiC35RQsy5vaXyj+0uoKcs6n22mPG7UT6BaxynXVP+/ah5TgOcp0umFPaunENORGvNK8H3ILkEFiXpHIkS91XyK7zXLm1nH+BFYPrjGJUdTLQsZIV+wpBI+K9uSsxpNY+3KFj7zjlUVG5Vl7coWPgHRqCyuxr647eFLMqOjfRzckRo3KsuqklFZVpSOypIYPyrLirYU5x2+/6gsK9tSsy5vMb7tfTlWtqU4//DpRz6ZzzZTHXe6UVkqcZxyTYyhFsdc7Nysnmq0KZoQjUhGZSmX70O2CxW1a2yUHsRzsdjK3KjlreBwifOek2AaFRvHPJPwyKZdx8K4/C8ljv3FJZa4xAGLe1QWpSBTpVk94/Q/roRGKk+jlGUu56vUmJfD86KacnkZhSiK8xTfzQmPjoycr0IsJApIVykpFyJu5BNqvqT5ihALSmvSoy0tlzwhFppUxsOXpFwsEvIpNR8Tmq8IIeKtLeXTmvLrHYYQYo5SaY8gkKRcLB6SWc6VJOVCLCjtaZ+WpCTlQiw0yZQikZCkXCwukl3Ohe9D51KULx/yQiwEHWmfZknKhVhwEklFMikjOovFRxLz2QoCV1MuSbkQC0I249OUkPNViIUmkVCkUpKUi8VJEvPZCALo7EJ58iEvRNwpoEOSciEWpFTKI5GU5iti8ZLEfCaJhBt9RZJyIWJPAdlMQCYhtW1CLCSektFXhABJzKeXSEK2U5JyIRYABXQ2BaQDScqFWEiCQJFKK5SSpFwIScynkkxBRyfKkw95IeJOknIhFqZkSjp5ClFKEvPJpDPQnpW7dyEWAE9BZyYgJUm5EAuGUpBOe/gyRrkQ49Q9MddadwK3A4cArwLaGNMzYZ0TgX8D2oAC8A/GmNujv90EnAHsjVa/xhizcd4BtbSiWtrmvbkQonY8BUuaApK+JOVCzEW5n73l8H1FOq1QniTlQkwUh0+z/wo8YIw5Angg+n2iAeCjxpg1wPnAP2utO0r+/l+MMSdGX/NLypVyteSSlAuxIPgKuiQpF2K+KvHZO2eJpCLT5ElSLsQU4vCJdglwc/TzzcD7J65gjHneGPNC9PMbwA5gacUiUMq1J880VWyXQojqCTzoak6QkKRciPmq6WevAtIZT8YnF2IGdW/KAiw3xmyLfn4TWD7dylrr04Ak8FLJ4n/QWv8t0V2/MWZ4im2vB64HMMaQzWZBKfwlS1GpdLnlACAIArq6uiqyr4VisZV5sZU3bpK+ojMT4EuNmxDlqMRn76x4niKdUXhyzgoxo5ok5lrrnwEHTPKnz5f+YoyxWms7zX4OBL4NXG2MCaPFn8NdVJLA14HPAl+YbHtjzNejdQBsz5490LEE1dsHvX1zKdKUurq66O7ursi+ForFVuZGLe+KFSvqHcKM0oFLyqVjthAzq/Jn78R19q/4ApJJj3TGj/U522iVLY1UnkYqy2zVJDE3xpwz1d+01tu11gcaY7ZFJ/+OKdZrA34EfN4Y85uSfRfv+Ie11t8EPjPrwDqWoFKpWa8uhKifTMIjm473B7wQcVLNz95JjjW+4qunh1TaI5dX9A+UUYgaaLTKlkYqT6OUZS4VX3FoynI3cDXwxej7XRNX0FongTuBbxlj1k/4W/HConBt5DbN9sCSlAtROdUe5aEzE4fLlRANo6zP3plkmmQWTyHmIw69ML4InKu1fgE4J/odrfUpWutvROto4F3ANVrrjdHXidHfvqO1fgp4CugC/q/ahi+EiNRllAchxLyU+9k7LUnKhZgfZe2UzcoanX3jjTcqvtNGeewyF4utzI1a3uhR27w/TbXWm4EzSx6N/8IYc9QM2zwJXF4c+WEaDX2+Shz7i0sscYkDxsdS7vlaZVU5X6slTv/jSmik8jRKWeZyvsqzYSFEpVR0lIeJncmq0QEoLh2LJI79xSWWuMQB8YpFCFEdkpgLIWatlqM8TOxMVo1ak7jUxkgc+4tLLHGJAyatMRdCNBhJzIUQs1bLUR6EEEKIxSYOnT+FEI2hOMoDVGGUByGEEKLRSY25EKJSvggYrfXHgNdwIzqgtT4F+BNjzMcZG+Vhidb6mmi7a4wxG+sQrxBCCBErkpgLISrCGLMLOHuS5Y8DH49+vgW4pcahCSGEEAuCNGURQgghhBAiBiQxF0IIIYQQIgYW9QRD9Q5AiBiK7YQl9Q5AiBiS81WIBUJr/YloGOBpLeYac1WNL631hmrtO65fi63MDV7euGro/6XEEd9Y4hLHFLHEVd1fq7l8xel/LOVp3LIQTZg3k8WcmAshhBBCCBEbkpgLIYQQQggRA5KYV96M7Yca0GIr82IrbyOLy/9S4thfXGKJSxwQr1gaSaO9ro1UnkVXlsXc+VMIIYQQQojYkBpzIYQQQgghYkBm/qwgrfX5wFcAH/iGMeaLdQ6parTWq4FvActxQ2N93RjzlfpGVX1aax94HNhqjLmw3vGI8mitPwj8HXAMcFo0S2mtY6j7dUNrfSNwIbDDGPPWWh+/JI7YXFe01mngISCF+6xcb4z57/WIJYpHrj0VFqf3W7ni9n6thEZ5z2utXwV6gQKQN8acMt36UmNeIdEb6KvABcCxwBVa62PrG1VV5YG/MsYcC7wd+LMGL2/RXwDP1jsIUTGbgMtwH2g1F6Prxk3A+XU47kRxuq4MA+82xpwAnAicr7V+e51iAbn2VEOc3m/litv7tRIa6T1/ljHmxJmScpDEvJJOA140xrxsjBkBbgMuqXNMVWOM2WaMeSL6uRd38qysb1TVpbVeBbwP+Ea9YxGVYYx51hizuY4hxOK6YYx5CNhd6+NOEkdsrivGGGuM6Yt+TURfdemUJdee6ojT+61ccXq/VsJifs9LU5bKWQm8XvL7FuBtdYqlprTWhwAnAb+tcyjV9s/AXwOt9Q5ENIxFe92YSRyuK9ETjQ3A4cBXjTH1ikWuPVUWh/dbuWL0fq2ERnrPW+CnWmsLfG2m2T+lxlyURWvdAtwBfNoYs6/e8VSL1rrY/nZDvWMRc6O1/pnWetMkXw37RGuhi8t1xRhTMMacCKwCTtNa17z9vVx7qi8u77dyxeH9WgkN+J5/hzHmZFyTxT/TWr9rupWlxrxytgKrS35fFS1rWFrrBO5i9h1jzPfrHU+VnQ5crLV+L5AG2rTWtxhjrqpzXGIGxphz6h3DNBbddWMmcbyuGGP2aK0fxLXD31Tjw8u1p4ri+H4rV53fr5XQUO95Y8zW6PsOrfWduCaMU/ZrksS8ch4DjtBaH4r7YF0HXFnfkKpHa62AG4BnjTFfrnc81WaM+RzwOQCt9ZnAZxbqRULEyqK6bswkTtcVrfVSIBclORngXOCfah2HXHuqJ07vt3LF5f1aCY30ntdaNwOeMaY3+vk9wBem20YS8woxxuS11p8E7sMNe3ajMebpOodVTacDHwGe0lpvjJb9jTHm3jrGJMScaK0vBf4XsBT4kdZ6ozHmvFodPy7XDa31d4EzgS6t9Rbgvxtjbqh1HMTrunIgcHPUbtcDjDHmnjrEIaonTu+3csn7NZ6WA3dqrcHl3LcaY34y3QYy86cQQgghhBAxIJ0/hRBCCCGEiAFJzIUQQgghhIgBScyFEEIIIYSIAUnMhRBCCCGEiAFJzIUQQgghhIgBScxFxWmtX9Vax3lSFyEWpXqdm1rrQ7TWVmstQ/QKMYVKnp9a63dqrTfPct0zo2FSRQzIRVIIIYQQooEYYx4GjqrEvrTWNwFbjDH/rRL7E9OTGnMhhBBCCCFiQGrMRdVorf8OeCtQAN4LvABca4x5Mvr7Z4FPAW3AG8CfGmMe0FqfBnwFOAYYBO4A/tIYM1LzQgjRgKLZAT8LfAxYBjwPvN8Y87rW+ivAZUA77pz9tDHmYa31CuAlYKUxZne0n5OA+3GzDhaAvwH+DyAD/AT4c2PM3poWToiFLa21HgRWG2O6tdafB/5PoNMYs09r/fdAqzHm01rrFPAPgAZSwJ3AfzbGDEZT2d9ijFkFoLU+GbgBOBx3bobAC6W14Frrv8JdFwq4GVC/qbW+HvgwYLXWnwYeNMZcVIsXYrGSGnNRbZcA3wM6gVuBH2itE1rro4BPAqcaY1qB84BXo20KwH8GuoA/As4G/rTGcQvRyP4SuAJ3w9wGXAcMRH97DDiRsXP2e1rrtDHmDeDXwAdK9nMlsN4YkwOuib7OAt4CtAD/Wu2CCNFghnDn4BnR72cArwGnl/z+y+jnLwJH4s7Xw4GVwN9O3KHWOolL2m/CndffBS6dsNoBuJvxlbgb9q9qrbPGmK8D3wG+ZIxpkaS8+qTGXFTbBmPMegCt9ZeBvwLeDmzD3eEfq7XeaYx5tbiBMWZDyfavaq2/hrsY/XPNohaisX0c+GtjTLFz2JPFPxhjbilZ7//VWv83XFvVJ3GJ+pXAv2utFbAOV5tG9P3LxpiXAbTWnwM2aa2vrWpJhGg8vwTO0FrfBRwP/N/R7w8CpwIPReff9cDxJU+w/hF3jn5uwv7ejsv3/sUYY4Hva60fnbBODviCMSYP3Ku17sOd97+pSgnFlCQxF9X2evEHY0wY9fxeET0a/zTwd8AarfV9uOYqb2itjwS+DJwCNOHepxv237UQYp5W45ql7Edr/RlcjdkKwOJq1LuiP98B/C+t9YG4mroQeDj62wpczV7Ra7hzd3mlgxeiwf0S9xl4MvAUrrnYDbgE+0VjzC6t9TLc5+MGrXVxOwX4k+xvBbA1SsqLXp+wzq4oKS8awD31EjUmTVlEta0u/qC19oBVuPbkGGNuNca8AzgYlwD8U7TqvwHPAUcYY9pw7VZVLYMWosG9Dhw2caHW+p3AX+ParGaNMR3AXqLzzxjTA/wU+BCu5vy2kg/7N3DnctFBQB7YXqUyCNGofoWrrb4U+KUx5hnc+fRexpqxdOP6YK0xxnREX+3GmMmS6W3AyqiWvWj1JOtNxc68iqgUqTEX1bZWa30ZcDeuo+cw8JuojflK4D9wbeoGGbvTbwX2AX1a66OB/wTsrHXgQjSwbwB/r7V+BngROA7Yijv38rjzLdBa/1dcjXmpW3EdxA4G3l2y/LvAZ7XWP462/0fgdmNMvqRGTwgxA2PMgNZ6A/BnwPuixb8C/gT3NKv4BPrfgf+ptf6kMWaH1nol8FZjzH0TdvlrXN+tT2qt/y3a52nAL2YZ0nZcvxFRA1JjLqrtLlztWg/wEeCyqKNYCtdxpRt4EzcyRLFd3GdwtXG9wL8Dt9c4ZiEa3ZcBg6v93od7TJ4B7sON2PA8rinKEPs/8r4bOAJ4szjCUuRG4NvAQ8Ar0bZ/Xr0iCNHQfgkkgEdLfm/FnV9Fn8XdWP9Ga70P+BmTjF0ejWh2GS6p3wNcBdyDqyibjRtw/cH2aK1/MPeiiLlQ1soTClEd0XCJhxtjrqp3LEIIIYRwtNa/Bf4/Y8w36x2LGE+asgghhBBCNDCt9RnAZtxT6g/jRnv5SV2DEpOSxFwIIYQQorEdhWu+1gy8DFxujNlW35DEZKQpixBCCCGEEDEgnT+FEEIIIYSIAUnMhRBCCCGEiAFJzIUQQgghhIgBScyFEEIIIYSIAUnMhRBCCCGEiAFJzIUQQgghhIiB/x8zydLgqS69gQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10b729fd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axes = plt.subplots(1,3,figsize=(12,4))\n",
    "sns.regplot(x='lpsa',y='Ip',data=treinoCP,ax=axes[0])\n",
    "sns.regplot(x='lcavol',y='Ip',data=treinoCP,ax=axes[1])\n",
    "sns.regplot(x='lweight',y='Ip',data=treinoCP,ax=axes[2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vemos que, apesar do modelo com todas as variáveis indique relação entre elas e o prognóstico, os modelos individuais acima deixam bem claro que essa relação não foi exatamente modelada. O modelo linear não consegue se ajustar aos dados. É percepível que o domínio do modelo não é restrito ao domínio da variável prognóstico ${0,1}$. Outra questão a ser levada em conta é que o modelo também não pode ser generalizado para mais de duas categorias. Assim, usamos como alternativo a regressão logística.\n",
    "\n",
    "No modelo logístico, ao invés de modelarmos o valor da variável dependente em função das variáveis independentes, modelamos a probabilidade condicional da variável dependente em função das variáveis independentes, $P(Y=1|x_1,x_2,\\dots,x_p)$. Por exemplo, poderíamos modelar $P(prognóstico=mau|lpsa)$. Nesse caso, poderíamos dizer que o prognóstico do paciente é mau se $P(prognóstico=mau|lpsa) > 0.5$, ou se quisermos ser menos alarmantes, $P(prognóstico=mau|lpsa) > 0.8$. Por conveniência, expressamos $P(Y=1|X) = p(X)$.\n",
    "\n",
    "Como discutimos, modelar $p(X) = \\beta_0 + \\beta_1X$ fornece um péssimo modelo. Na regressão logística, usamos a função logística como modelo:\n",
    "\n",
    "\\begin{equation*}\n",
    "p(X) = \\frac{e^{\\beta_0 + \\beta_1X}}{1+e^{\\beta_0 + \\beta_1X}}\n",
    "\\end{equation*}\n",
    "\n",
    "Após algumas manipulações algébricas, podemos verificar que\n",
    "\n",
    "\\begin{equation*}\n",
    "\\frac{p(X)}{1-p(X)} = e^{\\beta_0 + \\beta_1X}\n",
    "\\end{equation*}\n",
    "\n",
    "A proporção à esqueda é chamada de **chance**. Se tomarmos o logarítmo da chance (**logit**), temos o nosso modelo linear\n",
    "\n",
    "\\begin{equation*}\n",
    "\\log\\left(\\frac{p(X)}{1-p(X)}\\right) = \\beta_0 + \\beta_1X\n",
    "\\end{equation*}\n",
    "\n",
    "O modelo é ajustado pela máxima verossímilhança ao invés do método dos mínimos quadrados. Queremos estimar os valores de $\\beta_0$ e $\\beta_1$ tais que, quando substituídos na definição de probabilidade condicional do nosso modelo, o valor seja o mais próximo possível da probabilidade real. A função de verossimilhança é\n",
    "\n",
    "$$L(\\beta_0,\\beta_1) = \\prod p(x_i)^{y_i}(1-p(x_i))^{1-y_i} = \\prod_{i,y_i=1} p(x_i)\\prod_{i^\\prime,y_{i^\\prime}=0} (1-p(x_{i^\\prime}))$$\n",
    "\n",
    "Para encontrar o valor máximo, derivamos a função com respeito aos coeficientes e igualamos a zero. Como a derivada de produtos é mais complexa que a de somas, usamos o log da função de verossimilhança.\n",
    "\n",
    "$$\\log(L(\\beta_0,\\beta_1)) = l(\\beta_0,\\beta_1) = \\sum_{i,y_i=1} \\log(p(x_i)) + \\sum_{i^\\prime,y_{i^\\prime}=0} \\log(1-p(x_{i^\\prime}))$$\n",
    "\n",
    "O modelo pode ser resolvido por qualquer método iterativo para o cálculo do zero de funções, por exemplo, o método de Newton-Raphson. Não vamos explorar as propriedades matemáticas do modelo, tampouco resolvê-lo. O importante em nosso caso é que existem ferramentas para fazê-lo em Python."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimization terminated successfully.\n",
      "         Current function value: 0.527651\n",
      "         Iterations 6\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>Logit Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>        <td>Ip</td>        <th>  No. Observations:  </th>  <td>    67</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>               <td>Logit</td>      <th>  Df Residuals:      </th>  <td>    65</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>               <td>MLE</td>       <th>  Df Model:          </th>  <td>     1</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>          <td>Thu, 03 May 2018</td> <th>  Pseudo R-squ.:     </th>  <td>0.2013</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>              <td>10:32:26</td>     <th>  Log-Likelihood:    </th> <td> -35.353</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>converged:</th>           <td>True</td>       <th>  LL-Null:           </th> <td> -44.260</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th> </th>                      <td> </td>        <th>  LLR p-value:       </th> <td>2.434e-05</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>z</th>      <th>P>|z|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>   -1.9992</td> <td>    0.760</td> <td>   -2.631</td> <td> 0.009</td> <td>   -3.489</td> <td>   -0.510</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lpsa</th>      <td>    1.0833</td> <td>    0.311</td> <td>    3.484</td> <td> 0.000</td> <td>    0.474</td> <td>    1.693</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                           Logit Regression Results                           \n",
       "==============================================================================\n",
       "Dep. Variable:                     Ip   No. Observations:                   67\n",
       "Model:                          Logit   Df Residuals:                       65\n",
       "Method:                           MLE   Df Model:                            1\n",
       "Date:                Thu, 03 May 2018   Pseudo R-squ.:                  0.2013\n",
       "Time:                        10:32:26   Log-Likelihood:                -35.353\n",
       "converged:                       True   LL-Null:                       -44.260\n",
       "                                        LLR p-value:                 2.434e-05\n",
       "==============================================================================\n",
       "                 coef    std err          z      P>|z|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept     -1.9992      0.760     -2.631      0.009      -3.489      -0.510\n",
       "lpsa           1.0833      0.311      3.484      0.000       0.474       1.693\n",
       "==============================================================================\n",
       "\"\"\""
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "treinoCP = sm.add_constant(treino.copy())\n",
    "treinoCP['Ip'] = treinoCP.prognostico.map(lambda x: x=='mau' and 1 or 0)\n",
    "model = smf.logit(\"Ip ~ lpsa\",treinoCP).fit()\n",
    "model.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAENCAYAAAAVPvJNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xd4W+ed4Pvvi0KAnZSoYhXbsi0XufcaxyVucUsU63VJnNjxxJn0zM7cuXfuvfusN3ef3czd3cx4Ns4kiuOWxHFeOXbsxL0l7o4tJ44sF1muKlahRLGBIAng3T9ekIAokARJ4JwD4Pd5HorCwcE5v3MAnB/PW5W1FiGEEGKskN8BCCGECCZJEEIIIQqSBCGEEKIgSRBCCCEKkgQhhBCiIEkQQgghCpIEIYQQoiBJEEIIIQqSBCGEEKKgiN8BzJB0AxdCiOlRk61Q6QmCzZs3+x3CHjo6Oujs7PQ7jJKopmOB6joeOZbgCvrxLFiwoKj1pIhJCCFEQZIghBBCFCQJQgghREGSIIQQQhQkCUIIIURBnrRi0lrfAlwEbDPGHFbgeQXcCHwaSADXGGNe9SK2YmXWrMY+cg90boWOeajzlhM6/Fi/wwqcseeJgw6Ht9fscd6mcz4n3Ha8HpSCgcS03p/MmtXY39wGWzdDJg3ROojFobkVkgPQ3eVWnL8QtfxL4257NMZNH8Fg0m1LKZi/EI45ZcJ47QfvwGP3wUA/hEIQq4e998sdZ/42bYatoTBYO7r9ieLa4/zl73/k/7u6IJ3CtR5XEI1CJAr9vW6/4OKC3PnZazHqvOUA0/p+pH93Fzx2H1sHk25751xK+OIr9ox3mp+RaX0OAv499zJG5cWMclrr04E+4I5xEsSngW/hEsSJwI3GmBOL2LT1oplrZs1q7J0/gUgE6mIwNAipFOqqrxZ8Y4LexG0qpnIse5ynnl3up7UNmttGzxunnAXPP1n0+Zx02+EodG13K86aA+HIuNsrdDyZNauxt94IiV53wc1k9gwgHMY1G7fQ2IK65tt7bHs0xtSwSyg2u51QyG0XBW3thePdtTN7EbbZdbNi9TA8BPX1kEjktpkvFHIX+HHi2uP8pVK5/Tc2uwSQyYAKue1b67ZZ6DzkUyF3PKlsEmxoLPr9hGxy+P2v3WkNhbOJD7joctS+S6f0ndvjGKcQR6len69c14BSxZht5jppPwhPipiMMU8DOydY5VJc8rDGmBeBNq31Xl7EVgz7yD3uDYnF3RchFodIxC0Xo/Y4TwMJ9xFMJHY7bzx235TP54Tb7u12FysVgp7uKb8/9pF7YHAge4EcuTiP+e5kMi5JqBAMJApuezTGgQS5Ppwqt02bGT/ewYHdL/4qu//Bgdxxjtcv1NoJ49ottlh89/339rjf4C7Q+ducUDZZJhIuxoHE1L8fj93nNhOOuIQUjrjHj903re/cTL+nlfA99zrGoHSUWwhsyHu8Mbvs47Eraq2vB64HMMbQ0dFR9uC2d3WimlpQKnfRsOFGbFdnwf1HIhFP4vLCVI5l7HlKpVPuL8N0ikjEfdRsuJH0to8JNzQWfT4n2zaQK/rI7mu87RU6nu1dnWQyabe9VPbCmL3+5QLEfSFDIUinCRXY9kiM6c6tuQvs6HbsaHyF4k2NJpG8napscgmFIT1U8LyMxjZBXPmxKaWy5y67/9QwhKK54x6NY/zd5Y7L5oqkLKPvMUz+fgKuWCkcBqXc5pRyxzqYJDTF79zYY5xKHMW8fvbs2e54rWX05Fjy3i+727kLp1PMbm7KPs6dTJu/ft6vvBXyH+zxXOeObajG5t1jjMWxO7YxK1a350EVSvShSW8cRgUlQRTNGLMSWJl9aL0oysm0d7gig1g8t3AwCe2FbyNrtohp7HkKRyA1BJE6UqnshTFb1pxO9Bd9PifbNpC78IYjbl/jbK9gEVN7hyviyWRyF+Wx36uRC2ImA6EwmQLbHo0xHHHxjGxHZTdgrXuuULwj+1X5dxzZ35l09q98W/gLr5gwrj3O30h84Labf9yjv5k4SYwkzHAkWyxF7j2GSd9PwMUylASlXHKw1h1rLF7Ud85mssVh2WKxTOusvNdkz/1QEprb2P7eercsY3PP2fz/Z8g0trg7urpY7uCHBqG5lc61f53gZOypvb2drq6uKb2mGLa5Fds7EmNWNsYdH75f3EbCYRbOKu6PvqC0YtoELM57vCi7LBDUectdue1g0n2oBpOu3C9bOSecPc5TfYP7njU07HbeOOfSKZ/PCbfd3Jq9SGSgpXXK7486b7kr67eZXNHO2Ktj9i90bAbqGwpuezTG+gZyRVQ2t00VGj/eWH2uqAdyiSBWnzvO8YqMlZowrt1iG0zuvv/mlry6kvDu25xQNvM1NLgY6xum/v0451K3mXTKJanUsLuAn3UhnHkhDA25CvvUMCT6YXAQTj4Lu+1j7JZNsO1j2L4FOrfBju1w0pmuvmag3x1rMgHDw3DK2e71iYRrcJBMuhiHBt0+hofc+qed42IZyh7H0KB7z087d5Jz4aHTznUxDQ26c5YccMdw8JHYda9jX/sT9sU/YJ96kMxDd5O59w4yd/6YzC3/QuZH/5XMD/4jmf/yd0XvzpNKagCt9b7A78eppL4Q+Ca5Sup/M8acUMRmPamkhqm1HKjVOwgIfium8Y6nqFZMCpg3jVZMIQXzptGKKV4Pi8dvxcRIK6bs9oPeislmMu7ilkm7/aTTZB6+B/74sLvg1cXg1LMJnXmhi3fdWnj2UejqhPYOOO1cQgceOvFnZBqvKeXrR0z3DsIOJqGvJ++nF9vXA4k+6Ot1v7t2uHqk/DqjKVr8wCtQRCW1V62YfgWcAXQAW4H/BEQBjDE/zjZz/SFwPq6Z67XGmFeK2LRnCWIqajlBBF01HU9Qj8Wmht1f5OlU7nc67X7GUa4iGb+MPR6bybiLek+X+2OjZxd2pCVeb3fuZ2hw+jsNh6G+0d3N1TdAPO93vB4Vr3d/DNQ3sujq6yEoo7kaY66c5HkLfMOLWIQQpWGtzRbPDLuf4ZT77VGpRJC4BLALdnZCVyeJZILMxxtd3daunS4xTNZseKx4vWuG3NTifjc2QWMzqqHRPW5ock2LGxpdYqiL7VZ5Pa5wePJ1siqukloI4Q+bSrmEMJwtu8+vlK4B1lp3oe/cCp1bsTu2w45tsHO7K/ZJ587HwEQbitdDSzu0tLmf5lZUc6srzmxudfVCjS2oaLTsxzQZSRBCiD1Ya7NJYChXkTvVv4ArlM1kXB3E1s2w/WPsti3Q+TFs3+rOw2Ra2ojMmUeqqRXaZqHaZkPrLGhth9Z2VH7LrICTBCGEALJ3CENJ11poaLAmiopscgC2bIQtm7BbNsBI66iJEoFS0DYbOubCrLmo2XNcj/hZc6BtNioapbVK6lQkQQhRw+xI09fkQHF/HVcwO5iEjzfAxg+wmz+CzR+5IqLxRKIwZz7M2Qs1Z777f8c8mD0HFfG/+McLkiCEqDE2k872B8i2oa9CNpNxF/+P3sNueA82vA/bPx7/rqh9NsxfjJq/0DVJnr8Q2jtQoaB0FfOHJAghaoQdHnbt6JMDVVd8ZNNpd3fw4XrXn+TDd11/kkJmz4WF+6AW7A0L9nZ9OeL13gZcISRBCFHFrLVuML1Ef1XdLdhMxtUdvPc29v118OH6XGe+fPWNsHgJavF+sHhfWLAPqr7B83grlSQIIaqQTaVcL+mB/qppfWR7u+Gdtdj1b8K7b7m7obHaZ8O+S1H7HAB77+96dRfTN0AUJAlCiCpi02no74GByi9GspmMq1B+ew28s9YVIY3VNgv2Owi15CBYshTVOsv7QKuYJAghqoDNpN1YPQOJik4MdngY3n0T+9Zr8NYaNw5UvljcJYQDDoH9D4FZc+QOoYwkQQhRway1rn6hvyc7lHXlsUNDDL76Ipk/PQPrXt9zPKK5e8GBh6MOPBT23h81haEixMxIghCiQtnkgGvKWYFDXtjUMLzzBnbNK/D2Gvryk4JSrh7h4CPg4CNQs+b4F2g1CSnXtyMam3zdLEkQQlQYOzwEvd2kBxMVlRystbDhPexfXoLXV2enZs0KhWH/g1CHHgMHH4lqbPIv0GoRiUBdnUsI0TpUZOqXe0kQQlQIm067IaGT2aHgGhv9DahItmcX/PkF7KsvuIHtRijl6hMOO5a2U86ke2jYvyArncreHdTVubk66upQoZkXxUmCEKIC2EQ/9HVXTD2DTadh3evYV551LZDyK87nL0IdfRIcfpwbxRQINTbBUOWPXeSZkHLT7Y7eIUTL0utbEoQQAWZTKTfEdIV0crM9u2D1c9hXnnNxj2hsgiNPRB1zMmreQv8CrFSRSLb+IJsUIlFPWm9JghAioOxAwk1CE/C7Ble38D72xadg7au5jnlKwQGHoI49DQ46fFpl4DUpHHbTvUbr3F1Cme4OiiHvmBABYzMZNxVlcsJpZ3xn02lY+yr2+Sdg04e5Jxqb4JhTUcedKi2QJqOyc39ni4mI1gWqGa8kCCECxCYT0Nsz4fzNfrODSXj1eZcYdu3MPbFgMeqks+DwY2tmOOwpC4ehLjaaFIIwa9xEJEEIEQA2NQw9M5y0vszsQD+8+AfsC0/lRkpVCg45EnXKp2Dv/aRX81gjdQYjCSFAdwfFkAQhhM9sfx/09QR2iAzb34d97jH409O5EVMjUTjmZNQpZ6Nmz/U3wKAYKS6qixHumAvReMXPJyEJQgif2EwauncVHqY6AGyiD/vs4/DSH3J3NrE4nPhJ1MlnoZpafI3Pd0q5O4SY64hGtG70DkrF4qhQgdFmK4wkCCF8YAcHoXtnIIfitskB7HOPw/NP5BJDvN4VI510Rm3PpxCNQl18NClUe5GaJAghPGYTfa4iOmBFSnZ4GF5+GvvHh9wAgJBNDGfDyWfV7qxrdXUQr4dYfcXVIcyUJAghPGKtdf0aEonJV/aQzWRgzSvYx+5zdzXgikxOORt16qdq844hmk0K8dpLCvkkQQjhAZvJuCahAWulZD9cj33o7lw/hnAYjvsE6pPnjw6DUTOiUVfHEq+XZrpZkiCEKDObyUDXDhgOznAZtnsn9uF73KiqIw49BnXOpbXTKkkp1ychFqvJ4qNiSIIQooxsJp1NDsEYqdQOD8Gzj2OfeTgX06J9Uedfhtpnf3+D80Io5O4SYnGoi1V8M9RykwQhRJnYdBq6OgMzZ4Nd9zr29792MQE0taDOWw5HHF/dF8q8lkeqrvjJcoQkCCHKwg4Pw64dgRgyw3bvxD54N7zxZ7cgHIaTz0adcQEqFvc3uHJQKtvqKF6yeRFqlWcJQmt9PnAjEAZuNsZ8f8zzewO3A23Zdf4vY8yDXsUnRKnYwWS2j4O/zVhtJgN/ehr72G9zleNLDkRdfCVqznxfYyuLujqIN7hK5mq+I/KQJ2dRax0GbgIuAJYBV2qtl41Z7f8FjDHmaOAK4EdexCZEKdmBhGut5Hdy2LoZe/P/xD7wa5ccGptRl12Luva71ZUconXQ3Apz5qNmzUE1NEpyKCGv7iBOANYbY94D0FrfBVwKvJG3jgVG+u63Aps9ik2IkrCJfjdMt58xpNPw9CPYPz6YK9469hTUectR9ZUxRemkpI+CZ7xKEAuBDXmPNwInjlnnBuBRrfW3gEbgU96EJsTMBSI5bNmEvfd22Jz9qs2ag7r0KtR+B/saV0nU1UHM1SvIxEPeCdKZvhK4zRjzP7XWJwM/11ofZozZbbAarfX1wPUAxhg6Ojp8CHVikUgkkHFNRzUdC5TneDL9vWQGE9DeXtLtTiYcCdPe3o5Npxl4/H4GHrzb3TUoRfyMC2i4+PKKabUzciz5VCyOqm9AVeCdQrV8b7xKEJuAxXmPF2WX5bsOOB/AGPOC1joOdADb8lcyxqwEVmYf2s7OzrIEPBMdHR0EMa7pqKZjgdIfj+3vdeMq+aC9vZ2d767D/uY2+Og9t3DWHNTyLzK0zwEM9SegP1jDeoynvb2drq6uvIrmOIoQDCTdT4UJ+vdmwYIFRa3nVYJ4GViqtV6CSwxXAFeNWecj4GzgNq31IUAc2O5RfEJMme3rgb5ef/ZtLcnnn8TefXuuhdKJZ6DO/Syqrs6XmKYtGiXU0gaRyptQp9p5kiCMMSmt9TeBR3BNWG8xxqzVWn8PeMUYcz/w98BPtdZ/h6uwvsYYE6zhLoXIsr090O9Tchjox953J/1rX3ULmltRn/0iaunYhoEBFg5DfcPouEeh5hbUYHCGIhGOsgEbcniK7ObNwWvsFPTby6mopmOB0hyP7e+D3u4SRTTFfX/4LnbVLblRV5cd7SqiG5p8iWfKwmFobIb6ht3mUpDPmbeyRUyTTmYRpEpqIQLPDg666UG93m8m45qvPvk7N49ENErjimtJHHxUZUxaUxeDhsbanVOiQkmCEKJINpVyf7l7fNdt+3uxd98K6990C+YvQq34MvGDljHQ1eVpLFMSCrn+Cg2NMnx2hZIEIUQRbCbjyxSh9oN3sOZnuSKtE053I69GA3rBVWp0TgVi8cq4uxHjkgQhxCSstW74DA+H7LbWwvNPYB+91yWlWBx16edRhx/nWQxTopSrdG5slpZIVUQShBATsNa6+Rw8nAnODiaxv/1FbjKfeQtRV14fzIl8QiGob4TGRhk1tQpJghBiHG6a0B0w5F3zS9u5FXvnj2H7FrfgyBNRl1wVvL4N4TA0NkG8QQbHq2KSIIQowJfksO5114Q1OQDhMOqCFa7OIUjl+OM0UxXVSRKEEGO4OgfvkoO1Fp55FPv4fa6FVFOLK1LaO0BTgIZC7o6hoUkSQw2RBCFEnlydg0fJYXgIe+/PYc0rbsGifVFXfhXV0ubJ/icVCuXuGKQoqeZIghAiy+sKadvb7eobNn7gFhx9Euriq4LRhHUkMTQ0yh1DDZMEIcSI7i7vksPmDdhf/jv0dIFSqPOWwyln+38xVgoaGqGxSVolCUkQQkB28L3kgDf7evM1Vxk9PAR1MZS+DnXQ4Z7se1wh5YbZll7PIo8kCFHzbKLfs5FZ7QtPYh+621VGt81CfeHrqHkLPdl3QZGIK0qK1/t/9yICRxKEqGl2cNCTkVltJoN9aBW8+Ae3YNG+qM9/DdXUMuHryqq+AVraJDGIcUmCEDXLDibdEBplHnzPDg25IqW3XnMLlh2N+tw1/nV+U8rNIdHQ6M/+RcWQBCFqkk0moHtX+ZNDfx/2lz+CDe+7Bad+ys365leT0bo6aG4LRkspEXiSIETNsYl+6NlV/v3s3I6944ewY5trqfRpjTrpjLLvt6BwGJpbUPEGf/YvKpIkCFFTMj27vEkOmzdgf/5DN7lQJIq67FrUoUeXfb8FNTZDU7PUNYgpkwQhaoLNpKF7F5mG8s9oZt9723WAG0y6Hsif/zpqHx+GzYhGXSV0NGAD/YmKIQlCVD2bGnY9pNNpKHOCsGtfxa66FdIpaGlHfelbqLl7lXWfe1DKjZvUKHcNYmYkQYiqZoeHXHLwYCY4+8qz2PvvdBXfHfNR13wL1Tqr7PvdTTwOTa2oiHy1xczJp0hULTs06EZlzZR/Dmn7zKNu9jdwfRyu/gaqoans+x0VDrvipFjcu32KqicJQlQlz/o4WIt97D545hG3YP9D3FDdXl6o62LQ1i5jJ4mSkwQhqo4dSLiWSuVODpkM9vd3wcvPuAXLjkatuNbbsYyaW1CNzd7tT9QUSRCiqti+Hugr/7hKNp128zi89pJbcMwpbmrQsEd/xYfDhDvmoXq9GUNK1CZJEKIqWGvdXcNAovz7Sg27oTPe+ItbcPKZqAtWeNdiaKRIKRYDSRCijCRBiIpnhwZdckilyr+v4SHsr1bCO2vdgk9egDr7Yu+SQ2MzqtnHAf5ETZEEISqWzaSht8eTuwZwicj+8t/hvbcBUOd8BnX6eZ7sG6VcK6V6GSpDeEcShKhIu3V+82J/yQHsL26CD98FQF2oUSed6cm+CYXc3BF1MW/2J0SWJAhRcbzs/AauVZS943+5uaOVQl18Jer4T3iyb6JRaJ0lHd+EL+RTJyqKV/0bRvc30I+9/X/Bpg9dcvjM1ahjTvZk38TjbrgOv4YGFzXPswShtT4fuBEIAzcbY75fYB0N3ABY4DVjzFVexSeCz/b3eTL72+j+Ev3Y22+EzRtccvjcNagjT/Bm59K/QQSAJ3+aaK3DwE3ABcAy4Eqt9bIx6ywF/gk41RhzKPBdL2ITwWetxfZ0eZwc+rC3/qtLDqEQasWXvUkOSrn6BkkOIgC8unc9AVhvjHnPGDME3AVcOmadrwA3GWO6AIwx2zyKTQSYzaRdfUPCm5ZKkJ0F7tYbYcvGbHK4DnX4ceXfcSgE7bNR8fIPSS5EMbwqYloIbMh7vBE4ccw6BwJorZ/DFUPdYIx5eOyGtNbXA9cDGGPo6OgoS8AzEYlEAhnXdPh5LHZokPTOTmhscD8lEI6EaW9vH/f5TF8PPT//IektGyEUpunL3yF2lAd3DuEw4dlzpzQVqHzOgqtajidIldQRYClwBrAIeFprfbgxZrfpv4wxK4GV2Ye2s7PT0yCL0dHRQRDjmg6/jsX297nZ2EpcGd3e3k5XV9e4+7S35d05XP43JPZZSmKc9UsmHnfzRHdPrQhNPmfBFfTjWbBgQVHreZUgNgGL8x4vyi7LtxF4yRgzDLyvtV6HSxgvexOiCAKbSbte0cmkt/tNjE0OX0EtO6q8Ow2FoLlVOr+JwPIqQbwMLNVaL8ElhiuAsS2UfgtcCdyqte7AFTm951F8IgDsYBK6uzzr3zC630T/7nUOXiSHujrXv8Grwf2EmAZPKqmNMSngm8AjwJtukVmrtf6e1vqS7GqPADu01m8ATwH/hzFmhxfxCX/ZTMa1UvKw89vovgf6sbf/Wy456L/xIDnEoG22JAcReMp61OGoTOzmzZv9jmEPQS9/nIpyH4tNDrjmqx4NmZFfB2GTA65YadOHuaashx1b3gBicdeMtQSD+8nnLLiCfjzZOohJP4RBqqQWNcSm09DrfV3D6P6TA274jJEe0p+7pvzJIR53xUpejfwqxAxJghCes4l+6Ov2ZK7ogvsfTLqB9za875LD8i+hjji+vDsdGVNJkoOoIJIghGdsati1UBoa8i+GkSG7R0Zl/czVqKPGdskpsVDI1TlIchAVRhKE8ITt73VTgfpY52WHh+ld+T/g/XUAborQcg+8NzJ0hlRIiwokCUKUlU0NQ/cuGPbvrsHFkcLetZLhda8DoC663Jshu5tbZR4HUbEkQYiyKVdv6CnHkU5jzc0wkhwuuAx14hnl33FjE6qhsfz7EaJMJEGIkrOpFPR0+VrXMBpLOo29+1Z48zUAGi65kuTxp5d/x7E4qrm1/PsRooxkJhJRUnYgATu3BSM5ZDLYe38Or692C868kPpzxw4iXAaRCLSOPyCgEJVC7iBESQShhVI+ay32d7+C115yCz5xLurMC8u/45EWSzILnKgCkiDEjFhrXeukRJ/vdQ0jrLXYB1fBK8+6BSefiTrnM+VvZjrSYknmjxZVQj7JYtrs8LCraxge9juUUdZa7GO/hRefcguO/wTqghXe9EFobZMWS6KqSIIQ0xKEfg0FPfUAPPOo+//RJ6EuusKb5NDciorLsN2iukwpQWit24ALgQXAZuDBkSlCRW2w6TR07wxMXUM++8yj2KcecA8OO9b1kvaiLqCxCdXYVP79COGxor89WuuzgA+AbwPHA9/CTexzdnlCE0FjB5OwIxgtlMayLz6FffRe9+DgI1GXXetNcqhvkOasompN5Q7ih8D1xhgzskBrvQK4CTi41IGJ4HAV0T3Q3+d3KAXZV57FPpD9WB6wDHX5dd4MbRGLQ0tb+fcjhE+m8ifWAuA3Y5bdC8wvXTgiaGxyADq3Bjc5/OUl7P13ugdLDkRd+VVUJFr+HdfVlWxeByGCaioJ4ufAN8Ys+xpwR+nCEUFh02nsrp2kd2z3bDKfqbKvr8bec7urKN97P9Tnv4aqqyv/jqN1MjqrqAlTKWI6GvhbrfU/4uaVXgjMBV7SWj89spIxxoNxDES52EzG3S2M9Guoj/sdUkH2rb9iV93iYly4D+rqb6JiHsQarYN26QgnasNUEsRPsz+iStlkAnp7AnvHMMK+8wb2rp+6+avnLUR98VuoeH35d1xXJ72kRU2ZNEFkWy8BbBhvHWPMkyWLSHjOZjJumIzkgN+hTMq+vw57548hnYKO+ahrvu3NiKmSHEQNKuYO4meTPG+B/UoQi/CBHRqE7q7A3zUA2I/exf7iR5AahllzUNd+B9XUUv4dRyU5iNo0aYIwxizxIhDhLZvJuKariX6/QymK3fgB9o4fwtCgaz107XdRXjQxjUalzkHULBlqowZVSl3DCLt5A/b2f4PBJLS0ueTQNqv8O5bkIGqcJIgaYjPpbF1D0u9Qima3bMLefqOrH2lqcclh1pzy73h02G6ZS1rULkkQNaKS6hpG2O1bsLfd6IrBGptccuiY583OW1q96Y0tRIBJgqhy1loYGXm1gtjOrdhb/sXFXt+IuuY7qLl7ebPz+gYZmVUIJEFUtSDO11AMu3M79tZ/dZXo8XrXlHX+Im92Hg6DDL4nBCAJoipZa11P6CDO1zAJ29Xp7hx6dkEsjvrSt1EL9vYugJZ2qZQWIksSRJWxqWHo3gXDwRuSezJ21w6XHLq7ssnhW6hF+3oXQGMzKiYzwgkxQhJEFbH9fa5YpsLuGgDsrp3Yn/0L7NoJdTHU1d9ALfaw/2Usjmr2oNOdEBXEswShtT4fuBEIAzcbY74/znqfA+4GjjfGvOJVfJXMDg9D765ATuRTDNu909057NoB0TqXHPY5wLsAwmFolXkdhBjLk8JWrXUYN7HQBcAy4Eqt9bIC6zUD3wFe8iKuSmczaWxPV2BneSuG7e7C3vKv0NUJ0ahLDvsu9S4ApaS/gxDj8Ko27gRgvTHmPWPMEHAXcGmB9f4/4J+ByunJ5QNrLba/103kk0j4Hc60jd457NwOkSjqC19HLTnQ2yBa21BRDyYYEqICeZUgFrL7aLAbs8tGaa2PARYbYx7wKKaKZAcSLjH09kCm8uoaRozWOezcnr1z+DpqP49nrm0HJeGAAAAYRklEQVRqlv4OQkwgEJXUWusQ8APgmiLWvR64HsAYQ0dHR3mDm4ZIJFLSuKy12EQ/tq8HGwJavKtMDUfCtLe3l3Sb6Z2d9Nx+I7arE6J1tPztPxI96LCS7mM8I8ejGhoJt8/2ZJ/lUurPmZ+q6Vigeo7HqwSxCVic93hRdtmIZuAw4A9aa3DzXN+vtb5kbEW1MWYlsDL70HZ2dpYt6Onq6OigFHG5/gz9rk+DT0NktLe309XVVbLt2a4deRXSrlipb+5CKOE+JtLe3k5Xfz/U1aMC+NmZilJ9zoKgmo4Fgn88CxYsKGo9rxLEy8BSrfUSXGK4Arhq5EljTDcwmm611n8A/qGWWzHZoUHo7a64XtATsTu35/o51MV8qXNQkYgbLlzmkxZiUp7UQRhjUsA3gUeAN90is1Zr/T2t9SVexFApRlsm7eysruSwYxv2Zz/IJYcvftP7CmmlCLVLiyUhiqVsBXaqymM3b97sdwx7mO7tpU30u45umUwZopqeUhQx2W0f58ZWisXdHNJ7+zAJYUsrc/beN9C3/lMR9GKMqaimY4HgH0+2iGnS2+hAVFLXOjeoXmUOjzEZ+/EG7G3/5upR4vVubCUvh88YEa9HNTR5v18hKpgkCB/ZTDo77Wfl9mWYiN34gZsJLjkADY2oL30HtWDx5C8stXAYvJieVIgqIwnCBzaTcRfNgBUnlZL94B3sL37kpgltanHzOcwrruVEyckIrUJMiyQID9nBQUgmYHCgoju5Tca+sxb7q5+4SvaWdtS13/FuJrixGhpkhFYhpkkShAfsYNLdLVRRq6Tx2LWvYlfd4vptzJrj7hz86pAWCkGTjNAqxHRJgigjOzToEkOFDqQ3VfbPL2Dv/bkbbnzuXi45+Dk7W0urNGkVYgYkQZSBHUxid26vmcQAYJ97HPvwb9yDhfu4fg5+thqKxWWcJSFmSBJECdlkAvr7SQ8maiY5WGuxj98PTz/sFiw5EHXV36Li9f4FFQpJqyUhSkASxAzZTAYGEmPGS2r0NSav2EwG+/u74OVn3IJDjkStuM7/4bNbWlFhKVoSYqYkQUyTG0ivD/r7qrap6kTs8BB21a3w5l/cgmNORl3yef8vzPF6KVoSokQkQUyRtRYG+l1i8GmEVb/ZgX7sL38MH653C047B3XuZ/0fAE86xAlRUpIgihSEobeDwHZ3Ye/4IWxzY2CpCy5DnXK2z1FlSYc4IUpKEsQkCtcx1Cb78UbsL25y40aFw6jlX0IdcbzfYTkNjdIhTogSkwQxDptOuzuGgb6q7vVcLLv+TexdK93QGbE46sqvovb3eIrQ8YTD0iFOiDKQBJHHplLuApgcqMqRVafLvvoC9r5fuMr4ljbU1d9EzV84+Qu90twqRUtClIEkCLLFSP297o6hsufHKCmbyZC4/y7so791C+YtRF39DVRraeeonpF4vb99LoSoYjWfIGwyAb09NV+/MJYdHsL+5nYG1r7qFhxwCOryrwTrYhwKQYuPQ3kIUeVqNkG4AfR6pSipANvbjb3zx7DxA7fghNNRn9b+93EYq6VNxloSooxqLkFIYpiY3fShSw49u0ApGpZfzcCRJ/nfx2GslrZg3c0IUYVqIkG4zm3ZpqqplN/hBJZdsxp77+1uWPK6GEpfR/1Jp5Oc4ZzUJdfcgmqojeFMhPBTVScIm8k2VU301+RwGMWymQz2yd/BH7MD7rV3oD7/Nf9mgJtIYxOqsdnvKISoCVWZIGxq2N0tDAxIq6RJ2IEE9u5bYd3rbsGSA1FXfMXfobrHUxfzd34JIWpMVSUIO5h0YyQNDfodSkWwWze7+oad292Ck85Enf+54FVGg+sM1xag5rVC1ICKTxCjQ2EM9Ev9whTY1/6Eve+XrrI+EnEjsR59kt9hFaYUtM2SFktCeKziEwSdW2QojCmwqZSb+e2lP7gFre1u2IyF+/ga14SaWlDROr+jEKLmVH6CkORQNLtrB/bXN+f6NxywDHXZtajGANY3jIjXBzs+IapY5ScIURT7xl+w997hxplSCnXGp+GMTwd7DKNoHQRpWA8haowkiCpnU8PYR++FF55yCxqb3F3DAcv8DWwy4bCrdwhaBz0haogkiCpmt2/BrroFPt7gFiw50CWHoM+6NlIpHcTWVELUEEkQVchaC688i31olesVrZQrTgp6kRK4WFvbpVJaiACQBFFlbH+va7765mtuQess1GXXoPZd6mtcRQmF3J1DncwMJ0QQeJYgtNbnAzcCYeBmY8z3xzz/H4C/AVLAduDLxpgPvYqvGtg3X3MT+/T3uQWHH4e6+EpUfYO/gRUjHIb22ahI1O9IhBBZnpQ3aK3DwE3ABcAy4Eqt9dha0j8DxxljjgDuBv5/L2KrBnYgQeaeO1yv6P4+1zT0c9egVny5MpJDJAKz50hyECJgvLqDOAFYb4x5D0BrfRdwKfDGyArGmKfy1n8R+IJHsVU0+/YaV6TU2+0W7HcQavkXUa2z/A2sWOGwGxxQekkLETheJYiFwIa8xxuBEydY/zrgoUJPaK2vB64HMMbQ3h68dvLhSLjscWX6++i/+zaGXn7WLaiL0XDJlcRPP7ekFdFlPZZQiHDHXE8rpCORCB0dHZ7tr5zkWIKrWo4ncJXUWusvAMcBnyz0vDFmJbAy+9B2BW2uAqC9vZ1yxWWthb++7FoojdQ1LDkQ9ZkvkJw1h2R3d0n3V7ZjUcrVOXT3lH7bE+jo6KCzs9PTfZaLHEtwBf14Fiwobih/rxLEJmBx3uNF2WW70Vp/Cvh/gE8aY2RI1jFsVyf2/l/B+mzJXCyOOvezcNxpwW++OlZrm7RWEiLgvEoQLwNLtdZLcInhCuCq/BW01kcDPwHON8Zs8yiuimBTw/Ds49inH3L9GgAOOQp1oUZV4lAULa2oeAVUngtR4zxJEMaYlNb6m8AjuGautxhj1mqtvwe8Yoy5H/jvQBOwSmsN8JEx5hIv4gsy+85a7O9/nZuzoaUNdeHlqGVH+RvYdDU2BXMyIiHEHpSt7BnX7KZXX/Y7hj2Uotze7tjmhuV+669uQSgEJ5+FOvNCVCxegiiLU9I6iHg9qs3f1lVBLxueCjmW4Ar68WTrICYd6CxwldS1ziYHsH94EF58CtJpt3DfpaiLrgjmHNHFisdlZFYhKowkiICw6TS8/Az2qQfcfNrghp04bzkcekxlj2oai7shPyr5GISoQZIgfGathTf+jH3sPtiRrZuvi6E+cR6cenblD1pXF5Nhu4WoUJIgfGKthfVvYp+4HzZlh5xSCo49FXXWRajmVn8DLIV4vRuZVZKDEBVJEoQP7PvrsE/8Dj5cn1t40OGocz+LmruXf4GVilLQ3CKtlYSocJIgPGKthffXuTqGD97JPbHkQNSnLkXtvZ9/wZVSJOLqG6Iy8J4QlU4SRJlZa+GdtdinH4YP3809sXgJ6uyLYb+Dq6cIZqS+odJ6dQshCpIEUSY2nYbXV2OfeRS25o0qsvf+qDMvhP2rKDEANDRAc1t1HZMQNU4SRInZ5AADTzyHffJB6N6Ze2L/Q1CnnwtLDqqui6hS0NSCapT6BiGqjSSIErGdW7F/+iO8+gKJwaRbqJQbM+n081AL9/E3wHKIRqGlXeobhKhSkiBmwGYysO517Et/zI2wCq4s/piTUSefhZo1x78Ay0UpaGyCxubquhsSQuxGEsQ02O6dsPp57OrnoSdvnKK22agTT6ft7AvpHhz2L8BykrsGIWqGJIgi2eFheOuv2D+/4O4W8gc5POAQ1IlnwIGHoUIhQg1NMBi8iYxmRO4ahKg5kiAmYK2Fj97D/vVPsOYVGEjknmxshmNOQR13anUWI+WL1rlhxuWuQYiaIgliDGstbN2Mff0V+OvL0LUj92QoBEsPRR19suv5HKny0xdShNraUTGZ3EeIWlTlV7jijCaFN/4Ma1ZD55bdV5i/CHXUiXDE8dUxRlIx6hugqYVQYzMMyOyvQtSimk0QNpOBje9j33wN3vhLbsa2Ee2z4fDjUEecUNnzMExVXR00t1b+KLJCiBmrqQRhBxLw7lvYdWtg3evQ37f7Cm2z4bBjUIcdCwv2rq3K2HDYJYZ4vd+RCCECoqoThM1kYPNHbljt9Wthw/uQyey+0twFsOwoN8fz/EW1lRTA1as0NUN9Y+0duxBiQlWVIKy1sH2LGzX1vbfgvbchObD7SuEILFmKOugIOOgwVHuHP8H6LRRyzVbrG2VwPSFEQRWfIOzmDfDReuwH690w2v29e640aw4sXYZaeqgbXrsu5n2gQREOQ0MTNMgdgxBiYpWfIP79v+65sKnFJYL9D3bDabfP9j6woIlEXN+NeL0kBiFEUSo+QQDQ0uaG0V5yICw5EDrmyUVwRCzu7hZicb8jEUJUmIpPEOrv/4ubwUwSQk4o5Pox1DdWf2c+IUTZVPzVQ7VJ8RHgxkqKxSDeALG4JEwhxIxVfIKoaUq5ocXj9RCLoUJhvyMSQlQRSRCVJhRy9QqxGNTFpYmqEKJsJEEEXTjsRlOtq4NoTEZUFUJ4RhJEkEQiEIm6hBCJQjQixUZCCN9IgvBDKJRLBpGom6UtEpWKZSFEoHiWILTW5wM3AmHgZmPM98c8HwPuAI4FdgCXG2M+8Cq+kgqF3EU/HIZQePffEbkrEEJUBk8ShNY6DNwEnANsBF7WWt9vjHkjb7XrgC5jzAFa6yuAfwYuL3UsmXVr4dlHoasT2jvgtHMJHXjo5K976gF47gkYSrpK4tPPJ3T+cnfhD4WyCcD9jsydi4p2zizONauxj9wDnVtdx7/zlhM6/NiSbg8o6T4A0r+7Cx67DwYHIFYP51xK+OIrCsey6SNIp1wy3WsxHHQ4vL0mtxwLKPd8c6trtTWQmDDWyba9dctGGBoaXV7sMZfq/RhvO7stj9cXPNaxr02uuAb2WTrlGEqpVOclufoF0qtuK+lnUcycsvlzK5eJ1vpk4AZjzHnZx/8EYIz5b3nrPJJd5wWtdQTYAswxxkwUoN306st7LlUq9zP6OETm7TXw21+64p26OnehSKfhsmsIHXoUqBAo3O9QKPe6Bww8YNxzoTBk0u7addHle1z8ADo6OujsnH6CyKxZjb3zJ9k4YzA0CKkU6qqvTv+iNHZ7/X3u+BoaJ9zHVI4l/bu74Pe/nvA8jcaSGobebtwKuIviwADU10MymX2tde+Dtbn3c9YcN+BigViL2vbgINjsiL4tbRCOTnpeS/V+jLcdTjkLnn/SLU+loCs7N0nese62Tva1YZshc/lXfLuQlvK8hH79U9IqVJLPexDM9BpQbgsWLAD3TZ2QV20kFwIb8h5vzC4ruI4xJgV0A5P3gps7H+buBfMWoOYvdD/zFqDm7oWaM9/9dMxDzZ4Dzz/hmoc2NLqy/4ZGlyj++BCqoQlV34CKN6BicVS0DhWJosJhePx+dyrDEXeRCkfc48fuK8W52YN95B73pYvFsx3g4hCJuOWl2t7ggPsLtUT7ANz5mOQ8jcYykMhbJwSJhFs3kcgldnDJAdww7SoEPd3jxlr0tvOXF3HMpXo/xtsOj92XW97bnf1DZfdj3W2d0ddGZ/Z+zVBpz0u0tJ9FURIVV0mttb4euB7AGMOcufOKfu32rk5UU8tulcE23Ijt6qSjY/xhv7cOJl39Qf6FKxSGwWTB10UikQm3V644p7K9VPav+0jeUByF9jGVYynmPI3Eku7cmrtLC4XcX/2hOkgPQSgKqWxisKP/uPXSKSKRSMFYi9r2yJ1IdlvhhsnPa6nej/G2k972MeHs6LqpdMrFBrsda/46oyJhQtP8TJRCSc9LcwsRSvN5D4KZXgOCwqsEsQlYnPd4UXZZoXU2ZouYWnGV1bsxxqwEVmYf2qncxmXaO6C7y/2FMmIwCe2T3A7G4q7uIf/LmUlDLF7wdTMuYppunFPZXigMFlKp1IT7mNKxFHGeRmMJR1wdQSiUuzvIpLO/M24b1mZvgrP/z2QgHHExF4i1qG3n35GEI6QT/ZOe11K9H+NuJxZ3ccTiudhh92PNXycrnBomM83PRCmU8ryovh7Skbw+PjP4vAdBhRQxTcqrIqaXgaVa6yVa6zrgCuD+MevcD3wp+//LgCcnqX+YMnXecleeO5h0F4rBpCvrzFbYjuucS90fsemUe1065R6fc2kpw5t5nFPZXqzeDehXon0ARZ2n0VjqG/LWyUBDg1u3oSF3EYdcsgmF3HotrePGWvS285cXccylej/G2w7nXJpb3tzqYhtzrLutM/ra4Zm9XzNU2vMyXNrPoiiJ8A033FD2naxYsSKzatWqd4BfAt8CfmGM+Y3W+nurVq1qXrFixdurVq1aA3x+1apV/w04CvjbFStWdE2y6Rt6ewtMEDQONW8BzFsAGz+A7p0wqwP12asnrQgLHXSYK+T48F0YHnIX1ws+V7CCGqChoYFEIlF0XKWKc0rbu+xa1FEnTrqPqRxLMedpNJYtG13FcSjkKpEX7uMqYoeH3PJweLR/CPUNMHuuu3gOD40ba1HbHkoCyi2fv6io81qq92O87YQ/cW5ueaLPzY0+5lh3Wyf72parv87wgYdNKYZSKuV5adr/IAbffbskn/cgmOk1oNyam5sB/vNk63nSiqmM7ObNm/2OYQ9Bv72cimo6Fqiu45FjCa6gH0/QWjEJIYSoMJIghBBCFCQJQgghREGSIIQQQhQkCUIIIURBFd+Kye8AhBCiQlV9KyYVxB+t9Wq/Y5Bjqf7jkWMJ7k+FHM+kKj1BCCGEKBNJEEIIIQqSBFEeKydfpWJU07FAdR2PHEtwVcXxVHoltRBCiDKROwghhBAFVdyEQUGmtT4fuBEIAzcbY77vc0jTprW+BbgI2GaM8W/I0BLQWi8G7gDm4ZpGrzTG3OhvVNOntY4DTwMx3Hf4bmPMf/I3qpnJzlv/CrDJGHOR3/HMhNb6A6AXSAMpY8xx/kY0fXIHUSLZD/hNwAXAMuBKrfUyf6OakduA8/0OokRSwN8bY5YBJwHfqPD3ZhA4yxhzJG5o/PO11if5HNNMfQd40+8gSuhMY8xRlZwcQBJEKZ0ArDfGvGeMGQLuAsozo5AHjDFPAzv9jqMUjDEfG2Nezf6/F3chGjsnesUwxlhjTF/2YTT7U7GViVrrRcCFwM1+xyJ2J0VMpbMQ2JD3eCNwok+xiHForfcFjgZe8jmUGcnesa4GDgBuMsZU8vH8K/CPQLPfgZSIBR7VWlvgJ9lpkiuS3EGImqG1bgJ+A3zXGNPjdzwzYYxJG2OOws3vfoLWuiLribTWI/Vcq/2OpYROM8Ycgytu/obW+nS/A5ouSRClswlYnPd4UXaZCACtdRSXHH5pjLnH73hKxRizC3iKyq0vOhW4JFuxexdwltb6F/6GNDPGmE3Z39uAe3HFzxVJiphK52VgqdZ6CS4xXAFc5W9IAkBrrYCfAW8aY37gdzwzpbWeAwwbY3ZpreuBc4B/9jmsaTHG/BPwTwBa6zOAfzDGfMHXoGZAa90IhIwxvdn/nwt8z+ewpk3uIErEGJMCvgk8gqsENcaYtf5GNX1a618BLwAHaa03aq2v8zumGTgVuBr31+lfsj+f9juoGdgLeEpr/VfcHyaPGWN+73NMwpkHPKu1fg34E/CAMeZhn2OaNulJLYQQoiC5gxBCCFGQJAghhBAFSYIQQghRkCQIIYQQBUmCEEIIUZAkCCGKpLX+QGv9Kb/jEMIrkiCEEEIUJAlCCCFEQTLUhhBTpLW+ATgMNyHMp4F3gGuNMa9ln/8/gW8DLcBm4OvGmCe01ifgJpQ6BBjAjQ31H7LDwwsROHIHIcT0XAqsAmYBdwK/1VpHtdYH4YZcOd4Y0wycB3yQfU0a+DugAzgZOBv4usdxC1E0uYMQYnpWG2PuBtBa/wD4e9xsdR/jpgJdprXeboz5YOQFY4a0/kBr/RPgk7j5EIQIHEkQQkzP6ORQxpiM1nojsMAY84zW+rvADcChWutHcMVIm7XWBwI/AI4DGnDfv2qaB0FUGSliEmJ6Ruf+0FqHcPN/bAYwxtxpjDkN2Ac3u9jIUNz/DrwFLDXGtAD/N6C8DFqIqZA7CCGm51it9XLgflyF9CDwYrYOYiHwHJDEVUaHs69pBnqAPq31wcDXgO1eBy5EseQOQojpuQ+4HOjCzTWx3BgzjKt/+D7QCWwB5pKdEAf4B9wkUr3AT4FfexyzEFMi80EIMUXZZq4HVPLMZ0IUQ+4ghBBCFCQJQgghREFSxCSEEKIguYMQQghRkCQIIYQQBUmCEEIIUZAkCCGEEAVJghBCCFGQJAghhBAF/W9Q8MOJfO2w7QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10b7ef390>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "sns.regplot(x='lpsa',y='Ip',logistic=True,data=treinoCP);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimization terminated successfully.\n",
      "         Current function value: 0.423807\n",
      "         Iterations 7\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<table class=\"simpletable\">\n",
       "<caption>Logit Regression Results</caption>\n",
       "<tr>\n",
       "  <th>Dep. Variable:</th>        <td>Ip</td>        <th>  No. Observations:  </th>  <td>    67</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Model:</th>               <td>Logit</td>      <th>  Df Residuals:      </th>  <td>    61</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Method:</th>               <td>MLE</td>       <th>  Df Model:          </th>  <td>     5</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Date:</th>          <td>Thu, 03 May 2018</td> <th>  Pseudo R-squ.:     </th>  <td>0.3585</td>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Time:</th>              <td>10:32:59</td>     <th>  Log-Likelihood:    </th> <td> -28.395</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>converged:</th>           <td>True</td>       <th>  LL-Null:           </th> <td> -44.260</td> \n",
       "</tr>\n",
       "<tr>\n",
       "  <th> </th>                      <td> </td>        <th>  LLR p-value:       </th> <td>6.717e-06</td>\n",
       "</tr>\n",
       "</table>\n",
       "<table class=\"simpletable\">\n",
       "<tr>\n",
       "      <td></td>         <th>coef</th>     <th>std err</th>      <th>z</th>      <th>P>|z|</th>  <th>[0.025</th>    <th>0.975]</th>  \n",
       "</tr>\n",
       "<tr>\n",
       "  <th>Intercept</th> <td>    4.5834</td> <td>    2.935</td> <td>    1.562</td> <td> 0.118</td> <td>   -1.169</td> <td>   10.336</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lpsa</th>      <td>    1.0838</td> <td>    0.539</td> <td>    2.011</td> <td> 0.044</td> <td>    0.027</td> <td>    2.140</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lcavol</th>    <td>    0.0818</td> <td>    0.433</td> <td>    0.189</td> <td> 0.850</td> <td>   -0.767</td> <td>    0.930</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lweight</th>   <td>   -1.6921</td> <td>    0.890</td> <td>   -1.901</td> <td> 0.057</td> <td>   -3.437</td> <td>    0.053</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lcp</th>       <td>    0.9065</td> <td>    0.412</td> <td>    2.198</td> <td> 0.028</td> <td>    0.098</td> <td>    1.715</td>\n",
       "</tr>\n",
       "<tr>\n",
       "  <th>lbph</th>      <td>    0.2431</td> <td>    0.265</td> <td>    0.917</td> <td> 0.359</td> <td>   -0.276</td> <td>    0.762</td>\n",
       "</tr>\n",
       "</table>"
      ],
      "text/plain": [
       "<class 'statsmodels.iolib.summary.Summary'>\n",
       "\"\"\"\n",
       "                           Logit Regression Results                           \n",
       "==============================================================================\n",
       "Dep. Variable:                     Ip   No. Observations:                   67\n",
       "Model:                          Logit   Df Residuals:                       61\n",
       "Method:                           MLE   Df Model:                            5\n",
       "Date:                Thu, 03 May 2018   Pseudo R-squ.:                  0.3585\n",
       "Time:                        10:32:59   Log-Likelihood:                -28.395\n",
       "converged:                       True   LL-Null:                       -44.260\n",
       "                                        LLR p-value:                 6.717e-06\n",
       "==============================================================================\n",
       "                 coef    std err          z      P>|z|      [0.025      0.975]\n",
       "------------------------------------------------------------------------------\n",
       "Intercept      4.5834      2.935      1.562      0.118      -1.169      10.336\n",
       "lpsa           1.0838      0.539      2.011      0.044       0.027       2.140\n",
       "lcavol         0.0818      0.433      0.189      0.850      -0.767       0.930\n",
       "lweight       -1.6921      0.890     -1.901      0.057      -3.437       0.053\n",
       "lcp            0.9065      0.412      2.198      0.028       0.098       1.715\n",
       "lbph           0.2431      0.265      0.917      0.359      -0.276       0.762\n",
       "==============================================================================\n",
       "\"\"\""
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model = smf.logit(\"Ip ~ lpsa + lcavol + lweight + lcp + lbph\",treinoCP,).fit()\n",
    "model.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora observamos que o modelo se ajustou melhor aos dados. Da mesma forma que nos modelos de regressão linear simples e múltipla, os p-values estão relacionados com a hipótese nula de que os coeficientes não estão relacionados à variável dependente. Por outro lado, a interpretação dos coeficientes não é tão direta. Podemos interpretá-los como: valores positivos indicam que o aumento de $X$ multiplica a chance por $e^{\\beta_1}$ unidades, enquanto valores negativos a dividem."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Predizendo os valores do conjunto de teste\n",
    "\n",
    "Uma vez que obtivemos o modelo, agora podemos fazer previsões para o conjunto de teste e avaliar sua precisão."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "7     0.180741\n",
       "9     0.126993\n",
       "10    0.206372\n",
       "15    0.412494\n",
       "22    0.900738\n",
       "25    0.358521\n",
       "26    0.536305\n",
       "28    0.302081\n",
       "32    0.376385\n",
       "34    0.413963\n",
       "36    0.321005\n",
       "42    0.349197\n",
       "44    0.751478\n",
       "48    0.659462\n",
       "49    0.505584\n",
       "50    0.743541\n",
       "53    0.857350\n",
       "54    0.919921\n",
       "55    0.557184\n",
       "57    0.950451\n",
       "62    0.975161\n",
       "64    0.974580\n",
       "65    0.538045\n",
       "66    0.781062\n",
       "73    0.800235\n",
       "74    0.977351\n",
       "80    0.895569\n",
       "84    0.982613\n",
       "95    0.998565\n",
       "97    0.999024\n",
       "dtype: float64"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.predict(teste)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Como podemos ver, o modelo retorna as probabilidades associadas a cada objeto no conjunto de teste de terem um mau prognóstico. Podemos converter tais probabilidades em classes. Podemos arbitrariamente dizer que se a probabilidade for maior que 0.5 então a classe é mau, e bom caso contrário."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "7     bom\n",
       "9     bom\n",
       "10    bom\n",
       "15    bom\n",
       "22    mau\n",
       "25    bom\n",
       "26    mau\n",
       "28    bom\n",
       "32    bom\n",
       "34    bom\n",
       "36    bom\n",
       "42    bom\n",
       "44    mau\n",
       "48    mau\n",
       "49    mau\n",
       "50    mau\n",
       "53    mau\n",
       "54    mau\n",
       "55    mau\n",
       "57    mau\n",
       "62    mau\n",
       "64    mau\n",
       "65    mau\n",
       "66    mau\n",
       "73    mau\n",
       "74    mau\n",
       "80    mau\n",
       "84    mau\n",
       "95    mau\n",
       "97    mau\n",
       "dtype: object"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.predict(teste).map(lambda x: x > 0.5 and 'mau' or 'bom')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Podemos agora comparar tais valores com os valores originais para avaliarmos a precisão do modelo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.DataFrame(dict(prob=model.predict(teste),\n",
    "                  pred=model.predict(teste).map(lambda x: x > 0.5 and 'mau' or 'bom'),\n",
    "                  real=teste.prognostico))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "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>pred</th>\n",
       "      <th>prob</th>\n",
       "      <th>real</th>\n",
       "      <th>gleason</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.180741</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.126993</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.206372</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.412494</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.900738</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.358521</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.536305</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.302081</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>32</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.376385</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>34</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.413963</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>36</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.321005</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>42</th>\n",
       "      <td>bom</td>\n",
       "      <td>0.349197</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>44</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.751478</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>48</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.659462</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>49</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.505584</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.743541</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>53</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.857350</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>54</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.919921</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>55</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.557184</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>57</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.950451</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>62</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.975161</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>64</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.974580</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>65</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.538045</td>\n",
       "      <td>bom</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>66</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.781062</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>73</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.800235</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>74</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.977351</td>\n",
       "      <td>mau</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>80</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.895569</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>84</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.982613</td>\n",
       "      <td>mau</td>\n",
       "      <td>9</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>95</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.998565</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>97</th>\n",
       "      <td>mau</td>\n",
       "      <td>0.999024</td>\n",
       "      <td>mau</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   pred      prob real  gleason\n",
       "7   bom  0.180741  bom        6\n",
       "9   bom  0.126993  bom        6\n",
       "10  bom  0.206372  bom        6\n",
       "15  bom  0.412494  mau        7\n",
       "22  mau  0.900738  mau        7\n",
       "25  bom  0.358521  bom        6\n",
       "26  mau  0.536305  bom        6\n",
       "28  bom  0.302081  mau        7\n",
       "32  bom  0.376385  bom        6\n",
       "34  bom  0.413963  bom        6\n",
       "36  bom  0.321005  mau        7\n",
       "42  bom  0.349197  mau        7\n",
       "44  mau  0.751478  mau        7\n",
       "48  mau  0.659462  mau        7\n",
       "49  mau  0.505584  bom        6\n",
       "50  mau  0.743541  bom        6\n",
       "53  mau  0.857350  mau        7\n",
       "54  mau  0.919921  mau        7\n",
       "55  mau  0.557184  mau        7\n",
       "57  mau  0.950451  mau        7\n",
       "62  mau  0.975161  mau        7\n",
       "64  mau  0.974580  mau        7\n",
       "65  mau  0.538045  bom        6\n",
       "66  mau  0.781062  mau        7\n",
       "73  mau  0.800235  mau        7\n",
       "74  mau  0.977351  mau        9\n",
       "80  mau  0.895569  mau        7\n",
       "84  mau  0.982613  mau        9\n",
       "95  mau  0.998565  mau        7\n",
       "97  mau  0.999024  mau        7"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df['gleason'] = teste.gleason\n",
    "df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Podemos ainda sumarizar os resultados numa tabela, mostrando os acertos e erros do modelo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "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>real</th>\n",
       "      <th>bom</th>\n",
       "      <th>mau</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>pred</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>bom</th>\n",
       "      <td>6</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mau</th>\n",
       "      <td>4</td>\n",
       "      <td>16</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "real  bom  mau\n",
       "pred          \n",
       "bom     6    4\n",
       "mau     4   16"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.crosstab(df.pred,df.real)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Essa tabela é conhecida como **tabela de contingência**. A diagonal principal mostra os acertos do modelo, enquanto a secundária mostra os erros. Os exemplos que eram realmente maus e foram classificados como tal são chamados de *verdadeiros positivos* (TP), já os bons classificados corretamente são chamados de *verdadeiros negativos* (TN). Os bons classificados de maus são chamados de *falsos positivos* (FP) e os maus classificados de bons são *falsos negativos*.\n",
    "\n",
    "A **acurácia** do modelo é a proporção de acertos que o modelo teve. Além dessa métrica, outras duas métricas de qualidade são **sensitividade** e **especificidade** que revelam a capacidade do modelo de classificar corretamente os exemplos que eram positivos (maus) e os que eram negativos (bons).\n",
    "\n",
    "\\begin{align}\n",
    "ACC & = \\frac{TP+TN}{TP+TN+FP+FN}\\\\\n",
    "SEN & = \\frac{TP}{TP+FN}\\\\\n",
    "ESP & = \\frac{TN}{TN+FP}\\\\\n",
    "\\end{align}\n",
    "\n",
    "Dessa forma, podemos avaliar a qualidade do nosso modelo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "acc= 0.733, sen=0.800, esp=0.600\n"
     ]
    }
   ],
   "source": [
    "tabcont = pd.crosstab(df.pred,df.real).values\n",
    "acc = np.sum(tabcont.diagonal())/teste.shape[0]\n",
    "sen = tabcont[1,1]/(tabcont[1,1]+tabcont[0,1])\n",
    "esp = tabcont[0,0]/(tabcont[0,0]+tabcont[1,0])\n",
    "\n",
    "print(\"acc= {:.3f}, sen={:.3f}, esp={:.3f}\".format(acc,sen,esp))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Assim observamos que nosso modelo acertou 73% das classificações, sendo mais preciso com exemplos positivos (mau prognóstico), acertando 80%, que com exemplos negativos (bom prognóstico), acertando 60%. "
   ]
  }
 ],
 "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
}