{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Week11\n",
    "\n",
    "Up until this point in lecture we discussed how to model the conditional distribution of a random variable with a normal distribution---linear regression. \n",
    "Linear regression can be used to understand how a set of covariates $(x_{1},x_{2},\\cdots,x_{p})$ relate linearly to a random variable $Y$. \n",
    "We assume that $Y$ is a normally distributed variable.\n",
    "Because we assume $Y$ can be modeled as a normally distributed random variable we expect that our $Y$ data is some set of negative, positive decimal data (i.e. -2.12, 3.14, 78.2, etc).\n",
    "\n",
    "But what if our $Y$ data is instead a set of $0$s and $1$s that represent the absence and presence of some phenomena? Our plan is to explore logistic regression.\n",
    "\n",
    "The goal of logistic regression is to model a set of $Y$ data with binary responses (Yes/No), (Presence/Absence)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data setup\n",
    "\n",
    "Suppose we have a dataset of $p$ covariates and a single target variable $y$.\n",
    "Then our dataset can be represented as \n",
    "\n",
    "\\begin{align}\n",
    "\\mathcal{D} = [  (x^{1}_{1},x^{2}_{1},\\cdots,x^{p}_{1},y^{1}),(x^{1}_{2},x^{2}_{2},\\cdots,x^{p}_{2},y^{2}), \\cdots, (x^{1}_{N},x^{2}_{N},\\cdots,x^{p}_{N},y^{N}) ]  \n",
    "\\end{align}\n",
    "\n",
    "where $x_{a}^{b}$ corresponds to the $a$th covariate from the $b$th datapoint.\n",
    "Note that our setup above is identical to our setup for multivairate linear regression (MLR). \n",
    "\n",
    "The difference between MLR and logistic regression (LR) is that in LR each $y_{i}$ is either the value 0 or 1. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## A model for $Y$\n",
    "\n",
    "Lets start to model $Y$ by first searching or a random variable that generates the value 0 or 1.\n",
    "We know that a Bernoulli distributed random variable with parameter $\\theta$ will generate the value 1 with probability $\\theta$ and the value 0 with probability $1-\\theta$. \n",
    "\n",
    "That is, if $Y$ is a Bernoulli distributed random variable then \n",
    "\n",
    "\\begin{align}\n",
    "    Y &\\sim \\text{Bern}(\\theta) \\\\ \n",
    "    p(Y=1) &= \\theta \\\\ \n",
    "    p(Y=0) &= 1 - \\theta\\\\\n",
    "    p(Y=y) &= \\theta^{y}(1-\\theta)^{1-y}\n",
    "\\end{align}\n",
    "\n",
    "The expected value of $Y$ is $\\mathbb{E}(Y) = \\theta$ and the variance of $Y$ is $\\mathbb{V}(Y) = \\theta(1-\\theta)$.\n",
    "\n",
    "When we explore MLR, we started by assuming our $Y$ had a normal distriubution $\\mathcal{N}(\\mu,\\sigma^{2})$ and then modified $\\mu$ so that it was a function that depended on parameters $\\beta$ and $x$ data.\n",
    "We chose to modify $\\mu$ because the expected value of $Y$ was $\\mu$. \n",
    "\n",
    "Let us take the same approach with our $Y$ above and model the conditional distribtuion of $Y$ given parameters $\\beta$ and $x$ data as \n",
    "\n",
    "\\begin{align}\n",
    "    Y_{i}|\\beta_{0}, \\beta_{1},x_{i} \\sim \\text{Bern}(\\theta(x))\\\\\n",
    "    \\theta(x) = \\beta_{0} + \\beta_{1} x\\\\\n",
    "\\end{align}\n",
    "\n",
    "But we have a problem. \n",
    "The parameter $\\theta$ is constrained to be a value between 0 and 1, yet the quantity $\\beta_{0} + \\beta_{1} x$ can take any value from negative to positive infinity. \n",
    "We need a way to constrain our values for $\\theta$ to be between 0 and 1. \n",
    "\n",
    "### The logistic function\n",
    "\n",
    "One method to constrain $\\theta$ to be between 0 and 1 is to use the logistic function. \n",
    "The logistic function $f$ is\n",
    "\n",
    "\\begin{align}\n",
    "    f(x) = \\frac{e^{x}}{1+e^{x}}\n",
    "\\end{align}\n",
    "\n",
    "For $x$ values that approach positive infinity the logistic function approaches the value 1. \n",
    "For $x$ values that approach negative infinity the logistic function approaches the value 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAApa0lEQVR4nO3deXhU5fnG8e9D9gUSSNi3sO+iEEDUautSse5Wi1rXurYuVX/VSrVWbdVa69LWutZdK6JiCxVxr1bFStj3HUzYISH7nvf3xww2kAGGmJOTmdyf6zrXzNlm7rkY5sk57znva845RERE9tbG7wAiItIyqUCIiEhIKhAiIhKSCoSIiISkAiEiIiHF+h2gqWRmZrqsrCy/Y4g0sHZ7KQB9O6b4nESkoTlz5uxwznUMtS5qCkRWVhY5OTl+xxBpYOKTswB47arxPicRacjMNuxrnU4xiYhISFFzBCHSUl137AC/I4g0igqEiMeOGpDpdwSRRtEpJhGPLdlUyJJNhX7HEDlonhUIM3vWzLaZ2eJ9rDcz+7OZrTazhWY2qt66i81sVXC62KuMIs3h7ulLuXv6Ur9jiBw0L48gngcm7Gf9ScCA4HQl8DiAmXUAfgOMA8YCvzGz9h7mFBGREDwrEM65T4H8/WxyOvCiC/gSSDezrsCJwPvOuXznXAHwPvsvNCIi4gE/G6m7A7n15vOCy/a1vAEzu5LA0Qe9evXyJqWIiE/q6hzFFTUUVVRTWF5NaWUNZVW1lFTWUFZVQ0llLWWVNWSkJnD+uKb/DYzoq5icc08BTwFkZ2drYAsRabHKq2rZUVLJztIqdpZUsrOkih2lgcfC8mqKyqspqqimqLwmMF9RTUllDeEM2XNYr/SoKxAbgZ715nsEl20EvrvX8n83WyqRJnbLhEF+RxAP1dY5thdXsrmwnM2FFWzaFXjcUljB5sJytgeLQVlVbcj9k+JiSE+OIy0pjnaJcXRLT2Rwl7a0S4oLTImxwcc42ibGkhwfQ0pCbGCKjyE5Ppb4WG9aC/wsENOAa81sMoEG6ULn3GYzexe4t17D9PeBSX6FFPm2Rvfu4HcE+ZYqqmvJKyhj3Y4yNuwsZd2OUjbsLGP9zlK2FFZQU7fnn/mJcW3olpZE1/REsnt3ICMlnozUBDJS48lMjScjJfA8IyWBpPgYnz7VgXlWIMzsVQJHAplmlkfgyqQ4AOfcE8AM4AfAaqAMuDS4Lt/MfgvMDr7U3c65/TV2i7RoczYEvr4qFC1fdW0da7eXsmJrMSu2FLFiSzErthaTV1C+x6metKQ4sjJTGNWrPT07JNE1LYmuaYl0TUuiW3oiaUlxmJl/H6SJWLSMSZ2dne3UWZ+0ROqsr2Wqrq1jxZZi5ufuYkHuLhZtLGTN9hKqawO/ibFtjL4dUxjYuS39O6WSlZFCVmYKWRnJpCfH+5y+6ZjZHOdcdqh1Ed1ILSISrl1lVfx3XT7/XZvP/NwClmwqorKmDoAOKfEc0iON7w7qxOAubRnUpS19O6aQENtyT/80BxUIEYlKRRXVfLlmJ1+uzWfW2p0s31KEc5AQ24ZDeqRx4eG9GdkznUN7ptOjfVJUnBJqaioQIhIVnHOs2V7CR8u38dHybeSsL6CmzpEQ24bsrPbcdPxADu+XwSE90lr9kUG4VCBEJGI555iXu4u3F27m/aVb+Tq/DIDBXdpy5dF9OWZgRw7tla6C0EgqECIeu+PUoX5HiCrOOZZsKmL6wk28vXAzeQXlxMe04agBmVx5dF++N7gT3dOT/I4ZFVQgRDw2rFua3xGiQn5pFVPn5vHa7FxWbSshto1x1IBMbjx+ICcM60y7xDi/I0YdFQgRj322ageggYMao67O8fmaHUyenct7S7ZQXesY1Sude88cwUnDu9A+JXouN22JVCBEPPaXj1YBKhAHo6K6lqlzN/LMZ2tZs72U9OQ4Ljw8i4ljejKoS1u/47UaKhAi0mLsKKnkhS/W8/KXGygoq2Z493Y8PHEkPxjRVQ3NPlCBEBHf7Sip5KlP1/LSrA1U1NRy/JDOXH5UH8b26aD7E3ykAiEivtldGF6ctZ6qmjrOOLQ71xzbn34dU/2OJqhAiIgPKqpreeazdTz28WrKq2s547DuXPu9/vRVYWhRVCBEPHbvWSP8jtBi1NU5/rlgIw/MXMGmwgpOHNaZX04YrMLQQqlAiHhMp0sCFuUVcvs/FrEgr5AR3dN4eOKhjOub4Xcs2Q8VCBGPfbB0KwDHD+3scxJ/lFTW8NB7K3n+i3VkpCbw4DkjOfOw7rRpo8bnlk4FQsRjT/9nLdA6C8T7S7dyxz8Xs6Wogh+P68XNJw4mLUl3PEcKFQgRaXKF5dXcNW0JU+dtZHCXtjx6/ihG925/4B2lRVGBEJEm9fnqHdz8+gK2Fldy/XEDuO7Y/sTFtPE7ljSCCoSINInKmlruf2cFz36+jr6ZKbz50yM4tGe637HkW1CBEJFvLa+gjGtemcuCvEIuHt+bW08aQlK8usaIdCoQIh57eOKhfkfw1IfLtnLTlAXU1TmeuGA0E4Z38TuSNBEVCBGPdYvSwWtq6xwPvreCx/69hmHd2vHYj0fROyPF71jShFQgRDw2fcEmAE4d2c3nJE2npLKGGybP44Nl2zh3TE/uPG0YiXE6pRRtVCBEPPbylxuA6CkQeQVlXP5CDqu2lXD36cO4aHyW35HEIyoQIhK2ORvyueqlOVTW1PHcJWM4emBHvyOJh1QgRCQs7y3ZwrWvzqNbWiKTrxxD/07qYyraqUCIyAFNmZ3LrVMXMqJHOs9dMoYOGgu6VVCBEJF9cs7xxCdruX/mco4e2JHHfzyKlAT9bLQW+pcW8djjF4z2O0KjOOe4d8Yynv7POk4b2Y0/njOS+Fh1mdGaqECIeCwST8c457hr+lKe/2I9F4/vzW9OHabuuVshFQgRj72ekwvAOdk9fU4SnvrF4bKj+nD7yUMwU3FojVQgRDz2xpw8IDIKhHOO30xbwouzNnDFd/rwqx+oOLRmKhAiAgSKw53B4nDl0X2ZdNJgFYdWTi1OIgLAg++t5IXgkYOKg4AKhIgAT3+6lkc/Xs25Y3rqtJJ8QwVCpJWbMjuXe2Ys4+QRXbnnzBEqDvINtUGIeOz5S8f6HWGfZi7ezK1TF/KdAZk8NHEkMbqUVepRgRDxWEsdWW3OhnyunzyfkT3TefLC0STEtsyc4h+dYhLx2Euz1vPSrPV+x9jDhp2lXPHiHLqlJfLMxWNIjtffitKQCoSIx/61cDP/WrjZ7xjfKCit4tLnZuOc47lLx0bknd7SPPRng0grUllTy1UvzSGvoJxXrhhHn0wNESr75ukRhJlNMLMVZrbazG4Nsf5hM5sfnFaa2a5662rrrZvmZU6R1sA5x6Q3F/HV+nz++KORjMnq4HckaeE8O4Iwsxjgr8AJQB4w28ymOeeW7t7GOXdjve2vAw6r9xLlzrlDvcon0to889k6ps7byE0nDOS0KBn+VLzl5RHEWGC1c26tc64KmAycvp/tzwNe9TCPSKv12aod3DtjGROGdeHa7/X3O45ECC/bILoDufXm84BxoTY0s95AH+CjeosTzSwHqAF+75z7R4j9rgSuBOjVq1fTpBZpYq9dNd7X98/NL+PaV+fSv1Mqf/zRSHXbLWFrKVcxnQu84Zyrrbest3MuGzgfeMTM+u29k3PuKedctnMuu2NHDZ4usreyqhqueDGHujrHUxdmk6rR4OQgeFkgNgL1+zfuEVwWyrnsdXrJObcx+LgW+Dd7tk+IRIynPl3DU5+uafb3dc5xyxsLWbm1mL+cP4osXbEkB8nLAjEbGGBmfcwsnkARaHA1kpkNBtoDs+ota29mCcHnmcCRwNK99xWJBB8u28aHy7Y1+/u+8MV6/rVwM784cRDHDNQRthw8z443nXM1ZnYt8C4QAzzrnFtiZncDOc653cXiXGCyc87V230I8KSZ1REoYr+vf/WTiOzfgtxd3DNjGccN7sTVRzc4OysSFk9PSDrnZgAz9lp2x17zd4bY7wtghJfZRKJVYXk11/x9Lp3aJvKgGqXlW1CLlUgUCbQ7LGBLYQVTrh5PerK60ZDGU4EQ8VhiXPP1kvrc5+t5d8lWbj95CKN6tW+295XopAIh4rEXftI840HMz93Ffe8s44ShnbnsqD7N8p4S3VrKfRAi8i2UVNZw/avz6NQ2kT+ePVKjwkmT0BGEiMf+/OEqAK4/boBn73HXtCXkFZTx2lXjSUuO8+x9pHUJq0CY2RFAVv3tnXMvepRJJKp8vnoH4F2BeGfRZl6fk8e13+uvHlqlSR2wQJjZS0A/YD6wuysMB6hAiPhsS2EFk95axCE90vj58d4doUjrFM4RRDYwdK8b2UTEZ3V1jl+8voDK6joemXgocTFqUpSmFc43ajHQxesgInJwnv18HZ+t3sGvTxlK346pfseRKBTOEUQmsNTMvgIqdy90zp3mWSqRKNLeg5vVlm8p4g8zV3D8kM6cN7bngXcQaYRwCsSdXocQiWZPXDi6SV+vuraO/5uygHZJsdz/wxG6pFU8c8AC4Zz7xMw6A2OCi75yzjV/15QiAsBjH69hyaYinrhgNBmpCX7HkSh2wDYIM/sR8BVwDvAj4L9mdrbXwUSixf0zl3P/zOVN8lpLNxXxl49WcerIbkwYrqZB8VY4p5huA8bsPmows47AB8AbXgYTiRZzNxQ0yetU19bxi9cXkJ4cx12nDWuS1xTZn3AKRJu9TintRF10iDS7xz5ew9LNgVNLHVLUS6t4L5wCMdPM3uV/Q4JOZK8xHkTEW7tPLZ2mU0vSjMJppL7ZzH5IYNhPgKecc295G0tEdvvfqaV4nVqSZhVWX0zOuTeBNz3OIhKVuqYlfqv9d59aevLC0bTXqSVpRvssEGb2mXPuKDMrJtD30jerAOeca+d5OpEo8Mi5hzV639Xbinn048CppROH6dSSNK99Fgjn3FHBx7bNF0dEdnPOcdtbi0mOj+WOU4f6HUdaoXDug3gpnGUiEtpd05dw1/QlB73f1Lkb+e+6fG49aTCZuiFOfBBOG8QerWJmFgs0bd8BIlFs6aaig96noLSKe2YsY1SvdCZmq68l8cc+jyDMbFKw/eEQMysKTsXAVuCfzZZQpBW6f+ZyCsuruefMEbRpo76WxB/7LBDOufuC7Q8POOfaBae2zrkM59ykZswo0qrkrM9n8uxcLjuqD0O66loQ8U84d0R/ZWZpu2fMLN3MzvAukkjrVV1bx21vLaZ7ehI3aIQ48Vk4BeI3zrnC3TPOuV3AbzxLJBJl+nZMoW/HlLC2feazdazYWsydpw0jOT6s25REPBNWX0yN3E9EgPvOOiSs7XLzy3jkg5WcMLQzJwzt7HEqkQML5wgix8weMrN+wekhYI7XwURaE+ccd05bQhsz7lR3GtJChFMgrgOqgNeCUyVwjZehRKLJpKkLmTR14X63eW/pVj5cvo0bjx9I9/SkZkomsn/hdNZXCtzaDFlEotLa7aX7XV9SWcOd05YwuEtbLjkyq3lCiYThgAXCzAYCvwCy6m/vnDvWu1girccj769kc2EFj54/irgYDbUiLUc4jc2vA08AfwNqvY0j0ros2VTIc1+s57yxvRjdu73fcUT2EE6BqHHOPe55EpFWpq4u0BlfelIct04Y7HcckQbCKRDTzexnwFsEGqgBcM7le5ZKJIoM7Rb6bui/f/U183N38fDEkaQlxzVzKpEDC6dAXBx8vLneMgf0bfo4ItHnN6c2vGx1e3El989czvi+GZxxaHcfUokcWDhXMfVpjiAirck9by+lsrqO3505HDN1xictUzhXMV0Uarlz7sWmjyMSfW6YPA/438hyn6/ewT/mb+L6Y/vTr2Oqn9FE9iucU0xj6j1PBI4D5gIqECJh2FxY8c3ziupabv/HYnpnJPOz7/X3MZXIgYVzium6+vNmlg5M9iqQSDR74pM1rNtRyos/GUtiXIzfcUT2qzF35ZQCapcQOUjrdpTy2MdrOHVkN44e2NHvOCIHFE4bxHQCVy1BoKAMBaZ4GUokGv36H4tJiG3Dr08e4ncUkbDsb8jR3aOk/xF4MDjdBxztnAurbyYzm2BmK8xstZk12MfMLjGz7WY2PzhdXm/dxWa2KjhdvPe+IpFiVO/2tE2M5bPVO7hlwiA6tUv0O5JIWPZ3BDELGAVc7py78GBf2MxigL8CJwB5wGwzm+acW7rXpq85567da98OBAYlyiZw9DInuG/BweYQ8dvVx/TjuAc/YWSPNM4f19vvOCJh21+BiDez84EjzOysvVc656Ye4LXHAqudc2sBzGwycDqwd4EI5UTg/d13a5vZ+8AE4NV97bB2eykTn5y1x7JTDunKheOzKK+q5ZLnvmqwz9mje3BOdk/yS6v46csNh7i44PDenDqyG5t2lXPja/MbrL/iO305fmhn1mwv4VdTFzVYf92xAzhqQCZLNhVy9/SGH/uWCYMY3bsDczbk84eZKxqsv+PUoQzrlsZnq3bwl49WNVh/71kj6NcxlQ+WbuXp/6xtsP7hiYfSLT2J6Qs28fKXGxqsf/yC0XRIief1nFzemJPXYP3zl44lKT6Gl2at518LNzdY/9pV4wF46tM1fLhs2x7rEuNieOEnYwH484er+Hz1jj3Wt0+O54kLRwNw/8zlzN2wZ+3vmpb4zWWhd01fwtJNRXus79sx5ZuBeCZNXdigx9Sh3dp9c4PaDZPn7XElEQT+qv9lsHuLq1+aQ0FZ1R7rj+yfyfXHBYb8vPjZr6io3rMbsuOGdOLKo/sBNPjewZ7fvRMe+oQdJZV0SUvg/Ke/BPTd03eveb57jfndq29/BeJq4MdAOnDqXusccKAC0R3IrTefB4wLsd0PzexoYCVwo3Mudx/7Nrjd1MyuBK4ESO3a7wBxRJrfgtwCthVXkhDbhhQNISoRxpxz+9/A7DLn3DMH/cJmZwMTnHOXB+cvBMbVP51kZhlAiXOu0syuAiY65441s18Aic653wW3+zVQ7pz7477eLzs72+Xk5BxsTBHP1NTWcdqjn7NqWzEje6bzxtVH+B1JpAEzm+Ocyw617oCXuTamOARtBHrWm+8RXFb/tXc653Z3APg3YHS4+4q0dM9/sZ6lm4vIykghRt1pSATycnSS2cAAM+tjZvHAucC0+huYWdd6s6cBy4LP3wW+b2btzaw98P3gMpGIsGlXOQ+9v5LvDepIh5R4v+OINIpnBcI5VwNcS+CHfRkwxTm3xMzuNrPTgptdb2ZLzGwBcD1wSXDffOC3BIrMbOBudS8ukeTu6Uupc467Tx/udxSRRgvnRrkzgY+cc4XB+XTgu865fxxoX+fcDGDGXsvuqPd8EjBpH/s+Czx7oPcQaWk+XLaVmUu2cMuEQfTskMyR/TP9jiTSKOE0Us93zh2617J5zrnDvAx2sNRILS1BWVUNJzz0KcnxMbx9/XeIj9UY09Ky7a+ROpzr7kJ9w3W9nkgIf/5wNRt3lfPalYerOEjEC+cbnGNmD5lZv+D0ELD/uytEWqEVW4r523/Wcs7oHozrm/HN8ouf/YqLn214w5JISxdOgbgOqAJeC06VwDVehhKJNHV1jtv/sYi2ibFM+sGenfFVVNc2uBtWJBKEMx5EKRBW53wirdXrc3KZvb6AP5x9iC5rlaixzwJhZo84527Yq7vvbzjnTguxm0irs6OkkntnLGdsVgfOHtXD7zgiTWZ/RxAvBR/32b2FiMC9by+jrKqGe84cTps2umNaosc+C4RzbndD9KHOuT/VX2dmPwc+8TKYSCT4fPUOps7byHXH9mdA57YhtzluSKdmTiXSNMJppA41WM8lTZxDJOJUVNdy21uLyMpI5prv9d/ndlce3e+b7plFIsn+2iDOA84H+phZ/T6U2gHq9kJavb9+vJr1O8t4+bJxJMbF+B1HpMntrw3iC2AzkElguNHdioGFXoYSaelWbS3miU/WcOZh3TlqwP670tg9qMvuQW5EIsX+2iA2ABvM7HgCYzHUmdlAYDDQcAgrkVairs7xq7cWkZIQy20nDznwDiIRKpw2iE+BRDPrDrwHXAg872UokZZs9z0PvzppCJmpCX7HEfFMOAXCnHNlwFnAY865c4Bh3sYSaZm+ueehTwfOydY9DxLdwioQZjaewPjUbweXqUVOWqV7gvc83HvmcEyjxEmUC6dX1hsIjNnwVnDAn77Ax56mEmmB/rNqO2/N28j1x/anf6fQ9zyEcsohXQ+8kUgLdMDxICKFxoMQL5VW1nDiI58SF9OGd37+HV3WKlGjUeNBqC8mkf954N0V5BWUM+Wq8QddHMqrAj25JsWrqEhkUV9MIgeQsz6fF2at56LxvRnbp8NB73/Jc4GxIHQfhESaA/bF5JxTn0vSalVU13LLmwvplpbELRMG+x1HpFkdsJHazBbR8BRTIZAD/M45t9OLYCItwZ8+XMXa7aW8+JOxpCZopF1pXcL5xr8D1AJ/D86fCyQDWwjcMHeqJ8lEfLZ4YyFPfRoYQvTogR39jiPS7MIpEMc750bVm19kZnOdc6PM7AKvgon4qbq2jpvfWEiHlHhuP3mo33FEfBFOgYgxs7HOua8AzGwM/7tRrsazZCI+euLfa1i2uYgnLxxNWnLct3qts0frjmuJTOEUiMuBZ80sFTCgCLjMzFKA+7wMJ+KHZZuL+PNHqzj5kK6cOKzLt369c7J7NkEqkeZ3wALhnJsNjDCztOB8Yb3VU7wKJuKHyppabnxtPmlJ8fz29OFN8pr5pVUAdEiJb5LXE2kuB+yLyczSzOwh4EPgQzN7cHexEIk2f/pgFcu3FPP7s0Y02Q/6T1+ew09fnnPgDUVamHA663uWwCBBPwpORcBzXoYS8cOcDfk88ckaJmb35Pihnf2OI+K7cNog+jnnflhv/i4zm+9RHhFflFbWcNOUBXRNS+L2UzQIkAiEdwRRbmZH7Z4xsyOBcu8iiTS/+95Zxtf5ZTz4o5G0Tfx2Vy2JRItwjiCuBl6s1+5QAFzsXSSR5vXJyu28/OXXXH5UHw7vm+F3HJEWI5yrmBYAI82sXXC+yMxuABZ6nE3Ec/mlVdz8+gL6d0rlFycO8uQ9Lji8tyevK+K1sDuXcc4V1Zu9CXikydOINCPnHDe/voBdZdU8d+kYz8Z4OHVkN09eV8Rr4bRBhKKxFiXiPf/Fej5cvo1JPxjMsG7eXbm9aVc5m3ap2U4iT2O7p4yOYeik1VqyqZD7ZiznuMGduOSILE/f68bX5gMaD0Iiz/5GlCsmdCEwIMmzRCIeK6uq4bpX55GeHMcD54zETAfEIqHsb8Cg8EdlF4kgd01byrodpbxy2Th1fyGyH41tgxCJSP+cv5HXcnL52Xf7cUT/TL/jiLRoKhDSaqzcWsytby5iTFZ7bjh+oN9xRFo8jaEorUJxRTVXvzyHlIRYHj1/FHExzfe30RXf6dts7yXSlFQgJOo557jljYVs2FnGK5ePo3O7xGZ9f3X8J5HK0z+jzGyCma0ws9VmdmuI9TeZ2VIzW2hmH5pZ73rras1sfnCa5mVOiW7PfLaOdxZv4ZcTBvnSlcaa7SWs2V7S7O8r8m15dgRhZjHAX4ETgDxgtplNc84trbfZPCDbOVdmZj8F/gBMDK4rd84d6lU+aR3+u3Yn972znAnDuvh2qudXUxcBug9CIo+XRxBjgdXOubXOuSpgMnB6/Q2ccx8758qCs18CGrxXmkxeQRk/e2UuvTsk88A5h+h+B5GD5GWB6A7k1pvPCy7bl8uAd+rNJ5pZjpl9aWZnhNrBzK4MbpOzffv2bx1YokdpZQ2Xv5BDVW0dT12UrS68RRqhRTRSm9kFQDZwTL3FvZ1zG82sL/CRmS1yzq2pv59z7ingKYDs7Gx1/yEA1NU5bnxtPiu3FvPcpWPp3ynV70giEcnLI4iNQM968z2Cy/ZgZscDtwGnOecqdy93zm0MPq4F/g0c5mFWiSIPvb+S95Zu5faTh3LMwI5+xxGJWF4eQcwGBphZHwKF4Vzg/PobmNlhwJPABOfctnrL2wNlzrlKM8sEjiTQgC2yX/+cv5FHP17NuWN6cumRWX7HAeC6Ywf4HUGkUTwrEM65GjO7FngXiAGedc4tMbO7gRzn3DTgASAVeD3YgPi1c+40YAjwpJnVETjK+f1eVz+JNPDFmh3c/PpCxvbpwN2nD28xjdJHDVCXHhKZzLnoOHWfnZ3tcnJy/I4hPlm+pYhzHp9Fl7RE3rj6CNKSW06j9JJNhQCejjkh0lhmNsc5lx1qnfpikoi3aVc5lzw7m+SEGJ7/ydgWVRwA7p6+lLun6wBYIo8KhES0wvJqLnnuK0ora3j+0rF0T9dQJSJNpUVc5irSGGVVNVz2/GzW7SjlhUvHMqRrO78jiUQVHUFIRKqoruXyF3KY+3UBj0w8TGM7iHhARxAScapq6vjpy3OYtXYnD54zkpMP6ep3JJGopAIhEaWmto6fT57Hxyu2c8+ZwzlrVMvvvuuWCYP8jiDSKCoQEjGqa+u4acoC3lm8hdtPHsKPx/U+8E4twOjeHfyOINIoKhASESprarnu7/N4b+lWbj1pMJdH0ChtczbkAyoUEnlUIKTFq6iu5aqX5vDJyu3ceepQLjmyj9+RDsofZq4ANB6ERB4VCGnRdnfb/eW6ndx31gjOG9vL70girYYKhLRY24srueyF2SzeWMhDPxrJmYe1/AZpkWiiAiEt0trtJVzy3Gy2FVfw5IXZnDC0s9+RRFodFQhpceZ+XcBlz8/GzHj1isM5rFd7vyOJtEoqENKizFi0mRtfm0/XtESev3QsWZkpfkf61u44dajfEUQaRQVCWoTaOsdD76/grx+vYVSvdJ6+KJuM1AS/YzUJdfMtkUoFQnxXVFHNDZPn89HybZw3tid3njaMhNgYv2M1mc9W7QA0cJBEHhUI8dXKrcVc/dIcvs4v47dnDOeCcb1azEhwTeUvH60CVCAk8qhAiC+cc0yenctd05eQmhDLK5ePY1zfDL9jiUg9KhDS7Ioqqpk0dRFvL9zMUf0zeWjiSDq1TfQ7lojsRQVCmlXO+nxunDKfTbsquGXCIK4+uh9t2kTXKSWRaKECIc2ivKqWB95dwXNfrKN7ehJTrhrP6N66v0GkJVOBEM99tS6fm99YwIadZVw0vje/nDCYlITW89W796wRfkcQaZTW879Uml1+aRUPvLucybNz6dE+iVevOJzx/VpfQ3S/jql+RxBpFBUIaXK1dY7Js7/mgXdXUFxRw2VH9uGm7w8kOb51ft0+WLoVgOPVn5REmNb5P1Y889W6fH739lIW5hUyrk8H7j59OIO6tPU7lq+e/s9aQAVCIo8KhDSJFVuK+cPM5Xy4fBud2yXwp3MP5bSR3aLupjeR1kQFQr6V3PwyHvlgFVPn5ZGaEMstEwZx6RF9SIqPnq4yRForFQhplNXbSnjs36v55/xNxLQxrvhOX356TD/ap8T7HU1EmogKhByUhXm7ePzfa5i5ZAsJsW24eHwWVxzdh65pSX5HE5EmpgIhB1RVU8c7izfzwhfrmfv1LtomxnLNd/tz6ZFZUdMlt5cennio3xFEGkUFQvYpr6CMKbNz+ftXuewoqaRPZgp3nDKUs7N70C4xzu94EaNbuo6uJDKpQMgeSitrmLFoM2/OzePLtfmYwbGDOnHREVl8p3+m+k1qhOkLNgFw6shuPicROTgqEEJFdS2frtzOO4u3MHPxFsqra+mdkcxNJwzkzMO607NDst8RI9rLX24AVCAk8qhAtFIllTV8vHwbMxdv4eMV2yirqiUtKY4zDuvGD0f1YHTv9rqHQaSVU4FoJZxzLN1cxKcrd/DJym3M2VBAda0jMzWeMw7rzknDu3B43wziYtr4HVVEWggViCjlnGPN9hJmry9g9rp8/rN6B9uLKwEY0rUdPzmqD8cO6kR2Vgdi1K4gIiGoQESJsqoalm0uYu6GXXy1Pp+c9fkUlFUDkJESzxH9MzlmYEeOHpBJp3YavU1EDkwFIgIVV1SzdFMRizcVsWRjIYs2FrJmewl1LrA+KyOZ44Z0ZmxWB7Kz2tMnM0XtCT56/ILRfkcQaRQViBZsZ0klq7eVsGpbCavrTVuKKr7ZpnO7BIZ3S+MHI7oyvHsaI3umaXznFqaDuh+RCKUC4aPyqlryCsrILSgjN7+c3Px6zwvKKK6o+Wbb5PgY+ndK5Yh+GfTrlMrQru0Y1r2dikEEeD0nF4Bzsnv6nETk4KhANLG6OkdRRTX5pVXsLK1ia1EF24oq2Vpcwfbg47aiSrYWVVBUrwAAJMa1oUf7ZHq2TyI7qz29M1Lo3ymV/p1S6douUTepRag35uQBKhASeTwtEGY2AfgTEAP8zTn3+73WJwAvAqOBncBE59z64LpJwGVALXC9c+5dL7Pu5pyjorqO4spqSipqKK2spbiymtLKWkoqqymprKWovJpdZVUUlAUe80ur2FVWTUFZFYXl1d+0BdQXH9OGjm0T6NwugX4dA0cCndol0qN9UqAodEiiY2qC2gpEpMXwrECYWQzwV+AEIA+YbWbTnHNL6212GVDgnOtvZucC9wMTzWwocC4wDOgGfGBmA51ztU2dc2dJJec9/SUlFTWUVAamUD/we0uMa0P75PjAlBJH1/Qk2ifH0T45nvTkeNonx5GRGigIndsmkp4cpx9/EYkoXh5BjAVWO+fWApjZZOB0oH6BOB24M/j8DeBRC/yKng5Mds5VAuvMbHXw9WY1dcik+Bj6dUwlJSGW1N1TYiwpCbG0TYhtsDw1IZa2ibEkxmlAHBGJbl4WiO5Abr35PGDcvrZxztWYWSGQEVz+5V77dt/7DczsSuBKgF69ejUqZHJ8rC5DFBEJIaIbqZ1zTwFPAWRnZ4dxYkik+T1/6Vi/I4g0ipcd72wE6l+20SO4LOQ2ZhYLpBForA5nX5GIkBQfozG6JSJ5WSBmAwPMrI+ZxRNodJ621zbTgIuDz88GPnLOueDyc80swcz6AAOArzzMKuKZl2at56VZ6/2OIXLQPDvFFGxTuBZ4l8Blrs8655aY2d1AjnNuGvAM8FKwETqfQBEhuN0UAg3aNcA1XlzBJNIc/rVwMwAXjs/yN4jIQfK0DcI5NwOYsdeyO+o9rwDO2ce+9wD3eJlPRET2TZ3/i4hISCoQIiISkgqEiIiEZIGLhiKfmW0HNvidoxEygR1+h2hm+sytgz5zZOjtnOsYakXUFIhIZWY5zrlsv3M0J33m1kGfOfLpFJOIiISkAiEiIiGpQPjvKb8D+ECfuXXQZ45waoMQEZGQdAQhIiIhqUCIiEhIKhAtiJn9n5k5M8v0O4vXzOwBM1tuZgvN7C0zS/c7kxfMbIKZrTCz1WZ2q995vGZmPc3sYzNbamZLzOznfmdqLmYWY2bzzOxffmdpKioQLYSZ9QS+D3ztd5Zm8j4w3Dl3CLASmORzniZXb1z2k4ChwHnB8dajWQ3wf865ocDhwDWt4DPv9nNgmd8hmpIKRMvxMHAL0CquGnDOveecqwnOfklgUKho88247M65KmD3uOxRyzm32Tk3N/i8mMAPZoPhgqONmfUATgb+5neWpqQC0QKY2enARufcAr+z+OQnwDt+h/BAqHHZo/7HcjczywIOA/7rc5Tm8AiBP/DqfM7RpCJ6TOpIYmYfAF1CrLoN+BWB00tRZX+f2Tn3z+A2txE4LfFKc2YTb5lZKvAmcINzrsjvPF4ys1OAbc65OWb2XZ/jNCkViGbinDs+1HIzGwH0ARaYGQROtcw1s7HOuS3NGLHJ7esz72ZmlwCnAMe56Lwhp1WOrW5mcQSKwyvOual+52kGRwKnmdkPgESgnZm97Jy7wOdc35pulGthzGw9kO2ci7QeIQ+KmU0AHgKOcc5t9zuPF8wslkAD/HEECsNs4Hzn3BJfg3nIAn/lvADkO+du8DlOswseQfzCOXeKz1GahNogxC+PAm2B981svpk94XegphZshN89LvsyYEo0F4egI4ELgWOD/67zg39ZSwTSEYSIiISkIwgREQlJBUJEREJSgRARkZBUIEREJCQVCBERCUkFQkREQlKBEBGRkFQgRDxiZmOC410kmllKcHyE4X7nEgmXbpQT8ZCZ/Y5A/zxJQJ5z7j6fI4mETQVCxENmFk+gD6YK4AjnXK3PkUTCplNMIt7KAFIJ9DuV6HMWkYOiIwgRD5nZNAIjyfUBujrnrvU5kkjYNB6EiEfM7CKg2jn39+D41F+Y2bHOuY/8ziYSDh1BiIhISGqDEBGRkFQgREQkJBUIEREJSQVCRERCUoEQEZGQVCBERCQkFQgREQnp/wELcPW7vIDv9AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def logistic(x):\n",
    "    import numpy as np\n",
    "    e = np.exp(x)\n",
    "    return e/(1+e)\n",
    "\n",
    "fig,ax = plt.subplots()\n",
    "domain = np.linspace(-5,5,10**3)\n",
    "\n",
    "ax.plot(domain,logistic(domain))\n",
    "\n",
    "ax.set_xlabel(\"x\")\n",
    "ax.set_ylabel(\"Logistic function\")\n",
    "\n",
    "ax.axvline(0,ls=\"--\")\n",
    "ax.axhline(0.5,ls=\"--\")\n",
    "\n",
    "ax.set_yticks([0,0.25,0.50,0.75,1.0])\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Because logistic function maps values from negative infinity to positive infintiy to numbers between 0 and 1, we can use this function to give us valid $\\theta$ values. \n",
    "\n",
    "\\begin{align}\n",
    "    Y_{i}|\\beta_{0}, \\beta_{1},x_{i} \\sim \\text{Bern}(\\theta(x))\\\\\n",
    "    \\theta(x) = L(\\beta_{0} + \\beta_{1} x)\\\\\n",
    "\\end{align}\n",
    "\n",
    "where $L$ is the logistic functon. In other words,\n",
    "\n",
    "\\begin{align}\n",
    "    Y_{i}|\\beta_{0}, \\beta_{1},x_{i} \\sim \\text{Bern}(\\theta(x))\\\\\n",
    "    \\theta(x) = \\frac{e^{\\beta_{0} + \\beta_{1} x}}{1+e^{\\beta_{0} + \\beta_{1} x}}\\\\\n",
    "\\end{align}\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## A way to write logistic regression that looks more familiar \n",
    "\n",
    "At times it can be convienant to rewrite our logistic regression model above to look more similar to the more familar multivariate linear regression. For multivariate linear regression, the expedcted value of our target variable was equal to $\\beta_{0} + \\beta_{1} x$. \n",
    "\n",
    "or logistic regression our expected value is \n",
    "\n",
    "\\begin{align}\n",
    "    \\mathbb{E}(Y) = \\frac{e^{\\beta_{0} + \\beta_{1} x}}{1+e^{\\beta_{0} + \\beta_{1} x}}\n",
    "\\end{align}\n",
    "\n",
    "Let's rearrange the above equation to isolate $\\beta_{0} + \\beta_{1} x$. \n",
    "\n",
    "\\begin{align}\n",
    "    \\theta = \\mathbb{E}(Y) &= \\frac{e^{\\beta_{0} + \\beta_{1} x}}{1+e^{\\beta_{0} + \\beta_{1} x}} \\\\ \n",
    "    1+e^{\\beta_{0} + \\beta_{1} x}(\\theta) &= e^{\\beta_{0} + \\beta_{1} x} \\\\ \n",
    "    \\theta + \\theta e^{\\beta_{0} + \\beta_{1} x} &= e^{\\beta_{0} + \\beta_{1} x} \\\\ \n",
    "    \\theta  &= e^{\\beta_{0} + \\beta_{1} x} - \\theta e^{\\beta_{0} + \\beta_{1} x} \\\\ \n",
    "    \\theta  &= e^{\\beta_{0} + \\beta_{1} x}(1-\\theta) \\\\ \n",
    "    \\frac{\\theta}{1-\\theta}  &= e^{\\beta_{0} + \\beta_{1} x} \\\\ \n",
    "    \\log \\left(\\frac{\\theta}{1-\\theta}\\right) &= \\beta_{0} + \\beta_{1} x\n",
    "\\end{align}\n",
    "\n",
    "The ratio $\\frac{\\theta}{1-\\theta}$ is called the **odds** and so the quantity  $\\log \\left(\\frac{\\theta}{1-\\theta}\\right)$ is called the **log odds**. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## A 1-unit change in $X$\n",
    "\n",
    "Just like in MLR, we can look at how a one unit change in an x covariate will chaneg our expected value. \n",
    "Lets take the difference between \n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left(\\frac{\\theta}{1-\\theta}\\right) &= \\beta_{0} + \\beta_{1} x\n",
    "\\end{align}\n",
    "\n",
    "and \n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left(\\frac{\\theta^{*}}{1-\\theta^{*}}\\right) &= \\beta_{0} + \\beta_{1} (x+1)\n",
    "\\end{align}\n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left(\\frac{\\theta^{*}}{1-\\theta^{*}}\\right) - \\log \\left(\\frac{\\theta}{1-\\theta}\\right) &= \\beta_{0} + \\beta_{1} (x+1) - (\\beta_{0} + \\beta_{1} x) \\\\ \n",
    "    \\log \\left(\\frac{\\theta^{*}}{1-\\theta^{*}}\\right) - \\log \\left(\\frac{\\theta}{1-\\theta}\\right) &= \\beta_{0} + \\beta_{1} (x+1) - \\beta_{0} - \\beta_{1} x  \\\\\n",
    "    \\log \\left(\\frac{\\theta^{*}}{1-\\theta^{*}}\\right) - \\log \\left(\\frac{\\theta}{1-\\theta}\\right) &= \\beta_{1} \\\\\n",
    "\\end{align}\n",
    "\n",
    "We can simplify the left hand side using the logarithm rule\n",
    "\\begin{align}\n",
    "    \\log\\left(\\frac{x}{y}\\right) = \\log(x) - \\log(y)\n",
    "\\end{align}\n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left(\\frac{\\theta^{*}}{1-\\theta^{*}}\\right) - \\log \\left(\\frac{\\theta}{1-\\theta}\\right) &= \\beta_{1} \\\\\n",
    "    \\log \\left(\\frac{ \\frac{\\theta^{*}}{1-\\theta^{*}} }{ \\frac{\\theta}{1-\\theta}}\\right) &= \\beta_{1} \\\\\n",
    "\\end{align}\n",
    "\n",
    "The expression in side the logarithm on the left hand side of the equals is called the odds ratio. \n",
    "The logarithm of the odds ration is often called the log odds ratio. \n",
    "To be clear, the log odds ratio is related to a one unit change in x.\n",
    "\n",
    "We can also estimate the odds ratio by exponentiating both sides of the above equation\n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left(\\frac{ \\frac{\\theta^{*}}{1-\\theta^{*}} }{ \\frac{\\theta}{1-\\theta}}\\right) &= \\beta_{1} \\\\\n",
    "    \\exp \\left(  \\log \\left(\\frac{ \\frac{\\theta^{*}}{1-\\theta^{*}} }{ \\frac{\\theta}{1-\\theta}}\\right) \\right) = e^{\\beta_{1}}\\\\\n",
    "    \\frac{ \\frac{\\theta^{*}}{1-\\theta^{*}} }{ \\frac{\\theta}{1-\\theta}} =  e^{\\beta_{1}}\n",
    "\\end{align}\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### An example \n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 148,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/usr/local/lib/python3.9/site-packages/statsmodels/tsa/tsatools.py:142: FutureWarning: In a future version of pandas all arguments of concat except for the argument 'objs' will be keyword-only\n",
      "  x = pd.concat(x[::order], 1)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "const                         1.000000\n",
      "Pregnancies                   3.845052\n",
      "Glucose                     120.894531\n",
      "BloodPressure                69.105469\n",
      "SkinThickness                20.536458\n",
      "Insulin                      79.799479\n",
      "BMI                          31.992578\n",
      "DiabetesPedigreeFunction      0.471876\n",
      "Age                          33.240885\n",
      "dtype: float64\n",
      "Optimization terminated successfully.\n",
      "         Current function value: 0.470993\n",
      "         Iterations 6\n",
      "                           Logit Regression Results                           \n",
      "==============================================================================\n",
      "Dep. Variable:                Outcome   No. Observations:                  768\n",
      "Model:                          Logit   Df Residuals:                      759\n",
      "Method:                           MLE   Df Model:                            8\n",
      "Date:                Wed, 10 Nov 2021   Pseudo R-squ.:                  0.2718\n",
      "Time:                        22:56:30   Log-Likelihood:                -361.72\n",
      "converged:                       True   LL-Null:                       -496.74\n",
      "Covariance Type:            nonrobust   LLR p-value:                 9.652e-54\n",
      "============================================================================================\n",
      "                               coef    std err          z      P>|z|      [0.025      0.975]\n",
      "--------------------------------------------------------------------------------------------\n",
      "const                       -8.4047      0.717    -11.728      0.000      -9.809      -7.000\n",
      "Pregnancies                  0.1232      0.032      3.840      0.000       0.060       0.186\n",
      "Glucose                      0.0352      0.004      9.481      0.000       0.028       0.042\n",
      "BloodPressure               -0.0133      0.005     -2.540      0.011      -0.024      -0.003\n",
      "SkinThickness                0.0006      0.007      0.090      0.929      -0.013       0.014\n",
      "Insulin                     -0.0012      0.001     -1.322      0.186      -0.003       0.001\n",
      "BMI                          0.0897      0.015      5.945      0.000       0.060       0.119\n",
      "DiabetesPedigreeFunction     0.9452      0.299      3.160      0.002       0.359       1.531\n",
      "Age                          0.0149      0.009      1.593      0.111      -0.003       0.033\n",
      "============================================================================================\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1gAAAFgCAYAAACmKdhBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABgG0lEQVR4nO3dd3xUZdr/8c+VDgkklNAJoVfp1QZ2dF2xrGvZta+4Rbf67LrN3cftq9ueXV1F17WLZVdFRbGBYqV3CISWQiAJISGQPnP//sjgL8sGCTAzZ8r3/XrlxeTMIfOdyeTMuc65z3Wbcw4RERERERE5cQleBxAREREREYkVKrBERERERESCRAWWiIiIiIhIkKjAEhERERERCRIVWCIiIiIiIkGS5NUDd+3a1eXm5nr18CISIsuXLy93zmV7neNItO0RiT3a7oiIF4607fGswMrNzWXZsmVePbyIhIiZ7fQ6w2fRtkck9mi7IyJeONK2R0MERUREREREgkQFloiIiIiISJCowBIREREREQkSFVgiIiIiIiJBogJLREREREQkSFRgiYiIiIiIBIkKLBERERERkSBRgSUiIiIiIhIkKrBERERERESCRAWWiIiIiIhIkKjAEhEREQHM7GEzKzWzdUe438zs/8ws38zWmNn4cGcUkcinAktEPpNzzusIIiLh8ggw8zPuPx8YHPiaDfw9DJlEJMokeR1ARCLXwk2l/N87W3j0xsl0TEv2Oo6IxDjnHB9t28vTSwoZ3TuTm08fEO7Hf8/Mcj9jlVnAY675yNPHZpZlZj2dcyXhSSgSv5xz1DT4OFDfRHVdEwfrmzjY0ERtg4/aRh+1DT7qmvzUN/qob/JT3+SnoclPo6/ll6PJ56fR7/D5HE1+h8/vx+fA73c0+f38+Ypx9MhMO6GsKrBE5L/4/Y77FuXzhzc3M7xHRw7WN6nAEpGQKT9Qz7+WF/H0kgJ27K2hY1oSw3p08DpWa3oDhS2+Lwos+68Cy8xm03yWi5ycnLCEE4k2fr+j/EA9RZW17K6qY3dVHWUH6imvrmfvwQYqDjZQWdNAZW0j+2sb8R/joJqUxASSE43kpASSEhJISTSSEhNISjSSEozEhASSEoyEhMD3ZviDMHJHBZaI/IeD9U1899lVLFi/h4vH9uI3l46mXUqi17FEJMb4/c1nq55aUsAb63fT6HNMyu3EN88azAUn9SQtObq3O865OcAcgIkTJ2qstcS1ukYfeburydtTzZY91WwvP8iOvTUUVNTQ0OT/j3WTEoyuGal0Tk+hS0YKfTu3J6tdMpntkumQlkRGWhIZqc1f7VOSaJ+SSLuURNolJ5KanEBaciIpiQmkJiVgZp48XxVYIvKp3VV13PToUjaW7OenF47gxlNyPds4iUhsKq2u4/nlRcxdUkhBRQ1Z7ZO5dlouV07qy+DuEXnWqqVioG+L7/sElolIgHOOgooalmyvYEXBPlYWVLKl9AC+wOmnlKQE+ndJZ2B2OmcO60bfTu3oldWOHplp9OiYRqf2KSQkRPe+hwosEQFg/a4qbnpkGdV1jfzj+kmcMbSb15FEJEb4/Y7388t5ekkBb27YQ5PfMaV/Z7537hDOG9kjms5WzQNuNbO5wBSgStdfiUBtg48P8st5e9Me3ttcTnFlLQAd05IYl9OJs4d3Z2Svjgzr2ZG+ndqRlBjbffZUYIkIb2/cw21PrySzXTLPffVkRvTq6HUkEYkBtQ0+nl9RxMPvb2d7+UE6tU/mhlNyuWJSDoO6ZXgd77+Y2dPADKCrmRUBPwOSAZxz9wPzgQuAfKAGuMGbpCLea/T5eTevjJdW7+LNDbupa/STkZrEKYO6cMv0AUwb0IWB2RlRfzbqeKjAEolzzy4t5I5/r2Fkr0z+cd1EunU8sc45x8rMZgJ/ARKBh5xzvz3s/n7Aw0A2UAF82TlXFLjv98DnaJ5y4k3gW0595UU8V1pdx2Mf7uSJT3ZSWdPImD6Z/OXKscwc1YPUpMg9W+Wcu+oo9zvgG2GKIxKRdlfV8dQnO5m7tJDS6no6tU/msvF9OH9UTyb370xKUmyfnWoLFVgicWzOe1v59fxNnDa4Kw9cM4H2KeHdJJhZInAvcA7N3biWmtk859yGFqvdQ3Nb5EfN7EzgN8A1ZnYycAowOrDe+8B0YFG48ovIf8rbXc1Di7fx0qpdNPr9nD28OzefNoBJuZ10PadIlNtadoAH3t3KCyuLafI7zhjajasn5zB9aDbJMT7k71ipwBKJQ8457l6Qx32LtvK50T350xfHenXEaTKQ75zbBhC4rmEW0LLAGgF8N3B7IfBi4LYD0oAUwGgexrMn9JFFpCXnmq+venDxdt7bXEZacgJXTOrLjaf2p3/XdK/jicgJKt1fx5/e2swzSwtJTkzg6sk5fOW0AfTt3N7raBFLBZZInPH7HT95aR1PfVLAVZNz+OXFo0j0bnx0a3PKTDlsndXApTQPI7wE6GBmXZxzH5nZQprnnzHgb865ja09iOajEQm++iYf81bt4h/vb2fT7mq6ZqRy+7lD+NKUfnRKT/E6noicoCafn4c/2M6f3txCk9/PtdNyufXMQXTNSPU6WsRTgSUSR/x+xw//vZZnlhXy1ekD+cHModEwbOd24G9mdj3wHs0tkX1mNggYTnObZIA3zew059ziw3+A5qMRCZ7aBh9PfLyTOYu3UVZdz9DuHfj9F0Yza2yviL6+SkTabv2uKv7nuTVsKNnP2cO78dMLR9Cvi85It5UKLJE44fc7fvRCc3F125mD+O45QyKhuDrqnDLOuV00n8HCzDKAy5xzlWZ2M/Cxc+5A4L7XgGnAfxVYInLi6hqbC6v7391G+YF6ThnUhXsuH8Ppg7tGwrZERILAOcfDH+zgd69tIqt9Mvd/eTznjeyhv/FjpAJLJA74/Y4fv7iWuUsLufWMiCmuAJYCg82sP82F1ZXA1S1XMLOuQIVzzg/8kOaOggAFwM1m9huahwhOB/4cptwicaOu0cfTSwq4b9FWyqrrOXlgF/7+5fFMyu3sdTQRCaID9U18e+4q3tq4h7OHd+f3XxhNZw33PS4qsERinHPN11w9vaSQb5wxkO+dGzHFFc65JjO7FVhAc5v2h51z683sLmCZc24ezXPS/MbMHM1DBA+1SH4eOBNYS3PDi9edcy+H+zmIxKq6Rh/PLC3kvkX57Nlfz5T+nfnrVeOYOqCL19FEJMgKK2r4yqPLyC87wM8+P4LrT86NmH2FaKQCSyTG/fa1TTz1SQFfmzGQ28+NvGuunHPzaZ68s+WyO1vcfp7mYurw/+cDbgl5QJE4U9/k49mlhdy7cCu799cxObczf7piLCcP7Op1NBEJgY0l+7nmH5/Q0OTnkRsmcdrgbK8jRT0VWCIx7P53t/LAe9u4dlo/vn9e5BVXIhI5Gpr8PLe8kHvfyWdXVR0T+nXinsvHcMqgLtp2iMSoVYWVXPfwEtqnJDJ39jQGdcvwOlJMUIElEqPmLingt69t4qIxvfj550dqB0lEWuX3O+at3sXdC/IorqxlXE4Wv71sNKepeYVITFtVWMmXH/qETunJPPWVqZrXKohUYInEoNfWlvCjF9YyfUg291w+hgTv5rkSkQi2fOc+fvHKBlYVVjKyV0d+eckoZgzJVmElEuPySw9wwz+X0Dk9hWdumUrPzHZeR4opKrBEYszSHRV8a+4qxuV04v4vTyAlKcHrSCISYYr21fC71/N4efUusjukcvcXRnPZ+D46GCMSB0qqarn2H5+QmJDA4zdNVnEVAiqwRGLIjvKDzH5sGX06teMf102kXYom/RSR/+9gfRN/X7SVBxdvA+C2Mwfx1ekDSU/V7oBIPKhr9DH7seXsr2ti7uypmjw4RLRFFYkRlTUN3PjIUgAevn4SWe01d4WINPP7Hc+vKOKeBXmUVtdz0Zhe/OD8YfTO0pFrkXjhnONHL6xlbXEVD107kVG9M72OFLNUYInEgIYmP7c8vpyifbU8efMUcrvqiJSINPtk215+8eoG1hXvZ2zfLO6/ZgLjczp5HUtEwuyxj3by7xXFfPvswZw9orvXcWKaCiyRKOec44f/Xssn2yv48xVjmZTb2etIIhIBdu49yG/mb+L19bvplZnGX64cy0VjeqmBhUgc2liyn1+9upGzhnXjm2cO9jpOzFOBJRLlHv5gB/9aUcS3zhrMxeN6ex1HRDzW0ORnzntb+b938klKML53zhC+ctoAXZMpEqfqGn18e+4qMtsnc7c6C4eFCiyRKPbh1nJ+PX8j547ozrfO0hEpkXi3smAfd/xrLXl7qvnc6J7ceeEIundM8zqWiHjongV55O2p5p83TKJzuq7PDgcVWCJRqmhfDbc+tZL+XdP54xVjdURKJI4dqG/ingV5PPrRDnp0TOOhayfqGgsRYfnOfTz0/naumdqPM4Z28zpO3FCBJRKFaht83PL4chp9fuZcM4EMtVgWiVtvb9zDT19cR8n+Oq6d2o/bzxtKh7Rkr2OJiMeafH5+/MJaemam8YPzh3kdJ65or0wkyhxqs7qhZD//uG4iA7IzvI4kIh4oq67nf19ezytrShjcLYPnv3oyE/qpO6CINHvkwx1s2l3N/V/Wgdhw06stEmXmLi3khZXFfOfsIZw5TEOAROKNc47nlhXxq/kbqW3w8Z2zh/C1GQNJSUrwOpqIRIiSqlr++OZmzhrWjfNGal8h3FRgiUSRTbv38/N56zl1UFduPXOQ13FEJMx2lB/kRy+s5cOte5mU24nfXHoSg7p18DqWiESYexZspsnv+PlFIzU1gwdUYIlEiYP1TXzjyRV0bJfMn64YS6KaWojEDZ/f8dDibfzxzc2kJCbwy4tHcfXkHDW3EZH/smHXfv69sojZpw2gb+f2XseJSyqwRKLET19ax7bygzx50xSyO6R6HUdEwqSkqpbvPLOKj7dVcM6I7vxi1ih6ZKr1uoi07revb6JjWjJfn6GRLl5p04BtM5tpZnlmlm9md7Ryf46ZLTSzlWa2xswuCH5Ukfj13LJC/r2imG+eOZiTB3X1Oo6IhMnr60qY+efFrC6s4veXjWbONRNUXInIEb2/pZz3Npdx6xmDyGyvbqJeOeoZLDNLBO4FzgGKgKVmNs85t6HFaj8BnnXO/d3MRgDzgdwQ5BWJO9vLD3LnS+uZNqAL39RkwiJxoaahiV+8soGnlxRyUu9M/nLlWHUMFZHP5JzjD2/m0SszjWum9fM6TlxryxDByUC+c24bgJnNBWYBLQssB3QM3M4EdgUzpEi8avT5+fYzq0hJStB1VyJxYl1xFd+cu5Lt5Qf56vSBfPecIeoQKCJH9dHWvawsqOQXF48iLTnR6zhxrS0FVm+gsMX3RcCUw9b5OfCGmd0GpANnt/aDzGw2MBsgJyfnWLOKxJ2/vZPP6sJK/nb1OA0LEolxfr/jofe3cfeCPDqnp/DkTVM0JFhE2uyv7+TTrUMql0/o43WUuBesQ2JXAY845/oAFwCPm9l//Wzn3Bzn3ETn3MTs7OwgPbRIbFpRsI+/Lczn0nG9uXB0L6/jiEgI7dlfx7UPL+HX8zdx5rBuvP6t01VciUibLdtRwUfb9jL79AE6exUB2nIGqxjo2+L7PoFlLd0EzARwzn1kZmlAV6A0GCFF4s3B+ia+88wqenRM4+ezRnodR0RC6M0Ne/j+86upbfTx60tO4qrJfTVvjYgck3sX5tM5PYWrp2iEWCRoyxmspcBgM+tvZinAlcC8w9YpAM4CMLPhQBpQFsygIvHkl69upKCihj9+cQwd09QFSCQW1Tb4+MmLa7n5sWX0ymrHK7edxtVTclRcicgxyS89wMK8Mq4/OZf2KZqBKRIc9bfgnGsys1uBBUAi8LBzbr2Z3QUsc87NA74HPGhm36G54cX1zjkXyuAiseq9zWU8vaSAW04fwJQBXbyOIyIhUFhRw82PLWPT7mpmnz6A7507hNQkDesRkWP36Ic7SElM0NmrCNKmMtc5N5/m1ustl93Z4vYG4JTgRhOJPwfqm/jhv9cyIDud75wzxOs4IhIC728p59anV+D3O/55wyTOGNrN60giEqWqahv514oiLhrbi64ZqV7HkQCdRxSJIL99bSO7qmp5/qvTdJGqSIxxzvGP97fz6/kbGdQtgznXTCS3a7rXsUQkij23rJCaBh/Xn5zrdRRpQQWWSIT4cGs5T3xcwE2n9mdCv85exxGRIKpt8PHDf6/hxVW7mDmyB3/44hjSU/URLCLHz+d3PPrRDibndmZU70yv40gL2rqLRICahibu+Ndacru05/Zzh3odR0SCqGhfDbc8vpwNJfu5/dwhfOOMQWpkISInbPGWMgoravnBzGFeR5HDqMASiQB3L8ijoKKGZ2ZPpV2KhgaKxIqPtu7lG0+toNHn5x/XTeTMYd29jiQiMeLZZYV0Tk/h3BE9vI4ih1GBJeKx1YWVPPLhDq6Z2k9dA0VihHOORz7cwS9f3Uj/runMuWYCA7IzvI4lIjFi74F63tywh2un5ZKS1JZZlyScVGCJeKjJ5+dHL6ylW4dU/memhgaKxIK6Rh8/fmEd/1pRxDkjuvPHL46hg+azE5EgemFlMY0+xxWT+nodRVqhAkvEQ498uIP1u/Zz35fGa0JhkRhQUlXLVx9fzuqiKr599mC+eeZgEhJ0vZWIBI9zjmeWFjIuJ4sh3Tt4HUdaoQJLxCO7Kmv545ubOWNoNueP0vhpkWi3omAfsx9bRl2jnwevncg5I3S9lYgE34qCSraUHuB3l53kdRQ5AhVYIh75+bz1+J3jrlmj1FFMJMq9vXEP33hqBT06pjF39iQGddP1ViISGi+uLCYtOYHPje7ldRQ5Al0VJ+KBN9bv5o0Ne/jWWUPo27m913E8ZWYzzSzPzPLN7I5W7u9nZm+b2RozW2RmfVrcl2Nmb5jZRjPbYGa5YQ0vAjyztIDZjy9naPcO/OtrJ6u4EpGQafT5mb+2hLOGdydDc+lFLBVYImFW2+Djf1/ewNDuHfjKaf29juMpM0sE7gXOB0YAV5nZiMNWuwd4zDk3GrgL+E2L+x4D7nbODQcmA6WhTy3SzDnHX9/ewg/+tZZTB3XlqZun0iUj1etYIhLDPsgvZ+/BBmaN0dmrSKYCSyTM7luUT3FlLXfNGklyYtz/CU4G8p1z25xzDcBcYNZh64wA3gncXnjo/kAhluScexPAOXfAOVcTntgS73x+x50vrecPb27m0vG9eei6iaTraLKIhNi8VbvomJbE9KHZXkeRzxD3e3ci4VSwt4YH3tvGRWN6ac6rZr2BwhbfFwWWtbQauDRw+xKgg5l1AYYAlWb2bzNbaWZ3B86IiYRUXaOPbzy5gsc/3skt0wfwh8vH6GCJiIRcXaOPBet3c/6onqQm6eMukukTQSSM7nplA0kJxo8uGO51lGhyOzDdzFYC04FiwEdzk57TAvdPAgYA17f2A8xstpktM7NlZWVlYQktsamqtpFrH17C6+t389MLR/DD84erSY2IhMXbG0s52OBj1lgND4x0KrBEwmRRXilvbdzDbWcOpkdmmtdxIkUx0HKWxD6BZZ9yzu1yzl3qnBsH/DiwrJLms12rAsMLm4AXgfGtPYhzbo5zbqJzbmJ2toZVyPHZXVXHFQ98xMqCffzfVeO46dT4voZSRMLrlTW7yO6QqhEwUUAFlkgYNDT5uevlDfTvms6Np+Z6HSeSLAUGm1l/M0sBrgTmtVzBzLqa2aFt1Q+Bh1v83ywzO1QxnQlsCENmiUP5pdVc9vcPKdpXyyM3TOYiXWAes9rQ2TTHzBYGhiavMbMLvMgp8aWu0ceivDLOG9mdRE1eHvFUYImEwcMfbGdb+UHu/PwIjZtuIXDm6VZgAbAReNY5t97M7jKziwKrzQDyzGwz0B34VeD/+mgeHvi2ma0FDHgwzE9B4sDynfv4wv0fUd/kZ+7sqZwyqKvXkSRE2tjZ9Cc0b6vG0XxQ6L7wppR49N7mMmobfcwc2dPrKNIGankkEmJl1fX89e0tnDWsG2cM7eZ1nIjjnJsPzD9s2Z0tbj8PPH+E//smMDqkASWutZxA+LEbp5DTJb7nrYsDn3Y2BTCzQ51NW54dd0DHwO1MYFdYE0pcen39bjLbJTNlQGevo0gbqMASCbE/v7WZ+iY/P/qcGluIRJN5q3fxnWdWMbJXRx6+fhJdNcdVPGits+mUw9b5OfCGmd0GpANnt/aDzGw2MBsgJycn6EElfjT6/Ly9sZSzhnVTx9Iood+SSAht2VPN00sK+NKUHAZmZ3gdR0Ta6OXVu/j23JVM6NeJp2+equJKWroKeMQ51we4AHi8xXWin1JzHQmWJdsrqKpt5LxRPbyOIm2kAkskhH49fyPpqUl86+whXkcRkTZ6dU0J335mFRP7deaf10/SBMLx5aidTYGbgGcBnHMfAWmALsyTkHl93W7SkhM4fbAK9WihAkskRN7fUs7CvDJuPWMQndNTvI4jIm3w2toSvjl3JeNzsvjnDSqu4tBRO5sCBcBZAGY2nOYCSxPsSUg453hzwx6mD8mmXYqaZEULFVgiIeDzO3756gb6dGrHdSfneh1HRNrg9XUl3Pb0Ssb2zeKfN0xWcRWH2tjZ9HvAzWa2GngauN4557xJLLFuQ8l+du+v46zh3b2OIsdAnx4iIfCv5UVs2l3NX68aR1qyjjiJRLoF63dz61MrGd0nk0dumESGiqu41YbOphuAU8KdS+LTorzmk6Mzhmp4YDTRGSyRIDtY38Q9b+QxLieLC0drvgqRSPfG+t1848kVnNQnk0dvnEyHtGSvI4mIAPDOplJG98mkW4c0r6PIMVCBJRJkDy7eRml1PT/53HDMNNu6SCR7a0PzPFcje6u4EpHIsu9gAysL9jFDc2hGHRVYIkFUur+OB97dxudO6smEfpoMUCSSvb1xD197cjkjenbksRsn01HFlYhEkPe2lOF3cOYwFVjRRgWWSBD99Z18Gn1+vj9zqNdRROQzLNxUyteeWMHwnh157KYpZLZTcSUikeWdTaV0SU9hdO9Mr6PIMVKBJRIkBXtreHpJAVdM6ku/LulexxGRI1iYV8otjy9naI8OPH6jiisRiTw+v+PdzWVMH5pNQoIuN4g2KrBEguRPb20mKdH45lmDvY4iIkewKFBcDemRwRM3TSGzvYorEYk8q4sqqaxp1PVXUUoFlkgQbNq9nxdXFXPdybl076hOPyKRaNmOCm55fDmDslVciUhke39LOWZw6qCuXkeR46ACSyQI7lmwmYzUJL42faDXUUSkFfmlB7jp0WX0ymrHE1+ZQlb7FK8jiYgc0ftbyhnZqyOd07WtikYqsERO0PKd+3hr4x5uOX2AdtpEIlDp/jque3gJyYnGozdM1g6LiES0A/VNrCjYx6mDNLlwtNJU9SInwDnH3Qs20TUjhRtO6e91HBE5THVdI9f/cyn7ahp4ZvY0crq09zqSiMhnWrJ9L01+x2mDNTwwWukMlsgJWLylnI+3VXDrGYNIT9XxCpFI0tDk5+tPriBvTzX3fWk8J/VRq2MRiXyLt5STmpTAhH6dvI4ix0kFlshxaj57lUfvrHZcNSXH6zgi0oJzjjv+vYbFW8r57aUnqROXiESN97eUM7l/Z9KSE72OIsdJBZbIcXp93W7WFlfxnXOGkJqkjaBIJLnnjTz+vaKY754zhMsn9vU6johIm+yuqmNL6QEND4xyKrBEjoPf7/jzW1sYmJ3OJeN6ex1HRFp4/OOd3LtwK1dNzuG2Mwd5HUdEpM0+yC8HUIOLKKcCS+Q4zF9XQt6ear551mASNcO6SMR4Y/1ufvbSOs4e3o1fzBqJmf4+RSR6fLRtL53TUxjWo4PXUeQEqMASOUZ+v+Mvb21hULcMLhzdy+s4IhKwfOc+bnt6JSf1yeL/rhpHUqI+4kQkuny8bS9T+ncmQQdvo5o+fUSO0atrS9hSekBnr0QiyNayA9z06FJ6Zqbx8HUTaZ+irp4iEl0KK2oo2lfL1AFdvI4iJ0gFlsgx8Pkdf3l7C4O7ZfC5k3p6HUdEgNLq5omEkxKMR2+cTJeMVK8jiYgcs0+2VwCowIoBKrBEjsGra0vILz3At87W2SuRSHCgvokbH1nK3gMN/OO6SfTrku51JBGR4/LR1ubrrwZ3y/A6ipwgjaEQaSOf3/GXtzYzpHsGF4zS2SsRr/n8jlufWsHGkmoeunYiY/pmeR1JROS46fqr2KEzWCJt9MqaXWwtO8i3zhqijZ9IBPjjm3ksyivjrlkjOWOYJhIWkehVWFFDcaWuv4oVKrBE2uDQtVfDenTg/FE9vI4jEvdeW1sSmOuqL1+a0s/rOCIiJ+TjbXsBXX8VK1RgibTBy6t3sa3sIN86a7DOXol4bMueam5/bjVj+2bx84tGeh1HROSEfbytgk7tk3X9VYxQgSVyFD6/46/vNJ+9Om+kzl6JeGl/XSOzH19Ou5Qk7v/yBFKTEr2OJCJywpbuqGCyrr+KGSqwRI7i9XW72Vp2kFvPHKQNn4iH/H7Hd+auorCihr9/eTw9MtO8jiQicsL27K+joKKGSbmdvY4iQaICS+QzOOf428J8BmSnc746B4p46i9vb+HtTaXc+fkR2hERkZixdEfz/FcTtV2LGSqwRD7DwrxSNpbs52vTB2reKxEPvblhD395ewtfmNCHa6aqqYWIxI5lO/bRLjmRkb06eh1FgkQFlsgROOf46zv59M5qx8XjensdRyRubS07wHefWcVJvTP55cWjMNPBDhGJHUt3VDAuJ4vkRO2Wxwr9JkWO4KOte1lZUMlXZwzURk/EI9V1jcx+bBnJSQncf80E0pLV1EJEYkd1XSMbS/ZreGCMSfI6gEik+tvCfLp1SOXyCX28jiISl/x+x+3PrWbH3hoev2kyvbPaeR1JRCSoVhRU4ncwKbeT11EkiHRYXqQVKwr28eHWvdx82gAdMQ8xM5tpZnlmlm9md7Ryfz8ze9vM1pjZIjPrc9j9Hc2syMz+Fr7UEg5/f3crC9bv4YfnD+PkgV29jiMiEnTLdlSQmGCMy1GBFUtUYIm04t538slqn8zVU3K8jhLTzCwRuBc4HxgBXGVmIw5b7R7gMefcaOAu4DeH3f8L4L1QZ5XwWphXyj1v5DFrbC9uOrW/13FEREJi6Y4KRvTsSEaqBpXFEhVYIodZv6uKtzeVcuMp/UnXBi/UJgP5zrltzrkGYC4w67B1RgDvBG4vbHm/mU0AugNvhCGrhMmO8oN86+mVDOvRkd9eOlpNLUQkJjX6/KwqrGRCP529ijVtKrCONoQnsM4XzWyDma03s6eCG1MkfO5buJWM1CSum5brdZR40BsobPF9UWBZS6uBSwO3LwE6mFkXM0sA/gDcfrQHMbPZZrbMzJaVlZUFIbaEysH6Jm55fDkJCcacaybQLkVDdEUkNm0s2U9do5+Juv4q5hy1wGrLEB4zGwz8EDjFOTcS+Hbwo4qEXn7pAeavK+Haaf3IbJ/sdRxpdjsw3cxWAtOBYsAHfB2Y75wrOtoPcM7Ncc5NdM5NzM7ODm1aOSE/eXEdW0qr+etV4+jbub3XcUREQmb5zn0AjNf1VzGnLeOfPh3CA2Bmh4bwbGixzs3Avc65fQDOudJgBxUJh78v2kpqUoKu+QifYqBvi+/7BJZ9yjm3i8AZLDPLAC5zzlWa2TTgNDP7OpABpJjZAedcq2fZJfK9uLKYF1YW852zh3DaYBXCIhLbVhRU0qNjGr3UITXmtKXAam0Iz5TD1hkCYGYfAInAz51zrwcloUiYFFbU8OKqYq6d1o8uGalex4kXS4HBZtaf5sLqSuDqliuYWVegwjnnp/lM+cMAzrkvtVjnemCiiqvoVVhRw09fXMfEfp34xhkDvY4jIhJyK3buY3y/LK9jSAgEq8lFEjAYmAFcBTxoZlmHr6TrICSSPbR4GwkGs08f4HWUuOGcawJuBRYAG4FnnXPrzewuM7sosNoMIM/MNtPc0OJXnoSVkGny+fn2M6sA+NMVY0nSxN4iEuNK99dRXFmr4YExqi1nsI46hIfms1qfOOcage2BHaHBNB+d/pRzbg4wB2DixInueEOLBFvFwQaeWVbIRWN60zNTp+rDyTk3H5h/2LI7W9x+Hnj+KD/jEeCREMSTMPjbwnyW79zHX64cq+uuRCQurCgIXH+lDoIxqS2HCT8dwmNmKTQP4Zl32Dov0nyU+dBwniHAtuDFFAmtxz/aSV2jX2evRMJs+c4K/u/tLVwyrjezxh7eQFJEJDatKKgkJTGBkb06eh1FQuCoBVYbh/AsAPaa2Qaa56n5H+fc3lCFFgmmukYfj320gzOGZjO0Rwev44jEjf11jXxr7ip6d2rHXbNGeh1HRCRslu/cx6jeHUlN0lQUsahNs6i2YQiPA74b+BKJKs8vL2LvwQZmn64L60XC6Wcvraekqo5nb5lGhzRNiyAi8aGhyc/a4iqundrP6ygSIrqSWOKaz+94aPE2xvTJZOqAzl7HEYkbL61qbsn+zTMHM0HXIIhIHNlQsp+GJr+uv4phKrAkrr2xfjc79tYw+/SBmJnXcUTiQmFFDT95QS3ZRSQ+rQw0uBiXk+VtEAkZFVgSt5xzPPDeNnI6t2fmqB5exxGJC2rJLiLxbnVhJd07pqprcQzTJ5vEraU79rGqsJKbT+tPYoLOXomEw6GW7L+8ZJRasotIXFpVWMnYvllex5AQUoElceuBd7fSOT2FL0zoe/SVReSEqSW7iMS7fQcb2LG3hrF9df1VLFOBJXFpy55q3t5UyjVT+9EuRS1SRUKtOtCSvVdWO/5XLdlFJE6tKqoE0BmsGNemNu0isebBxdtIS07g2mlqkSoSDnd+2pJ9Kh3Vkl1E4tSqgkoSDEb3yfQ6ioSQzmBJ3Nmzv44XV+7i8gl96ZKR6nUckZh3qCX7bWcOYkI/TYcgIvFrVWElg7t1ID1V5zhimQosiTv//GAHTX4/Xzmtv9dRRGJey5bst54xyOs4IkdlZjPNLM/M8s3sjiOs80Uz22Bm683sqXBnlOjknGN1kRpcxAOVzxJXDtQ38eQnO5k5qgf9uqR7HUckpjnn+NELa/E7p5bsEhXMLBG4FzgHKAKWmtk859yGFusMBn4InOKc22dm3bxJK9Fmx94aKmsaGav5r2KePu0krsxdUkB1XRO3nK7JTUVC7d8rilm8pZwfnD9MLdklWkwG8p1z25xzDcBcYNZh69wM3Ouc2wfgnCsNc0aJUqsKmycY1hms2KcCS+JGo8/PP97fzpT+nRmjjZtISJUfqOcXr25gQr9OfHmKmslI1OgNFLb4viiwrKUhwBAz+8DMPjazma39IDObbWbLzGxZWVlZiOJKNFldWEW75EQGd8vwOoqEmAosiRsvr95FSVUdt0wf4HUUkZj3vy9voKbex+8uO4kETeQtsSUJGAzMAK4CHjSzrMNXcs7Ncc5NdM5NzM7ODm9CiUiriyo5qXemhkvHAf2GJS4455jz3jaGdM9gxhANlxcJpbc37uHl1bu49cxBDOrWwes4IseiGGg5+3yfwLKWioB5zrlG59x2YDPNBZfIETX6/KzftV/t2eOECiyJC+9tKWfT7mpuPm2AjqaLhFB1XSM/eXEdQ7pn8NXputZRos5SYLCZ9TezFOBKYN5h67xI89krzKwrzUMGt4Uxo0ShvN3VNDT5dYlCnFCBJXHhgXe30r1jKrPGHj6UXkSC6e4FeezeX8dvLxtNSpI+YiS6OOeagFuBBcBG4Fnn3Hozu8vMLgqstgDYa2YbgIXA/zjn9nqTWKLF6qJKAMb0yfI0h4SH2rRLzFtXXMWHW/fyw/OHaYdPJISW7ajg8Y93cv3JuYzP6eR1HJHj4pybD8w/bNmdLW474LuBL5E2WVNYRaf2yfTt3M7rKBIG2tuUmPfAe9vISE3iqik5XkcRiVn1TT7u+PdaemW24/Zzh3odR0QkoqwuqmR0nyzMdJlCPFCBJTGtsKKGV9fs4uopOXRMS/Y6jkjMunfhVvJLD/CrS0aRnqrBESIih9Q0NLF5TzVj1OAibqjAkpj2j/e3k5hg3HBKrtdRRGJW3u5q/r4on0vG9WbGUHXpFBFpaf2u/fgdanARR1RgSczad7CBZ5YWctGY3vTM1JhnkVDw+R0/+NcaOqQl89MLR3gdR0Qk4qwurARgtBpcxA0VWBKznvh4J7WNPmafromFRULlsY92sKqwkjsvHEHn9BSv44iIRJzVRVX0ykwju0Oq11EkTFRgSUyqa/TxyIc7mDE0m6E9NNGpSCgU7avh7gV5zBiazayxvbyOIyISkdYEGlxI/FCBJTHpXyuK2HuwgVtO10SnIqHgnOPHL6wD4JcXj1JnLBGRVlTVNLJzbw2j+6rBRTxRgSUxx+d3PLR4O6P7ZDJ1QGev44jEpJdW7eLdzWX8z3lD6dOpvddxREQi0triKgBG987yNoiElQosiTlvbtjN9vKDzD59gI6qi4TA3gP1/O/L6xmXk8W103K9jiMiErHWFFcCcFJvncGKJyqwJKY453jgvW3kdG7PzJE9vI4jEpN+8coGDtQ38bvLRpOYoIMYIiJHsraoin5d2pPZXnNxxhMVWBJTlu3cx8qCSr5yWn+SEvX2Fgm2D/LLeXHVLr42YxBDuquBjIjIZ1lTVKUGF3FIe6ASUx54dyud2idz+YS+XkcRiTlNPj//+/J6+nZux9dnqIGMiMhnKT9QT3FlLaM1PDDuqMCSmJFfWs1bG0u5dlou7VISvY4jEnOeWlLA5j0H+PEFI0hL1t+YiMhnOdTg4qQ+KrDijQosiRkPvred1KQErp3Wz+soIjFn38EG/vDGZk4e2IXzRnb3Oo6ISMRbW1SFGYzSGay4owJLYkLp/jpeWFnM5RP70CVDM6WLBNuf39pMdV0jd35+hLpzioi0wZqiSgZmZ5CRmuR1FAkzFVgSE/754Q4a/X6+cuoAr6PIMTKzmWaWZ2b5ZnZHK/f3M7O3zWyNmS0ysz6B5WPN7CMzWx+474rwp48PebureeKTAr40pR/DenT0Oo6ISFRYU1Sl66/ilAosiXoH6pt44uOdnD+qB7ld072OI8fAzBKBe4HzgRHAVWY24rDV7gEec86NBu4CfhNYXgNc65wbCcwE/mxmWWEJHkecc9z1ynoyUpP47jlDvI4jIhIV9uyvo7S6XsMD45QKLIl6c5cUUF3XxOzT1dUsCk0G8p1z25xzDcBcYNZh64wA3gncXnjofufcZufclsDtXUApkB2W1HHkzQ17+CB/L989Zwid0lO8jiMiEhXWFDU3uBjTVwVWPFKBJVGt0efn4fe3M7l/Z8b2zfI6jhy73kBhi++LAstaWg1cGrh9CdDBzLq0XMHMJgMpwNYQ5YxL9U0+fvnqRoZ0z+BLU3K8jiMiEjXWFlWSYDCipwqseKQCS6LaK2t2sauqjq9O17VXMex2YLqZrQSmA8WA79CdZtYTeBy4wTnnb+0HmNlsM1tmZsvKysrCkTkmPPz+DgoqarjzwpGauFtE5BisKa5iSPcOmjYmTukTU6KWc44H3t3G4G4ZzBjSzes4cnyKgZazQvcJLPuUc26Xc+5S59w44MeBZZUAZtYReBX4sXPu4yM9iHNujnNuonNuYna2RhG2Ren+Ov72zhbOGdGdUwd39TqOiEjUcM6xtqiKk3T9VdxSgSVRa/GWcjbtrubm0weQkKC20VFqKTDYzPqbWQpwJTCv5Qpm1tXMDm2rfgg8HFieArxAcwOM58OYOS78fkEejT7Hjy8Y7nUUEZGosquqjr0HGxitCYbjlgosiVoPvLeV7h1TmTW2l9dR5Dg555qAW4EFwEbgWefcejO7y8wuCqw2A8gzs81Ad+BXgeVfBE4HrjezVYGvsWF9AjFqdWElzy8v4sZT+6szp4jIMVobaHChDoLxSzOfSVRaV1zFB/l7ueP8YaQmaXxzNHPOzQfmH7bszha3nwf+6wyVc+4J4ImQB4wzfr/j5y+vJ7tDKreeOcjrOCIiUWdtcSVJCcbwnpo3MF7pDJZEpTnvbSMjNYmr1dlMJKheWl3MyoJKvn/eUDJSdQxORORYrSlqbnCRlqwDwPFKBZZEncKKGl5Zs4urp+TQMS3Z6zgiMeNgfRO/fW0To/tkctn4Pl7HERGJOs451hZX6fqrOKcCS6LOQ4u3kZhg3HhKf6+jiMSUvy/ayp799fzs8yPVOEZE5DgU7aulsqaRk1RgxTUVWBJV9h6o55llhVw8tjc9MtO8jiMSMworapizeBuXjOvNhH6dvI4jIhKV1gQaXKhFe3xTgSVR5bGPdlLX6OcWTSwsElS/nr+RRDN+MHOY11FERKLWmuJKUhITGNqjg9dRxEMqsCRq1DQ08dhHOzh7eHcGddOGSyRYPtxazmvrdvONMwbqzLCIyAlYU1jFsJ4d1OE4zqnAkqjx7NJC9tU08lWdvRIJGp/fcdfLG+jTqR1fOU1/WyIix8vvd6wrrtLwQFGBJdGhyefnwcXbmdivExNzO3sdRyRmvLx6F5t2V/ODmcPUUlhE5ATs2HuQ6vomdRAUFVgSHV5dW0JxZS23TB/odRSRmNHo8/PHNzczomdHPndST6/jiIhEtbXFzQ0uRvfJ8jaIeE4FlkQ85xz3v7uNQd0yOGtYN6/jiMSMZ5cVUlBRw/+cN1Rt2UVETtCaoipSkxIY3C3D6yjiMRVYEvEWbylnY8l+Zp8+QDuBIkFS1+jjr2/nM6FfJ2YMzfY6johI1FtbVMXIXh1JStTudbzTO0Ai3v3vbqV7x1Rmje3ldRSRmPHExzvZvb+O/zlvKGY6cCEiciJ8fse6XVUaHiiACiyJcGuKKvlw615uPKW/Wp6KBMmB+ibuW7SV0wZ3ZeqALl7HERGJetvKDlDT4FMHQQFUYEmEe+C9bXRITeLqKTleRxGJGQ+/v52Kgw3cfu5Qr6OIiMSENUWHGlyowBIVWBLBdu49yGtrS/jS1H50SEv2Oo5ITNh3sIEH39vGeSO7M6ZvltdxRERiwpqiStqnJDIgWw0uRAWWRLAHF28jKSGBG0/J9TqKSMy4/72tHGho4ns6eyUiEjSri6oY1TuTRDXjElRgSYQqP1DPc8uKuHR8b7p1TPM6jkhMKN1fx6Mf7uDisb0Z0r2D13FERGJCQ5OfDSX7GaPhgRLQpgLLzGaaWZ6Z5ZvZHZ+x3mVm5sxsYvAiSjx65IMdNPj83Hz6AK+jiMSMvy3Mp8nn+PbZg72OIiISMzbvqaahya8OgvKpoxZYZpYI3AucD4wArjKzEa2s1wH4FvBJsENKfNlf18ijH+3gvBE9GKixzCJBUVhRw9NLCrhiUl/6dUn3Oo6ISMxYXVQJwBgVWBLQljNYk4F859w251wDMBeY1cp6vwB+B9QFMZ/Eocc/2kl1XRPfOGOQ11FEYsaf39pCghm3namzVyIiwbSmsIpO7ZPp27md11EkQrSlwOoNFLb4viiw7FNmNh7o65x79bN+kJnNNrNlZrasrKzsmMNK7Ktt8PHw+9uZPiSbkzSWWSQo8kureWFlEddO60ePTF3TKCISTKuLKjmpT5YmbZdPnXCTCzNLAP4IfO9o6zrn5jjnJjrnJmZnZ5/oQ0sMmru0gL0HG3T2SiSI/vjmZtolJ/K1Gfq7EhEJptoGH1tKD6jBhfyHthRYxUDfFt/3CSw7pAMwClhkZjuAqcA8NbqQY9XQ5GfOe9uYnNuZyf07ex1HJCasLapi/trd3HTaADqnp3gdR0QkpqzfVYXP7ziptwos+f/aUmAtBQabWX8zSwGuBOYdutM5V+Wc6+qcy3XO5QIfAxc555aFJLHErBdWFlFSVcfXzxjodRSRmHHPG3lktU/mK6f19zqKiEjMWV1UBaCJ2+U/HLXAcs41AbcCC4CNwLPOufVmdpeZXRTqgBIffH7H3xdt5aTemUwfouGjIsGwZHsF724u42vTB9IxLdnrOCIiMWdNUSXdO6bSXXN2SgttugbLOTffOTfEOTfQOferwLI7nXPzWll3hs5eybF6dW0JO/bW8I0zBuoiUZEgcM5x94JNZHdI5dppuV7HEYkamvtTjsWaoirNfyX/5YSbXIicKOcc9y3MZ1C3DM4d0cPrOCIx4d3NZSzdsY9vnjmIdimJXscRiQqa+1OORWVNA9vLDzJWwwPlMCqwxHNvbyxl0+5qvj5jIAkJOnslcqKcc9zzRh59OrXjikk5XscRiSaa+1Pa7ND1Vyqw5HAqsMRTzjn+tjCfPp3a8fkxvbyOIxITXl+3m3XF+/n22UNISdJmXuQYaO5PabPVhZWYoXk75b/ok1c89dHWvawqrOSr0weSnKi3o8iJ8vsdf3prMwOz07lkXO+j/wcRaTPN/SktrS6sZGB2hpoIyX/RHq146m8L8+nWIZUvTOjjdRSRmPDWxj1s3nOA284cTKKG3IocK839KW3inGNVYaWGB0qrVGCJZ5bv3MeHW/dy82kDSEvWRfgiJ8o5x32LttK3czsuHN3T6zgi0Uhzf0qbFO2rZe/BBs1/Ja1SgSWe+cvbW+icnsKXpuoifJFg+Ghb85Db2acPJElDbkWOmeb+lLZaVVgJwDgVWNKKJK8DSHxavnMf720u44fnD6N9it6GIsHw90Vb6ZqRyuUacity3Jxz84H5hy278wjrzghHJok8qwsrSU1KYGiPDl5HkQikQ5ziiUNnr66Z1s/rKCIxYW1RFYu3lHPTqf015FZEJMRWF1UyqnemGnRJq/SukLBbUdB89mr26QN09kokSO5blE+HtCS+rCG3IiIh1eTzs7a4ijF9sryOIhFKBZaE3V/eCpy9mqqzVwJmNtPM8sws38zuaOX+fmb2tpmtMbNFZtanxX3XmdmWwNd14U0eOfJLD/D6+t1cO60fHdQuWEQkpDbtrqau0c/YnCyvo0iEUoElYbWiYB/vBs5epafq7FW8M7NE4F7gfGAEcJWZjThstXuAx5xzo4G7gN8E/m9n4GfAFGAy8DMz6xSu7JHkgXe3kpKYwA2n9Pc6iohIzFtRsA+A8Sqw5AhUYElY6eyVHGYykO+c2+acawDmArMOW2cE8E7g9sIW958HvOmcq3DO7QPeBGaGIXNE2VVZywsri7lyUl+6ZqR6HUdEJOatLKikW4dUeme18zqKRCgVWBI2KwNnr24+TWev5FO9gcIW3xcFlrW0Grg0cPsSoIOZdWnj/wXAzGab2TIzW1ZWVhaU4JHiwcXbALj59AEeJxERiQ8rCvYxLicLM03mLq1TgSVh85e3t9CpfTLXqnOgHJvbgelmthKYDhQDvmP5Ac65Oc65ic65idnZ2aHI6ImKgw3MXVLIRWN70adTe6/jiIjEvL0H6tm5t4bxOXE5Il3aSAWWhMXKgn0syitj9ukDdfZKWioG+rb4vk9g2aecc7ucc5c658YBPw4sq2zL/411j3ywndpGH1+bPtDrKCIicWFlQSUA41RgyWdQgSVh8cc3N9M5PUVnr+RwS4HBZtbfzFKAK4F5LVcws65mdmhb9UPg4cDtBcC5ZtYp0Nzi3MCyuHCgvolHPtzBuSO6M7i7JroUEQmHlYX7SEowRvfJ9DqKRDAVWBJyH24tZ/GWcr4+Q2ev5D8555qAW2kujDYCzzrn1pvZXWZ2UWC1GUCemW0GugO/CvzfCuAXNBdpS4G7AsviwlOf7GR/XRNfP2OQ11FEROLGip2VjOjVURO6y2fS3q6ElHOOexbk0aNjGl9W50BphXNuPjD/sGV3trj9PPD8Ef7vw/z/M1pxo77Jx0OLt3PywC6M7ZvldRwRkbjg8ztWF1Vy+YQ+R19Z4prOYElIvbOplBUFlXzzrME62iMSJP9aXkxpdT1fn6GzVyIi4ZK3u5qaBh/j++n6K/lsKrAkZPx+x90L8ujXpT2XT9TRHpFgaPL5eeC9rYzuk8kpg7p4HUdEJG4s/3SCYRVY8tlUYEnIvLq2hE27q/nO2UNITtRbTSQY5q/bzc69NXx9xkDNwSIiEkbLd1TQrUMqfTppgmH5bNrrlZBo8vn545ubGdq9A58f08vrOCIxwTnH3xdtZWB2OueO6OF1HBGRuLJ0xz4m5XbWwS05KhVYEhL/WlHE9vKDfO/cISQmaEMkEgyL8srYWLKfr04fSIL+rkREwqakqpbiylom6PoraQMVWBJ0dY0+/vLWFsb0zeKcEd29jiMSM+5blE+vzDRmje3tdRQRkbiybEfz9VcTc1VgydGpwJKge+qTAnZV1fH984bqNLpIkCzdUcHSHfuYffoAUpK06RYRCaflO/fRPiWRET07eh1FooA+pSWo9tc18reF+Uwb0IVTBnX1Oo5IzLhvYT5d0lO4YlKO11FEROLO0h0VjO2bRZKadkkb6F0iQfX3RVupONjAjy4Y7nUUkZixYdd+FuaVccMpubRL0XxyIiLhdKC+iY0l+5mY29nrKBIlVGBJ0OyqrOXh97dz8dhenNQn0+s4IjHjwcXbSE9J5JppuV5HERGJO6sKKvE7mKgGF9JGKrAkaO55Iw8H3H7eUK+jiMSM0uo6Xlmzi8sn9iWzXbLXcURE4s7SHRUkGIzLyfI6ikQJFVgSFOuKq3hhZTE3nJJLn07tvY4jEjOe+qSARp/jupNzvY4iIhKXlmyvYHjPjnRI00EuaRsVWHLCnHP8ev5Gstol8/UZg7yOIxIzGpr8PPFxAWcMzaZ/13Sv44iIxJ36Jh8rCvYxpX8Xr6NIFFGBJSdsUV4ZH27dyzfPGqwhTCJBNH9tCeUH6rn+lP5eRxERiUtriqqob/IzZYAaXEjbqcCSE9Lk8/Pr+RvJ7dKeL03p53UckZjyzw93MCA7ndM05YGIiCc+2bYXgMnqICjHQAWWnJDnlhexpfQAP5g5TJOfigTRyoJ9rC6s5PqTc0lI0ITdIiJe+GR7BcN6dKBTeorXUSSKaI9YjtuB+ib++OZmJvTrxMxRPbyOIxJTHv1wBx1Sk7h0fB+vo4iIxKVGn59lO/Yxpb/OXsmxUYElx+3ehfmUVdfz488Nx0xH2EWCpXR/Ha+uLeELE/uQkZrkdRwRkbi0pqiK2kYfUwaowYUcGxVYclx2lB/kH4u3c+n43ozP0cR7IsH05CcFNPkd12liYRERz3yyPXD9lc5gyTFSgSXH5ZevbiA50bhj5jCvo4jElPomH09+UsAZQ7uRq9bsIiKe+WRbBYO6ZdA1I9XrKBJlVGDJMVuUV8pbG0u57azBdOuY5nUckZjyaWt2TSwsIuKZRp+fpTsqmKr27HIcVGDJMWlo8nPXKxvo3zWdG07J9TqOSExxzvHPD3YwMDud0warNbuIiFdWFVZS0+DjVE2TIcdBBZYck0c/3MG2soP89MLhpCYleh1HJKasLKxkTVEV15+cq8YxIiIe+iC/HDOYqgYXchxUYEmble6v4//e3sIZQ7M5c1h3r+OIxJxHPlBrdhGRSPBBfjkn9c4kq73mv5JjpwJL2uwXr26k3ufnzs+P9DqKSMzZs7+O+WtL+OKkvqSrNbuIiGcO1jexsqCSkwdqeKAcHxVY0iaLt5Tx8updfH3GQPqrs5lI0D35SQE+57h2Wj+vo4iIxLUl2yto8jtdfyXHTQWWHFVdo4+fvriO3C7t+er0gV7HEYk59U0+nvpkJ2cO7Ua/LjqAISLipffzy0lJSmBirub5lOOjcShyVPe/u5Ude2t4/KbJpCWrsYVIsL26poTyAw1cr86cIiKe+yC/nIn9OmmfR46bzmDJZ9pefpD7Fm7l82N6cdrgbK/jiMScQ63ZB3XL0HAUERGPlVXXs2l3NadoeywnQAWWHJFzjjtfWkdqUgI//dxwr+OIxKQVBZWsLa7iOrVmFxHx3HubywA4XQeV5QSowJIj+veKYhZvKef284bSrWOa13FEYtIjH+6gQ1oSl47r7XUUEZG4t2hzGV0zUhjZq6PXUSSKqcCSVpVV13PXKxuY0K8T10xVVzORUNizv47X1pZwxUS1ZhcR8ZrP71i8pYzTh2STkKARBXL8VGBJq342bx21DT5+d9lobWREQuTJj3cGWrPneh1FRCTurSqspLKmkRlDu3kdRaKcCiz5L6+vK2H+2t186+zBDOqW4XUckZhU3+TjyU8KOGtYd3K6tPc6johI3Hs3r5QEg9MHq8GFnBgVWPIfqmoa+elL6xnRsyOzTx/gdRyRmPXK6hL2Hmzg+pNzvY4iIiI0X381tm8WWe1TvI4iUU4FlvyHX766gYqDDfz+C6NJTtTbQ0LPzGaaWZ6Z5ZvZHa3cn2NmC81spZmtMbMLAsuTzexRM1trZhvN7IfhT398nHM88mFza/ZTBnXxOo6ISNwrP1DPmqIqDQ+UoNAetHzq7Y17eG55EbecPoBRvTO9jiNxwMwSgXuB84ERwFVmNuKw1X4CPOucGwdcCdwXWH45kOqcOwmYANxiZrlhCX6CVhTsY21xFderNbuISER4N6+5PfuMoWrPLidOBZYAsPdAPT/411qG9ejAt84e7HUciR+TgXzn3DbnXAMwF5h12DoOONQvNxPY1WJ5upklAe2ABmB/6COfuKc+KSQjNYlL1JpdRCQivLlhD907pjKqlw4wy4lrU4HVhiE83zWzDYHhO2+bmfp6RxHnHD9+YR37axv50xVjSU1K9DqSxI/eQGGL74sCy1r6OfBlMysC5gO3BZY/DxwESoAC4B7nXEVI0wbB/rpGXl27i4vG9lJrdhGRCFDX6OO9LWWcM6K7OidLUBy1wGrjEJ6VwETn3Giad3p+H+ygEjr/XlHM6+t3891zhzC8pybWk4hzFfCIc64PcAHwuJkl0Hz2ywf0AvoD3zOzVjuzmNlsM1tmZsvKysrClbtV81btoq7Rz5WT+nqaQ0Rap4PK8eeD/HJqGnycM6KH11EkRrTlDNZRh/A45xY652oC334M9AluTAmV4spafj5vPZNzO3PzaeoaKGFXDLSsNPoElrV0E/AsgHPuIyAN6ApcDbzunGt0zpUCHwATW3sQ59wc59xE59zE7Gxvx9c/s7SQ4T07cpKucxSJODqoHJ/e3LCHDqlJTBugpkMSHG0psNoyhKelm4DXTiSUhIfP77j92dX4neOey8eQqNPiEn5LgcFm1t/MUmhuYjHvsHUKgLMAzGw4zQVWWWD5mYHl6cBUYFOYch+XdcVVrC2u4spJfdXcQiQy6aBynPH5HW9t3MP0odmkJKk1gQRHUN9JZvZlmo8g332E+yNmmI7A/e9u5aNte7nz8yM00al4wjnXBNwKLAA20twtcL2Z3WVmFwVW+x5ws5mtBp4GrnfOOZqPMmeY2XqaC7V/OufWhP9ZtN2zywpJSUrg4rFqbiESoYJ2UFn7PNFhVeE+yg80cO5IDQ+U4GnLFdZtGcKDmZ0N/BiY7pyrb+0HOefmAHMAJk6c6I45rQTNsh0V/PHNzXx+TC++OFHXgoh3nHPzaW5e0XLZnS1ubwBOaeX/HaC5VXtUqGv08cLKYi4Y1YPM9slexxGRE9TioPL01u7XPk90eGP9HpITTe3ZJajacgbrqEN4zGwc8ABwUeBaCIlglTUNfPPplfTOasevLxmloUoiYfDauhKq65q4YlKO11FE5MiO9aDyRUc6qCyRzznHq2tLOGVQVzqm6cCXBM9RC6w2DuG5G8gAnjOzVWZ2+DUUEiGcc3z/+TWUHajnb1ePo4M2KCJhMXdJIbld2jN1QGevo4jIkemgchxZVVhJ0b5aLhzdy+soEmPaNAlLG4bwnB3kXBIij3+8kzc27OEnnxvO6D5ZXscRiQvbyg7wyfYKvj9zqM4Yi0Qw51yTmR06qJwIPHzooDKwzDk3j/88qAxQ4Jy76Ig/VCLWy6tLSElM4NyR3b2OIjFGs1zGkVWFlfzylY2cOawbN53a3+s4InHj2WVFJCYYXxivZmMikU4HleOD3++Yv7aE6UOzNTxQgk79KONE+YF6vvbEcrp1TOWPXxyjo+giYdLo8/P88iLOHNaNbh3TvI4jIiLAsp372L2/jgtH9/Q6isQgncGKA00+P7c9tZKKgw3862snk9U+xetIInHjnU2llB+o58pJ6tYpIhIpXlmzi7TkBM4eruGBEnwqsOLA7xfk8dG2vfzh8jGM6p3pdRyRuPLM0kK6d0xl+hC1ABYRiQSNPj+vrinhrGHdSU/VrrAEn4YIxrhX15Qw571tXDO1H5dN0PUfIuFUUlXLorxSLp/Ql6REbW5FRCLBorwy9h5s4NLxmvRdQkOf+DFsXXEVtz+3mvE5Wfz0whFexxGJO88vK8Lv0GTeIiIR5LllhXTN0MgCCR0VWDFqz/46vvLoMjq1T+b+ayaQkqRftUg4+f2OZ5YVcsqgLuR0ae91HBERAfYeqOedTaVcOr63RhZIyOidFYNqG3zc/Ngy9tc18tB1k+jWQZ3LRMLtw617KdpXyxWTcryOIiIiAS+u2kWT3/EFXTYhIaQr+2KM3+/43nOrWFtcxYPXTGREr45eRxKJS3OXFpDVPplzR6hDlYhIJHDO8dyyQsb0yWRI9w5ex5EYpjNYMeYPb+Yxf+1ufnzBcM7Wjp2IJyoONvDG+j1cMq43acmJXscRERFgbXEVm3ZX8wVdFyshpgIrhjz+8U7uXbiVqyb35aZT+3sdRyRuvbCymAafnys095WISMR47KOdpKckcvHYXl5HkRinAitGzF9bwp0vrePs4d34xaxRmJnXkUTiknOOZ5YWMLZvFsN6aIiuiEgkqDjYwLzVu7h0fB86pCV7HUdinAqsGPBhfjnfnruKCTmd+OtV49UVR8RDKwsr2bznAFfq7JWISMR4ZmkhDU1+rpnWz+soEge0Jx7l1hVXMfvx5eR2bc8/rptEuxRd7yHipWeWFNI+JZELx2gIiohIJPD5HU98vJNpA7qouYWEhQqsKLZlTzXXPbyEzHbJPHbjFDLb65S3iJcO1Dfx8ppdfH50LzJS1aRVRCQSvL1xD8WVtVyrs1cSJiqwotTWsgNc9eAnJCQYj980mR6ZmutKxGuvrN5FTYOPKyZreKCISCRwznH/u1vpndWOc9RdWcJEBVYU2lF+kKsf/BhwPH3zFAZkZ3gdSUSAuUsLGdI9g3F9s7yOIiIiwNId+1hRUMkt0wfoGnUJG73TokxhRQ1XP/gxjT7Hk1+ZyqBuGkssEgk27d7PqsJKrpiUoy6eIiIR4r5F+XRJT+HyCRpZIOGjAiuKbC8/yJVzPqam0ccTN01haA8VVyKR4pmlhaQkJnDJuN5eRxEREWDDrv0syivjhlNy1QRMwkpXYUeJjSX7ueYfS/A7xxM3TWFEL82vIxIp6hp9vLCymHNHdqdzeorXcUREBLh3UT7pKYlcMzXX6ygSZ3QGKwqsLNjHlXM+JinBePaWaYzqnel1JBFpYeGmUiprGvniRA1BERGJBOuKq3h1TQk3nNJfXZYl7HQGK8J9uLWcmx9dRpeMVJ78yhT6dm7vdSQROcwLK4vJ7pDKKYO6eh1FRESAuxfkkdkumZtPH+B1FIlDOoMVwV5aVcz1Dy+ld6d2PPfVaSquRCJQZU0Di/LKuGhMLxIT1NxCRMRrn2zby7uby/j6jIFkttPZKwk/ncGKQM457lu0lbsX5DG5f2fmXDOBrPa6rkMkEs1fu5sGn5+Lx6q5hYiI15xz/H5BHt07pnLdyblex5E4pQIrwjT5/Pz0pXU8vaSQWWN78fsvjCY1SZ1vRCLVi6uKGZCdzqjeajwjIuK1F1cVs3znPn576UmkJWv/SbyhAiuCVNU0cuvTK1i8pZxvnDGQ750zlAQNORKJWMWVtSzZXsH3zhmiua9ERDxWXdfIr+dvYkyfTDUdEk+pwIoQm3bvZ/ZjyympquW3l57ElZNzvI4kIkfx0qpiAGZpeKCIiOf+8tYWyg/U89C1E3WAWjylAisCzF9bwu3PrSYjNYm5s6cxoV8nryOJyFE453hxZTET+nUip4sa0IiIeGljyX7++eEOrpzUlzF9s7yOI3FOXQQ91Ojz89vXNvH1J1cwrEcHXr7tVBVXIlFiY0k1m/cc4OKxvbyOIiIS1xqa/Hz32dV0ap/C988b5nUcEZ3B8krRvhq++fRKVhRUctXkHH5+0Qg1sxCJIi+tKiYpwfjcaBVYIiJe+us7W9hYsp8Hr51Ip3R1XRbvqcDywGtrS/jBv9bgHPz1qnF8fox20ESiid/veGnVLqYPyaazPsxFRDyzqrCS+xZt5dLxvTlnRHev44gAKrDC6mB9E7+av5GnPilgTJ9M/nrVeF27IRKFPt6+l9376/jR54Z7HUVEJG5V1jTwjSdX0L1DKj/7/Eiv44h8SgVWmHy8bS//8/xqivbVMvv0Adx+7lBSknQJnEg0emnlLtJTEjlnuI6Wioh4we93fOeZVZRW1/HsLdPIbJfsdSSRT2kPP8RqG3z878vruXLOxySY8ewt0/jRBcNVXIkEmNlMM8szs3wzu6OV+3PMbKGZrTSzNWZ2QYv7RpvZR2a23szWmllaqPPWNfqYv7aE80b1oF2KrpsUEfHCX9/JZ2FeGT+9cATjctQgTCKLzmCF0Luby/jpi+soqKjh+pNz+f7MobRP0UsucoiZJQL3AucARcBSM5vnnNvQYrWfAM865/5uZiOA+UCumSUBTwDXOOdWm1kXoDHUmRduKqW6vomLNfeViIgnXlhZxJ/e2swl43pzzdR+XscR+S/a2w+BPfvruOuVDby6poQB2ek8ffNUpg3s4nUskUg0Gch3zm0DMLO5wCygZYHlgI6B25nArsDtc4E1zrnVAM65veEI/OKqYrI7pHKy/qZFRMLuw/xyvv/8GqYO6MxvLzsJM00oLJFHBVYQNTT5efzjnfzpzc00+Px875whzJ4+QO3XRY6sN1DY4vsiYMph6/wceMPMbgPSgbMDy4cAzswWANnAXOfc71t7EDObDcwGyMnJOe6wVTWNLNxUxpen9iMpUcN8RUTCaUXBPm55fDn9u6bzwDUTtX8lEUsFVhA453hzwx5+89omtpcf5PQh2dx10Uhyu6Z7HU0kFlwFPOKc+4OZTQMeN7NRNG+/TgUmATXA22a23Dn39uE/wDk3B5gDMHHiRHe8QeavK6HB5+eScRoeKCISTst3VnDdw0vpkpHCozdOVlMLiWgqsE7Q2qIqfj1/Ix9t28vA7HT+ef0kZgzN1ilrkbYpBvq2+L5PYFlLNwEzAZxzHwUaWXSl+WzXe865cgAzmw+MB/6rwAqWF1YWMyA7nVG9Ox59ZRERCYoP8suZ/dgyunVM4+mbp9IjM+T9jEROiAqs47R+VxV/fmsLb27YQ+f0FH4xayRXTs4hWcOGRI7FUmCwmfWnubC6Erj6sHUKgLOAR8xsOJAGlAELgO+bWXugAZgO/ClUQYsra1myvYLvnjNEB1BERMLk2WWF/OjfaxmQnc7jN02he0cVVxL5VGAdo7zd1fz5rc28tm43HdKS+O45Q7j+lFw6pulUtcixcs41mdmtNBdLicDDzrn1ZnYXsMw5Nw/4HvCgmX2H5oYX1zvnHLDPzP5Ic5HmgPnOuVdDlXXequbeGuoeKCISeo0+P3cvyGPOe9s4dVBX7vvyeO1rSdRQgdVGqwsrmfPeNuavKyEjJYlvnTWYG0/trzHAIifIOTef5tbrLZfd2eL2BuCUI/zfJ2hu1R5yL64sZnxOFjld2ofj4URE4tauylpue3oly3fu48tTc/jZ50dqhJBEFRVYn8Hnd7y1cQ8PLd7G0h376JCaxNdnDOTm0waQ1T7F63giEiYbS/aTt6eaX8wa6XUUEZGY5ZzjmaWF/Gr+Rvx+x/9dNY6LxvTyOpbIMVOB1YrqukZeWFnMw+9vZ8feGnpnteOnF47gikl9yUjVSyYSb15cWUxSgvG50fqgFxEJhc17qvnfl9fzQf5epvTvzO8uG61uzBK1VC0EOOdYWVjJ3CUFvLy6hNpGH2P7ZnHvecM4b2R3zXkjEqf8fse81bs4fUg2ndN15lpEJJhK99fxl7e38PSSAjJSk/jlxaO4enIOCQlqJiTRK+4LrH0HG3hxVTFzlxSSt6ea9imJXDSmF1dO7su4nE5exxMRj32yvYKSqjp+eMFwr6OIiMSMwooaHly8jblLC/H5HddM7ce3zh6iA1kSE+KywDpQ38Qb63fz8updLN5STpPfMaZPJr+59CQ+P6aXhgGKyKdeXFlMekoi5wzv7nUUEZGo5vM73s8v54mPd/L2xj0kmPGFCX34+oxBaiAkMSVuKokD9U28t7mMV9bs4u2NpdQ3+emd1Y6bTuvPrDG9GdFLE4eKyH+qa/Qxf10J543qQbuURK/jiIhEHb+/+RKM19eV8NKqXZRW19M5PYWvzRjIl6b0o1dWO68jigRdTBdYuypreXvjHt7aWMpHW/fS4PPTNSOFKyf15aKxvRjXt5PG+IrIES3cVEp1XZPmvhIROQal1XV8tHUvi7eU8+7mMsqq60lKMGYMzebS8X04c1g30pJ10EpiV0wVWLUNPpbuqOCDreUs3lzOhpL9AOR2ac91J/fj7OHdmdCvkxpWiEibvLiqmK4ZqZw8sIvXUUREIlJ1XSN5u6tZV1zFmuIqVuzcx469NQBktkvm1EFdOWdEd84Y1k1zh0rciOoCq8nnZ3VRFR/ml/PB1nJW7KykwecnOdEY17cTd5w/jLOHd2dgdjpmOlMlIm1XVdPIwk1lfHlqPx2UEYkjZjYT+AuQCDzknPvtYfenAo8BE4C9wBXOuR3hzhkuzjn21zaxq6qWon21FFbUUFBRw7byg2wtPUBxZe2n63bNSGVs3yyunpLD5P5dOKl3JokaKSRxKKoKrIqDDazYuY8VBc1fa4qqqGnwATCiZ0euO7kfJw/qyuTczqSrUYWInID560po8Pm5eJzmvhKJF2aWCNwLnAMUAUvNbJ5zbkOL1W4C9jnnBpnZlcDvgCvCn7btnHM0+hy1DT5qGps4WO/jYH0TB+qbqK5rZH9tE1W1jVTWNlBxsJGKg/XsPdBA2YF6SvfXU9vo+4+fl56SSP/sdCbmduLq7jkM7d6BUb0z6d4xVQe0RYiCAuuD/HL+tbyIFQX//5RzYoIxomdHLp/Qh0n9OzNtQBe6ZKR6nFREYsmLK4sZkJ3OSb0zvY4iIuEzGch3zm0DMLO5wCygZYE1C/h54PbzwN/MzJxz7kQffMn2Cl5bV4Lf7/A5h981N4nwHfpyjia/w+dzNPn9NPkdTT5Hg89Pk89Pg89PY5OjvslHQ5OfuiY/9Y0+aht9+NuQLjHB6NQ+mc7pKXROT2FMnyyyO6TSMzONHplp9M5qR9/O7emSnqJCSuQzRHyBlV96gPe2lDEupxNXTMphfE4Wo/tkqaOXiIRMbYOPXVW1XD6hr3YiROJLb6CwxfdFwJQjreOcazKzKqALUN5yJTObDcwGyMnJadOD55ce4PnlRSSYkZhgJJiRYJCUYCQkNC9L+vTfBJISjeTEBJISjPTUJLISE0hONFKTEklJSiAtOYHUpETaJSeSlpxAu5Qk0lMSaZeSSEZqEhmpSXRIS6ZDWhJZ7ZPJSE3SNk8kCCK+wLp6Sg7XTuunP3gRCZt2KYm8e/sZNPj8XkcRkSjlnJsDzAGYOHFim85uXT0lh6untK0YE5HIFfFXbicnJqi4EpGwS0gwtREWiT/FQN8W3/cJLGt1HTNLAjJpbnYhIgJEQYElIiIiEiZLgcFm1t/MUoArgXmHrTMPuC5w+wvAO8G4/kpEYkfEDxEUERERCYfANVW3AgtobtP+sHNuvZndBSxzzs0D/gE8bmb5QAXNRZiIyKfaVGBpTggRERGJB865+cD8w5bd2eJ2HXB5uHOJSPQ46hDBFnNCnA+MAK4ysxGHrfbpnBDAn2ieE0JERERERCSutOUarE/nhHDONQCH5oRoaRbwaOD288BZps4UIiIiIiISZ9pSYLU2J0TvI63jnGsCDs0J8R/MbLaZLTOzZWVlZceXWEREREREJEKFtYugc26Oc26ic25idnZ2OB9aREREREQk5NpSYGlOCBERERERkTZoS4GlOSFERERERETa4Kht2jUnhIiIiIiISNu0aR4szQkhIiIiIiJydGFtciEiIiIiIhLLzKtLpcysDNgZ4ofpCpSH+DGUQRmU4T/1c85FbJvQY9z2RMLvK1z0XGNTvDzXaNvuRNLvRVlapyytU5b/1Oq2x7MCKxzMbJlzbqIyKIMyRFaGaBFPr5Wea2yKp+caTSLp96IsrVOW1ilL22iIoIiIiIiISJCowBIREREREQmSWC+w5ngdAGU4RBmaKUN0iafXSs81NsXTc40mkfR7UZbWKUvrlKUNYvoaLBERERERkXCK9TNYIiIiIiIiYaMCS0REREREJEhipsAysx1mttbMVpnZssCyzmb2ppltCfzbKYSPPzTw2Ie+9pvZt83s52ZW3GL5BSF47IfNrNTM1rVY1upzt2b/Z2b5ZrbGzMaHMMPdZrYp8DgvmFlWYHmumdW2eE3uD2GGI77+ZvbDwOuQZ2bnhTDDMy0ef4eZrQosD9Xr0NfMFprZBjNbb2bfCiwP63simpnZzMD7It/M7vA6Tyi1tu2MJceyfYxmx7r9k9Azs8sD22C/mU1ssfwcM1se+LtbbmZnHuH/B+33d6QsgfuO+lloZv3N7JPAes+YWcrxZjns57b6+djKeiHfTrX19Q7H58OR9p9aWS9kr8vRnqeZpQZ+f/mB90ZuMB+/xeO0uk9z2DozzKyqxe/uzlBkOSbOuZj4AnYAXQ9b9nvgjsDtO4DfhSlLIrAb6Af8HLg9xI93OjAeWHe05w5cALwGGDAV+CSEGc4FkgK3f9ciQ27L9UL8OrT6+gMjgNVAKtAf2AokhiLDYff/AbgzxK9DT2B84HYHYHPg+Yb1PRGtX4G/363AACAl8D4Z4XWuED7f/9p2xtLXsWwfo/nrWLZ/+grb72Q4MBRYBExssXwc0CtwexRQfIT/H7Tf32dkadNnIfAscGXg9v3A10Lwen36+djKfSHfTrXl9Q7X5wNH2H8K1+vSlucJfB24P3D7SuCZEP1eWt2nOWydGcAroXx/HOtXzJzBOoJZwKOB248CF4fpcc8Ctjrndh51zSBwzr0HVBy2+EjPfRbwmGv2MZBlZj1DkcE594Zzrinw7cdAnxN9nGPN8BlmAXOdc/XOue1APjA5lBnMzIAvAk+f6OMcJUOJc25F4HY1sBHoTZjfE1FsMpDvnNvmnGsA5tL8GkkUOsbtY9Q6xu2fhIFzbqNzLq+V5Sudc7sC364H2plZqhdZaMNnYeCz60zg+cCioP/NhOvzMQjC8vkQ7v2nVrTlebbcjj4PnBX4PQbVZ+zTRLRYKrAc8EbgdPvswLLuzrmSwO3dQPcwZbmS/9xI3Bo4zftwGIeiHOm59wYKW6xXRHjeqDfSfJbkkP5mttLM3jWz00L82K29/l68DqcBe5xzW1osC+nrEDhlPw74hMh7T0SqeHs9Wtt2xjqvPhu84MXnj7TdZcAK51z9Ee4P9e+vLdu7LkBlix3+UGwTW/t8bClc26mjvd5efD4cvv/UUqhel7Y8z0/XCbw3qmh+r4TMYfs0h5tmZqvN7DUzGxnKHG0RSwXWqc658cD5wDfM7PSWd7rmc4gh70kfGJd8EfBcYNHfgYHAWKCE5lPgYRWu534kZvZjoAl4MrCoBMhxzo0Dvgs8ZWYdQ/Twnr/+LVzFfxbeIX0dzCwD+Bfwbefc/pb3ef2ekIjymdvOWBfjfwuRtP2LSWb2lpmta+XrqGc1AjuBvwNuOcIqx/T7O5EsodTGXId/Ph4uKNupo2QJ699LW16XVvafDhc32+/P2qcBVgD9nHNjgL8CL4Y53n9J8jpAsDjnigP/lprZCzSf3txjZj2dcyWBIU+lYYhyPs1Ho/YE8uw5dIeZPQi8EoYMcOTnXgz0bbFen8CykDCz64ELgbMCOzIEjtTVB24vN7OtwBAg6BeufsbrH+7XIQm4FJjQIlvIXgczS6Z5Q/Skc+7fgcUR8Z6IAnH1ehxh2/met6lCzovPhrDz8PMnbjjnzj6e/2dmfYAXgGudc1uP8LOP6fd3nFnasr3bS/PQ8aTAmYpj2iYeLVdrn4+t/IygbKfa+hp9xusdtM+HNrwu13PY/lMrPyNU2++2PM9D6xQFfoeZNL9Xgu4I+zSfallwOefmm9l9ZtbVOVceijxtERNnsMws3cw6HLpN88WB64B5wHWB1a4DXgpDnP84CnPYtSyXBHKFw5Ge+zzgWms2FahqMVQmqMxsJvB94CLnXE2L5dlmlhi4PQAYDGwLUYYjvf7zgCutuQtO/0CGJaHIEHA2sMk5V9QiW0heBzMz4B/ARufcH1vc5fl7IkosBQZbc9esFJqH/M7zOFNIfMa2M9Z58dkQdh5+/shnsOaOcK/S3Gjlg89YLxy/v6N+FgZ27hcCXwgsCvbfzH99PrYUru1UG1/vsHw+HGn/6bB1Qvm6tOV5ttyOfgF450iF4In4jH2aluv0CKyHmU2mub4JSbHXZp/VASNavmjucrI68LUe+HFgeRfgbWAL8BbQOcQ50mn+hWa2WPY4sBZYQ/ObsWcIHvdpmk9nN9I8TvamIz13mjvF3Utzd5i1tOgmFIIM+TSPz10V+DrUbeaywO9pFc2ndT8fwgxHfP2BHwdehzzg/FBlCCx/BPjqYeuG6nU4leYhT2tavPYXhPs9Ec1fgddrc+A1+bHXeUL4PFvddsbS17FsH6P561i3f/oKy+/kksDvoh7YAywILP8JcLDF9nkV0C1w30OHtsHB/P0dKUvgvlY/C4H5/P9uhwNoLrzyab4EIjWIr1Nrn4+9gPktHjvk26kjvd4tswS+D/nnA0fefwrb69La8wTuornoA0gLvBfyA++NASF6LY60T/PVQ+8b4NbAa7Ca5qYgJ4ciy7F8WSCYiIiIiIiInKCYGCIoIiIiIiISCVRgiYiIiIiIBIkKLBERERERkSBRgSUiIiIiIhIkKrBERERERESCRAWWiIiIiIhIkKjAEhERERERCZL/B8/nxbjyJrqtAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 864x360 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import statsmodels.api as sm\n",
    "import pandas as pd\n",
    "\n",
    "d = pd.read_csv(\"diabetes.csv\")\n",
    "d\n",
    "\n",
    "\n",
    "X = d.drop(columns=[\"Outcome\"])\n",
    "X = sm.add_constant(X)\n",
    "\n",
    "averages = X.mean()\n",
    "print(averages)\n",
    "\n",
    "model = sm.Logit( d.Outcome, X)\n",
    "model = model.fit()\n",
    "\n",
    "print(model.summary())\n",
    "\n",
    "def fromLogOdds2prob(vs,var,model):\n",
    "    def invlogit(x):\n",
    "        import numpy as np\n",
    "        e = np.exp(x)\n",
    "        return e/(1+e)\n",
    "\n",
    "    probs = []\n",
    "    for v in vs:\n",
    "        x = averages\n",
    "        x[\"{:s}\".format(var)] = v\n",
    "\n",
    "        logodds = model.params.dot(x)\n",
    "        p = invlogit(logodds)\n",
    "        probs.append(p)\n",
    "    return probs\n",
    "\n",
    "\n",
    "glucoseLevels = np.arange(40,210)\n",
    "probsGlucose = fromLogOdds2prob( glucoseLevels, \"Glucose\", model )\n",
    "\n",
    "pregs = np.arange(0,20,1)\n",
    "probsPregs = fromLogOdds2prob( pregs, \"Pregnancies\", model )\n",
    "\n",
    "dpfs = np.linspace(-12,4,100)\n",
    "probsdpfs = fromLogOdds2prob( dpfs, \"DiabetesPedigreeFunction\", model )\n",
    "\n",
    "fig,axs = plt.subplots(1,3)\n",
    "\n",
    "axs[0].plot(glucoseLevels, probsGlucose )\n",
    "axs[1].plot(pregs, probsPregs )\n",
    "axs[2].plot(dpfs, probsdpfs )\n",
    "\n",
    "fig.set_size_inches(12,5)\n",
    "fig.set_tight_layout(True)\n",
    "plt.show()\n",
    "\n",
    "predictions = model.predict(X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 146,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimization terminated successfully.\n",
      "         Current function value: 0.502215\n",
      "         Iterations 6\n",
      "                           Logit Regression Results                           \n",
      "==============================================================================\n",
      "Dep. Variable:                Outcome   No. Observations:                  768\n",
      "Model:                          Logit   Df Residuals:                      765\n",
      "Method:                           MLE   Df Model:                            2\n",
      "Date:                Wed, 10 Nov 2021   Pseudo R-squ.:                  0.2235\n",
      "Time:                        22:56:02   Log-Likelihood:                -385.70\n",
      "converged:                       True   LL-Null:                       -496.74\n",
      "Covariance Type:            nonrobust   LLR p-value:                 5.967e-49\n",
      "==============================================================================\n",
      "                 coef    std err          z      P>|z|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "const         -7.5156      0.605    -12.416      0.000      -8.702      -6.329\n",
      "Glucose        0.0352      0.003     10.693      0.000       0.029       0.042\n",
      "BMI            0.0763      0.013      5.723      0.000       0.050       0.102\n",
      "==============================================================================\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/usr/local/lib/python3.9/site-packages/statsmodels/tsa/tsatools.py:142: FutureWarning: In a future version of pandas all arguments of concat except for the argument 'objs' will be keyword-only\n",
      "  x = pd.concat(x[::order], 1)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEKCAYAAAAVaT4rAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAACIXElEQVR4nO2dd5hTxdfHv5Nk07NLL9JRBLGhIh0bUkTFgthFsWDvYkVF/dnFhmJ/RUUURUQERREQBSmColIVlN7r9pZ83z9Osmn3ZrM1C8zneeZhubnl5ObeMzNnTlEkodFoNJqDB0uqBdBoNBpN9aIVv0aj0RxkaMWv0Wg0Bxla8Ws0Gs1Bhlb8Go1Gc5ChFb9Go9EcZFSZ4ldKtVVKLYlomUqpO5RSdZRS05VS/wT/rV1VMmg0Go0mHlUdfvxKKSuATQA6A7gZwG6Szyil7gdQm+R9VS6ERqPRaABUn6mnF4A1JNcBOAfAB8HtHwA4t5pk0Gg0Gg0AWzVd52IAnwT/bkhyS/DvrQAaGh2glBoKYCgAeDyeE9q1a1flQlaEoqIirF69Grm5uWjatCkaNjT8WhqNRlNtLF68eCfJ+rHbq9zUo5SyA9gM4EiS25RSe0nWivh8D8mEdv6OHTty0aJFVSpnZZCTk4PBgwdj4sSJuP766zFq1CikpaWlWiyNRnOQopRaTLJj7PbqMPWcAeA3ktuC/9+mlGocFKoxgO3VIEO14PF48Pnnn+P+++/HW2+9hf79+2Pv3r2pFkuj0WiiqA7FfwnCZh4AmAzgyuDfVwL4qhpkqDYsFguefvppvP/++5g9eza6du2K1atXp1osjUajKaFKFb9SygOgN4CJEZufAdBbKfUPgNOD/z/guOqqq/DDDz9g+/bt6Ny5M3766adUi6TRaDQAqljxk8whWZfkvohtu0j2ItmG5Okkd1elDKnkpJNOwoIFC1C/fn2cfvrpGDNmTKpF0mg0Gh25W9UcdthhmDdvHk466SQMGTIEDzzwAAKBQKrF0mg0BzFa8VcDtWvXxrfffovrr78ezzzzDAYNGoScnJxUi6WpShYuBG64ARg8GPj6a0B39poahFb81URaWhreeOMNvPzyy5g0aRJOOukkbNq0KdViaaqCZ54BTj0VeOcd4KOPgEsuAS64QCt/TY1BK/5qRCmF22+/HZMnT8bff/+NTp064bfffku1WJrKZPNm4LHHgNzcsKLPyQGmTwe+/z61smk0QbTiTwFnnnkm5s6dC5vNhp49e+LLL79MtUiayuKHHwCrNX57djYwcWL8do0mBWjFnyKOOeYYLFy4EEcffTTOP/98PPvss9CF7w8A3G7AYvBaWa2Az1f98mg0BmjFn0IaNmyIWbNm4eKLL8b999+Pa665BoWFhakWS1MR+vc33m63A1ddVa2iaDRmaMWfYlwuF8aNG4dHH30U77//Pnr37o1du3alWixNeXG7gSlTgPR0aT4f4HQCI0cCRx+dauk0GgDVlI+/ouwvSdoqyrhx43D11VejadOmmDJlCmp6RlJNAvLzZTE3Px/o1QuoWzfVEmkOQlKZpE2TJJdeeilmzZqFrKwsdO3aFTNmzEi1SJry4nQCAwYAF16olb6mxqEVfw2ja9euWLBgAZo2bYq+ffvi7bffTrVIGo3mAEMr/hpIy5YtMXfuXPTp0wfXX3897rrrLvj9/lSLpdFoDhC04q+hpKenY/Lkybj99tvx0ksv4ZxzzkFWVlaqxdJoNAcAWvHXYGw2G15++WWMHj0a06ZNQ/fu3bF+/fpUi6XRaPZztOLfD7jxxhvxzTffYN26dejUqRMWLFiQapE0Gs1+jFb8+wl9+vTBvHnz4PF4cMopp2D8+PGpFkmj0eynaMW/H9G+fXssWLAAHTt2xMUXX4zHH39cp3nQaDRlRiv+/Yx69erhhx9+wODBg/Hoo4/i8ssvR35+fqrF0mg0+xG2VAugKTsOhwNjxoxBu3bt8OCDD+K///7DpEmT0KBBg1SLptFo9gP0iH8/RSmFBx54AJ9//jmWLFmCzp07Y+nSpakWS6PR7Adoxb+fc8EFF+Cnn35CQUEBunXrhm+//TbVImk0mhqOVvwHAB07dsTChQtx6KGH4qyzzsKoUaNSLZJGo6nBaMV/gNC0aVP8/PPPOPvss3Hbbbfh5ptvRnFxcarF0mg0NRCt+A8gvF4vJk6ciHvvvRejR4/GmWeeib1796ZaLI1GU8PQiv8Aw2Kx4Nlnn8W7776LmTNnolu3bvj3339TLZZGo6lBaMV/gHLNNddg+vTp2Lp1Kzp37ow5c+akWiSNRlND0Ir/AOaUU07BggULUKdOHfTq1QsfffRRqkXSaDQ1gCpV/EqpWkqpCUqplUqpFUqprkqpOkqp6Uqpf4L/1q5KGQ522rRpg3nz5qF79+4YPHgwhg8fjkAgkGqxNBpNCqnqEf8rAKaRbAfgWAArANwPYAbJNgBmBP+vqULq1KmD7777Dtdeey2efPJJXHTRRcjNzU21WBqNJkVUmeJXSmUAOAnAewBAspDkXgDnAPgguNsHAM6tKhk0YdLS0vD2229j5MiR+OKLL3DyySdjy5YtqRZLo6ke1q8H7rgD6N4duP564O+/Uy1RSlFVld1RKdUBwNsAlkNG+4sB3A5gE8lawX0UgD2h/5vRsWNHLlq0qErkPBj5+uuvcckll6B27dr4+uuv0aFDh1SLpNFUHcuWAd26AXl5QFERYLUCTicwfTrQtWuqpatSlFKLSXaM3V6Vph4bgOMBvEHyOAA5iDHrUHodw55HKTVUKbVIKbVox44dVSjmwcfZZ5+NuXPnAgB69OiByZMnp1gijaYKufNOICtLlD4A+P1ATg5w442plSuFVKXi3whgI8lQuagJkI5gm1KqMQAE/91udDDJt0l2JNmxfv36VSjmwcmxxx6LhQsXon379jj33HPxwgsv6Nz+mgOTOXMAo2f7r7+AgoLql6cGUGWKn+RWABuUUm2Dm3pBzD6TAVwZ3HYlgK+qSgZNYho3bowff/wRF1xwAYYNG4brrrsOhYWFqRZLo6lcfD7j7Q4HkJZWvbLUEKraq+dWAB8rpf4E0AHAUwCeAdBbKfUPgNOD/9ekCLfbjU8//RTDhw/He++9h759+2L37t2pFkujqTxuuQVwu6O3OZ3AVVcBloMzlKnKFncrE724Wz2MHTsW11xzDVq0aIEpU6bg8MMPT7VIGk3FKS4GrrkG+OwzGeUXFAB9+gCffgq4XKmWrkpJxeKuZj/j8ssvx8yZM7Fnzx506dIFs2bNSrVINYqcHGDLFkDHv+1n2GzABx8Aa9YAX3wBrFwJfPXVAa/0E6EVvyaK7t27Y+HChWjcuDH69OmDd999N9UipZzcXGDwYKBuXaB1a6BJE2DixFRLpSkzhxwC9OoFtGiRaklSjlb8mjhatWqFX375Bb169cJ1112HYcOGwe/3p1qslHHFFcDnn4uFID8f2LpVts2bl2rJNJryoRW/xpCMjAxMmTIFN998M1544QWcf/75yM7OTrVY1c7WrcA334jCjyQvD3hGuyVo9lO04teYYrPZ8Nprr2HUqFGYMmUKevTogQ0bNqRarGpl0ybAbo/fTgKrV1e/PBpNZaAVv6ZUbrnlFkydOhX//fcfOnXqhF9//TXVIlUbhx8eDviMxGYDevSofnk0mspAK35NUvTr1w+//PILnE4nTjrpJEyYMCHVIlULPh9w772AxxPeZrGIW/j9Oq+sZj9FK35N0hx55JFYsGABjj/+eAwaNAhPPfXUQZHm4dFHgdGjgfbtgfr1gfPOA379FWjVKtWSaTTlQwdwacpMfn4+rr32Wnz88ce44oor8M4778DhcKRaLI1GE4MO4NIkz9y5wNlnA0cdBdx8MxCzoOt0OvHRRx/hiSeewEcffYRevXpBZ1DVaPYf9IhfE82nn0p4e6hCV1qaGLh/+83QtvHZZ5/hyiuvROPGjTFlyhS0b9++mgXWaDRm6BG/pnT8fuDWW8NKHxCXlsxMMXQbcOGFF2L27NnIzc1F165d8d1331WTsBqNprxoxa8Js2FDtNIPEQgAM2eaHtapUycsXLgQLVu2xJlnnonRo0dXoZAajaaiaMWvCVO7toz6jWjYMOGhzZs3x5w5c3DGGWfg5ptvxm233Ybi4uIqEFKj0VQUrfg1YTIyZFE31kPH7RZn9lLw+XyYNGkS7rrrLowaNQoDBgxAZmZmFQmr0WjKi1b8mmjefx84/XQpVJGeLqlrH3gAuPDCpA63Wq0YOXIk3n77bUyfPh3dunXDf//9V8VCa2okpDgFzJsH6MpuNQrt1aMxZvNmae3aAV5vuU4xc+ZMDBw4EGlpaZg0aRK6detWyUKWkZwcWcdo0sS8HJ+mcvjzT5k97t4NKAVYrcBHHwFnnZVqyQ4qtFePpmwccgjQsWO5lT4AnHbaaZg/fz4yMjJw2mmnYdy4cZUoYBkggYcekrDbTp2ABg2A224zX8/QVIyCAuC004D164HsbCArC9i7F7joImDt2lRLp4FW/Joqpm3btpg/fz66dOmCyy67DI888ggC1V3C6pVXgJdfllzKWVmSY/m994DHH69eOQ4Wpk0zNu0UFwP/93/VL48mDq34NVVO3bp18f333+Pqq6/GE088gUsuuQR5eXnVJ8Dzz8e7qebmSmdQE0yd69aJWcRul4X0q68G9u1LtVTlZ+dO4/qUhYXY/l8OLrpIlpCcTpkEbN9e/SIe7Ggbv6baIIkXXngB9913H0488UR89dVXaNSoUdVf2OEwHoEqJWaJtLSql8GMrCzgsMOAXbvCpie7HTjySGDxYpFxf2P1auDoo+Oq1xR5aqGtZwM27PYi5OlrswHNmgGrVqX2ZzhQ0TZ+TcpRSmHYsGGYOHEili5dik6dOuHPP/+s+gt36GC8/fDDU69txo6VRefI9YbCQuCff4Cff06dXBXhsMOAq66KzmXtdmNy05uwM8+DyPCO4mKZIHz1VbVLeVCjFb+m2jn33HMxZ84cBAIBdO/eHVOmTKnaC770kphQIkfPbjcwalTVXjcZliwRxR+L3w8sX17t4lQao0cDH3wgrsHduwPPP48VFz+G7Oz4GUx2NrBiRQpkPIjRil+TEo477jgsXLgQbdu2xYABA/DSSy9VXW7/bt2AOXOAAQOAli2Bfv2AGTOA3r2r5npl4dhjpROKxWoFjjii+uWpLJQCBg4Epk+Xe3/TTWh3lC1qEhDC692/v+r+iLbxa1JKbm4urrjiCkycOBFDhw7Fa6+9hrRUm1+qk8zMsI0/tCBqt4sm/P33/dPGb0JhIdC2LbBxI6Js/E2bio3fqLaxpmJoG7+mRuJ2u/H555/jgQcewNtvv40zzjgDe/bsSbVY1Ud6OrBwocxCbDZxdbn4YuDHHw8opQ+IYp83Dzj3XPnbbpe/583TSr+60SN+TY3hww8/xLXXXotWrVph6tSpOOyww1ItUvVCHnDK3oyQ2jlIvm7K0CN+TY1n8ODBmDFjBnbt2oXOnTtj9uzZqRapekm1Fly3TpLxnX22xD5U4cxLqdR/3YOZKh3xK6XWAsgC4AdQTLKjUqoOgPEAWgJYC+BCkgmfMD3iP7hYs2YNzjrrLKxZswZvvfUWhgwZkjphtm4V75QNG4CTTxbbxIG4BvHLL0DfvhLXUFQkyfnS0yWWoEmT8p2TFLfU/HyJS7BaK1dmTamYjfhBssoaRLHXi9n2HID7g3/fD+DZ0s5zwgknUHNwsWfPHvbu3ZsAeN9999Hv91e/EHPmkB4P6XSSAOn1kh06kNnZ1S9LRSksJDdsIPPy4j8LBMjDD5fvGNlsNvKqq8p3vZUryXbtSLdb7luDBuQPP1TsO2jKDIBFNNLNRhsrq5ko/lUAGgf/bgxgVWnn0Yr/4KSwsJA33HADAfDcc89ldnUq3ECAbNYsXhk6neQTT1SfHBUlECBfeIH0+UQJu93kvfeSkR3pzp2k3R7/XQGyXr2yX7OwkGzYkFQq+lweD7l+feV9N02pmCn+qrbxE8D3SqnFSqmhwW0NSW4J/r0VgGFpJ6XUUKXUIqXUoh07dlSxmJqaSFpaGkaPHo1XXnkFkydPRs+ePbFp06bqufg//4iLZSz5+cDHH1f99TMzgW+/lejdimQRff994JFHJDVEbq60114DnngivI/TaZ6zqDzZWb/7Tq4Te06dpK3GUNWKvwfJ4wGcAeBmpdRJkR8GeyTDJ47k2yQ7kuxYv379KhZTU1NRSuG2227D5MmT8c8//6BTp05YvHhx1V/YbjdXhrEVyiqb994DGjUSt84zz5RkNn/8Ub5zPfmkcYK6F18Mxw14POJOGrt24XYDN91U9mtu3WrcWRUUiBO/JuVUqeInuSn473YAXwLoBGCbUqoxAAT/1bn5NKVy5pln4pdffkFaWhpOOukkTJw4sWov2LIlcOih8a4nbjcwdKjhIZXCkiXArbdKCunMTBmpb9kC9OkDlKeG8ZYtxttzcqKTqI0ZAxxzjHQC6ekyCzjnHODOO8P77NsHvPOOdCY//WTeMXbvbvyZ1yspHDSpx8j+UxkNgAeAL+LvXwD0A/A8ohd3nyvtXNrGrwmxdetWdunShQD49NNPMxAIVN3FVq0SW7XPR7pcYh8//3yyqKjqrnnjjaTVGm9r9/nI774r+/k6dza23TdvLvb/SAIBctEicsIEcvXq6M8WLCDT0+UeWCxir+/XT+z5RlxxhewTup7LRR53XMn+hYXmh9ZY/H7jxfEaDKp7cRdAawB/BNsyAA8Ft9cFMAPAPwB+AFCntHNpxa+JJDc3lxdffDEB8KqrrmJBQUHVXayggPzyS3L0aPL336vuOiEGDTJW1D4f+dlnZT/fzz+Lso48l9styj1ZAgGyadN4mdxu8s03jY/x+8kxY6TjOeYY8plnyJwcrl9P9u0rDkM2m/xd49d7CwvJu+6SjsxiIdu0Ib//PtVSJUW1K/7KbFrxa2IJBAIcMWIEAfCkk07ijh07Ui1S5fDxx9Ej5Uhvoq1by3fO+fPJ008Xl8quXcnp08t2/J9/ikumUYd04olJnyY/n2zSJHpCY7WShxwin9VYhgyRGUtsp/frr6mWrFTMFL+O3NXslyil8Oijj+KTTz7BggUL0KVLF6xcuTLVYlWcQYMkY2cojaVSsq4wfDjQ0NABrnQ6d5Ysmdu2SaBWiuzsX34pywSR675+vyxlfPllSkQqnV27gHHjZM0lkrw84KmnUiNTJaAVvyYl5OfHFWgqFxdffDF+/PFHZGVloUuXLvjh66/37yLqaWnArFnA668DZ5whtQmnTpVi8aniyCOBjIz47aEykUmyerVx6YGcHPGerZFs2GDsxUXu1/UStOI/gNm6VfTFaaeJo8jq1amWCFi7VuTx+aT17g2sX1+xc3bp0gULH38czXJz0W/AALzldgP33Vc+L5iagN0OXHkl8M03wCefAKecklp5LBbgiy/E28fjkf97PEDPnsA11yR9mqOPhmE+fo9HPquRtG5tXLbTagU6xmdC2G8wsv/UtKZt/GVn9Wqydu1wQKbVKqbjOXNSJ1NurjjJWCyVbOOdNYt0u7kPYH+JC+EdNhuLb7qpskTXkOTeveTbb0vk8uzZ8V5BpVBUJFkcIoOE7XaybdvkHaWWLCFPPplMSyPr1iVHjKhaJyuS5D33xC+Qezzk8uVVfOGKA724e3Bx7rnGa3FHHJE6mT76yHiN0Ocj33lHPBkbNSJbtCD/9z9xqEmKU04pOVkxwDuCyv9Mi4WZmzdX5VeqGHv2kLffTjZuLF4zDz8svWN5KS4mJ08mH3qIfOstct++Cov35pvk8OHklCly+oqyaxd53XVkRoa0666TbcmwenX88+N2lz+dUNL4/eRLL8lv5HLJ87Z4cRVftHLQiv8gw+EwVvxKkZmZqZFpxAhjmSwWslYtGcVFun2feWaSJzZwNXwDoBXg0W3bcu3atebHFheTU6eSr75K/vhjmUex5aagIH7463SSPXuWT4bsbPL448Oa0e2Wm/rHH+US77ffRDGHBrperzjw5OSU63TmZGcnfdIbbhAX0Njnx+Egt2ypZLkOEMwUv7bxH6CYmbdJMdGmgmOOCad+scCPyzAW36EP7uRIFGYVoKgovG9enqxxJpWp4Ljj4iJsbwDwrdOJ9Vu3olOnTpg/f378cVu3Si3Aiy8G7rlHvF3S04G775bPqpJJkyR9QaT9OD9fyi3OnVv28z37rCw2ZmfL/3Nzgb17gUsvLfOpSLkl+/aFsz1kZwN//SWZHiqF1auBHj2AWrVk4bhXL1lITcCiRcbPtdMJ/P13Jcl1kKAVfzVCijK75hpps2aZR71XFJvNeLtS4RQt1c3ZZ0t9VXsaMQnn4E3cgD6Yjp2sg1x/vOeEUpIOHpAU8Z98Alx2GXDXXcCKFRE7Pv645I+PxONB74cfxvz58+H1enHKKafg008/jdplw2X3Y9Z/LbEry4bNhXXxWPGDuDj7HYx6sQiZ7bsA//5byXcggvnzw0o6kqIi0XBlZexYYzepNWuAzZvLdKr16411cH4+8OGHZRctjpwcoFs3qblYXCxt9mzZFtn7A9KZ9eoFpKXh2D8+gFXFe2wVFEh2DU0ZMJoG1LR2oJh6brklHJujlPx9yy1Vc63+/Y3NKoceWjXXS5bdu8mXzp7BLISDlJ7C/XQix9D2P2OGRMl36hS+dzabmII+/TTixAsWyKqfx0O2bk2++26JyWTHjh3s2bMnAXDEiBHMywtw4LnFdCKXGdjNL3EOfdhHB3LFSoJsHoKN3Hz20Kq7Ea+9Fr9gGPrSX35ZplMVF5PZjQ41/sHLYQdZvz5cgiC2tWtXplMZ8957xkFqPh/5xRfh/TZtkjQRwfTOK3E4PciOOsTlIi+8sBJkOkCBtvGnliVL4oP/Qg9uOc2wCVmxQt6ZkE3UYhE9M2NG5V8rkpwcWRRMyP33R92ErWhAH/ZF3Rebxc82bWRdbfRoYx3p9Sa/Fpqfn88rr7ySAHj44ZfQ6RAlXxu7eCT+iju3DYW8yj62orfDnD17xAYfeVGrVUJby5DEJieH7NiRfCztCeYg5gFTSvLjlIOjjopPp+9ykc8+W67TRfPAA8a9is1GPv989H4xi1Xz0JnHq9+oVIAej2RSqNFRvylGK/4U89RTxrm3bDby6aer5ppr18qMomNHyZn1119Vcx1SankMGCBrlWlp5JFHyiDckOeei3uhF+F4tsdftCOfdhSwV9sNDDnknHyysZ5IT5f12GQJBAJ86qmnCYBAVwJbmYY8KhQbnr+uJUl3k/Ly119S0St003r0kB+tDPzvfzI6dyCPP6InM+FhAWzMUl4G6teXRHPlYPlyqcHi88kz6vGQp55aSUr288+N3bu83ugcOH36GP/wGRks/vBjBp54guzWTfIbzZ1bCYIdeGjFn2JeecV4+ux0ikPJ/kwgIPor0isn9B5v2GBwwKZNxtMfgDtQl3sdDcj//ivZ/cwzjd9/r7f0dCm7d0s/c/755COPiBkD+IKAi0ALAn8RCBiev2mtanJ/2rkziWmSMW3bRsoc4En4kQ/hcX5ouZJZp59DPvgguW5duc6dm0uOGyf376efKtHhqbBQSj1GPjAOh3glRVYGu+8+48pgTqe4wIZeKKVkSvh//2d8vc2bJT/RmjWV9AX2H7TiTzGbN5ubevZ3V7T5841Ntg6HuJQb8s030baoyBc45qApU4zPb5RZOJJ168j69cP33eGQEawoy0UEGhPwEZga/7tYC/jE4ymo81sW/vqLsz1ncA8y+DcO47V4iw2xmRvQhNkI2sbsdrl5v/ySammj2b1b/DPr1JEf6Y47ov2MN2wgzzvPWOkfeqixv7LPF502ubiYvOYaOSYjQx6Evn33z5rJ5UQr/hrAhAmi13y+cGrzsmTHTQa/X5I4Vqfd85NPzJM3DhyY4MCCAnLmTEnZe9FFkgXxp5/idgsEZPDndMq98/kk0WRppqsLL4yOEg61tm1DawYbCBxHwELglYiRf4Bnn12GALJU8PffpNdLP8KG+Cy4uRgdWAADZ/fDD6++GIWKsnu3/MCxtlGlJOLr6KONH7b0dHLhwvB5nnsufnHI4SCvvDJ5WULTnhdekLD3/eUeBtGKv4awb5+YOCdMqPxAqg8+ELus0ymDm9tuq9piF99/T/buLe+b0XvodpMjRyZ/vu3b5R2bODEY07N5Mzl2rESj5udzwwbJWvztt8mF6ft8xnJZrRJ4KXolm8C5FLv/jQQK2aBBGW7Cd9+JLapLF1E0O3aQkyaJoLFplAOBylMcV15puGgUMPrCIYW3bVvlXLuqee45U1MgXS7yhBPMP4s05xjVEAg9AO++K4tQiX6PZcskL4TXK2Ypj0ce+Bo9IohGK/4DnKlTjett3Hxz1VzvmWeMPW0i361GjSS9SzK89lr0iH6E/UkWpznlpfP5JPGQ6WqxMQ0aGMuWliYd4rXXhszEfgL3EgAtlt589NE9yV3gySejb4LdLlMMn0/kdjplpLh9uyxApqXJjTnzzMR296IimUYNHEhefTU5b178Pu3aGX65QKwrTqRsFUzhUG2Y5RuJfLBjHz6bTXx+IzEbkQBhN7cuXcxHYEceGe/aVNbRTIrRiv8Ax6zCnstVOWH2mzZJAaq8PJmJm/l5h1r79nJMMvzxR/QArzt+DtuoI1vdulFTmLw8kWnjRuPzPvhg/MDRbg/7fefkkGeesJlOlccM7KUNb1LBRrv9CLrda9izp7HOJSkJZkq7CSFF0axZ9EJmqFc0+mGKisR9JrSoEVJQscpmwIB4pRQ6d8yXDths3HBUP958s3SwcZ3xhAmSxMnrlQepLK5SVcEDDxgv6oZaenrYdh+ymR53HBmbl+mcc4xtfbEzoauvjpdh7VrzWUdlJbzauVNMV7Vry7N9++1kVlblnDuIVvwHOA0bmuudcjp1kBRnkz59wqNxr1dmEYkGU4AEkCXLnXdGWy3ex+Ao23XUCx+sHvXGG+HJgNMpBaZ2744+b36+DK5dLtnP45F8MyUONGPGkG43/0MLzsCp3IxG/BZ22pBOoC6Bn+h2R5uNS/jmG1kwLE3xWyzx7k6ACD9mTPx5x483r8D155/k0KGSxa5Nm6jz+qG4ydmaORcMltlF8Ev7PV4uTzuaLT3bCcip69aVJQKS4gkTO3p2uVKr/NevN180AuTHnDZNerCZM8mlS43P8/ffyf1GTme8yefff80Vf9u2Ff+OBQXkYYeV7tlUQbTiP8A56yzjAWBGRsXs/H36xA++nE5jXRapN8oSm3DdddHHT8D5xidOTycnTeL06fG6ym4ne/UyPv+yZbJ28OuvEe93IGDaWy5AOpvDQSssrI2X2Lt3+JCStYWFCxMrp8hmZn4ZNkzOlZ9PjholZof69Y339XjivaCC7T+0ZEv8S5e1gE5ngFdfTRb8tYocN46P9vmFNmsgTpyTTgp+ITN7WOfOpf9w06aRZ5whgSJPPJG8XS8ZFiwQty2zezpgQHLX27jROLNbbOdsVHj+UINoaJdLgicqyqefmscy/PBDxc8fRCv+A5wlS+IHim63TO3Ly6ZNyVkzYt+hBg2ST7VLymJtpOwX4ZOolA5RL93evaZxPU6nSdyAEVlZpgohAHA3wFOgCICNnDdx2DA/vV5RmkceSf44KyAjttJMCU6n8U30eiVPdVGR1MFNtGAS6tlMets8ODgZZ0bdpiuukK+ZaIE7f9tec6Xo9Sa+f888E/2jOZ2SKqMy1xH8frH3G90/u10ymSbD+eeb/05KkaedZnzcokVhU1LonnTqVDm202HDzH/nSlxD0Ir/IOC338RNuW5d8Xj7/POKne/3380Vh5muuPLKMijfIH6/vJshPZKmiviDpTcL7J6wlnK5JMc8RfEaXT89Xd7VpC9aihmgEOA1kNz+FstAIiKfkNtNLpm6Uey9odF4Wlp4gTd0Q848U/aJnDbZbOJxkpcnLkzJzByMwr4jWi6cbIINJZscDjF91a5tfEhaGlmYV2yev7tZM/N7t2ePsTJ2ucQjpzIJBMynsy4XuXJl6eeYPdv4eLtdblCi6Obdu8nXX5d1h6+/rpyiBKR4FZnlKypjrqZEaMWvSUhOjqwFRJqF8vKSt2aEzCHlxe8Xz6SrrpI1hO+/Lebb/b/k5+4r+GGd2/jJg3+WmD7vuMN47c/jKWMdk1ivHIMWAPgYvAQUgY4ENhEQ3X7hhRTF9Ndf5M8/y8V/+03yZAweLG6ofr8oj6uukms5nWKDDy1E3nij+fVDi5d160qPnmB2sQcZ7IJfou7F8uXkrbfG63abTdY96fdHzSKmoQ+742cego08xz7VPIfU9OnmneZJJ5X3ETCnZ0/ja2VkiDttaZh1HE2bygJrKsjKkuC12HJ0zZtXqg+2VvwaQ4qKxN/f5RK9lJ4uxYZCvP668cAktvXoUXky7d0r72SkZcPtlvguUnRmvXrxn48aVcYLBQKi/DMySJvN1Ad+uaU9gckEPASaEPiNgHgumZ533DhRWCecIFN3sx7psceMezGvl7z3XgmWKCyUqUyCTioXzrgMp0OGiF7r1CnsXerzSSzXtm2UDim481hcQndE5kuFYno8MuuLw8iuCIhyvfjiMv4ISfDQQ8Yzk2Qzj5ot0losqfHJX76cPPZY+d1DHZLNJgtqZi5q5UQrfo0hRuVE3W7RWyG++04WTtu2NTYzu92SQTOKQEAOvPBC8uyzxS89yWnyyJHG76rTGc5htnmzdFhHHCGyTZtWgZtQXCy9zaBBcQom4PbwWtv/Bf+7hEAzAm5aLJN4+eUm57vuumjFGAo6MhrJrV9vrNDr1o0Pv/7iCxklxuybDQ+vwJi4U7hc4vUYCEhA9GuvyX0qcRoJKn4/FBtgq6Fu7NfP4PsFAmJvizU/ud3xqSECAVnEGTRIFmTHjy+7uWTbNunpI6/ndpPJ1lSuW9dY8TscZHExly0jv/pKSjtWOdnZ8hvGzkDq1q1Y2U0TyqX4AdyVqCU6tjKbVvxVQ2Gh+SDSLO/6F1+IQgl1AB6POKPEpYi4555o5efxiI+n30+/X0aSS5aIDpg5U97hO+8US8lZZxnLlJ4ena690snMFE0Xyu3idJL33cdbbwlE3KctBDoRULz77ucYiPUGWbXKeHTq9cYUEAgT+OZbFvrqMDfNx/w0D/Mbt4jK1b1ypTjPuFxk/foBPjr4Xxaec4Gsordvz41Pvs/atY0TzTkckorbMKaiuJh0OrkN9elAnuHxdrvJQv2GDTJqDU0TvV6xW8dyxx3xz8GAAWWPYF63TsxnDRvKCGT06OTdHh94IH4k4XAw6/IbeOqp8lF6uvzcAwdW8STg/ffNvXk++qjSL1dexf9oopbo2MpsWvFXDbt2mcfJpKebH/fPP5JSf8gQ8rPPDNInrFlj6smybOS3bNxYnnOPR146l0sGQKFYpe7djWcWHo+JT31ls26d5GUJBgb4/eLE0rCh3K8uXXJ5+ukXEgCvvvpqFoQ0xd9/S2CW0Q0N2V4i2LRJ2uDBZLq7iCdiATuoJXQ5AyVp6TdulD4ocoDocsVbVA45xPyyPp90AOefbzCofOop5ts8dBkUwgk1m03c2g1ZuVJG+Uaj1VWrjKduleyyWCoFBbKoERnwdfLJvPLSwrg+2uWSLK5VRqLC008+WemX06YeTRyBgLmeOvVU42OKi5PwZnv7bdOpxFtpN5kqmMhRamy/YbOJp1JNyZHl9/v58MMPEwBPOeUU7tq+nWza1DxXjt1ODh9OUuKwjjwyHA9htGbrcIgpuFs349M5ndGBeYk8FiOPufbamC8SCJBPP81b7W/QFVPdKrKVq57L66+b29fvvLPc9z5pAgExrYRmBqtWSZTy+++z+PU36UjzG4pWplxNZeWbb5KrRVBJlHfE/2qilujYiHNYAfwOYErw/60ALACwGsB4APbSzqEVf9nYvp189FFxT77ppsTeauPGRevoUGbk2JF1fr6cy+USU+vhh4uJxpDx4w39QIutaXzK9nCpit/tlspKzZrJ9ULBWbE5z6qCQEDWUb/4IqokgCljx46l3W5nmyZNuMyZYBXc6ST/+4/79sUX3jK7Bw0amMd+ZWSIDgmxbBlL4gxKU/5GpozCfD8vOK+QZrUJrNZy3MyxY42VnN1OPv54OU5YBj75RKqZ2WzyLD72mMzgjj+e9HqZ78ygxaQAj8dTybIEAuIKeu65YrNr1Sp6ZON0SsBcJUbshiiv4i8E8BuA+wEMBnBlZEt0bMQ57gIwLkLxfwbg4uDfbwK4sbRzaMWfPOvWhTN0hkbKbnd0BP6uXZI7rHv3sEeZyyUuzWecIZkrY7nwwvjBm9ttUjYyJ8cwp0OhzcVDsZoWFPN4LOLxWESFYvbCdL6JoXwVt/AE/Eqvl/zwQ3lf/v2XfOcdcYcfMkRKSpbGv//KYLOsKeh37JCRbcgt3+kUM0xpa5Fz5sxhfV86awGcYaBJAqGbtWMH3347OS8phyNxuhqXS0xukSxdKinsGzUyH/2npQXiUluE+Ptv8+vZbNH7rlolRa+yFycw9WRlGQeCuFzJ9arlxSxjYciTJrjtOCyOE81iEV+ESuXGG6N/dLdbUm40by7pNx5+uHKCwgwor+KvC+AGALMATAdwLYBaiY6JOb4pgBkATgMwBYACsBOALfh5VwDflXYerfiT57LLjGN9Dj1UFOmqVaLsjWzoLlf0KDLE5s3GJnuLhbz0UhNBfvlFLpSeXmJXXfrwJzzd+TO3oCH3wcd98DAHTubAST/AIliYDTcfsvyPf/4pi8+xSSiVIl9+2fiSfn98srp69cyTxWVni7ILvXNnnBF/X5xO80BKv19cSNu0Idt4fuYRULQBfMdsqP3ccxw+vHSlHxphm31mtxSy3xH/ydTOBDPX9Vb4l4GBF0QXLInALAdTKPfS5s3ioOR2+pluyaQb2RztuMN8cffnn2VEEXoOPJ6KRxaWhlna5pj2K06gF5m0I7/kJ6pdO75DrRDLlxubuzyeqi+AzXIq/qgdRYnfA2AzgCuSPGYCgBMAnBJU/PUArI74vBmApSbHDgWwCMCi5s2bV/kNOlAw81yz22VEe+qpic0BxxwTf85ffjGP1zn22ATCFBbKw/3NN2R2NgM7dzHXWnpEWC6cPAx/87y+xjZnpYyTGA4ZYnzK2CDU4mIxMbtcoq/cbomgT5TSxSg9y3XXRQ8sm2MJPehBALwFdhbHnsTr5ZyR802D4kLJ5Nxu6YSilX+ACn46kcvr8BZznHUSKo/ly4P1ci1izrCgmG5kcxr6yIVM8nUvXBh/Hxo0CHeOJ5xA2mzR5iA3sjkLJxu7cxYXiwdOmzZi4vjf/6rEbTEKA5dXs7YWzXkPnmefPgE+8kgVmBNffdU878k991TyxeKpkOIHcDyA5wEsAfAegPZJHHMWgNHBv8us+CObHvEnT6tW5oo/M7P0BUCXK/6cO3caP7s2mySLTJacl9/mJPtA3oNn2RvTOBHnGC6GZsPNa/E2zezNgHGd4kSmkchR/4gRxuncS9MT77wTPsfGjWbZDopoxY0EwP6wMitmh4DXy9OP2RZ1P10uCbJ6/30xTWdmygJwpIyDMYa7UCs+a2nduqZVaf79l7ze/QGPxe+8EJ/yN3SIvqiJDauggHzqKTFzjR9PLh+7mD+3voKLfKfwIcuTrIXdMd/Zz3PwpXEA1wUXRH+R0JdNppJOeTn5ZOMf0GaLfwEsFpNghUriww/N1zmeeqrqrhukvKaexwEsBjA2qMhtifaPOfZpABsBrAWwFUAugI+1qadqeekl48yVAweKqccsNUuomaUav+WW+PP6fMnXr/7mG9KTlk8bCkoU+hN4yLBMYCa8vAwfUQqkGMtplBIm0UwmtG4RCJSeUtqs1asXvtaUKebn8SCTh+E2AlYqHMPT8Q7/Q4sSxZf/5PN85BHppFu3lsSWRpaXzz+XhWCfj/zRcorxxXy+8Ch71iypeh/K9z90qPmXsVqTGnnPvfMz5sDNYojCzIGT69GUdbEj6nSdMF/+iEzZ8Pvvxt5dXm/589EUFEh1s8sukxGzUa6euXONbfwjRsj0JWRv93gkG6qpr2olkJlpvM7hdkvwXhVTXsUfALAGwF/B9mew/QXgz0THxpznlIjF3c9jFndvKu14rfiTp7hY6ko4HGKecbslnUIog+2ll5orSLdbcoZFEQiQ779Pf9sjONL9EJs6d9Dj8rNPHxmVcssWefkMRo+rV4tr8t13G88Y2mIFcxBv/8yGm+nYy7SITiK2GVURNKu0p5R4CV13nUT7lkfph3RliL/+Mps9BXgk/gwGRE0jkE6gEWvhW+5FsKcoQ1m0wkIJasvpepq54p8/XyrGxBZgidm3CBG2oySKiRTmFHK3is/ylg87n8Wwkk1O5PJxDI9P0jZqlLmZ47bboi+2Y4cUX69bl2zcWJR0bFRgbq7YmkKK22aTaxpF9f30kyz4hBZSQ8FRmZlSzOH668UDoDqqkv38s3yv0DqHzyflOauB8ir+FolaomNjzhOp+FsDWAhx5/wcgKO047XiLzsbNohzw7Jl0dv79TN+D51OGUjF8cgjpMfDAqRxN2oxACUP7oIFsmDgcIRHTpMmMTdXOpm33gr7qScaiQ/FG8yFg/vg5V74mAUP++JbAqQLOazlK4o75oYbjL+zVOwzNw9FdgTlUfx16oSvtWaN8Xna4y+6kRWxbSmBlgScHIJ+CSN4E/Lxx8buQA0aSKfbu3dSXyIf9qQLrfw98S9mwjg96wq0LVH6zbGWux2N4tMyf/65sZnD6WRJhBopCwiNG8f/SLFpl1991XgGkZ5uEDoeJBCQmdDLL4tLZXlMTAUF4gpakSCSwkKRY/p004X1qqDCi7slB4idXpX1uIo0rfgrh0DAXOlZrTKKHTJEgobuv5/csjqb+c4M3oDRdCKXduSzOdbya8sAGcHEuMDkWd083vYHrdayKde62M5b1aschE/pQSaBAC0o5rvHv859+2SAXL++pL+PnZHs3Cku2t27ky5n6Uq/Ii0yCdwzzxivC3THT1HJzqRtJ9CdAPhEgwYMmCmpRPj95EUXieKz2aQT8PnErEEmDt2NaMVQDBx1lLgYlpKcZvPCDcyF8Yh9vrUrO7bZyxGHfcTdHU6VRdvIwihbt9K0cILHE15FLSqSdQGzHjqyzrJZNFt6evg+RJKVFc5Q53DI/WrZMvmaoAUFYuMM5Shp0qT8OUNyciQh0IQJlVuwphTKO+LvAuBHABMBHAdgadBevx1Av0THVmbTir9yKCpKrBfc7rD5wm4n66QXcYB1Slw4/4mYz2IV72/oh2IBbHwA/6MV8SP1RC0tLcDI0bpNFbNd2wAfeyy6CHuDBpLjhxRzT+PGYWtCLewqszJv1kwighMVnLFaw8WySOkAzAJSj8If9ESN+KW5sIsn+E4kAF5++eXMN1D+q1fLjCxkwbj22ug64IsXk12OzqJV+ZnuKuB9d+SF87716pXUFy4x/9hsohBLyYHxh6973DpMNjz85a4ELpl+v0T4GfmkNmwoZpjQfr17J/Y4uPXW8HnNZjVmaUTvvDN+UctqlRTXyXDNNcbBK7NnJ3d8iO++k3sdMvW4XFWSl8eI8ir+RQD6ABgEYA+ALsHt7QD8nujYymxa8ZdOXp4E8CRw7SZZuldPdAsYRjf2w1RmK+MopLVoFvSLrvjo2+k09tRp3FisG3fcEZ50dMBvfBL3G4y2zZvNFg6WHD9eYga8XompOf988UJcuDA6tiYz01zpAzLrOM65jPaIpGcWFLOO2sVda3bzf//7HwGwe/fu3B7xY+3eTdat7Y+63w5rIbt08jMQkE4h1mricpGXXBI8wc8/l17Fy6h17Jjwedn+11Yudx3HbHi4F+nMhZMzuz3IgD+B2eOHH4xlsdmi1wC+/bb0gg933BHe/8sv481dSskquZEZxqyMpc0Wb27Zt0/sk3ffLeHs27aZjwZOPz3hPYtizx7je+FyJe8ZUQHKq/iXRPy9Iuaz3xMdW5lNK/7EvPJKuPB4KBlXdnb8fosXl62iVkj5x26rj23MR7x7UCGsfAD/Mz1XKUWkkm4WS4AjRoRLoh6FP5kFDwth5fmYQBdy6EUmbSigMvEMstkkyrU0pk4VU3Pr1rJoPm6cuTePzycDz73r9nJw65/pQB6tKGKv2ov4zw9rS8752Wef0el0slWrVlwaLBQ+8ql8ulV8ojSPNZfz54tlxsi05HQGU7j7/bJ2kMin1agplVSqgFUT/uSip77jzpU7Sr9pr79ufr1Bg8L73Xln6bJFhoYHArJCH7IlWiziwbB8ubEcZuXHbLbo3nzlSjFdhjoVr1eias06pZYtS78HId5/P6qz2oTGHI0b+JrlVq6/++Xkz1NOyqv4fzP62+j/Vdm04jfnq6/iBxQOR/T7RcrsOtFItSyK32rx85sjo9MuF8LKjWhs4OMd/b4Znc9MOSdqDksBmzeSiMsvcB6LI/zbl+EIfojL+Q36sgN+o8cSrVA9HlHkpQXrvP569L21WkmfL2B6H88/L8DAV5PFXtO1KwPXXMPAmA/EYyWGhQsXslGjRkxPT+e0adN4eaeVhuf0IpPvj1jLzp2MZ1AZGeTMR2aJ+6bDIVOgZIISQs3trvzMd889Z369UAgwKS5fZv7FSkmejki2bRNbX+T383jEA8iI664zDlFv0CA690eXLvGLUmb30WKJf7kS8dprJTOHjWjMn9CDN2A0fdhLp7XQMB6lMimv4vcDyASQBaA4+Hfo/0WJjq3MphW/OV26mChGB6Nyshx/fNkUa+TzHzm7VipAr5f8e1VAPE1OOIGFjZryJdzG+thWhk5EbPpO5PBePG2aDz6h8lf5dLnI9TDx4wTo96bz+zfX8LnnxEHpmWcks8A334SLuhiRn288O7KhkBdYJrCn+jlqu9VKPm+5N75IvMUiP8brr8ddY/369Tz22GNpsVh4XrMr48xUR+Iv/oLODAB82vIgbSiMk8dpL+ZG56ExQtpEILudAYABq40BiyV+yuVyxbtVVgbvvmtuU7zqqvB+GzYYm0EsFvLFF+NdhB980LijcDrFpBLLjh1iBoq9Rij/99y5sgBs1lE6HPEdgs0W9GNOkr//Ji2WqMC7LHj4PU6ngt8w51JlUmlePaloWvGb07y58TPr9cozF6K8ZpaLLxbT52GHyeiyf3/x/olk9GiyPDb9NOTzCTzEffDyYowz3c9sRpCOPbwPT/EvHGF8oM0Wlbs4P1/MYJF1Vs4/38ATMBDg0v99Sa+KX6QFyEPxD3OUmyc6/2RGhuiRFlZzDxgCDLhchgojKyuLZ599NgHQiWupgnljDsFG5sBZshjrh+J/aMHDEZ4ZuFzkRQ1nGl/TbpeaAuvWyQ+2dasEcUR++fPOM3eDrAirVxvbxz2e+MCtqVNFnpB/e8OG0Z48kXTqZPxdMzLCC8ax5OdLDmyj45o2lVmE0awg9PzEdmBOp7wQpbFqlYw0Bg82fPky4WUfTGNaGvn002W4t2VEK/4ayK5dNM2UmAyvvWY+WKlVK7rSn9k6l1mz28VvPTKJYmYmOWaMzOTnzw9bCKS2RLzib4gtvBGvszt+Nu0Y7sNTfBBPlHjkpKXJLCZykHY0ltAoiteJbK5AW+bBXqIgV6AtL8I4NlfreFLjVVEpzu+5J97c5XKFU6bk50sSsuKnnuE2V/OS5F2x7WTMYhEsHG+9mEOGBAujYBz3mfi8E2CxsjJw112Gv2NxcTHvvuEGAmBdnEALdvIDXB4XgBUAuAUNaVXFTE+XkrwFLQ83vqbXa5zKdMMG8SWPme6EvnuZ3dyzsqRTiTUX3XJLfEbKnj2N00QUFIinzLx5iVOhDhpk7Cdc2rA5kY2zfn2JcIxVzna7+WgplKAqO9v4u7/zjlyzFJPbc7ibVmvVZqjWir8GsXy5mF7sdlF0XbuWHjW+bh35wQcSg1JQIKYKMycOt5t8773o4595JnmnD6XIE0+U9A+hALCFC2Vg5fWGzT/nnCMDuGHDDN4N/MYdqMNsuPg1zqQXmfG6CfvYHP/RgdySbQ5HOGV56L27FS8bVIiSpGWXYGyJeWUpjqAXmbREuJK63eEUz2bpkH0+ceZwu8kMZz6z4GUAoAvZcbMNN7L5LfoGr9eedrscdzq+T6j4CXDTmddF/yiFheSbb0oQwskn851zz6UNYHuA/yb6gSIrNV10kbFZxeNJKiVDcTF53z3FdLtFV9WqZZwHKY69eyUPiN0uP1qLFtGFRAIBiU494wzJnfP22xWvafjNN8b3IzYLXywNG5rfS4uFPOWU+BX7Qw4xf2Fat5Y85Q6HtGbNxEOJNE9sFdNy4eAwPEuXS7zxqgqt+GsImZnibBA5cLFYZG3O6L0IBESxOp1hV+C6dcULz+w5NgoM9ftFuSWj+EPKP1Qa8fzzJXbFaL+QXCFlHNruQg6Pxe/Mh51+KHbEwig7vgN5bIp1dBt0CG63FHm57bbgvcFm/h+uZAZ2M3bm4EAu++NrEuA5mEil4mcGdevKoC7Rdw2944fiH2bCy1VoQxeymYYCOpBHH/YxHXv5Jq4jISmkx+LSkntuRRE3oXF8ErVgy4SXr5zxbfQPcuqp0crF4+GMU05hLY+H9S0WzjUTOLKy/fLl8RVY3O7S7QeBADl6NDM9DRkAuA7NeAnGlhxuGMUdycknx3sQud0V02I5OTIb+ekn45H/ffcZj8IdDmMbP0UPzzr9CeZZE4x6lIpX1qEczUbXatkyfq3B7ZYAk3HjSndRhaQlae3cVLVlHkmt+GsKZoU4fD7jNOVTpxrvbzYL9fnM1578/lKfR4MW4GW2cYZBSaU1D7L4Pq4kIQta9+A5pmMvvdjH6/AmL8WHpp3JE0/IgNLjkf83cuxmF/xCZRBX4EQu16MpG2CL4flstsRejpEDZg+ymAsn/0XLklnGs7iHi3AcCxC2BWfBzfZYSkDO7XSS7bCca9CK+UHTU8hUkwkPJ+B83n5rhNtkAh/2VW43D1OKdoBjjYTduTP6h/3003Aq4rQ0Cb+OcNH84YdgDn23pOj54guSr73GQMyINhtuno8JBKSeuSlmtXStVgl6Kg/jx4cr4IQi9WKDy0480fgHTE83tPGvXSuJ9dyOYr6HIcyDg3uQbto5x7XmzUWm0MPj8YjST1Sc4vzzzcPW09JY7PYyz57Or/q/yT//qPo6olrxVyF790oUeDJecfffb/pMGGacPPts8/2NlL/Pl3i97qijknvmQ+0Z3Ms56Eof9pWj0yAHYXzJf4ph4UjcUeoxoZiESIXscpn7z2dgD+egG4/CH6bnPA6LOAo30oOsEu8YhWJ60gpojZkljMJNzIabR+IvKkhQ1YN4gttRl0Wwcj46sTPmlezfsKEkwhT9EGAH/MZTMIN98A1H4Uaeia+p4OdLLwV/hO3bpY5lgpuwE+ApAAHwYYD+kHKJNQj/+Wf8yMDjkdqblAF0XPCpK8A8bz3D64Zy8Hi9CR7i774zL9DQo0fpL0Es//xj3JH4fDICeP11uWcXXmhu4w+ln1i+XFJBN2nCC+tOpyXit22ILTwZs9gb06JfJDNF7XTKAtCll4qXw7vvSodt9t1r1TI387hc8ps0aiQ9sMcji8vz5pX9fpUBrfirgJ07xcslNOJr0aL0ojoTJ5rXWjaqYXvqqcbPkc8nijBy/cjtlsSDJCWl47vvynAvOPpbuVJczCOf80Q5dXzYx1w46YdiQ2wus9K3oZC34eWSDdlwsxPml3pcyHQau10pkyAm5HI3arG/msr4ReQArShkNtz0A/wPLXg7XmIPzOZNGMVVOCzOjm9FEZ/D3VyCo1kfW+nFvoQpKI46KlFEtJ+D8CnfxHVcX/vocK6dJBIaFQC8Oqj8LwSYa7dLVslIBg40V4Y5OTzuuPiP7MgvSbMcf800Xo13edWxv8U/jDNnygN/1FHmJpeHHirze8QHHzT3rLFa5bu43VIKzSjn+Mkny3lChYeDP0YG9pg8lwXMsgcL2bRvb+75EHKB9XolT1BurrzARj92osVgq1W+o1GH4fOJl0cVoRV/FXDiifEmBI8ncXHzoiLxLotUbE6nLPAazRjeecd4jcnhECeJ9u1l9tmrV3BtLT9fcpG43dJ8PrJNG/72zRZ6PNHPplISff7hh8bvXQf8VpJK+Cf0oBeZQfOH30DBMm6bCzn8HcewCFZmw83ncVepSj/0viX6PNLc40YWh+EZEmBtk1w9FhRzKxqUbAgAJQFfO1CXVhTEfQ8XcjgfHZkPO8fiYvbEjzTyLHI4Agn1t0IxX8KthgVnkmkBgM8DVAA7AdzicETncTepvFPs9pHLl5s4tAS4CY2Nj4NiNtwsdrrlOcrPF7t9z57RP0zsj2S1ymJKeUpYJaobEPtyffSR2G88HnkJ+vcPu8bFmFkOwca4U1yNd7kFDRlQShTx6NFyztI8H+x2ecnM4g48HnPvgb59xTPDaMTnchnGeFQWWvFXMkuWmKcjueWWxMfu2ycLrY0by2zvoYfMay0XFMjsOTIFuVLhDic0ILrlFglUnN3rMfqdrjih5mf0NnwmTzjBvKZuXeyI8k3fiTochZs5HI/zgcM/i/JYcyGbTuTQgyymYy/TsZc3WN7gfXiaD2MEj8Xvcedv0cJYyaelBdha/ct+mGLiUhlgfWxjOyznO7i6RKkaBTiF9t8Qoej8AIfjsYicQtGBZSdjJn9FuG7rVfg/umKCq2wo4HScxmdxD0uLYSgpUmLWnE5pCXqQSQDdAJspxSUvvhh+QPr2Ndw/F06e1GFvSVqL2HYV3mM2oh/guM4plAzIzB3SYpEF0Pr1ycsvTxwRF4HfL5kYliwJTkYnT05qQZTp6ZIhr3dvuV8NGkgVq337JFdRjC3wMQyP8ga7HB/EfWe63bLwNmWKeEzUqVP2vN3t24sHk5HidzikmPqzz5q7dz78cFL3rTxoxV/JfPWVuc25T5/KvVZRkSz8XnmlzDjNZsUAucEkijUfdkOXSqWk0zEzTX6oBscXS3G7yYULuWKFmECvvJIcc96XXGNvxyU4ir+gC3ehNp/CfQmVYqgTi93uQRbfxjUciM9Mjg3wBdwVtdEPRWVijrEjL6oIST7SuAN1eRj+jjonQLbG6qjo2x2oaxhVrODnIIwP5iZKrPiPwNLEiuO66yRo6aSTzEeNAH8D2EQpelwuTp48WR6OOXPiRiDZcPMNXM+0tPhi9aHvWhs7eQnGcg1asQBpsoZg1EqbfhkVaU7A/PniIeb1SmvShJw/t1hG0wm+uzwYnvhOKLTYZSBnIWwciM/pRC6BAP+DSbRj48ZhATdvLr1MXWzr0EGm60YmAJ9PFgDnzjX+fl5vtBtsJaMVfyWzfr2xsnS5yl9KMztb/O2PO05G+Z98Em3+mTSp9DQs22BsryyA3dAUkp4u577xRuOBXRoK+BJuZzbcLISNe+u15tTbppXEBm3fLgMwu5102YtZDzv4NO7lELxHq+kIPPq9jd3mRC63oT7PwFST4wK8EJ8wT7lIl4tZysdtqE+bScDVaFwfN5othoU/4LS489ZVu6IS0C3C8ab2/eZYy2Vox8SKP8BOmMetqMdCgzKT9HrFdYsUv/4xY5jV6xw+eeh77KB+Z0/8yPEYxABAvzWNG49oz44dO1IpxZEjRzIQCJBff82NaS1ZCBuz4OFzuCe8gK2MOlc/j8SfdAfjFGpjlxRoKYuyC7VkUxxTPC6N0mD4fOSe7YVic+zbV+yeiXL4lEG+HLj4P9xPgNEVyGLPGfmiHX102e7BGWeEv+BFF4Xt/SecIGttpJz/zDOjO2m3W2IIkkiSV1604i8ngYCMUiZNiq/fcO210b+jzSYeHuVZq8nPl+ctUvl6PKKQSXLRouQCsN7A9VFuh6G2pcHRcce7XOJlRIpJ6brrpDMzuo4VReIRYw2UBPtcdZUEMZYlJ1hkS1OFrOfLo0WJ50wa8ulELt/DEBLgh7g8oRtpR+df3PjwG/zg3In02gt4CNYb7meWSqEYlpi005KAbW2HASWK5y+0N1XsjbGRBHgqfjDdB5DOU6GYc9E5WpZQAMGnn5YontxcsRxEDircyOLNeJUzLb34xK3bmJOTw4EDBxIAhw4dysLCQjZtEqAHWYZptI2aFUWcjxN5gXUi27cn19U6Oj5SOC2NebCbLgTT4ZCRep8+YisvJWDsjZfzeWbad+yHb6JMMB4P+dboYnLsWFGOAwdKB1C+rIIRv6/i/6khJebCVWhjvG+LFtGCLl2afGFmt1scKCIpKDC23RYViffFCSfI6O7VVyse1FYKWvGXg40bZaocci92OCSTbGhw4PdL2oS2bWW96eijxST69ddl78Q/+MB4Juh0SlTvRRclN9iph+38Dy3C5opgFZOihb/x6qvDqVocDjHRhNI6BAJi6mza1DhFiVFzOss+K45uYWVpRz67YC7XRkzHs+FiRywoSV6WhoIoxeb1Bjjj1aXctejfoMNErL2etCPXdDRbgDTDPEDHttonIzGAi3Ac00xmEodhVVDBWPgEHmQt7En4G6Uhn4/hYf6H5lyPptyFWuEPg/bBd981s3gEeBT+5Ex1GgMOB/316vHBnj0JgKeddhqHDNkdNXuqhd3sirlsZtIZ2lAQnoH4/TIy9flKehy/28MttqbMSzQTiHSFdLvF28dksWrxcz9wLzK4F+nci3RmwcMB+JKApNl+4rD3o7+4xyPFXMpWQCKqBQAekbGJ7mAnMwifGpstjSIezaYnSsloPj1djo0sy1YD0Yq/HHTqFO+h5fFIcF4k774rg5PQvh6PDFwSpR2J5dJLjZ9fr1ecDswidY2aCzl8rNk7EsjzzDNRlcnXriU/+yw6Bw9JvvBC+ep4VFaOfZE7m0vRngT4KQYxHXvpw17akccM7OKVeI9NsI4f4xJuQUMuRXs+iP8F/ffj3TgdyOXbuMbQoyYfdn6MS0o2WVDEERjOrWjA7yx9SwKKVuIwGnnzNMJGTsQ5URtzjziO7w36li6Vm8T3DbAnfizpoAMAOXs2L7jAeP8W+I/74IuyxefAxestPQikMS3tcHq9f9PpCPBp3MtcOLkXGSywOPmNpX9M5k8psHMZPmTB0RHv1rZt8rwMGcJn27zDB9VT5qN9ox7ORBHu+Xc3s2OzlgblPwQb6bEXcLbNoJh8BZR+SMZsR22+gaEcgnf5HO7mNtTjFjRigXLQf3jbxKUUp0yR7xR6yF0uGRktWCD+92YeGTUIrfjLyLp15gueXbqE99uzx3hG6vHIM7VrlyQ2e+edxKU+hw0ztnf7fLL2c889Zauxccwx0Qm3AgFxJXa7pTNxOsnbb5d9CguTn9km06wojCuDmMw77EAeX8Zt/BNHxaUotqCYbmRzJ2qzKEIZZcFjGmPgQB63IzpQKQCJpP0DR5eseTiRy9G4vmQff0iJuVz8A0fTjryomcap+CEYF6DE9g5Fv1Jku3acYj8v6WA3Kwp5HBazGBbmwsnA+efz9tuNTWev4BZDE14unKyDLwnUJVCHDza6goUq+kHJVw6OU5fQyN122CXr457FrVvJu60vxXvARP6YZh4Gp5wSd76frnw3Pl01wDw4OEw9x3741tjdNS0tsSdDac3koSuqVTf5KfnSpeJuevrpkmlz5Mjo4jA1HK34y8hff5l7mLVrF95v4kRzpdmtm3QKXm849cArrxhf759/4kfcSonXQ3GxOBvUrVu25/6QQ8Ieds8/bzyid7tlMTqRObU0E1PkqD8NBWyEzdyMhnwMwwmIIkvm/fUik+/jSt5gedNwQdWJPM7AyVEbp6C/6eKrG9n8Fy3jPhiN6+kqWTvwc0uEj39Ua9OGxS4ve2A2J2IAa2E3bcjnTtQxvRGFsLEhNiddXMaLTE7FGcxHGv/uOYQrVhgr/lEwjvTdgwyegpkE/mET+JgG8D2D/fLgMCxL6fXGx4+sXpbPTJOEc35AHmaz6eEFF8Q927MGjDQ0GfmhOA4XGS96h4SrVSs+SjGZwvJut/lIKS2tbAXPly8X11GfLxxMduGFZZvSpwit+MtIUZGxonU4ZOQc4ptvzBW/kRnE5TKvFDdlirgRezwyWFFKnrU775R1s3XrxHMsGYUCyPHdusm5G5jottLeESCcRdT8OoGgoguwKdZxDrqSELfC050/87XXklun8yKT++DjAJuxN48HWZyIc0s2/IIuCWvsNsFGw5FkIWwcjDF0IpetsNrcjl2rFgOffc69SKcf4DbU5/V4vdQsnCvQlodjBd3Ipg/76Lbk0YpiGi8A+zkCjzAHLr57/Xz+9JPxuokXmcwxWKTOhZNNsZ525PFrdGNvSKTvsJCSjtivoUEuI6PKi/4Vq5iljEc9+Ta3FHowChxzOMhZs6JPtncv/+51fbxtHVK0/Q8cY34va9USm+QNN4iJpX17yWY6cWL8QojFIvsfcojYRb/80tw7p1at5JV2ICCLeLGjH49HpvE1HK34y0GorGFIgbvd8rxH5tDPz5fnyEhZGik7my1xvMbmzaLsI58zpzMcG7BvH9mmTfL2eLtdChGVZot3uSrsRCGKBH5mYA+3oCH9UJze4mqSkblsjJpEyv7oPoOsU4cTL5vAJ22PcBieZTOsK9nPimIujlAUZ2AKI+3vIXOMDQV0I5vfoXfwQOMvvx31uBAnmEfVHnkkl3y6Is7kYaTEiOgAqACkA1iMDiy0ufg/y3DDRWIvMvkGhvI967X85BNyyNnbDWcL6djLr3FmnBwTcL7oXORxA5qwCOCNEOV/LsDs4L7brI1p1PF06GDwEO7dy+I041X7gq4nyT6NGnELGvJF3MHheJw/4iT5/r/8Ej5PICAXsNv5I3pGmXuy4OFcX18GlIkN0Os1L8hCykvkcMioy+uVafj6GLOVUV1St1uCqZLl77/NX4w2barUFbMy0Iq/nCxdKi6V/fuLmSYzM36fuXPDSQW9XnkeBwwwfl6UEnu+GU88YR4fEMp4m50t3kT9+5OXXUZm+MxGkyLL5s3idlma0jbqwMrTHMjlwxhBAlzU7lKSogNOPtnsmAAv6bGOOROnkUOGMOB2sxiK+bAzBy5eVFKdK0AHcvglzuYeZLAV/ok6z/UYzf6YwlvwClciXKAkB27m2X30w8LdqFWyYFkMC3fb6jNweFtDBbH66c/YyfUnMxEe/W7AIdyB+KlgFjzcjQzDL1gIK3PhZH1si0o3oeCnF5nsi2/ZoJ6fRdfdwEutnxjeIxdy+AXOYzZc9EMxBy6+gltLXBX7YWrJzgGArwK0ADwO4AaA/6S1ow/72NsynbfiFZ5hmUav28+5c42fw79OuDKug8uGm2vf+Z6cP5/foTfdyA4GR/npQRYH4EsWd4h4V2fMKLGXBgD+iuM5F134C7pwfv/H6F/0m/FLYrcnXhALsX27TJN//dU8Q+Inn0i+fKVkCj9yZNlqDC9bZh5YZrGw6ORe/HhMIUeNMq57w127pDjG6NFR1eCqC634q5jcXFnMHTNG3EDNArzc7sQJ+QYMMH7GfD7xxonj33+5Ob0tT8CvhiPFI46Q3WbNKn2WYLNV3iLvGZjKLHi58Y3JJaLeeKPZ/oGgH38B++FbbkajGIXjilkwDZRkzYwc8W81sdUXwspXcAvrYCftyGc69rAvpjIDe2hHAetnFPD3U+8I229r1WJg1Gts1UqU8xaEC3n8gNPYD1OZCS+z4GEBbMyGm+NwMU/EvKiFZwIsgK3Ehr0KbXgiFgT9+sU0lpYm3mMbP5pJejycjLPoMYiwlrQTYtv+CmexMTYxzSbf3WolN0fIGGrfAPQBbAzwV4gnU65yMl85mGvzMb9Vu/gUzyTz8sharny+jWuYCydz4OJ21OPl+JAXXEAWPP6MYQI0D7I41hWRlvnVV839fRs3ltF0aKQTynTodIpvc2VTUFC+ovJ+v8hq8qBnw83b7KPpdMrjc8MNEZf58kvZGIo4djrFa6oa0Yo/BbzwQtjNMxSrc/PNiY955BHjd8XtNnEmuPJK0mrlHmSwLZaXpGVwIYc+X4CLF4d3XbCg9JF/olxTyTY78ngPnuOm7oOipsJjxpgr/pLOB4U8FP9EuRHuRToH4nPjzgqFJcr/JxhXW/k/xBcyj50hud3kh2/l8t3HNrJJw6KSNRZASi1mwcNcOLkWzelELjOwh9fgHd6LZ3gCfiUgM51xuKgkF38RFOejI/ch2l6+A3W5w9OC218ZF9a7Q4eSStEPxUEYHwxc8zMNBbSiiE7kRK1v+J1ubpy5ivte+5A5Ay6i/xDjVB1/AWwB0AVwQuznaWmySBnDqlVhhxgXcngINpaY0WrVImeP+oPp2Gv4e/TGd3Jwkybi/paoGpXDISP2VavkZXn1VRk1VQcffRROFtW6tcmoKsjs2QntoPPQOer9+eormrv7Wa3hqmUXXEBu2VKlX1Mr/hTx558SHXvXXdHmTzM2b44fdTsc4cyzcbRoUbJjyDf9Zozi07bhnPnhOr7zjriCvvOOmKmWL0+s2K1W8v/+z9z90ucrza3UT68tj+vG/Rw3wnrvPSMPoXgTlQ/7SsobhhT/uZhoes00FNDpJM/yzGSuin7ZsuBmC/yXVIdVq5Z0ALWxi+2xNGjGkM8aYxPvw1N8FTezJf4t+ewQbOQ9eJZP4T5+j9OYG7FQ7IfiTtTmP2jNgkjPFatV8rJH+oHfdFPJzQkAnIFTeTte5L14mv0wlenYy68QLM4QWvRp3bpkGue32FiANH6Fs7gL0ZWjtgLsArH7P4WYZGxpaXG/0/Ll5vfI7Q7mQ1NGsxKZ6UVtKM0lzG4vPathZZKbKyHqsa5Tbjf5/vtiTzXyz//hB1MPhznoFrXp7LMpHUtpiedsNin2UhUF74NUu+IH4ASwEMAfAJYBeCy4vRWABQBWAxgPwF7aufZnxV8e/vornAXX6ZTUEJmZYi7MyZF/8/KCO3fuHPUwBQA+iQekXqwKK1WPR6KLTzxR3kULikrMDbHPY6Jo3JC51OzzQw81d3OOL8puvC7hRA5H44YI5e1J6L2TgT1saBef/P6YwlVowwAkm+gwPJswl350C3AsLmUeHNwHH7Pg4R14IWoflzPAEZYRrIttPAtfMQcu5sFBf6xCDbYcuPgYhvMLywX0W4N5+Pv3j7f3zp0bVxGLEFOCD/voseQwB8FcGaefbvojbUc9ZmA3FyIc8eeH4kx0ogdnEQDPQD3mR3ZCfj8LC8UJoLhYHGnMnAEOOSTo8VYn3qzoQSY/x8D4g0pT/i6XpOokxc1y+HCJ2j3uOBmxxC6gBgIxL0GSbNok3kEJi1D4pBMYMSK6QwwEDPP2Z8HDq/B/UZv79qWMnpKZOnu98RGhlUgqFL8C4A3+nRZU9l0AfAbg4uD2NwHcWNq5DjbFH8Lvl+dt5kxRqKFBSigt8xVXkNmffh1lvP8Ql5kqSXneA7wPT7Enfkw6r0tZWu/e5t/n43eyDZR9vPL3IIszcTKz4GYOXOyPKQmv6UQeT4wr8BLp7ROr+I07nGZYH7WguQxHsC52RKV7rlePPL7herqDrqexJzFS/hMsg6Q+SSBg6gWSl0e+WeeBYEdiZ3bwu5+rvqTbLW7D9PulxGCCxZpMeHkU/mRrrI6SJR82XoV3eAiuJQAeidrcoRQDffty+PCwU0KdOpJOxshbMy2NJTViR48O3cdwOxR/m2f5TORWZrNJMElurrhORnZqbreYM0N8/70UoAiZS664QrwdkmHgwOSjgd1uCckPUVwcVYM3AOmUJ+C8qPfI45GUQ9yyJami6wSko6siUmrqAeAG8BuAzgB2ArAFt3cF8F1pxx+sip+UabfZe+50yuCRL75YUq/wCLU84TPWDT/zLHxlWLu2Mtr555t/l1cummOidKNnAe2wjHfieV6HN1kXO0q2q6AvfOQithvZvBGj6EBpaRLC1zgLX0UlCQPElj05xl3yCCyLWzB3OqVeQC9MLylSE9siF3gLbS7uuKv0dK1jxojSaINVvAMvcijeZD1sZ1pajDOA2ep/sOXAxVZYQzeyuQZh7R0IdgpZ8PBpNKcNVrYGOPbYIWzu2h6n8558UuQJDTYcDnGLz8qSEX89g8qNHmTxC5wX/4HVKqNlM+Ufimw0GyU7nRLh+Mcf8S9DqBhLMpQ1sdRhh4WPDeUyivj8J3SnG9m0Kynk4/GIKCUR86+/zpKiFWazDK83icr25Sclih+AFcASANkAngVQD8DqiM+bAVhqcuxQAIsALGrevHmV3ZiaznXXJR4sOZ1Bq0FuLvnHH6xTK3HEaNhunTiHfHmax5M4tfhprdYkdV0Liug0nLX4+TJu5tmYxHTsYUv8y1dtd/IKjCn9nBaylq+I6e4i7rHW4RT0ZwcsZjr2shPmcwZOjUqJEFls3ehcp2CmqeJfHVK4Sskocft285sSCJArVnDwucZpHrxe6RRKSBDBVwwLl+AYAtKRrUMzw/0CAH+GYgOAGQAnwMs2WBW1W7NmorNCA2S7XQbau3ZJCQCj/GWAmNriNrrdctB11xkf5HLJ4laihFUffijFXoxG7E5nckVgyhqoEspZTkqnY2CzX4+mfLL5m7zzTqlvHDehW7lSpkl33y0L3pFrCzabmJ4OJBt/1EWAWgBmAeiRrOKPbAfziL+HsaNKScvIkNl/iLPOSmTCLN3MUt5mtYbNAEY8+SRpsZiVbCxNzvB2K4r4RYR3Cy+5hF5PaecM0IE8OpDHIbYPWQgbi2DlVzibz+EefoN+LIYlyjTyNw4zNJl5kMWr8H/8Pwzm7sjsmsGWBQ8fVo/LDenVi1y5knv3GnpNymp/s2ak283htqcMq42FcjX9/ruYnd/qNZ5708LD7SJYmQkv98HH/9CcLfEvFfwlxV8SlXxcC/BogFaAt6B91MdGz5DdLo4CP/5o7vbbG9+FHwijPCVffBEuCRqyp4cSpT34oLHngM8ni6snnpjcS2DG4MHG5zd7YSLtloGALMTG7hNaFE6GbdvEi8puF6V/3nnJxStUgJR79QB4BMAwbeopG/ffn3iG6nTKolyIUI6h5OtVBJiG/JKiHJYkiqeYtcMOMzZhr16dvLkzseKX5kUmN+AQPpz2FDsdlVOmGbwL2RyK0WyN1fRhH20ooBf72B5L+QeOZHbQxh8A2AvfRx3bEQu5FxnMh4PzcSL74ltmwstMeJkHB7Ph4hu4npdfFiADAW7YQPbrmc0LbRN4mfUT9my/M7zwvX171AhyHZqV1B4ILbpbLNIv9O0b/j1VMFjqR8vJJCRG4Sncx96YVmKWUijmXXiBv6CLeR6cYMsEeCbE4we4nQiaAM1M4a1bi0u8Ud1wj6OQH43aIwsWkydLuuPIhzNEZiY5YYK0yIjI//6LN+VYLOK55veTd9xh7FnjcBhfJ5ZduyTtg9crytfnk4d2+PDoDiFUQze04Bxi8WLO9vbnudav2Ekt5AjbE9x1zlVlz9kTCJQvpqAcpGJxtz6AWsG/XQB+BnAWgM9jFndvKu1cB7Pi37JFFtyMXkS3m7zttvhjVq6UWXOrVvKcJ8qzY0ceR+JOTkNvbkM9/ooO5Vb8AHnvvWE5fvxR1t4SF2spz6wjQCdy6VDxI2QFP4fgHU7CAD6He9gS/8btY0ERbTEF1u3IZ29M43o0YRGs9EOxEDZ+gXPpshVSwc9NKpwcLADwUPzDDOzmYIzhbXiZ7bGUAPncc2LnvbzBd8yEl3vh4z74mAsnb3G+LYV6Xnwxrjech078Gd1ZDAvz4ODkOoP51H3G/vIeZHMF2nIOupXktTdqn+GCuKCy2FYM8GakBZV/fzqd+0yfmVClxW++EctJqNONs2+Xlx9+kICp0GzhhBPCZpz166XHiXwZ3G7y1luTP7/fT377rUTwjh8vpjOvN2zXcjplEdggodZbb5Fud7jmg9PuZ9Om5Su8VF2kQvEfA+B3AH8CWArgkeD21kE3z9XBTsBR2rkOFsW/caMUcsnIkKRqDz4o5r81a+RZDKVTdrlEqb/2WukDh169zN95B3J5Ld7iKhxWYq/ugl8qpPgbNJDrPvSQvJNlrVtdFuVvtF2hiB3wGwmprZsFD3vgp6SOnYWT4srz5Vvd/GnAc3xp8G8scsoIfT2a8jd04HK0ZVusCBaY31OyfnLWWeTgc/YapiLOgYvHOFayjjOHN2B0nLkoyjRjt3NJWkdTeUtrPTGbi9GBfkg1KkJMQ3FlKNOcHFfvFrrdbxCwsnXro3j88Wvj1pbcblF+ITZtkk7u3ntFX1faINbvl9FLbO4dUqJ9Bw6UtZPWreUlKG++nJDZJfJLOp1iz4ohN9fYLd/hqFKnnAqTclNPRdrBoPj37pWyjZEvm9Mp0/yPP5bOIFQgqUuX5AL+lixJ5L0WYH1spQO5TMdeOpHL6/BGUnVygQAH4jMOw7PsiIVRn9WrJzP2ZEw7LpfkywqZPCurk3AiN2phczVaM9ZzKPaYWthtWqkrs35rjrnlV25xteDJmEknctkUa3kzXuUCdOQSHMPvcDrvw1MlyuAK9VFcxC4hmUEfx0Oi15HPdlhuXgsWsmZQls7YZpMAogsb/xSXXK4IVk53nMFAly5hO7vbTZ56aknQ0vTp05mRkcF69RqwWbN5IWcxOp1SajNZHVtURL78suROa9WKvO++smVCrnKKisynwnXrxu2+YIH5usZxx6VA/iTRir+G88orxm6bDke8ErXZTLIqRvDrr4lNPBZLICoyFZAMj8l53RTzLHxVkqfmU1xYYl++9VYp4Wjmgmq1ymdOp/iLk6IQzjvPzBwUX06xtOZBFlegbYTCs/BPHMWPcBmPwh8licUij2mC9YaFTghws2pMu7WYnTCPaShgD/xU4ho5CWexD6bxOCyO8kS6Dm8ZFjLxQ/FF3F6yyYtMfhlTySuy5Vi8vBLvJ/W9HQ6p/EaS2R26Ge6U7wsWIfn9d0lgZhBtt2LFCh566KF0OBx89NFx/PRTKf9pxpQpEkjcsSP59NNitj/vvOhnwOEQF/0qdGApGwUF5qMijydu99WrzZ2CylBvvtrRir+Gc/HF5orSaCTs8Ug6CDPatTNXECHX4rIo09im4OdQvMkAZFR6CT5m48byPn3yifG0OBS0+tZb0c4Mv/ySOIGcNabWbmmtMTbRj/BNC5k2imBhNtw8BTPi7kfLQwpY5I73UcxHGkfhZgJ+OpFDK4pKqno9juF0RxSDj4wQboZ1hkXes+DhyZgVsSnA4fVGS0isQU+dDU/crMqsnXNOxFqpmZayWqPzipuwY8cO9gzW9B0xYgQDJnacESOiXe+dTlmLNbp8qIxojaFLl/iXy2o1LCZDSjI9o0wP331XzXKXAa34azhm6ZjNzB8ZGYl95hOZTc48s+LlTAEZWX+KC0mA6w8/rcQMkJlpfn6vV2S/+uqwi+PNNyeW14NMw1F6eCYg/6Yhnx5k8QeE67cauTOuxOFRmw5tkst8d23SaqUfKDG9ZMPN/9CCdbGDR0ESk3XFXO6FjztRp9SgsYfwBLPhZjEs9AeV/oe4jLEzGLebnDlxD4tq16M/Ij99Lhz8Gd3j9je8Rx4ZxJeQyNUpN9f4odmxQ5KlDR1KjhnD/D17eNVVVxEAL7nkEubGHLd9u/Fl7HbzfE7XXGN86ZSwdKk8jKFeyu2WRSqjtQWKebVjR9ktPV0OGzmyekUuK1rx13C2bo23ISZShk5nYm+CRMfu3Ck1JIw+K6ud/STMkj8iaq3OmFF6rExamshQWBiVn8z4uyKX36BvsIJUvBK8H0/xMnzIh/A41yLsa23mw14Ea5SPvgXFJMBFOJ5tsIKP4yGOxaW8Aa+V7HcpPqALOeyGOdyLdE7FGYapiWPbiVjA13Aj38HV7Gf5LqESb+/6l5NxdkmuoFG4KWGOosjm9UpwaQlmP4DFYvzgLFkS1mahE7ZqxcD27XzmsssIgF3S0rj13HMlipbk11+Xrfqc0ynxHKREAD/0kASFHXqoDHzKmnqnUti+XexTl14q9tYkFiJWrpSwAaPaHDUNrfj3A5YsEe81WzCXl1nEbiiHVCIaxqdnJyAjtECAfPRR48+bNCldyUS2jlgow80PPyQpLs2jRiVZatErrtyJ0s+4kcUb8Dq/Rn+mxbhghtpnuCAqjXNY8Rv3JjlwRpllamEPCfCIoDumceeTwydxHzOwmztRhwvRsSQFdqLmRjYfxgj2tM+vsHktUatbV0xoU6cGXSrNqt40bGi4Qhs49tj4jjItTcwhbje/gKR2bg7wT6+X/Pdfzp9vbNJTKlw+NPb33rxZnpHjjoue4bpckpiwmtzbq4WiIjF9phKt+PcjsrPFFm42kkpUujE/XwIOzUbQtWrJex+RzTmqHXus7JNMYkE78viE7VHO6X4vB13gZ+PG4VrByZqShg+Xl/2228K1C0LrGofW28Ox1sH8HcckHPm2h1TJKo5Q9DnKzZ29L6LfGVs20cnXcUPEpgCvVB9wKxoEF7fNZb0VL3MKzuBdeIHZcLENVhlk/gyP6NNQwGZYVxLk1cSyudIVvtMp98vhCDvqNGlCrpm4xLjsYIyRPRAgXx2xiwVmdYcj2iKAhwD0ApzauzcDAZm1xf7WLpcs+HbqJHKF3I9Dacm/+sq4w/B6JfZjf6a4WGLXevaUftNikbrXy5YFd8jPFxvttGnmJrdKRCv+/YxBg4zfv/R0ydZpxogRiUfbaWliqzSbTbhcMoV95x1JwdC1q9m5ArQoP21WP2228vmZh1KwhPj9d/Lxx8nnnw/mHzr7bK7E4TwEG5gojQPg55H4i19iALeiARejAy9QE9ipU0AWEJxOMiOD+crJz3BBnILvmPE3d6U1NEyZEHXvvcUMHHec3MS0NP5q6chGanNw4VnkO7L5XrZWa9gYm3gDXuc2SCrfXDj5AP5X5nsU6sCtKGR7/MUWwYC0Bg3k97nwwvjf22IJuhj++qu42zRoIOm7p06Ne16efpps4Mo09WiKHUFsBHg8pKzjnXe+zBUrAjzqqLDd2+eL7lu2bJH4q8iR/IMPmj2bAT7/fEXfnNTx++9yq2M7QqVkKWHHhB/lj/T0cK3gr7+uUpm04t/P+O4741F33brm0ZHFxcZZEyObxyPTT7MkW61aRZ/T75cXNTRCK+uicFpaKKOlseJPNOjZ1PvKoB29fB2LwxGMd9ixg5wzh03UJsP9LCpAdu7MnpafDVI4R3+X998nj2yTz4Z1Cuh2xcvlthfyJ2dvwxO8h6vK/B2sVvK8w/+kF5n0YR9dyGE7rOCbj27mY4+Zm+acTtM1yqjnJZR64Xv0YmFsPEFoOhGxbQsa8gTMpAUDCIA22w18/fVCLl8uWUTzs4tkgfjQQ8VT6ZZb4tIpvDmqkG4VnwDPZ8nip+NqdvFyMwoLDdP1lzSXM8Bn0x4y+CCYoK6KOCAVfyAgZpGypsqoKeTlmfs1BwISQOh0hqfwGRnm9Xqzs2V9IJEScbkkkIYMR9ZGKS23jPRjKSiQ8qEXX1x2xd+ypYzgzT4//HCZ3r/7roz+9+wJX/f+81bQXor5pbTvG5m00SyozGKRtDLXdFlKjyXHtKNxOEo3gSkV4MXW8XEfZEYU7LDZkr+PLkdx3EzEGixPCQRMTXpud8karCl794Y75MbYxNVozX3wMRsuiUHo1Ys899yoG9cRCyPKXd5HALRYenPKlOAPN3Bg9IOVliZ2xYic+Xvfm8Ba2MNILy2FYjZQ25j3VQ32jUzAd98ZD6YysIeX40Neg3d4s+X1+B2cTknfUUUccIp/2jSJ2LbZ5AW//fbUL6Qky5o1svZmtYr8ffualxr9919Rxp99ZlwRLsRDDyWOllVK7Oih2UJxsZSDDNWC9nqlFkYof9TChZI08YcfZPqauNyieUtPT07Jud0ig9tNTpokMvY6rWLZQ0Pf+f77JTi1Ii6sZVm36FN/MQvtYeWXBzv/waElAXO9e8tvnsy5rCWV0qK3e5DFw7HC9LhmzUpfKPX7o2eIFhSzN77jdXiLg49YKDvl5IjHi8PBlY5jeCp+4LmYyEYIrVe8TyCNXm87rp42zdjO6HaHo/VI8oEH+AeO5pH4S3IuIY/HYxFX2dpLDoj9kM8+i1f8Z+JrZsPNfcFgP8MIbaUSp7WtIAeU4l+4MH606nJFF+qpqWRnx9sBrVYZFBUWlv+8Zou1gCh2s3q/OTnSuYRmHtu2SSIuj0ce5KrLtWPeXC6JMTJLxljWlqieQTLNbpeMvGaui7E67u03irjj/uf5Dw7jejTlSNzB2thV8lt8/z05a1bioLVwM+780rGXxweLvEe20Kzk55+Te24k8Vj8/Z8xI2bHP/9kbsOWwaRz6cyFk89gWFC+2bRa67Cu18vZZgtMl1wSPteYMSW2w01ozC0IuqD5fDK13A/Zti164JWOvYaR24Yvp9k0vhI4oBT/OecYKySn0yTveTnZs0dyQlXmTMKsyFBFn3kzxW+1RngUxDBtmqz51asnI+JffpHZfWUEd1VU0b7wAjlxYvlnGpXZDjtMnoPS8g+5XOTxx4f90e++O/q39ngkgVuopOb550d/7nAk30k5kcsM7I7aZrHIYnxZU7yPHSvf0emUReHp0w12OuYYBmIejCx4eB6+oNNJ3nrrP2zbvDnTAI4xegivuy58rpwccSuNPJ/NJnbBCqf3TB1PPx3uRC/FWMNcTXFK/5JLqtSH9YBS/G3DaViiWnp6TPRiOcnJIS+6SF5Er1fO++abFT8vKWYHI9nT0io2y3344XjFpJTY/Y34/HPjkV4qRvhGLVRXONVyAOGMo337JpapV6/oIKRAQEb3F10kg5Xx46PXo/x+iWM45xzx4vr2W3LIELPfIDzyV/Az3SR4zOmMXiepFP7+23R68qM6mc2aSUzY7l272MvtJgA+AETX3w3ZY0OsXSs2r5C98+yzq3SRs7r4+WcpFPbyMe+x0GEwwlNKFr7PO09GeuXNLJokB5Tiv+wy88jAffsqfrMuuiheibrdht5wZWb8eHMf5ork/MjJkdG71yv3xuuVkfzKlfH7BgJiA061Qt1f2umny33LzBT/bLP92rQx/31yc8Vzb9KkiIjPoiLxfezTR/JoTJzIM/ubrWv4aUWxob0/sqWnx0TwVgaLFpm6gW1u2CFqll24bh2HNmxIADwfYHbsSzR3bvS5i4r261G+KZs3G08RPR5y9uxqE+OAUvwrVsQrT7c7ughIedm92zzNSY8eFT9/QYF0+JG2a4dDsm1WtPP3+2WE+dRTok/MXCVzcipu9z4QmlKizxwOMXUY3RO3W1Lyhvj8c/PzNWlifL9DHh8h9223mxz/aUCi9GJsQWc3/930/Mmsd5SWyqNcFBYal9xyOsn//S9u98AHH/Alu50K4AkAN0Xe8FtuMT7/3LmyeFfFI+BqY9EiURhWaziq0eORZEXVGJ58QCl+Ukw6p50m97J5c0kTUBn3c9Uq4xE5IHb0ymDnTvn9MzKkutYtt1TOTCVZ/P7kFiors6V63cCs+Xyy0Lp8uSxyP/aYxDK43RJxGTtAnTvXvNPs3j3+Xu/ebTzwc9mLuc4Vb7M83zLRtJMq7bu43eS111bRQxOyDYa+vNstvrhGD+7YsaTXy68hUb5NAC4OfYlIcw8pC021a8sP4fWK7/+iRVX0JQz45huxh9auLdO5ZFfFEzFypNyf0ENvs4kL4o8/VntOigNO8VcVBQXGStFikTKC1cnOnTKLOfxwCX8fN67ynpsRI4zNtjXFxl+dLTTyd7vFsy7yHk+bJvnnDjtMklauWSNBdLHncDolFiCSv/8mmzY1vqbdUsincW/cB6/bbqU7Lb4YTmjQaPZ7+XzkAw9UsdXkr7/IG2+UFerRo839i3fvLnHr/AOS38cN8Eu7XUb1ITZuNH4Ia9VK7Lscw6+/StGX8ePLmOjt00+Ng1lmzSrDSWLYutW4p/d6JY9FNaMVfxl4443o58Fqlc5g9erqk2HfPpnJRC4mejziKVIZhCJyPZ5wudFUK+Ca0CJdzt98M/o5sNlkYDhhQrijsNnk30svjbZSFBSQjRqZd6QKAT5geSbugyxPQ7aqnxX1u7vdYhWK1VFpaeSJJ4qyr3HJzT77TB4sl4tb7HZ2VooK4LPPPhvO7f/008Z2VZ9PijqUQlEROWCA3BeHQw6rW9fciy2KQMC8V+7YsfzfOzjbMTzvkCHlP2850Yq/jEyZIlP95s1lMbm0KMjK5oUXjGNhStIQVBL5+TKQK6u/fChwrrIUbp06iVPIR7eKBXaV1po1k/titJ5ps4mZ7rbbZGCnlMQ9xC6oTpqU2E/f4/ZzruPU+A98Pu5al8V77hGT09FHi6+93y/rNpElOHv2lKzCNZatW8nXXydHjmTuH3/woosuIgAOGTKEBQUF5J13Gt8cp1Nq6ZbC6NHx91gpmSGX2hHm5iZOWFVeJk40fnAsFnloqhmt+Pcz+vQxfibT0yX7X2Wye3fZXCc9HvFXr4zgKotFFFxBgZhX775bMgGbpTBOQy6diM/zkkwrixnrwgvNB26hWVLkNq9XzEAhRo40P7fVKjOEwNdTwiu+Pp+4YRlE2gUCEkF9zz2St/7778kNG8x/zwULxDyVkUG2by8mkJqA3+/nI488QgA8+eSTudOsVJvLJUVSSqFDB+P7m0y6CgYCxgvWgHhflJfcXGNbsctVOb7mZUQr/gqQlyeRedUxnV61SuyVPXsam188nqpZ++rWLTnFaLFIzp7kR+eJW1qa+MmPGBGOXP7vP6MgtwAVinkzXuUNeN00N3+i5nAkX3IylCK6LN8j0mHlf/8z37dr14hnKT9fbMo//2xooC8uFpN6SD/a7aLYzJI6GkW1u90yOq4O/H6ZoRx7rLi3PvhgfFzB2LFj6XA4eNhhh3Fl587RAns8SZtEjjrK+P56POL5VypPPGF8sz74oKxfO5rZs6M7dKeTfPXVip2znGjFXw7y8sRLwuGQ1rix5K+pKu67T56RUA7z2AfaZpOHvSo6oNmzk7PzN20qjhBmg6WKtKOPDn+3uXPDAycLitkfk7kLtUmAG3EI62BnMFlY5csRqQPK0sF16xa+n2+9ZX4/jTwazfj4Y+NI7/R04wR/vY0Tg7J27epxl7/iivgi623axLsWz507l/Xr12etWrX4w913y0jn9NNlepLkA/7ss8bvSTJ5ikhKL/Xoo9KrOp2yqPzKK2X+zobk5oq9b9y4lNrjtOIvB5deGv9gud2V4/EVy88/m9uEPR55Lrt1K39w4/r14qESaY6IpFu30hW/UhJg+e+/iTNdVkTZRo5k580LLjyjiN+iLzPhoR9gHhxc52zDof3WsVUr6ZCrorpVrVryfUMR3IlmRLEj/v/+M+403O7E9RRi6dfP+Hrp6Qb5dGjscQRUefZfkuZpLTweyb4ay3///ccjjzySNpuNb7/9dpmvl5sbDloE5Npeb7wLbqkUFsqU/gAMJNOKv4zs2GGu3M44o/KvN3SosWKx2USOhg2lXmlZn83CQol3CJ07lMYhcgSWrI3f7RbXOTIuWy8BeemMRqdl6QwuuCBa/pBnTYbPz7Od3/OdWvdw573PRiWk2blT0ryErl3W4DSj+26xhGXZtk2cVMxs/qHvvny5pM5o1EhmREcdFZ0Gw+uVtYOyzNjOPtv4ej6fcQCoWWputzt+hrBvn5iAbr1VckiVwYPSkI8+Mr9Hl15qfMy+ffvYr18/AuBdd93F4jLmWC8uloped98t2Y1r9GJ3CtCKv4yEak8bPcSHH17517v66uQUr9kLZMbllxufq1ev8D47d5or/lA64mOPjS6Ll58vL5vPJ5/37En+8Yd4wsWeoywLx0bm3exsGSX/9pu50szKEiU2aBA5bFji8pORLS1N7kVkbJJS0vr3D9uKf/3VvHiNzyeynXlm9AzRahUTyxVXyO8wZUrZA1O/+sq4M61Xz3gQ8PXXxmbre+6J3m/1ajlH6Nwej5jxKjIrmDHD+B7Z7WLrN6OoqIi33XYbAfCss85i5v5QxXw/QSv+MpKdbWw/tFqrJv3z9OnJ1bl1OsWMkCyJRr+RSui44+IVpcMhOfuT5fffza+X7Ci8IrEzsTz6aLy5RSlR9iHzTfv2kuJg2bL4BXWLRfZZtkxGlkYVltxumZUsXWqeir4idTYCAfL66+XcIVOGz5fYnPHBB+EaCh6PBAHGDqSN6hPYbNHZk8uK3y8eWrG/tccTXRDHjNdff51Wq5XHHHMM161bV35BNCVoxV8OHnssWhmHIjz//rvi587LExe9WbPEHBMISHRuaYqxrO6cic6VlRXeb9kyGZ2Gvm9IKe7dm/y1nnzS/FqllYQMKdpRo5K/XmkEAhKJG6pg5nSKSW3KFPKll2TNI9T5mZWtVErMWqT8VqH1ltA9Ov10+f0+/th8RnDhhRX/LsuWybrjhx9GJHlLgN8vMzmjGg9FReYdscdTMTnXrpWgslDluCZNjNcizJg2bRrT09PZsGFDLohMkqQpF9Wu+AE0AzALwHIAywDcHtxeB8B0AP8E/61d2rlSpfgDAamx2q6dBBgNGJBkVGApTJ0a9vZKTxeF++abyQVEud1lk8HsnFZrvNlk3z7xRrnvPvFeMisMs3atuF/edJOYIkKjyQceMJc72VwzofTX//4r+eq93nBh8fLWRcjJEft7ok5s40bzexVKy0yKvf+llyS99vffhzuOefOMZ2wOh9yrmkRxsflieEZG5Vxj40YZIJUn59qyZcvYqlUrOp1Ojq8pQQj7KalQ/I0BHB/82wfgbwDtATwH4P7g9vsBPFvauVLtx1+ZbN5srGDS0kpXjg4HedJJZbveM/FZAQjIgl55mDJFFHQoeMvjkTKSBQWyAGr2HZLxunG5RLHu2CHeKZGmCJcrPPKuCnbtMg9I69Ch9OMDATGXxZ7D4YheG6kpnH++saw33ZRqyYTt27ezR48eBMDHH388nOZBUyZSbuoB8BWA3gBWAWjMcOewqrRjDyTFP3KksZufmWK0WKTZ7bKwW54snk8/HfYusdtl8bM875FZAjuHg7z5ZvNCLmlpiaN8Q/bryy+XxeEWLYz3d7kqx8wWy5o14jVl9Bu43ZKYMhl27RKFGuqwlAqnthg2rPLlrgjbt4uTQsgE5vWKR1BNWlfNz8/nFVdcQQC87LLLmFemDGwaMsWKH0BLAOsBpAPYG7FdRf7frB1Iiv/BB42Vn81m7P0SWswti5lj506JHO3bV7Lgrl4tU+69e+MX+crCnDnmtmGzkX6DBmKXTqT4vV4J3irN1JWeLgnSKpsePYxdTu12MeuUhT17zBd5zeoepwq/X9Y5Ro2S9YuaOKgOBAJ88sknCYDdunXjtm3bUi3SfkXKFD8AL4DFAM4P/n9vzOd7TI4bCmARgEXNmzevwltTvfz0k7Et2OUST4uQH7TVKsqirJHe69eL90loATItTa73008Vl33ixMSK2ahdcYUklUtk6glFK5d2LrdbXEYrk8xM804p0rafLOPGGS/yKiWzIk35+Pzzz+lyudiyZUsuTSKPj0YwU/wWVCFKqTQAXwD4mOTE4OZtSqnGwc8bA9hudCzJt0l2JNmxfv36VSlmtdKjB9C3L+DxhLd5PMBllwE//ACMHw9cfTVw663AL7/Iv2XhoYeA3buB/Hz5f1ERkJMDXHONqKCKkJ4OKGX8mdl2AJg4EUhLM/+8oEBaIux24LjjgGOOKV3OslDRe2J0PrNzVva1DiYuuOACzJ49G/n5+ejatSumTZuWapH2b4x6g8poEDPOhwBejtn+PKIXd58r7VwHkqmHlCn2+PES8HPOOZLSo7Km2WZuk3a72HV/+kk8Z8pTDGjHDmNzlMViPKJ3uSTBWMeOiUf8aWmlf37ppWVzLS0LRukqHA7yjjvKfq6IGiRxs5U5cypf9oON9evXs0OHDrRYLBxVmb6/ByhIgVdPDwAE8CeAJcHWH0BdADMg7pw/AKhT2rkONMVflbRsaa48jztOTEkul/zboUN85sTSuPvu6MhQi0U8cN59N7xIa7UmnzvH4ZCo4ERup0nlV68A//wjZp2Qmc3rlXQL5S2HOXZs2HwVWtwtTyeiMSYrK4sDBgwgAN58880sOgBz7FQW1a74K7NpxZ88L7wQH7Jvt4unTOxo3W4veznJQECUfPv2kpNm8OBwVObmzZJS+p57Slf4oWC4hx8WP/upU81zI3m9lbNGkYicHIl4fewxCZCryCI4KfnyR46UoLbKXpfQkMXFxRw2bBgBsE+fPtxbVdPB/Ryt+A8SiotFmTud4gXjdksGQ7MRtcNR+aPpG24wV/gNG0oWSyOzx8CBxsf4fBIZq9HE8u6779Jms7F9+/ZcY5Z69iDGTPFX6eKupvqxWoEPPwRWrgQ++giYNw+YPx/w+433Ly6ufBny8sw/a9QIGDUK6N49/rNevQC3O357cTFw4omVJ5/mwOGaa67B999/jy1btqBz586YO3duqkXaL9CK/wClRQtgwICwF0zfvtIpRGKxAL17J/bIKQ+33Wb+2TXXmH82eDBQv7548IRwu4HzzgPatKk8+TQHFqeeeirmz5+P2rVr47TTTsPYsWNTLVKNRyv+g4TXXgPq1g27kbrd8v/Royv/WscfD5xzTvz2Nm2Am282P87jARYvBm68EWjWDGjXDnjmGZnBaDSJOPzwwzF//nx0794dV1xxBYYPH45AIJBqsWosSsxANZuOHTty0aJFqRZjvycrCxg7FliyBDj2WODyy8U3v6r44gvguefER/+aa0TpW/RQQ1OFFBYW4qabbsJ7772HQYMGYcyYMXAb2Q8PEpRSi0l2jNuuFb9GozmQIIkXX3wRw4YNQ8eOHfHVV1+hcePGqRYrJZgpfj3+0mg0BxRKKdx9992YNGkSli9fjk6dOmHJkiWpFqtGoRW/RqM5IBkwYADmzJkDAOjRowcmT56cYolqDlrxazSaA5YOHTpg4cKFOOKII3Duuedi5MiR2B/M21WNVvwajeaApnHjxpg9ezYGDhyIe+65B0OHDkVhYWGqxUopWvFrNJoDHrfbjfHjx2P48OF499130a9fP+zevTvVYqUMrfg1Gs1BgcViwRNPPIEPP/wQc+fORdeuXfHPP/+kWqyUoBW/RqM5qLjiiiswY8YM7N69G507d8aPP/6YapGqHa34NRrNQUePHj2wYMECNGrUCL1798Z7772XapGqFa34NRrNQUnr1q0xb948nHbaabj22mtx7733wm+WzfAAQyt+jUZz0JKRkYGpU6fipptuwvPPP4+BAwciOzs71WJVOVrxazSagxqbzYbXX38do0aNwtdff42ePXti48aNqRarStGKX6PRaADccsstmDp1KtasWYNOnTrhQM4PphW/RqPRBOnXrx9++eUXOBwOnHTSSZgwYUKqRaoStOLXaDSaCI466igsWLAAHTp0wKBBg/D0008fcGketOLXaDSaGBo0aICZM2fi0ksvxYMPPoirrroKBQUFqRar0rClWgCNRqOpiTidTowdOxbt2rXDI488gn///Rdffvkl6tWrl2rRKowe8Ws0Go0JSik8/PDD+PTTT7Fo0SJ07twZK1asSLVYFUYrfo1GoymFiy66CD/++CNycnLQtWtXfP/996kWqUJoxa/RaDRJ0LlzZyxcuBAtWrRA//798cYbb6RapHKjFb9Go9EkSfPmzTFnzhycccYZuOmmm3D77bejuLg41WKVGa34NRqNpgz4fD5MmjQJd955J1599VUMGDAAmZmZqRarTGjFr9FoNGXEarXixRdfxFtvvYXp06eje/fuWLt2barFSpoqU/xKqf9TSm1XSi2N2FZHKTVdKfVP8N/aVXV9jUajqWqGDh2KadOmYePGjejUqRPmzZuXapGSoipH/GMA9IvZdj+AGSTbAJgR/L9Go9Hst/Tq1Qvz5s1Deno6Tj31VIwbNy7VIpVKlSl+kj8BiC1qeQ6AD4J/fwDg3Kq6vkaj0VQX7dq1w4IFC9C5c2dcdtllePTRR2t0mgdVlcIppVoCmELyqOD/95KsFfxbAdgT+r/BsUMBDA3+9ygAS432qwHUA7Az1UIkoCbLp2UrPzVZvposG1Cz5ats2VqQrB+7MWUpG0hSKWXa65B8G8DbAKCUWkSyY7UJVwZqsmxAzZZPy1Z+arJ8NVk2oGbLV12yVbdXzzalVGMACP67vZqvr9FoNAc91a34JwO4Mvj3lQC+qubrazQazUFPVbpzfgJgHoC2SqmNSqlrADwDoLdS6h8Apwf/nwxvV5GYlUFNlg2o2fJp2cpPTZavJssG1Gz5qkW2Kl3c1Wg0Gk3NQ0fuajQazUGGVvwajUZzkFEjFb9SyqqU+l0pNSX4/1ZKqQVKqdVKqfFKKXsKZaullJqglFqplFqhlOpaU1JRKKXuVEotU0otVUp9opRypvLelSVthxJeDcr5p1Lq+BTI9nzwd/1TKfWlUqpWxGcPBGVbpZTqW92yRXx2t1KKSql6wf9X631LJJ9S6tbg/VumlHouYntK751SqoNSar5SaolSapFSqlNwe3U/c82UUrOUUsuD9+j24PbqfydI1rgG4C4A4yDBXwDwGYCLg3+/CeDGFMr2AYBrg3/bAdQC8ByA+4Pb7gfwbArkagLgPwCuiHt2VSrvHYCTABwPYGnENsN7BaA/gG8BKABdACxIgWx9ANiCfz8bIVt7AH8AcABoBWANAGt1yhbc3gzAdwDWAaiXivuW4N6dCuAHAI7g/xvUlHsH4HsAZ0Tcrx9T9Mw1BnB88G8fgL+D96fa34kaN+JXSjUFcCaAd4P/VwBOAzAhuEvKUj0opTIgD9Z7AECykORe1JxUFDYALqWUDYAbwBak8N6xbGk7zgHwIYX5AGoFYz2qTTaS35MMJVefD6BphGyfkiwg+R+A1QA6VadsQV4CcC+ASI+Mar1vCeS7EcAzJAuC+4RidGrCvSOA9ODfGQA2R8hWnc/cFpK/Bf/OArACMmCr9neixil+AC9DHu5A8P91AeyNeCE3Qm5WKmgFYAeA94OmqHeVUh4ADUluCe6zFUDD6haM5CYALwBYD1H4+wAsRs25dyHM7lUTABsi9ku1rFdDRltADZBNKXUOgE0k/4j5KOWyBTkcQM+gWXG2UurE4PaaIN8dAJ5XSm2AvCMPBLenTDYl6WyOA7AAKXgnapTiV0qdBWA7ycWplsUEG2Qa+QbJ4wDkICbDKGWOVu0+skG74DmQzukQAB7EZ0etUaTqXpWGUuohAMUAPk61LACglHIDeBDAI6mWJQE2AHUgJolhAD4LztZrAjcCuJNkMwB3IjhjTxVKKS+ALwDcQTKqgkt1vRM1SvED6A5ggFJqLYBPIWaKVyBTnFBeoaYANqVGPGwEsJHkguD/J0A6gpqQiuJ0AP+R3EGyCMBEyP2sKfcuhNm92gSxYYdIiaxKqasAnAXgsuBLCKRetkMhHfofwXejKYDflFKNaoBsITYCmBg0SyyEzNjr1RD5roS8DwDwOcKmpmqXTSmVBlH6H5MMyVTt70SNUvwkHyDZlGRLABcDmEnyMgCzAFwQ3C1lqR5IbgWwQSnVNripF4DlqBmpKNYD6KKUcgdHWiHZasS9i8DsXk0GMDjoydAFwL6I6W+1oJTqBzEzDiCZG/HRZAAXK6UcSqlWANoAWFhdcpH8i2QDki2D78ZGyCLhVtSA+xZkEmSBF0qpwyGODzuR4nsXZDOAk4N/nwbgn+Df1Xrvgu/lewBWkHwx4qPqfyeqchW7Ig3AKQh79bSGPCyrIT22I4VydQCwCMCfkIe9NmQdYgbkgfoBQJ0UyfYYgJWQFNYfQTwpUnbvAHwCWW8ogiira8zuFcRz4XWI18dfADqmQLbVEJvqkmB7M2L/h4KyrULQQ6Q6ZYv5fC3CXj3Vet8S3Ds7gLHBZ+83AKfVlHsHoAdkvesPiE39hBQ9cz0gZpw/I56x/ql4J3TKBo1GoznIqFGmHo1Go9FUPVrxazQazUGGVvwajUZzkKEVv0aj0RxkaMWv0Wg0Bxla8WsOSpRSDZVS45RS/yqlFiul5imlzlNKnaKCWWE1mgMVrfg1Bx3BQJpJAH4i2ZrkCZCAwaYJD9RoDhC04tccjJwGoJDkm6ENJNeRHBW5k1JqhFLqnoj/Lw0m14JSanAwR/ofSqmPgttaKqVmBrfPUEo1D24fFDz2D6XUT8FtViX5/38N7n991X9tjUawlb6LRnPAcSQkurRcKKWOBDAcQDeSO5VSdYIfjQLwAckPlFJXA3gVkmL3EQB9SW5S4eIu10BC8E9USjkAzFVKfU9JXazRVCl6xK856FFKvR4cjf+a5CGnAfic5E4AIBnK/94VUkAIkJQZPYJ/zwUwRil1HQBrcFsfSB6WJZA0AnUheWw0mipHj/g1ByPLAAwM/YfkzUpKGS6K2a8Y0YMjZ3kuRvIGpVRnSIGhxUqpEyB5WG4l+V15zqnRVAQ94tccjMwE4FRK3RixzW2w31pI2m0E6522ijh+kFKqbvCzkKnnF8giMQBcBuDn4OeHklxA8hFIIZ9QCcUbg2l6oZQ6PFjUR6OpcvSIX3PQQZJKqXMBvKSUuheijHMA3Bez6xcQc8wyiDnm7+Dxy5RSTwKYrZTyA/gdUt/4Vkh1tmHBcw4Jnud5pVQbyCh/BiRL5J8AWkLy6qvg/udWxffVaGLR2Tk1Go3mIEObejQajeYgQyt+jUajOcjQil+j0WgOMrTi12g0moMMrfg1Go3mIEMrfo1GoznI0Ipfo9FoDjL+H65XizzugYiRAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "model = sm.Logit( d.Outcome, sm.add_constant(X.loc[:,[\"Glucose\",\"BMI\"]]) )\n",
    "model = model.fit()\n",
    "print(model.summary())\n",
    "\n",
    "\n",
    "fig,ax = plt.subplots()\n",
    "ax.scatter( X.Glucose, X.BMI, color=[\"red\" if _ else \"blue\" for _ in d.Outcome.values]  )\n",
    "\n",
    "def f(x):\n",
    "    return (7.51-x*0.035)/(0.0763)\n",
    "ax.plot([40,200],[ f(40), f(200) ], color=\"black\")\n",
    "\n",
    "ax.set_xlim(40,210)\n",
    "ax.set_ylim(10,70)\n",
    "\n",
    "ax.set_xlabel(\"Glucose\")\n",
    "ax.set_ylabel(\"BMI\")\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Decision boundary\n",
    "\n",
    "We can use logistic regression to compute a **decision boundary**. A decision boundary is when choosing to label an observation as a \"one\" or a \"zero\" is ambiguous. \n",
    "Often for logistic regression choosing to label an observation as a one or zero is ambiguous when the probability that observation equals one is $1/2$. \n",
    "\n",
    "For example, suppose we fit the following logistic regression model\n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left( \\frac{\\theta }{1-\\theta} \\right) = \\beta_{0} + \\beta_{1}x^{1}_{i} + \\beta_{2}x^{2}_{i}\n",
    "\\end{align}\n",
    "\n",
    "and our fitted MLEs for our parameters are $\\beta_{0}=1$, $\\beta_{1}=-1$, and $\\beta_{2}=1/2$. \n",
    "\n",
    "The MLE for the above is \n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left( \\frac{\\theta }{1-\\theta} \\right) = 1 -x^{1}_{i} + (1/2)x^{2}_{i}\n",
    "\\end{align}\n",
    "\n",
    "Lets look at the observation $(x^{1}=1, x^{2}=0)$.\n",
    "To compute the log odds we can plug in the above values for $x^{1}$ and $x^{2}$\n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left( \\frac{\\theta }{1-\\theta} \\right) &= 1 -1 + (1/2)0 \\\\ \n",
    "    \\log \\left( \\frac{\\theta }{1-\\theta} \\right) &= 0 \\\\ \n",
    "\\end{align}\n",
    "\n",
    "We find the log odds is zero, but what does that mean for our estimated $\\theta$ (i.e. the estimated probability this observation is a one)?\n",
    "\\begin{align}\n",
    "    \\log \\left( \\frac{\\theta }{1-\\theta} \\right) &= 0 \\\\\n",
    "    \\frac{\\theta }{1-\\theta} &= e^{0} = 1\\\\\n",
    "    \\theta  &= 1-\\theta\\\\\n",
    "    2\\theta &= 1\\\\\n",
    "    \\theta &= 1/2\\\\\n",
    "\\end{align}\n",
    "\n",
    "Ah ha! When the log odds equals zero the probability we would assign to our observation equals $1/2$---we are at a decision boundary. We can build an algebraic equation for the $(x^{1},x^{2})$ pairs that are on the decision boundary by settign the log odds to zero. \n",
    "\n",
    "\n",
    "\\begin{align}\n",
    "    \\log \\left( \\frac{\\theta }{1-\\theta} \\right) &= \\beta_{0} + \\beta_{1}x^{1}_{i} + \\beta_{2}x^{2}_{i}\\\\\n",
    "    0 &= \\beta_{0} + \\beta_{1}x^{1}_{i} + \\beta_{2}x^{2}_{i} & \\text{when at a decision boundary}\n",
    "\\end{align}\n",
    "\n",
    "Finally we can plot this decision boundary on a $x^{1}$ and $x^{2}$ axis by putting the above algebraic equation in slope-intercept form.\n",
    "\n",
    "\\begin{align}\n",
    "    0 &= \\beta_{0} + \\beta_{1}x^{1}_{i} + \\beta_{2}x^{2}_{i} \\\\\n",
    "    -\\beta_{1}x^{1}_{i} &= \\beta_{0} + \\beta_{2}x^{2}_{i} \\\\\n",
    "    x^{1}_{i} &= \\frac{\\beta_{0}}{-\\beta_{1}} + \\frac{\\beta_{2}}{-\\beta_{1}}x^{2}_{i} \\\\\n",
    "\\end{align}\n",
    "\n",
    "The above is the equation of a line with intercept $\\frac{\\beta_{0}}{-\\beta_{1}}$ and slope $\\frac{\\beta_{2}}{-\\beta_{1}}$."
   ]
  }
 ],
 "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.9.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}