DSA-Python-Cap08-07-StatsModels.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# <font color='blue'>Data Science Academy - Python Fundamentos - Capítulo 8</font>\n",
    "\n",
    "## Download: http://github.com/dsacademybr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Versão da Linguagem Python Usada Neste Jupyter Notebook: 3.8.8\n"
     ]
    }
   ],
   "source": [
    "# Versão da Linguagem Python\n",
    "from platform import python_version\n",
    "print('Versão da Linguagem Python Usada Neste Jupyter Notebook:', python_version())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Statsmodels"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Linear Regression Models"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Para visualização de gráficos\n",
    "from pylab import *\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import statsmodels as st\n",
    "import sys\n",
    "import warnings\n",
    "if not sys.warnoptions:\n",
    "    warnings.simplefilter(\"ignore\")\n",
    "warnings.simplefilter(action='ignore', category=FutureWarning)\n",
    "warnings.filterwarnings(\"ignore\", category=FutureWarning)\n",
    "warnings.filterwarnings(\"ignore\")\n",
    "import statsmodels.api as sm\n",
    "from statsmodels.sandbox.regression.predstd import wls_prediction_std\n",
    "np.random.seed(9876789)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'1.21.0'"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.__version__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'1.3.0'"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pd.__version__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'0.12.2'"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "st.__version__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Criando dados artificiais\n",
    "nsample = 100\n",
    "x = np.linspace(0, 10, 100)\n",
    "X = np.column_stack((x, x**2))\n",
    "beta = np.array([1, 0.1, 10])\n",
    "e = np.random.normal(size=nsample)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = sm.add_constant(X)\n",
    "y = np.dot(X, beta) + e"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                            OLS Regression Results                            \n",
      "==============================================================================\n",
      "Dep. Variable:                      y   R-squared:                       1.000\n",
      "Model:                            OLS   Adj. R-squared:                  1.000\n",
      "Method:                 Least Squares   F-statistic:                 4.020e+06\n",
      "Date:                Fri, 09 Jul 2021   Prob (F-statistic):          2.83e-239\n",
      "Time:                        20:24:41   Log-Likelihood:                -146.51\n",
      "No. Observations:                 100   AIC:                             299.0\n",
      "Df Residuals:                      97   BIC:                             306.8\n",
      "Df Model:                           2                                         \n",
      "Covariance Type:            nonrobust                                         \n",
      "==============================================================================\n",
      "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "const          1.3423      0.313      4.292      0.000       0.722       1.963\n",
      "x1            -0.0402      0.145     -0.278      0.781      -0.327       0.247\n",
      "x2            10.0103      0.014    715.745      0.000       9.982      10.038\n",
      "==============================================================================\n",
      "Omnibus:                        2.042   Durbin-Watson:                   2.274\n",
      "Prob(Omnibus):                  0.360   Jarque-Bera (JB):                1.875\n",
      "Skew:                           0.234   Prob(JB):                        0.392\n",
      "Kurtosis:                       2.519   Cond. No.                         144.\n",
      "==============================================================================\n",
      "\n",
      "Notes:\n",
      "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n"
     ]
    }
   ],
   "source": [
    "model = sm.OLS(y, X)\n",
    "results = model.fit()\n",
    "print(results.summary())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Parameters:  [ 1.34233516 -0.04024948 10.01025357]\n",
      "R2:  0.9999879365025871\n"
     ]
    }
   ],
   "source": [
    "print('Parameters: ', results.params)\n",
    "print('R2: ', results.rsquared)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                            OLS Regression Results                            \n",
      "==============================================================================\n",
      "Dep. Variable:                      y   R-squared:                       0.933\n",
      "Model:                            OLS   Adj. R-squared:                  0.928\n",
      "Method:                 Least Squares   F-statistic:                     211.8\n",
      "Date:                Fri, 09 Jul 2021   Prob (F-statistic):           6.30e-27\n",
      "Time:                        20:24:41   Log-Likelihood:                -34.438\n",
      "No. Observations:                  50   AIC:                             76.88\n",
      "Df Residuals:                      46   BIC:                             84.52\n",
      "Df Model:                           3                                         \n",
      "Covariance Type:            nonrobust                                         \n",
      "==============================================================================\n",
      "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "x1             0.4687      0.026     17.751      0.000       0.416       0.522\n",
      "x2             0.4836      0.104      4.659      0.000       0.275       0.693\n",
      "x3            -0.0174      0.002     -7.507      0.000      -0.022      -0.013\n",
      "const          5.2058      0.171     30.405      0.000       4.861       5.550\n",
      "==============================================================================\n",
      "Omnibus:                        0.655   Durbin-Watson:                   2.896\n",
      "Prob(Omnibus):                  0.721   Jarque-Bera (JB):                0.360\n",
      "Skew:                           0.207   Prob(JB):                        0.835\n",
      "Kurtosis:                       3.026   Cond. No.                         221.\n",
      "==============================================================================\n",
      "\n",
      "Notes:\n",
      "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n"
     ]
    }
   ],
   "source": [
    "nsample = 50\n",
    "sig = 0.5\n",
    "x = np.linspace(0, 20, nsample)\n",
    "X = np.column_stack((x, np.sin(x), (x-5)**2, np.ones(nsample)))\n",
    "beta = [0.5, 0.5, -0.02, 5.]\n",
    "\n",
    "y_true = np.dot(X, beta)\n",
    "y = y_true + sig * np.random.normal(size=nsample)\n",
    "\n",
    "res = sm.OLS(y, X).fit()\n",
    "print(res.summary())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Parameters:  [ 0.46872448  0.48360119 -0.01740479  5.20584496]\n",
      "Standard errors:  [0.02640602 0.10380518 0.00231847 0.17121765]\n",
      "Predicted values:  [ 4.77072516  5.22213464  5.63620761  5.98658823  6.25643234  6.44117491\n",
      "  6.54928009  6.60085051  6.62432454  6.6518039   6.71377946  6.83412169\n",
      "  7.02615877  7.29048685  7.61487206  7.97626054  8.34456611  8.68761335\n",
      "  8.97642389  9.18997755  9.31866582  9.36587056  9.34740836  9.28893189\n",
      "  9.22171529  9.17751587  9.1833565   9.25708583  9.40444579  9.61812821\n",
      "  9.87897556 10.15912843 10.42660281 10.65054491 10.8063004  10.87946503\n",
      " 10.86825119 10.78378163 10.64826203 10.49133265 10.34519853 10.23933827\n",
      " 10.19566084 10.22490593 10.32487947 10.48081414 10.66779556 10.85485568\n",
      " 11.01006072 11.10575781]\n"
     ]
    }
   ],
   "source": [
    "print('Parameters: ', res.params)\n",
    "print('Standard errors: ', res.bse)\n",
    "print('Predicted values: ', res.predict())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7fd5c8d630d0>"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeAAAAFlCAYAAAAzqTv+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAABlE0lEQVR4nO3dd3gUVRfA4d9kU+lVhADSlC4goEREIkhRqtgAK/aKWBBQFBD9AMGGKDZEEUQE6Yj0UCT0Ih0pooReQktC2v3+OAkESELKbM15nycPyWYzcydL9sxt51jGGJRSSinlWn7uboBSSimVF2kAVkoppdxAA7BSSinlBhqAlVJKKTfQAKyUUkq5gQZgpZRSyg38XXmyEiVKmAoVKrjylEoppZTbrF279pgxpmR633NpAK5QoQJr1qxx5SmVUkopt7Esa19G39MhaKWUUsoNNAArpZRSbqABWCmllHIDl84BpychIYH9+/cTFxfn7qY4TXBwMGXLliUgIMDdTVFKKeUh3B6A9+/fT8GCBalQoQKWZbm7ObYzxnD8+HH2799PxYoV3d0cpZRSHsLtQ9BxcXEUL17cJ4MvgGVZFC9e3Kd7+EoppbLP7QEY8Nngm8rXr08ppVT2eUQA9iT9+/dn2LBhGX5/6tSpbN261YUtUkop5Yu8LgBPXR9F48ELqdh7Fo0HL2Tq+ijXnl8DsFJKKRt4VQCeuj6KPpM3ERUdiwGiomPpM3lTroPwBx98wA033MBtt93Gjh07APj2229p2LAhderU4d577yUmJobly5czffp0evbsSd26ddm9e3e6z1NKKaWuxqsC8NA5O4hNSLrksdiEJIbO2ZHjY65du5ZffvmFDRs28Pvvv7N69WoAOnXqxOrVq9m4cSPVq1dn1KhR3HrrrbRv356hQ4eyYcMGKleunO7zlFJKqatx+zak7DgQHZutx7Ni6dKl3HPPPeTLlw+A9u3bA7B582b69u1LdHQ0Z8+epVWrVun+fFafp5RSysOdOgWFC7vsdF7VAy5TJCRbj+fG448/zogRI9i0aRP9+vXLcBtRVp+nlFLKw/XtC1GuW1fkVQG4Z6uqhAQ4LnksJMBBz1ZVc3zM22+/nalTpxIbG8uZM2eYMWMGAGfOnKF06dIkJCQwbty4C88vWLAgZ86cufB1Rs9TSinl4dasgXvvhZUr5evevcGFGQu9KgB3rBfKoE61CS0SggWEFglhUKfadKwXmuNj3nTTTTz44IPUqVOHu+66i4YNGwIwcOBAbrnlFho3bky1atUuPL9z584MHTqUevXqsXv37gyfp5RSygMZAxER0LIlNGwICxbAnj3yvdBQuOYalzXFMsa47GQNGjQwl9cD3rZtG9WrV3dZG9wlr1ynUkp5tLZtYdYsKFUKXnsNnnsOChVy2uksy1prjGmQ3veuugjLsqzvgbbAEWNMrZTHhgLtgHhgN9DNGBNtW4uVUkopu5w/D0FB8nnr1nD33dCtG4TYv34oO7IyBP0D0Pqyx+YBtYwxNwI7gT42t0sppZTKvXXr4MYbYexY+fqll+CFF9wefCELAdgYswQ4cdljc40xiSlfrgDKOqFtSimlVM4kJ8OwYdCoEZw7B2U9L0zZsQ/4CWCCDcdRSimlcu/AAXjsMZg/Hzp1gm++geLF3d2qK+RqFbRlWW8DiUCG+28sy3rGsqw1lmWtOXr0aG5Op5RSSl3dypWwfLkE3kmTPDL4Qi56wJZlPY4szmpuMllKbYz5BvgGZBV0Ts+nlFJ2mro+iqFzdnAgOpYyRULo2apq1rc0GgObN8O0afDPP/DII9C0qWRSWrBA3vCLFZN/ixe/uABIOU9cHPz5JzRvDvfcA7t3w7XXurtVmcpRALYsqzXwJtDUGOPV1QeOHz9O8+bNATh06BAOh4OSJUsCsGrVKgIDA93ZPKWUE6QWdknNLZ9a2AXIPAifOQP9+sHUqbB3rzxWqpQEX4Dt2yWxQ1qBgTBuHNx3n81XoS44eFCC7oYN8rqULu3xwReytg1pPBAOlLAsaz/QD1n1HATMSyk2v8IY85wT2+k0xYsXZ8OGDYDUAi5QoABvvPHGhe8nJibi7+9VKbOVUleRWWGXSwLwuXMwZw7ExMDDD0O+fPDbb1CrlmRNatdO3uxT1a4N69fDiRNw/Lj8u3mzLAQCmDtXsi89/jiUKeP8C80L1q6FDh3g5En4+edLXw8Pd9XIYozpks7DPl3y5/HHHyc4OJj169fTuHFjChUqdElgrlWrFjNnzqRChQqMHTuW4cOHEx8fzy233MKXX36Jw+G4yhmUUu6UpcIuEyfC889LIK1bVwKwwyFDmxndlOfLJ8/NyMKFMGQIvPsutGkDTz0Fd92V8fFU5n79VW5mSpSQ4efMfvceyKNe9R49ZATBTnXrwqefZv/n9u/fz/Lly3E4HPTv3z/d52zbto0JEybw559/EhAQwAsvvMC4ceN49NFHc9NkpZSTlSkSQlQ6QbhMkRDptb70EowfL6kKJ06EJk0uPik3wXLwYHjySRg1Cn74AaZPh/BwWLQo58fMy/76C+rVg8mTZSrAy3hUAPYk999//1V7sgsWLGDt2rUX8kfHxsZyjQvziCqlcqZnq6qXzAFDmsIuO3fClCkwcKAMM9vdO73+egnEAwfCzJmyoAsgPl56IDffbO/5fM25czLPW6sWvPceJCR47SI3jwrAOempOkv+/PkvfO7v709ycvKFr1NLDhpjeOyxxxg0aJDL26eUyrnUed7UVdBVQgz/C9pLw3qtgVBZ2ezsHlVAgCwcSvX119C9u+xfHTLEK3t0Trdvn8z3HjkCu3bJkL+XBl/wsmpI7lKhQgXWrVsHwLp169ibsvqxefPmTJo0iSNHjgBw4sQJ9u3b57Z2KqWyrmO9UP7s3Yy9rfMz74eXadj35UtXNrtat27S4/75Z7jhBvjsM0hMvPrP5RXLlsnowN698P33Eny9nAbgLLj33ns5ceIENWvWZMSIEdxwww0A1KhRg/fff5+WLVty44030qJFCw4ePOjm1iqlsiQ+Hl59Fe64Q3qjS5dCxYrua0+BAjBoEGzaJKume/SQ3nBeZwx8+aW8ToULS5KN1peXJ/BOWo7QRfLKdSrlFYyRfbmTJ8uCq8GDIc20k9sZI3uNy5SBW26B06dlrtNDMzo5lTEyVJ+YKAUVihRxd4uyJVflCJVSyudYlrypN28ulXE8TWr7UvXqJQvDvvxSchvnBf/9J8G3fHkZlg8OBj/fGrT1ratRSqnMxMTIXCLIvl5PDL7pee456Q3fey888IAsQvJlixdD/foXh+Dz5fO54AsagJVSecXZs5L8omVLOHTI3a3Jnjp1ZO7zgw8k/3SNGpJz2tcYI9thmjeX4faRI2079NT1UTQevJCKvWfRePBCpq6Psu3YOaUBWCnl+06floU7S5bAd995RZ7gKwQEwFtvSarLunWhUiV3t8heMTFS1OLVVyXF58qVUK2aLYdOzf0dFR2L4WLub3cHYQ3ASinfFh0NrVrJG/ovv0DXru5uUe7UqCF1bitWlB5j164werQUoPdmycmS2WrgQMm3XaiQbYfOLPe3O2kAVkr5tlGjJGH/xIlw//3ubo29Tp+G/fvhiSegcWNYvdrdLcqemBgYMECqTBUoIO3v29f2+d4s5f52Aw3ASN7nDh06cP3111O5cmVeeeUV4uPjiYiIoG3btlc8f+bMmdSrV486depQo0YNvv76aze0WimVJa+9Jm/sHTu6uyX2K1wYIiIkr/Q//0iiim7dJJ+1p1u4UKpH9e8Pv/8ujzkpq1WZIiHZetxV8nwANsbQqVMnOnbsyN9//83OnTs5e/Ysb7/9drrPT0hI4JlnnmHGjBls3LiR9evXEx4e7tpGK6UyFxMDXbrAjh2ypadOHXe3yHn8/GS18I4d8OabEtgCAtzdqoxFR0sVqObNpe0REfDgg049Zc9WVQkJuDS3/4Xc327knQE4MlIyxkRG5vpQCxcuJDg4mG7dugHgcDj45JNP+P7774mJibni+WfOnCExMZHiKRvig4KCqFrVvS+iUiqNpCR46CGYMAG2b3d3a1ynUCHJIb1jBxQsKIk7OnS42Lv0FM8/Lz32Xr1kzrdpU6efsmO9UAZ1qk1okRAsILRICIM61b609rMbeF4ijvR6kw88IPv1YmJknuOvv2TC3s8PbrwRXnlFakIeOybZbdKKiMj0dFu2bKF+/fqXPFaoUCHKly/Prl27rnh+sWLFaN++Pddddx3Nmzenbdu2dOnSBT8f3KOmlNcxRlbRTp0quZQ7dHB3i1wvOFj+3b9fbkDatIG775aA16SJjAi42urVcoNQtapspXrjDdnn60Id64W6PeBezvuixqlTF1f7JSfL1y723XffsWDBAm6++WaGDRvGE0884fI2KKXS8emn8PnnEoS7d3d3a9yrYkXJKz1sGCxfLj3N6tXhwAHXnD82VlZnN2woc9MffiiPV6rk8uDrsYwxLvuoX7++udzWrVuveCxTy5cbExJijMMh/y5fnr2fv8y8efNMkyZNLnns1KlTplixYmbWrFmmTZs2mf780aNHTYECBa56nmxfp1IqexISjLntNmPuvdeYpCR3t8aznDtnzA8/GNOlizHJyfLYTz8Zs3jxxa/t9N57xhQtagwYU6OGMSNGGHPqlP3n8QLAGpNBTPS+HnBYmGSAGThQ/g0Ly9XhmjdvTkxMDGPGjAEgKSmJ119/nccff5x86ZS7Onv2LBFphrU3bNjAddddl6s2KKVs4O8Pc+fCTz/5ZNrCXMmXTxZq/fyzDEEnJ8O7717sFX/8MWzeLL3W7DJGpv9mz744OnnuHNx5JyxaJMd98UVb9/X6Cq2GBPz333+88MILbN++neTkZO6++26GDRtGZGQkd91114UFVwDjx49n0KBB7N69m5CQEPLnz89nn31GgwbpFru4wBOuUymftGsXvP02fPONbMtRWRMTI3ujv/764oLWd9+VfbnR0ZJ164Yb4Prr5d/gYKlEVLCg7KseMgR275aP1KnA+fNldbMx7plr9kBaDekqypUrx4wZM654PDw8nNh07gibNGniimYppa7m2DG46y44eRKOHtUAnB2pveLHHpPFWhs2SJYtkAVc48dLIE5rwgRZFBsXJ8+vXFlGIStXlp9NXdHshcF36voohs7ZwYHoWMoUCaFnq6pOX7SlAVgp5Z1iY6F9eylbt3AhVKni7hZ5r2rVLs27XKuWJPM4fhx27pSP8+dlQRXIbpSdO93TVidIzRWdmq4yNVc04NQgrAFYKeV9kpPh0UdhxQoZRr31Vne3yPdYFpQoIR8+/vvNLFe0BmCllErrwAEJvkOHSo1cN3DHkKVyDnflivaIAGyMwfLCOYOscuVCN6XyhLJlZY+rm+Z83TVkqZyjTJEQoqJjuSlqG43+3cSK8rVZF1rd6bmi3b5WPzg4mOPHj/tskDLGcPz4cYJTs9MopXJu4UIprpCUJCty3XTj7qnl7VTO9GxVldujNjNxXC9eWzqWcb+8TdjhnU7PFe32HnDZsmXZv38/R48edXdTnCY4OJiyZcu6uxlKebft22W4OTRUKui4cV+pp5a3U9m0ezesXEnHrl25vkg0lknGAZCUSN8CR6jp66ugAwICqFixoruboZTyZEePSj7jwECYOdPtSR1ShyzTe1x5sMhISeBkWVInYP58CAmBdu2o2bU9jPoM4uPxDwyUr53M7QFYKaUyFRcntXwPHpQ3zQoV3NwgGbJMOwcMnlHeTmUiMlKK/cTHy9elSsF778ETT0hykdQsixER8rxcZlnMCg3ASinPtn69JH346Se45RZ3twa4uNBKV0F7iYQECaxJKTdMfn7w8suSQS2tsDCXBN5UGoCVUp4tLAz27JEeiwfxxPJ26jKnT8Obb8Lff0v9gMBA6QEHBkKzZheelpwsHeQJE+TzESNc0zwNwEopzzR2rGRfevJJjwu+qU6ckJHxa66BYsXA4XB3i9QFv/8Ozz4re8Z79IAGDS4ZYjaNwli9SoLuxImSUC0oSErKuyqVtQZgpZTnmTtX5uZuuw0ef9xjIltSktSW/+MPmDMH/FYspymLWEQzVvmFUbKkBONrrpF7hoaJkTQ4G0HVZ8Mp2d51Q5t5WnS0DC+PHQs1a8KkSRemLkyjMDYEhzFhAvz6EOzdCwEB0Lo1DBokmU0LFnRdUzUAK6U8y/LlcM89ktz/t9/cHnwPHpRg+8cfMG8enD1xnjuI4H9FvyacqWBZJPp/wJw7hnDrn0M5F1OAM3sKkJSQRI34DRgg/vcg+jRdQLO3w2jeXKslOkVkpPRuGzaULGn9+klFp8BAQBY8v/aa5G9xOKRa4jvvyPq+okXd02QNwEopz7Fxo2w3Cg2VqOeud8aUpox8NJKif0UQQTinS1bh90IvcVPsbAJiz8DZAMCAMQQkx9P2hp1QpjnFzp6Fs2c5t3EzfgeTsQBDPHWXj2Bvy9E8EvoQ9bo34fEn/ChRwm2X51siIqBlS5nADQyU2sQplZn274fXX4dff4VKlaT6YqdOeMTvXgOwUspzLF4se3znz3fbvG9cnKzXWTp4GfOTm+FPIgQFw6Q5+D2/GR7uLGOVBQrIzULqop6uXS+soJ26PooJw3/l+7G9CUhKJNHhIL72GR7bNJVnor7l317l+OGtLpxo/RCd252j9vEIrDvCr1iBm1fyTefqOvftg4cflpXOIK/H8uXEhzXl009lp1FSkvzbs6eUNfYUlitTQDZo0MCsWbPGZedTSnmJtKteTp1yW47nJUvg6afB2rmdZflbUeLcv/INh0Oicp8+l/5A6rDnZftGGw9emG5u4cr5YEGVU5z5ahz5ls3hH66jtDlIEPFYwYH4LVxwSRBPb6/xoE61bQ3C7g7yubrOefOgSxcpTZmYKJE2MJDVQxbw6JdhbN8u90qffgruyvdkWdZaY0yD9L6nMxFKKfc6ckTK3UVGytduCL6nTsFzz0Hzpgk8dXQQWwLqUsLvpPRsHQ75Nzz8yh8MC5OgfFnPNTUl5brQ6nwZ9gDrQqsDsCcG6NqVgktm4Th8kHLP3E2wFY+DJIiLZeNrP5KYKMdwRb7p1OAXFR2L4WJRianro2w7x9Xk+DpHjZLVU6VLyz7xiAhOvTGQvmELuLl7GPHxMGMGTJvmvuB7NToErZRyn+hoaNUKduy4mCQhC+zstU2dCi++CIcOwaedV/PyL2/Bffcx+9m3mfLbUqpsXcOuGg24O7g8HbN4zCylqixZksDHusKYUZjz8VjJSdRZ8TVzrommzJghtuWbzux35a46uGnl+DobN5YV8sOHY/Ll58c/r+elEWEkJcGAAbL915OGm9OjAVgp5R7nzkHbtrBlC0yfLluOssCuUoCnT8sW4+mTzvNExQieWNGKhg1vhTfXMZVr5ByFKzE3rBIAS7NxjiynqkxJf2hFRGBuvoWtXy3m9klDsdpNpUfZp/j03rswlwWR7OSbvtrvyhOKSmQ5r3ZqpoyTJ+GHH6BaNRg1irNn4YXHJFFaeLh0jCtVcknTc02HoJVSrhcfL5WNIiPh559lKDGL7BiaPXZMEiEVm/wtx/OX58t9d9OwxF75Zr16uT5Hx3qhDOpUm9AiIVhAaJGQjOc0U4axrebNqDFxAPGbdrDx+vu4fv8e9v/QjJqrDvFC5K/cFLUt2/mmr3YdGQVzVxaV6NmqKiEBl241u+I6U/M4f/YZjBkjY8vISvX69WHcOOn1zp/vPcEXtAeslHIHhwPy5YNvv5XUQ9mQ215bVJTsWHl0x1u8mTwI6xwyx3vo0IXJQjt6hjlNVVm4Vjlu2TmW5YsTaNN1BeMXvUQwccT7B7Li24k0zcYxr3YdnlBU4qp5tRMTZTw5tYiCw4HZvIWvotrz6quSgWzBgvSn6D2dBmCllOvs2SPBrmxZyVCUg4wUuSkFuGePJGB4MOpjeiUNuviNpCRZzZyymMoTyg3e2jSASc8vw++deCwgMDGeuiM/h4fukvRNWXC16/CUohKZ3qzcfz8sWwb+/mAMJjCQvvPD+d8iGTgZMwZKlnRpc22jQ9BKKddYuRIaNYJHH5Wvc5gOKktDlunYskWmmU+dgudedMgYdEhIuqucc3oOu/k3D8cvJAjj5yAZB4VXzedwhZtJXrchSz+flevoWC+UP3s3Y+/gNvzZu5nn7TN++mmZ4F2yhP3PDuS+IgsYsiSMIUNg1izvDb6g+4CVUq7w22+SLKFMGUmSXzV3gSy7q6DXrIG7WiVTzfE3IxdWpVYtZO/xihUZ1n919/7YC1L2Gp++KZyvBxzikcjnWV22Ew1Xf8m11179xz3mOrLKGBg5Uoace/S48NCnn0KvXrLr6JdfXFo1MFcy2wesAVgp5TzGwMcfSwqiW26R1c4u7rIsWQId2iTyrXmKTvyG35bNcN11Lm2DXYyBnz47wRt9AqBgQaa8u57GIetkL7WLisg7VUwMPP+8jCt36ABTpnD8hEXre2JZszSEfNcfolaXnfTpWNmzbyLSyCwA6xywUh7C63oqWREXJ1tGOnWSYcQQ182hgnS2u3Q6z6TALrQ4O0XyEZYv79I22Mmy4NEexWjQEjp3BsfLz2NYCZaFFRwsq5G8NQhPmiRVjA4dkiXNffvy53KLjvclcvxoEEXv3ELBm/7h8HlytO3ME+kcsFIewBMyEtnq7FlJDxgSIkO8v/7q8uC7ZEgkq9sO4E9HE1qcmSJbWN55xzWFXp2sRg2ZUj/RoBUGC8sYTGysDN16o3nzZLHVoUMQGEhy8xYM/tCPpk3hTHwCpR7+k0L1/7nw0tmdEcxdNAAr5QFckXbQZWbMgBtvlCxFAMWLu7z+3l/fRNKgd3PeMQOoGbMa3n4bund3aRucLSQE7h7emuTAYBLxw2DBTz+RNPE3dzct644dk3/XrLlwY2SSkvixWwR9+sjASalHlxB07ekrfjTtFqup66NoPHghFXvPovHghV5z46oBWClniomBtWth61b5+tgxGSK87z549VX46CP49Vf89+5J98ddmZEo1/bsgXbtJPt9SIjkd3SDHTtgWo8IAomXsOTnB/nzu6UtThcWhn/EAmLfep+3Gy/mWb6ixfD27NmD/J9LrRDkaZKS4JNPZC5+3jyZvw4OJtnPQVxSID/8E87IkZL4qmyp9LdcpW6l8ubRI50DVspu/ftLcvjNmyUoGSO9wdGjpTdYtKjsifnjD0nHCHRt+QSDinbimjPHmTy2J6vK1eTP6+qy+8Zb3HklWTd1qlSl8feHYcOkt5nFvap2OnwY2rRO4n9JW/ELCoBEMi6k4CvCwigYFsb/DIwZ04RfukOj2uf4x785IWWLYY0cKa9FBqu9XW7rVskBumKFpCKtUYO44qH81HkBe3+IYHe5cIbPCKNOHXn61ZKFeEI+65zSVdBK2SEqSorIA9SuLYXBa9aEWrXko169K0uyGCPFCPbv548D53l12TEKnTjMOwtH0ejfvygRc0qeV7WqZIxq0sSll5Ql585J7/LAAXjrLfjgg4u/Bxc7exbCmxqe+utlnkv8QvasFC7sGUHHhf79F554AkIWzOC7kJcpFbtP9jobA0FB7luoFRkJ778Pc+fK6zJ8OHTpwrz5Fi++CH//DY88Al98AQULXvqjmS1QrNh7FulFMQvYO7iN0y/ranQVtFLOsnWrBJ558+QdpEwZSVCblTlPy5LecNGitK4NcddEMXROIC936EVooSDeq2xoFrVJ3jBTN3yOHQuffw4tWkgVoUaN3NLTZM8eeO01OHECFi+W6/7hB9e3I0ViIjzwALRYP4TnzBfwxhsweLDb2uNO5ctLjPvyy3bU7NmMKf5tuC1xMRbI3to0Gb9c4vx5Cb533y2r4v384McfOVCvDa91kWHm66+XNrdokf4hMsuU5QlZy3LMGOOyj/r16xulfMK//xrTrZsxfn7GFCpkzPvvG3P2rPPPO3myMbfeaozDYQwYU7CgMR06GBMb6/xz79tnzIcfGtOokZw7Xz5jBg82JiHB+efORHKyMU89Zcwj/Cjt6trVmKQkt7bJU+zYYcxTNZebGIJNEn4mMSjEmOXLjfn4Y2N++cW5r92+fcb06WNMyZLGtG9/4f9sssNhlrX5nylY0JigIGMGDMjdf98p6/aban1nm+t6zbzwUa3vbDNl3X77riUXgDUmg5ioQ9BKZdfRo1ChgnS7XnpJCrKXKOHaNkRHw6JFMGcO/POPzCcDvPKK9HJuuw2qV5fh65wuQDJGevPlysnc9fffy9xd/fqyPPXRRyWn82VcvZ954EAY8G4i/5ZqSJlaxWXzb2Cg087nSbLyu05MhOl9ItnxTQTTT4eTL/wWpu2/iQK7Ul7b7t1l1fratTkfrk/J1kV4uCw8HDFCkq6AzPO2aQM9emDi4zmfHMgdZgGFW4UxYgRUqZLb34Jn76HXTFhK2eHo0YtZnL77TsbLPC2j0qOPwpQpMiGaqksXKfkHMt5XpowEzpgYOHNGrqlyZRke/OYbeezQIUm0u3evvJm++KIU0D15MtNrvrz+LMiCmQxL8eXSDz9At25y2T98chLL4Sfzi3lAdn/XcXHw9dcwaBAcOZzMO/Vm8br1MYXWRcgT/PxkjnjGDHmNK1W6OJWSNsCmDdDGwOzZckOWmCg3PpUrw8GDksP5uefguuv46y+Y1TeSMzMi2FwinEdHhnHvvT6xJfuqNAArlVvTpsFDD0m2nmzUrnWL+HjYuRO2b4dt26SX8/jjsiUlf/4rt6b06CFbQmJiLvaWQ0LkzbZTJ0kJmMX0kY0HL0x3Pi60SAh/9m6Wq8u63Ny58OLde/nk2iG03PopgYWCr/5DPiSnv+tz5+DLL2HIEDh+HGZUeJk2/3yBhZHFWl27StayAgVkQWGpUnIzlpwsAfbdd2Xef98++YiJuXhwh0O21w0cyNEzwfz8s9wkbdggSxVeekmSXF2+yMqX6SIspXJj9Gh46ilo2FA+PF1g4MXV12n5+8tCsW3bpIeSPz8UKiQ9FpCge+yYvDvmcAjXjjq6WbF5M3zVYTareJQip85jHX4dCl1v6zlcITdDpzn9XefPL6m5n3tOFiIPH9yVZowikHiSrUA2lr+Xmp/dTr5df8HGjST9PhtHyk1bYtx5Ds9fSujJIzLF0bq19IK//BKSkjCBgSwu1olPHgzm99+lU1y//oUFzy6fqfF0GoCVyszQoVIMvGVLqehToIC7W5RzliVDixkNIVuWzPXmgitWpB45AgPvjGBiXFv8SMZKCJIbh+u9KwBfPoScmkACspbjOLe/64IFJUHYiy+GMbH3AuL+iGBcVDhLPwjDsqBOHShb/SzWrRP5delz+CcnkuDw59WSrbnj5fupX/paTp6UhfCBjgcI+DOCL7eGM+etMK69VgZWHnvsyvtAdZEOQSuVkQULUqq3PyjVWfLIwp7ccPYc8Pnz0LyZ4cvIutxo/pIHHQ5ZidWnT66P70q5Ha53xu86NlbyYyxZIqPMi5cmkZzooBGR3OFYwKKk5qwg/UVaQUEyW/HYY3K/6p+me+fJi6ScTYeglcqJZs1g4kS45x55k1dXlfqm6ow3W2NkXc+B5XupFrIXEgMuzkt6Yaar3A7XO+N3HRICd9whHwAV3phL3KHCbP+3OFtiHsIvKIGiIVtwBCcw+rm6FC0KxYpx4d/07lFz29P3ZVcNwJZlfQ+0BY4YY2qlPFYMmABUAP4BHjDGnHReM5Vykbg4eOEFeP11yWR1333ubpHXySxpQm4MHixrgwYMqETgMzslGcjixV6b6cqO4Xpn/a5ThZYIIsr/JMFlL317Dy0SQtu2WTuGN6eKdLasFGP4Abh82WdvYIEx5npgQcrXSnm3M2ckW8/o0TIOpzzG5Mnw+1tL+eXG//FOXyOZwW69VYadvTD4guQ4Dgm4dGQlbY5jT2BHG121MM8bXTUAG2OWACcue7gD8GPK5z8CHe1tllIulpwsyzSXLJF0j08+6e4WqRTr1kG/h3Yxw/8e7o8bg3Xu7NV/yAt0rBfKoE61CS0SgoX0Kp21Xzqn7GhjRj16j0sVeeKELLR0oSwtwrIsqwIwM80QdLQxpkjK5xZwMvXrdH72GeAZgPLly9fft2+fLQ1XylbvvisLeb74QoaglUc4cADurH+SGccaUbHQcfxWrrAndZJyGVcnZ8mWo0dlAtvfX/KHf/SR7KG3cUV9Zouwcl0POCXXZYZR3BjzjTGmgTGmQcksbuZXyqWSkmDVKikh8/zz7m6NShETA53axjPy6L1UtP7Bb+qULAdfby3Q7os8rqd/8KDsW27eXKYyliyRx198EVavdukNXk5XQR+2LKu0MeagZVmlgSN2Nkopl3I4JNNPYmLeyI3nBZKT4YO2kTy5/nua+C3B7/sfslyOUVfdeh5nLxbLkgMHZI/UggWypL5qValkVqmSfL9ixStLhjpZTnvA04HHUj5/DJhmT3OUcqHoaEkifPCgBOGgIHe3SCHvjcO7RPL2ouY8aY3GLzDgYrauLMhs1a3Kg5KT5d/ixWWj87vvSiq1bdtk2qlCBbc17aoB2LKs8UAkUNWyrP2WZT0JDAZaWJb1N3BnytdKeY/kZHj4YRg/HnbvdndrVBpDh8K1v35KMOfxM0mSuzoiIss/r6tuFSAjWp9/DnXrSgLsoCBYuhT695cthh4w2nXVIWhjTJcMvtXc5rYo5Tr9+8uw84gRUrpPZZkzsxqNGQOLes1mBr9h+QGWI9uJNlxWoP38+YujJmvWwPr1UnEgIEDaHBAgZfgCAuw9r7q6BQukNOeWLVK1LDpakmB7QNBNSzNhqbxnyhQZeurWTVc8Z5Mz51f/+AO+fSKSeX734nfjjVhDP5RFMdlMtNGzVdV0V93asr/2yBG5cZs+Xcoxbd0qubVnzID33rvy+adPSwAeNkxKAnXsKAUMvDmnuCc7d07meX/7TeZzp06F9u09LvCm0lzQKm9JTpaKRv7+kkUp2L4Sdnkh362zyg2uWgUvNN3CgoQmFLiuOI7ly6QMXg7Z/lr89ZeskI+MlEnqcuXkjb1nTwnAZ8/CqVNSCjIhQT7i46WigZ+fBOfhw6X+X3Cw9Mq6dJEPZR9joF07SdLy2mu2/n3nlNYDViqt6GhZjFG6tG2H9Oi9jmnkNjBV7D0r3T2HFrB3cJsctWnnTmjcGN5KfI/ugSNxrFju8tWo6frzT3lDv+02OHxYhpPbtZOKA3XqZL9XlZgIy5bJCMzUqXDTTfI5wNq18rWH9tQ8WnS0pI7t319ujIzxqN+jU/cBK+UVjJEUk3FxUKSIrcEXvGPlbepNQlR0LIaLw8fZ2SNrd1ajgwehVSt5v2y36h0cf21wf/A9eRKeeUYCb79+8lipUjLP26+fLOrJyRu8v78Mp3/2Gfzzj1SqB9i0SUZlmjSRoO9l3LrnevVquXEZM+bi786Dgu/VaABWecPYsZJoY8wYpxzeG1be2nGTYGf+4lOnoFPLs3z2XycWfr6FKtdbuRp2zjVjZFV8tWrw/feSGWn6dOecy7KgcGH5vHp1+OorKS5x220yT7x1q3POazM7bupyxBj45BMZOklKkmQanTs795xOoAFY+b79++Hll+XNzUk5nr0h360dNwl2ZTWKi4P72sfz3pZOtDXTqZV/b7Z+3ilmzICuXWVf6Jo1sh8qf37nn9ffX3rcf/8NH3wAixZBo0ZSHMTDuW3k55NPZI737rtl9bmXFuTQVdDKtxkjQTcxUYb8nFTX16krb21i1/ac3GY1OnUK3g7/k2EbnqcOm+D70WS5tp3d4uOlt1m3rrRh/Hi4/3731H/On18yMz3zjAytFiwo/3+/+kpuDFJ7zB7E5SM/iYlyw/LUU1CokPxte9GQ8+W0B6x829dfy3aRYcOylU0puzwu3206PKH83aFD8Fr9xXy6oakE34AASQnoDjt3yvxhs2ZyV+DnJ8OY7gi+aZUoAXfdJZ+vWydb5apXd95weC64bOTHGBmRaNRIFlAWKiRB2IuDL2gAVr6uaVMZqnr2WaefqmO9UP7s3Yy9g9vwZ+9mHhV8wf03Cbt2yZTddfuW4CBlpCA5OVtZrmyzerU05vBhWRfggb1LAOrXl7aWLCmrr7t2hWPH3N2qC1xyU3fqFHTqBG++KQv0EhPtO7ab6TYk5Zs8bCtCXrduHTzW8iCxyUFM+3AHNbs3l+HfwEDJWuTKObw5c+Dee+Gaa+RzG0vPOU18PAwaBO+/L2kU16/3mP/fTt3//tdf8lr984/0gF95xWOuO6t0H7DKez76SHoOP/zgEZvx87IFC6Bn+x1MO9+KYmHVyL/0D0loERGR7SxXtnjiCQlgs2dLOTpv8tdfkszjjjsk2cfJk3Ij4YuMkYQa+/bBr796bcpYDcAqb9m6Veb27roLJk/O0h1zXshi5Q6//gqfP7SC6aYthYo6cPzxuwyrusOpUzLUnJBwcR7Rm/3vf3KjOXy4DE17Wc8wQ3FxMsxcoADs3Qv58rl3e1ouaSIOlXckJEgu2IIFZfVoFoOvW/Yy+rgvvoCfHpzJvORmFCpfRDJcuSP4JidLpqQGDeDECVn45e3BF+Cee+CGG6SqV4cOst3O2/3zjyQkeeop+bpiRa8OvlejAVj5lsGDZQ/nyJFZ/sP1hixW3uTMGfjovkgOvvQ+3wc/R2DdmhJ8nbgKPUPx8VLz+eOPZUSkSBHXt8FZqleX1JbDhsH8+TI3/Ntv7m5Vzv2eMjry9995Jke27gNWvuPcOfjySxmOu+++LP+YM/cyJifDjh0yHb16NcQsiKTqwQi2XhPOoYphlCghC1xLlODC59deKyPo3jZ1bYykOP756YX8eLwtQVY8fgRgffije+Ypz56VBTxz58pwbe/evjNMm8rhkN59x45SLKJ8eXe3KPtOnZKdCt9/DzfeKDcRVaq4u1UuoQFY+Y78+WVxTWBgtn7MzvqxR45IVrzVq6XCz9q1FxMadQj6g4nxHXCYBJJOBTAidjARSU1YdPoGDpyVIdFGRBJOBH0CwnHcFkazZrLepmHDrF2Wu+ay9+2TZGPMmM5Yv0cJIQ7LGEhAfhHN3VA+/LXXZAXY999L6UlfVrmy3Gik6tFDRoDeeMPz6xHHxsLMmdCnD7z7rvfdeeaCLsJSvmH1ahm+8sv+rIodlYzOnIEhQ2DZ0EhujV/Iv45KVCqfRNOC6zjf5h6ue6gJ1T99Fr/vvrnyh2fM4HyLtsT+7xMKv/caAMmWP1+Vepf3Dz3NIa4lf35ZBJoakG+66cp8Ee6oyJSYKLUFvun7L8MSutMuaRqmYiWsA1HyTXdsM0p1/DisXCnpCvOSpCRJKDJpklRt+u47mf/2JCdOyBqN3r3lb/bMGVm34YN0FbTybVu2QL160Lev3EHnQE57jklJstOpb1+od2gWM632WCaZCwOdwcEy//j887LtpnVrCUwBARK5SpeGW26RIdqXXpIh9Mv+Jud/tJFpe29k25x/2fG3RVn20yZfBMlNw6neLYyWLWVxr7Nq9WZk5UrJb1Jh41TGOx4mKNDgN6C/9L7WrMnVNqMc9+RTq1499BAEBWX7vD5lyhR48UVJNvLaa/K34QlBbto0eO45OHpU5rAbNXJ3i5xKA7DyXUlJktFo927ZflSypMtOvXChvK9t3xjHTWHBTLjxA8p93Ve+6ecnAfWjjyR3barM9r9GRspQbWqCis8/l/mx7t3lGC+9BF98gcHCAAkE0IaZLPZvwe23w3q2ElL5CAHFzl1y2NzU6r1cXJxc95oRkZyfvZBNJZrx4nulaLmwF9awYVKcPpdy1ZPv3x8GDIBvvoGnn851W7xedLRkkBo7VqotXXut+5LUHDsm/5fHj5ee+ejRcuPs4zQAK9/12WfS4xo3ThZfucCOHdCzJyyfcYwhBT+gc+Bv5Nu7FWvzpksDaE6GXjML0Nu2yYlnzbrwUGL+wrzz0klmzrI4uXk/ByhD44KLaF5gDpHla7OpegWuq5JI5Nt35Ph6jx2TU06fDov/iOXxmBEMpo+kkwwJwcrmdV6td5vjnvzIkZI3uVs3GDXK9xZc5cbBgzLaYowMydetK3ePLrxhpWVLqfT0zjsy9JzNtRreSgOw8k1790KtWhKsZs50+htuXBx89VgkpybOpZzjAA85fiEw4SzWE0/I9qfixZ2f4SltL9nfX97M3n4bjOHsNWWJjz5LocRzWCSTQCAtmcOKwNupV9eifn0ufFSrJodLTJRBhMTESz8/dUrW9EybJnXOmyfP5dXgr2iWOIegxBgM0rPG4YCBA2UBTRZkpXdbsfcs0ntXyrQnP3EiPPigVDSaPPnSUQd1UUyMZAL79VcICZGh4DfekOBsJ2NkePmrr2DECChaVG5IS5aUlc55iAZg5ZvWrpVhxmnToFw5p57q7Fnoc0ckH65pRjBxEnxuv13eYKpXd+q5r5BekE9MhJ9/5tT7gyn097YLc9B7y9/MF/evZOPqeBqu/JyV5+sSQDz12EAE4axAft4imXzEEEg8QZwnnEU8zDhmVOxOqUda8dTpjyk78WOsDh0kd/Jbb+Wop5+V3m22e8Dnzskq4CpV5K4hX74stSVP275dtmb9/LPcrMyeLav7cismRo45YgRs3Cj7rqdOlaIoeZQGYOW7XDCfdeKEjNo1XzWIgdY7+CUnyRzv++9nuefnMml7yA6H3CB06wabN0Pt2heelvpX/2eL/qy+qx/l9y3l3s9uv/J4gYES7G+6ST5P/V3nsKefld5tjuaAt26VXlzRollui0LWTnz2mYzg5MsnQXnWLFkYFRYm/5Ytm/7PJiTI0HZ8vNz8nDghN2cnTkgv9+WXZVooj98QZRaAdZxGeZ9Dh2SB0ttvO/2P++BB6ND8LK/teIbKPdri91XgxZ5feLhTz50jYWHSI708ONaqJatOe/WC0aNlj65lcdv1R7jtVSCqEpQdKte1cCHMmCFZRJKS5FiXB9mwsBwNsWdlz3VqkL3qKuhdu2Tq4ZVXoEaNbLdFISMHw4df/LpkSbnJ+uILWb0PMl+xdas83r07rFghaS8PHZIb4GbN5P9csWKyHiM8XPbM6Rz8VWkPWHmf+++XALFxo1OLue/dCw82O8qX/7ahPmuxxoyBSpXcV8XHDpevtE5v+Dgrz8kh2/YqHzsmvbNTp6R378P5gt0iPl7+vlasgNOn5WYXJEfz/v0QGio947Jl5W/w9nRGTxSgQ9DKl0ydKknoP/hA5iGdZOtWePKOPYw91ooKAVE4Jk6Adu2cdj6XysrwsRMXk+U6W9f583DnnZJ8JSLC5/eRKu+mAVj5huhoGWq85hp583VSir3Vq+GlljuZeboJRQsm4j97pnf2dn2RMVLt6qef4JdfZOWzUh5MyxEq39Cnj2T1GTXKacE3IkKmtE4VLk9wu5b4r/hTg68nWbNG9ny/954GX+X1dBGW8h7du0vmHCfVlF39eSQHenxB57KP0P/PVhQM/ckp51G50LChDFHkgQxKyvfpELTyfImJTk+scHByJMXvbUogCRh/f6wlS7Tn60lWrpTiCnmtsILyejoErbzbiy9Kge7kZKccPjYWlj89mgASAGSLTkSEU86lcmDfPmjfXra4xMe7uzVK2UYDsPJsc+dKYv1y5XJUavBqjIG3u+6l+YkJsm/R4fDcPb550enTkl7y/HnJeJZH8gervEHngJXnOnVK9h1WqyaLbpzgiy+g5dTnCQr2w/rxF8kM5K17fH1NYqLUtd22Df74w/UpP5VyMg3AynO9/jpERcHy5VJX12ZLl8Krr0LXFqNp+c4eaNLY9nOoXJg4UXIUf/217PtVysdoAFae6ehRSbrx5ptSsN5mBw7AyPazqVyhJcMnlsavsM3VYFTude4s9WvtKBKglAfSAKw8U8mSsGWLVFOxWXw8fN30Z36OfoiDL3xO4cIv2X4OZ8l1FilvMGWKTDtUr67BV/k0XYSlPM8ff0gRgFKlICjI9sN/0mUVfXY9wdHqt1O63zO2H99ZUvMoR0XHYoCo6Fj6TN7E1PVR7m6afRYulAQbTkwzqpSn0ACsPMvMmXDXXfDtt045/ISPo3hkckdiCpem5JLfvGpV7dA5Oy4pYgAQm5DE0Dk73NQim23cCB07wg03wOjR7m6NUk6nQ9DKc5w4Ac88I3Vru3Wz/fDbRi3n5tcfophfNP4Rq6BECdvP4UwH0injl9njXmXfPrnxKlxYRkCcMPWglKfRAKw8xyuvyOKrmTNtH3qOXRhJhafvJJDz+Pn7Y8WesfX4rpCVWrpea+BAyYiybFnGBeCV8jE6BK08w5QpMHas1B296SbbD7+k10wCTDwOkrFSi8x7mZ6tqhIS4LjksZAABz1bOa8mssuMGAGLF0PNmu5uiVIuoz1g5RlKl4YOHZyy+GbF5AM0WvO5ZNKy8NpMV6mrnZ29CtplK62TkqTn+8orULQo3Hij/edQyoNpMQblXsZICkgniTlnWFPyLhrGLcEaPZrgA3s001UmUldap13sFRLgYFCn2vYGYWMkx/fIkfDjj/Doo/YdWykPklkxBu0BK/d69VXJvzxsmFMC8ay2I7k/dg47e3zJDY9p/diryWyltW0B2Bh44w0Jvm++qcFX5Vk6B6zc57ff4LPPZCjSCcF3zbgdtIl4g63lW3PDx8/Zfnxf5PSV1omJ8OST8PHH8NJLMGiQPcdVygtpAFbusWsXPPEE3HwzfPih7YePiYFve+7kuH8prlvwvVOHuX1JRiuqbVtpfeKEJNvo1w+GD3dKhSulvIUOQSvXi4uD+++Xoedff3VKMoy334ZvDraj67zWlKsSYMsx80IayJ6tqqY7B5zrldZnz0JICFxzjSTcKFw4ly1VyvtpAFaut3Yt7NwJEybAddfZcsi0wbHRzmMUnHKeF194mqZ32hd80wam1DSQgE8FYaestD52DO6+Gxo0gC+/1OCrVAoNwMr1GjeGvXulN2SDtMExKCaeAdPeJ9iKY2XHcOAGW87hksVJHqJjvVD7rmn/fmjZEvbsgXfeseeYSvkInYBRrrNtG/z0k3xuU/CFS4Njz19+o0rybnqE9+HTtfttO4dPp4F0lp075WZr/36YMwfatXN3i5TyKNoDVq5x7pzM+x45Im/ENub6TQ2Cj839g6eOjuPXou1Yd/N1WDYGR59OA+kM8fHQurWkl4yIcEp2M6W8nfaAlfMZA889B1u3wrhxtifaL1MkhFv2bKLf+hEYoN2ZOdwUtc3W4GhbGsiEBDh0CM6cke1XviYmRq4rMBC++QaWLtXgq1QGNAAr5zJGki2MHQvvvQctWth+ip6tqlJr1TEMflhAQFIit0VtuSQ4Tl0fRePBC6nYexaNBy/Mdg3djvVCGdSpNqFFQrCA0CIh2csONWkSNG0qC5BKl4ZChcDfXwIWSErGGjWgYUP5HfXuDdOne1eQnjMHatWSBBsAd94JVX0gT7VSTqJD0Mq5Vq+WLFcvvih7g5yguqMEH+3vxOvWFwQST6J/APUfvYemKcHRrhXMWVqctHcvzJ0LkZGwfDnMmgXXXy/bcOLipNzi9dfL0Oy5c7I1B6QCUM2a8tiRI5Ko4rvvpDoUwKhRksTi1lvleZ60f/bIEXjtNRndqFYN6tZ1d4uU8gqaC1o538KFkn/ZCUHDnD3Hf9fcxHBe4a2J9Sj2V8QVuZ4bD16Y7vxtaJEQ/uzdzJ6GHDsmNxijRkmvtWRJCZbvvy+9wuyKjZVgXqOGfN24sQR0kN7zHXdA587y4U5TpsBTT8mQ+ltvQZ8+tpeSVMqbaS5o5XqTJslK59tvh2Y2Bbl0bL+vL9Vjd3L7G7Up1iYM2lxZZMElK5iNkdSaL7wAL78MVapckX0rW4k8QkIuBl+QOrm7d0vPeulSGe4tUUICsDES6Js2lRuPAHv2PmcqIUHOU6SItPPrry9tr1LqqrQHrOyXuuUkPFw+d1IayFNzVlCw9a1MKfU89xz4IsMOtlN6wPHxEnTmzIEZM+Qaz56FAgXSfbrtVYaMkV5yvnyyx7ZqVRmiLlxY5pDvvBPatLGvuH1SEqxaJdc6c6YE+q+/vtgWTfWpVLoy6wF70ESS8gnLl0OnTjJP+euvzntjPn+ec52fZD9lqTp5UKaj27YWsk9OlgxeNWpA9+6yiOrkSfleBsEXMk/kkSOWJcEXoFIlOH5ceuD33SevwXPPwbp18v3162WOdsIEGdbO7k13796ycOzWWyVvd/HicMstl7ZFKZVtOgSt7PPXX9LrCg2FP/6wfbtRWpu//pOq0TsZ02kaT95aKNPn2pZe8b//JMCtWiXF42fPhlatshSAnD4MXqiQ3Ph06iQBNm2msS1bZGXyJ5/I1yVKyLz0b79BsWJSj/fnny9eh2XJPO6UKfK5MdKjbtdO9vYWLWpPm5XK4zQAK/uMGgX588O8eVCqVJZ/LLtFDhISoMu3zShQZjfzx5TP0jlsSa9YpIhc348/wkMPSTGJLHJpIg/Lkl5xqocfhgcfhE2b5OZh1SrJUpXq/Hk4dUo+N0Y+HA44eBDKlIEhQ+xvo1Iqd3PAlmW9CjwFGGAT0M0YE5fR83UO2EclJckbdnKyvGmHZj3QZXtuNCmJn19azkNfNWHaNGjf3o4LuIoFC6BRIwm+OZzvtH0OWCnlFZwyB2xZVijQHWhgjKkFOAA374lQLjdjBtSpA4cPyzajbARfyP7c6In+w+n61e30uj3SNcF3xAgpJjBwoHydw/nOXCfyUEr5nNwOQfsDIZZlJQD5gAO5b5LyCsnJsvWlXz9JNZiQkKPDZGdu1OzeQ/7/vc0sRzte/KlRjs6XZcnJ0LOnJMTo0MGWSj62VhlSSnm9HPeAjTFRwDDgX+AgcMoYM9euhikPdvq0LPbp1w8eeUT2qOZwu0tGc6BXPL58OecatyAx2SLqrS8pV96JK29jY+GBByT4du8ui5Xy53fe+ZRSeVJuhqCLAh2AikAZIL9lWQ+n87xnLMtaY1nWmqOpafWUd+vVS/aCfvqpLEgKyflCoixtEYqMxISHU+DwHgKtBJ5o8V+Oz5clhw/LTcWnn8Jnn2VrsZVSSmVVbvYB3wnsNcYcNcYkAJOBWy9/kjHmG2NMA2NMg5IlS+bidMrtEhPl3w8+kPSSr7yS6z2gWZobjYjApMwT+1vJ+C+LyNU5M5S6ErhCBVkl/MorzjmPUkqRuzngf4FGlmXlA2KB5oAucfZFqfO98+bJiuBixSTFpE2uNjc6Lu5G7iGIAOJJcvizsnRNmtp29hQnT8pK5/vukxuMQpnvLVZKqdzKzRzwSmASsA7ZguQHfGNTu5QnMEaSTTRoIPO9FStKMHahlR98yZ+f7KZl8O98fOsjdH3wfZ7bHZjtcoKZio+XOe1//pFEE0op5QKaC1ql78ABSfS/dKkE3oEDoWtX16YdPHiQU+WrsSmxNvd1+JjgahfXENhWycgYePJJGD0afvpJklYopZRNtBqSyrrTp2X4tWRJCU5ffCHl5gIDXdsOYzj32PMEJsbzwnVDCKp66QI+21I4Dh4swffddzX4KqVcSgOwEjt3yjDzkiXw99+S6H/pUrc1x4z/hfzzptHTMYQTd8fhf1nH27YUjtWqSQ+4f397jqeUUlmkATgvO3kSpk6VLUXTpkkC/ldfdfk87xXOnCH+2ZdZzy1E93iSgiGriU2T5yPHlYzSiomRm4x77pEPpZRyMS1HmJckJ0si/u3b5etdu+CJJ2DFCkk4sWePrHbOpKyeK5xIKMhjjrF8XGs0Xw0pbn8Kx3/+gRtugPHj7WqyUkplm/aAfd2qVbBjB8yfLyuajx6VWrEjR0L9+lIrtk4dz6npGhNDr175mHS2NWvHSg4MW1M4njoFbdvCuXNQr549x1RKqRzQAOzt/v5bar/u2QO7d8u/oaEwfLh8/8EHpcdXrJhssbn7bqlhC1I8oW5dd7X8SseOcb5GXTjaj9d6Pk2dOjYfPyFBUkzu2AFz5sj8r1JKuYkGYE83cyZERsKhQxc/ChWCRYvk+08/DYsXy+eBgVIHNm3GsfHjoXBhGXL15JSKkZEkP/s8fkcP8V+ZRnzazwnneP11mDsXvvsOmtmwhUkppXJBA7AnOHYMNm68+LF3rwRVy4Jff5UgWqqUfFx7LVSpcvFnBw2Snl3FitLz9btsWr+Rk6sG2SEyEsLD8YuPJxF/+r9x1v7aB8ZA8eLwxhuy6lkppdxMA7C7vfWWBNFUZcrInOy5c7IY6vPPZZ9qRr3XsDDXtNOZpk3DxMdjAQ7L0CguArD5uixLtlkppZSH0FXQrrZunQwbr18vX3fuDB9+KHmWjxyBqCj4/feLK5ELF/bsoWMbJAXnAyARB35BgRAebt/BjYFu3WToWSmlPIj2gF0hNlaGkkeOhJUrpXxfkyayCvfGG+UjD/vA712W0ohhD66lzivh9vbqR4yAH36Ahg2hZUv7jquUUrmkAdjZkpKgZk2Z161WTerLPvooFCni7pa537Jl7FwVzXvvteXBri2pM87mALl5M/TsCW3awPPP23tspZTKJQ3AzvLff1C2rAwfv/OO1JgND/ec/bbuduwYyQ92xv9ofsqXbsUXXwTYe/y4OCkeUbgwjBqlv3ellMfROWC7GQPffCPbfsaOlce6dYM77tAgkCplXjbp0FHuS/iFUWMC7B8QGDsWNm2SBWylStl8cKWUyj3tAdvp9Gl45hmYMAFatNA5x4x89hnMnMlrDKfZ6/W44w4nnOPJJ6FqVZlrV0opD6QB2C5r117MOvW//0GvXlfuyVWwezfmzTeZE9Sexde/xOoPbD7+sWNw5ozsi9bgq5TyYBqA7fLffxAfLwk0Gjd2d2s8lqlYic9qfceQzW2Y+7NFUJCdBzdSu3j5cknJ6eaiEkoplRkNwLkRHS0Bt0MH6NhRciyH2FSn1tcsXw4zZjDTas+r6x9l2DCoXdvmc3z3nZRV/OgjDb5KKY9nGWNcdrIGDRqYNWvWuOx8ThUbKwF37Vr4919Jc6jStfiH6dz61L34JyVyniCeqjaHMVua2jtCv2ePRPRbb5VCCzr8r5TyAJZlrTXGNEjve/oulRNJSfDQQ7Bsmayy1eCboanrozj82Uj8kxIl1SSJ1C06lukbo+w7iTGy+M3hgO+/1+CrlPIKOgSdXcbAiy/ClCmymveBB9zdIrebuj6KoXN2cCA6ljJFQujZquqF+r2//Tibr7ZEYLBIxI8EP3+WV67KlDk77Kvxe/68LLq67z4oV86eYyqllJNpAM6uOXPg669llXP37u5ujdtNXR9Fn8mbiE1IAiAqOpY+kzcB0PFaP4aO6sVpR0EeShhHnWtWs6llcdaFVseKjrWvEcHB8O239h1PKaVcQMfqsqtVK5gx49IKRnnY0Dk7LgTfVLEJSQydswNKlmRe5Va0SpzD7NLNGP3w7awLrQ5AmSI2LFYzBnr3Bl9ZV6CUylM0AGfV7NmSW9iyoG1bzWqV4kA6PdmQ+Dji9x9g/yF/ehz4iW35qnFNpzX4BSTL9wMc9GxVNfcn/+03GDIEFi7M/bGUUsrFNABnxfLl0KkTvP66u1vicS7vyfonJTJi+hAm/9yLe9vEkRDnYNi3pylf1g8LCC0SwqBOtXM//3v8uMzF168Pr72Wu2MppZQb6Bzw1WzbJj3ecuUu5nZWF/RsVfXiHLAx/G/OCJrvXs2gKsNZszmYWbOgdetSvILN+Zhfew1OnJA6v/7631gp5X30nSsz+/fLnG9QkCy+KlnS3S3yOKk92aFzdtB12lc8sGk+42v14a3NL/PFF9C6tRNOumgRjBkDfftCnTpOOIFSSjmfBuDMDBgg2a6WLJFtLipdHeuF0nHK17BiIntqd6Drpg/o0QNeeMFJJ2zcGIYPl72/SinlpTQTVmYSEmD7difkTPQxkZHQvDkmNo5Ygnm38QKGLA7D4XDCueLjITDQCQdWSin7aSas7Nq5U+YXAwI0+GbmzBno0wfmzcOcj8fCEEg8H9wZ4Zzgu2wZVK4MGzc64eBKKeVaGoAvFxcH99wDd90l+0xV+nbtgkaNYOhQjsYVJM4EkoADv+BAglqF23++uDipdORwSBBWSikvp3PAl+vbF7ZulUVXutc3fX/8AV26gMPBpmFzCB/YnJsLNWLUIxGU6RoOYWH2n/O992DHDnldtNKRUsoHaABOa/Fi+PhjWT3UsqW7W+OZvv9eeqI33sj0blO4/82KVKgAI34Po0xlJwRegA0b4MMP4fHH9XVRSvkMDcCpTp+WN/jKleXNXl0UGQkRERAeDrfdhnniSYaGfkqvHvm5/XapS1GsmBPPP2YMlCghdX6VUspHaABOdf481KoFb70F+fO7uzW2yaxSUZZERsIdd8jq4+BgEv5YwDNJ3/LDe/Dww/Ddd7JN2qk++gh69HBylFdKKdfSAJyqZEkpsuBDMq1UlJUgvGkTPP203JwAJj6eH7tF8MOeMPr3h3ffdfI0+Z49sujquuugfHknnkgppVxPV0EfPQr33w/79rm7JbbLtFJRZg4dkjrHN94Ie/eCvz/G4SAuOZAf94UzZgz06+fk4JucLFMCTZrIfmyllPIxebsHbAw8+yzMmiXduTRyPXTrAdKrVHT541PXR/H7N5OpsnUN+26oQ4sXHqTjDUVg/Xro2xfT41WWjdrBkvciWGyFM3BGGOHhLmj811/D0qWy6CsgwAUnVEop18rbAXjsWFlB9OGHlyTcyPXQrYcoUySEqHSCcGoFo6nro5j46Xi+H/cWAUkJsMTikfhk6NGFjtu3s2a9gzfuhcWLw6hePYzJk6FaNRc0/N9/4c034c47pReslFI+KO8OQR85At27w223XVHOLsdDt5eZuj6KxoMXUrH3LBoPXsjU9VG5bnZ29GxVlZCAS1NSXajFu2QJxR96gFHj3iIoKSHlP4Kh/p71DPxlHw896qBhQ9kS/eWXknzKJcHXGHjuORmC/uYb3YutlPJZebcHPHgwnD0L337L5XkTszJ0ezWe0IvuWC8UR2wMaz4bTdj6RRRLjufcgPe5o14oTF9LuYN7WVLxJsL3rMXPJJPgCOD34x1Z99EtbA2QBeG9ekGhQi5proiPl9KPgwZpAQyllE/LuwG4f39o1izdbt3Vhm6zIrNedGoAztI8c9o9uJdnmIqNlfnrefOgTBkICZH0kLffDrt3w8030+7ECdql/Zn/NgJ3Qbt2PNSrAFHRsdT7bzv1Vv7HnP86sXxbE0redIi1U0tTrlyWL9c+QUEy/6uUUj4u7wXg5GRISpJuXdu26T7lkiLzKS4M3WbR1XrRWeohz50L7dtLr9DhkBuGe++VMnxnz0LBglee4J13JACXKiUrmffsgfnz5bodjgv5rU+dtgijLt/OjGP67hZMjQsk+LpjVGixnE9eqHAh+Lp0Mdo770CHDtAg3cIhSinlU/LeHPC4cVLE/cCBDJ/SsV4ogzrVJrRICBYQWiSEQZ1qZyvwZNRbTn08tYd8U9Q2Xoj8lQb/baHS/p3M/WqiPDE5WYLv+fMSNBMTYdUqOHhQvl+gALRqdXGO1OGQvUHvvXfx+yNHSk8/KAgcDpIDAplwOJw770xJLNWnGESVoniN45R6YCX1nv2LT16ocEkPvc/kTURFx2K4eJPglLns6dPh/fdh9mz7j62UUh4ob9UDPnMGbrhBkjpERoKf8+4/Lu/hgvSiUwN5xd6zaLp7Dd9OHoh/sjzHAjaVqkztQ7vkB959V1ZoJyZKDdwFCy4dhk6pw3uhRm6a7584IaWMt2+Hs/MiSV4YwYQj4awgjOrVJba3aycj1hmVDmw8eGG6Q/GhRUL4s3czW35PABw+LKvQy5SRmwyt96uU8hGZ1QPOW0PQ778vSSamTXNq8IWLw8gZDd+WKRJCn4jvCUgJvsnAjGpNGNPuWSalHuS996QsYgZzwAkNwjg4egExv0ewvnA4i0aHsb2XBN2jRy8+LzAwjFtvDeOB3vBTO6hSJWvXYMditKsyBrp1k5ujn3/W4KuUyjPyTgDesQM++UTe7G++2SWn7Fgv9NJh6+XLodPLMHQoPVtV5dsVnXl/5if4JyeR4PBnfKN7eLhL+CXHSLo5jL+LhrFrF+z6FP7+W0rx7tolybuSksIACcwlSsiasg4d5N/UjwoVMu7lZsaOxWhXNX68DDt/8QXUqGHfcZVSysPlnQD8xReySnjQINeed9kymYvduBG2bIGiRWHbNjq2bQv9utO9dFmqbF3DrhoNePCZTnSoG8rOnbJuav58WLQIoqMvHq5wYbj+ermH6NpVerNVqkDVqhKA7WTHYrSruv9+WRT38MP2HVMppbxA3pkDTkyUAFinjuvOuXChZHMyRhZLvfqqDCtfVm3p8GGZvk0Nuv/9J49fdx20aCHpkKtWlUBbrJhrc1M4bRV0bCycO2f/XYNSSnmQvD0HfP48xMRIz9OVwRdg5cqLn/v5SbBJCb7GSNAdMEA6ySBNbNZMEmC0aAGVKmUebF2xReiKYXS79OoFkyfLTVHhwvYfXymlPJzvb0P69FNZ+XzokGvOFxMDTz4pY8fh4RAcLBOwgYGkVjFYulRK7LZoAf/8Ax98AKtXy8KpSZMkE2PlylcPvi7bImS333+Hzz+X4WcNvkqpPMq3e8AHDsDAgTIMfO21zj/f7t2SKGPjRllQ9Prr0s1NWcW8yhHGO60kv8a118Lw4VJuNzg4+6fKSqYtj3T4sCyEu/FG18/HK6WUB/HtANyrl8z9fvyx8881fTo8+qgMNc+aBXffLY+HhbEhJIx334UZM2QUeuhQeOEFyJcv56dzyRYhu6VuOTp9WubHc3LnoZRSPsJ3A3BkpJQbfPttmUx11jkiIqBIEYmo9evLGHKFCgCcPAnPPw8TJshT3n9fCjCll0Eyu1yyRchu587JuPqwYVCzprtbo5RSbuW7AXjGDBnn7d3bOce/PAtVz56ywjmlV7drl6Sa3rMH+vaV0egiRew7vUu2CNmtQAGYOdPdrVBKKY/gu4uw/vc/2LBB3vSdYd482UqTlCRBuGjRC8F36VJJ8Xj0qGwrGjjQ3uAL9uSrdpmTJ+Hxx2XFmWVpjV+llMIXe8AJCbB/v9SSLVXKOedISpLFVSBzvmlWOP/0Ezz1lIxCz5qV9bSPOeG0LUJ2On8e7rlHsoA9/viF4XmllMrrfK8HPGqUZK3YssU5xzcGevSAJUvk3/ffhwULSL4ljHfekXVYjRvDihXODb5eITlZFl0tXgw//HDhJkUppZSv9YDPnJGSfI0aOS+v8LBhMGKETOoOGwbISHS3rrLY6skn4csv7akp4NJavM7w9tuS63nQIMmbqZRS6gLfCsAffQRHjsiWIGfMMyYmyrjyAw9ImUBkW2vHjpL0asgQWYtlx6kvL2eYmmgD8I4gfPasvA7PPivbwZRSSl3CdwLwwYPSI73/frjlFuecw98f/vhDPvfzY9cuyfFx5Aj89ptMddrFaxNtpCpQQOZ98+fXRVdKKZWOXM0BW5ZVxLKsSZZlbbcsa5tlWWFX/yknWbxY5hz/9z/7j715s1SwP3FCVjoHB3PqlGwzOntWpoPtDL7gpYk2QHJqPv44xMVJmkl/37nHU0opO+X23fEz4A9jzH2WZQUCucjtlEudO0ty5eLF7T3u/v1w110S3M+ehWLFSEqCLl0k8+SCBdAg3ToXueOViTb27JG7knz54NQpzXSllFKZyHEP2LKswsDtwCgAY0y8MSbapnZlz9at8q/dwXfePImuJ05IAYHy5QGpVjR7NpRvu53Hfp9F48ELbS+C0LNVVUICHJc85tGJNo4flxuVhAT55ThrC5hSSvmI3AxBVwSOAqMty1pvWdZ3lmXlv/xJlmU9Y1nWGsuy1hw9ejQXp8vAsmWS1nDCBHuPu3w5tG4tq6ySkqTKEZLd8sMPoUj9fSRV3e20SkRelWjj9Gno0AH27YNp06BaNXe3SCmlPF5uhqD9gZuAl40xKy3L+gzoDbyT9knGmG+AbwAaNGhgcnG+Kxkjy47LlIF27Ww9NDNnyrAzyOrniAhWOcJ46ikoVOkkhe64dJ+xMxZIeUWiDZAFcNu2SRaSJk3c3RqllPIKuekB7wf2G2NSq85PQgKy60yeLBkv3nsvd6WF0tOuHYSEXKjle7RmOPfcA6VLQ+E2q7EcV95LePwCKbv99ZfcBFWtKvO/99/v7hYppZTXyHEANsYcAv6zLCt1UrI5sNWWVmVFQoIUWqhZEx57zL7jJifDV19BvXqywmrgQM7/voC2H4Rx6pRsbS1XJv2BA49eIGUnY6SmYr160usFWfGslFIqy3K7CvplYFzKCug9QLfcNymLNm+GY8dkUjYHW10yzDI1ciS89JIElC5dMI3CePoxWLVKOty1a0PPRC+sRGSXmBhJdj1+vCQkufded7dIKaW8Uq4CsDFmA+CETThZUK+eVNcpVCjbP5pRlqn8/+6hxZtvymrezp0BSa71008yyp261zd1Xtar00TmxL//StqvDRtkv3Xv3ppkQymlcsi7syTkcNgzvSxT58/HU+qVNyAoCL77DiyLP/6QLIr33y81fdPymgVSdtq8GfbulVrLbdq4uzVKKeXVvDsA51B6i6WeWj2VG/dtgXHjoEwZTp6UQj41a8Lo0Xm4o3fwICxcCA89BHffLQHY7uLGSimVB+XJAJxelqlFlRpQxorn8S5dAOn5HjkitRfyX7G7OQ/Yt082PI8aJYuuWraEkiU1+CqllE18rx5wFqTNMmUZ2eu7v0wlinw0BCyLxYvh22/h1VfhJtdurHK/AwekpmKVKvJLeOwx2eNbsqS7W6aUUj4lT/aA0y6iemDGd1Q/d5jYb0bRoV4ocXFSQa9CBRgwwL3tdKn4+ItFjCdPhhdegDfegHLl3NsupZTyUXkyAEPKIqp1f0DkL5Jy8uYKgCzu3bED5szx8aHn48cl3eby5ZLOMyBA5nrLlIGoKPsTmyillLpEnhyCBmDpUnj6aZnfjIiAyEi2bIHBg+Hhh2XK02ckJ8Pff1/8+sUXoUQJKbH40UeS1KRpU0m5CRp8lVLKBfJsD5ghQyT4AsTHk7wogqdnhlGoEHz8sXublmXJyXD0qMzbHjwIt98OBQpI9/277+DQIfk4cEASaBw+DNdcI2Uby5eHW2+Vak8heSSDl1JKeZC8GYATE2HtWvDzk/1FgYFMORFOZCT8+KMHrTdKToYtW2DjRhkWfvBBmZyeOVPmaA8evNhrBVizBurXl6C8eTNcey00bCj/1q59sT5vx47uuBqllFJp5M0A7O8vQW3ZMtixgyM1wun2SBh33gmPPOLuxiHDxb16wZIlMlebqnp1CcDXXgt33AGhofJRpoxUiUgtA/jww/KhlFLKY+W9AHzypGTQuuYa6NQJY+CZe6Qj+dVXbki4sXOnbDaOiJA52SeflGHkDRukIlN4ONxyi6xGTl0V1qCBdNWVUkp5rbwVgI2B++6TOc+ZMwHZcTNtmkwJV67swrZs3gwDB8LEidKuKlVkNTZIb3bPHhc2RimllKvlrQA8bZpstRkxAoDoaHj5ZahbF157zcVtefppCcK9e8t8btmyLm6AUkopd8o7Afj8eXj9dUnu/OyzgMS+w4elxm8OKhpmz/r10s3+/HNZ5fX991CqFBQr5uQTK6WU8kQ+G4Avr/f7zcH51NyzB+bOBX9/Nm6Eb76B7t1lStVp1q2TlFrTp8vc84YNsg2oenUnnlQppZSn88lEHKn1fqOiYzHAgZPn8JswgYNNW0rwA/r0kXj47rtOaoQx0uNt0EBWMw8YIPWLU86vlFIqb/PJHvDl9X6N5UfHh4dRJTiZWcCiRTB7thT7ceoI8Pr18MAD8PXXOa5drJRSyjf5ZABOW++3zOkjnAgpRFxAMFsTpWPaq5eseXrpJSecfN8+SEqCSpVkq1BgYB4uJqyUUiojPjkEXaZISmpFYxg+fSi/jO8DxlCmSAiTJsHq1fDee07IwLh4sQw5P/KIRPqgIA2+Siml0uWTATi13m/3P8fTIGoby8rXJSTQn1ebVeWtt2Qh9KOP2nhCY+CLL+DOO6XIwejRGniVUkplyieHoDvWC6XE6uU0/vNnDPDUuunUfu5htq4NZdcumDEDHA6bTnb+vFQXGjVKMleNHQuFCtl0cKWUUr7KJ3vAALdNHoUFWEBwciKN/tnCgAHQpAm0aWPjiZKSZGtR374wdaoGX6WUUlnikz1gAM6cuaTa0biocA4fhilTbBodjo2VbnS+fPDnnzLfq5RSSmWR7wbgZctgwQJYvZoTN4bzSucwOnWCsDAbjp2UBA89BKdPS+1dDb5KKaWyyfcC8LFjEiBLlZJFUXfeSf/u0mH93/9sOkfPntKV/uwzGyeTlVJK5SW+Nwfcr5/UxT19GoDdu6XM4JNPQtWqNhz/88/hk0/glVckj6VSSimVA74VgPfskQTPnTtfWAz1zjtSaKFfPxuOP2MG9OgBHTrARx/ZcECllFJ5lW8F4H79ICBAoi5SB2H8eCk1WKaMDcevWBHuuQfGjdOhZ6WUUrniOwF40yYJjC+/fCHa9uoFxYvLlG2unDolyTZq1YJJkyB//ty3VymlVJ7mOwF4/nwoUkSiLrBwoTz09tu5rIMQHS1Lp/v0saOVSimlFOBLAfjVV2HXLihWDGOgf3/pCD//fC6OGR8PnTrJcVu1squlSimllA8EYGMkQMKF2oIREbB0KfTuDcHBuTj2K69I7cJRo+COO3LdVKWUUiqV9wfguXPhhhskIUaKAQOgdGl4+ulcHPePP2T/0htvSHUjpZRSykbenYgjOVnmZq+77kIPdfFi+fj001z2fuPi4LbbYOBAW5qqlFJKpeXdAXjSJFi/HsaMkcL3SJ3fa6+FZ57J5bE7dpT9vlpWUCmllBN47xD00qVSBrBiRejaFZD0zwsXwptvQkhIDo87a5akmExO1uCrlFLKabwzAEdGQosWkvc5KgpWrQJk7rdUKXj22Rwe98QJeOopWXSVmGhfe5VSSqnLeOcQdETExQCZlAQRESw3YcyfD8OGSYXAHHn5ZQnqs2dfGNJWSimlnME7e8Dh4RIgHQ75NzycAQOgZEl47rkcHnPyZPj5Z0ljWbeujY1VSimlruSdPeCwMKn1GxEB4eGssMKYOxeGDMlhlsjYWHjhBahXTzNeKaWUcgnvDMAgQTgsDIABd0GJEhJDcyQkBH75RQ4SEGBfG5VSSqkMeG8ATrFqleTMGDQIChTIwQFOn5bSheHhdjdNKaWUypB3zgGnMWCAZKB88cUc/PDhw3D99TBypO3tUkoppTLj1QF49Wr4/Xd4/XUoWDAHB3jhBSk1qL1fpZRSLuaVAXjq+igaD15I0y6H8Q9JoEKTA9k/yNy5svK5f3+oXt32NiqllFKZ8boAPHV9FH0mb2LP9gBid5cif4M9DJz7F1PXR2X9IImJ0m2uVEnKGCqllFIu5nUBeOicHcQmJBGzrQx+QQkUqv8PsQlJDJ2zI+sHWbcO/v4bPvwQgoKc11illFIqA163CvpAdCwARcK3U6DePvyCEi95PEtuvhl274YyZZzRRKWUUuqqvK4HXKaIVFmwLAgoEnvF41e1cycYA6GhWmxBKaWU23hdAO7ZqiohAY5LHgsJcNCzVdWr//CePVC7Nnz8sZNap5RSSmWN1w1Bd6wXCshc8IHoWMoUCaFnq6oXHs9Ur17g7w9duji5lUoppVTmvC4AgwThLAXctJYuhUmTJHOHzv0qpZRyM68bgs6R5GTZbhQaCm+84e7WKKWUUt7ZA862PXvgv/9yWSxYKaWUsk/eCMBVqsCuXTmsVaiUUkrZz/eHoNeskcxXBQuCn+9frlJKKe/g2xEpKgqaNoU333R3S5RSSqlL+HYAfust6f2+/LK7W6KUUkpdwncD8Pr1MGaMrH6uWNHdrVFKKaUu4bsB+N13oWhR6NPH3S1RSimlrpDrAGxZlsOyrPWWZc20o0G2OHMG9u2TPb+FC7u7NUoppdQV7NiG9AqwDShkw7HsUbAgbNgg879KKaWUB8pVD9iyrLJAG+A7e5pjgz17IDpathwFBrq7NUoppVS6cjsE/SnwJpCc+6bY5Kmn4NZbpeSgUkop5aFyHIAty2oLHDHGrL3K856xLGuNZVlrjh49mtPTZc2iRfLx7LNa61cppZRHs0wOe4qWZQ0CHgESgWBkDniyMebhjH6mQYMGZs2aNTk631UZA7ffLkPQu3dDcLBzzqOUUkplkWVZa40xDdL7Xo57wMaYPsaYssaYCkBnYGFmwdfp5s2DZcvg7bc1+CqllPJ4vrMPeP58KF8ennzS3S1RSimlrsqWAGyMiTDGtLXjWDn24YeS/SooyK3NUEoppbLC+3vAxsC//8rnxYq5ty1KKaVUFnl/AJ4yBSpXhshId7dEKaWUyjLvDsDJydCvH1SqBA0burs1SimlVJbZkYrSfSZOhM2b4eefwd+7L0UppVTe4r094KQk6N8fataEBx5wd2uUUkqpbPHebuOaNbBrF4wfDw6Hu1ujlFJKZYv3BuBbbpEAXK6cu1uilFJKZZv3BmCA665zdwuUUkqpHPHeOWCllFLKi2kAVkoppdxAA7BSSinlBhqAlVJKKTfQAKyUUkq5gQZgpZRSyg00ACullFJuoAFYKaWUcgMNwEoppZQbaABWSiml3EADsFJKKeUGGoCVUkopN9AArJRSSrmBZYxx3cks6yiwz8ZDlgCO2Xg8d9Jr8Ty+ch2g1+KpfOVafOU6wP5ruc4YUzK9b7g0ANvNsqw1xpgG7m6HHfRaPI+vXAfotXgqX7kWX7kOcO216BC0Ukop5QYagJVSSik38PYA/I27G2AjvRbP4yvXAXotnspXrsVXrgNceC1ePQeslFJKeStv7wErpZRSXskrArBlWa0ty9phWdYuy7J6p/P9IMuyJqR8f6VlWRXc0MyrsiyrnGVZiyzL2mpZ1hbLsl5J5znhlmWdsixrQ8rHu+5oa1ZYlvWPZVmbUtq5Jp3vW5ZlDU95Xf6yLOsmd7QzM5ZlVU3zu95gWdZpy7J6XPYcj31NLMv63rKsI5ZlbU7zWDHLsuZZlvV3yr9FM/jZx1Ke87dlWY+5rtXpy+BahlqWtT3l/88Uy7KKZPCzmf5fdLUMrqW/ZVlRaf4f3Z3Bz2b6fudKGVzHhDTX8I9lWRsy+FlPe03Sff9169+LMcajPwAHsBuoBAQCG4Ealz3nBeCrlM87AxPc3e4MrqU0cFPK5wWBnelcSzgw091tzeL1/AOUyOT7dwOzAQtoBKx0d5uvcj0O4BCyb88rXhPgduAmYHOaxz4Eeqd83hsYks7PFQP2pPxbNOXzoh54LS0B/5TPh6R3LSnfy/T/oodcS3/gjav83FXf79x9HZd9/yPgXS95TdJ9/3Xn34s39IBvBnYZY/YYY+KBX4AOlz2nA/BjyueTgOaWZVkubGOWGGMOGmPWpXx+BtgGhLq3VU7VARhjxAqgiGVZpd3dqEw0B3YbY+xMFuNUxpglwInLHk779/Aj0DGdH20FzDPGnDDGnATmAa2d1c6sSO9ajDFzjTGJKV+uAMq6vGE5kMHrkhVZeb9zmcyuI+U99gFgvEsblUOZvP+67e/FGwJwKPBfmq/3c2XQuvCclD/WU0Bxl7Quh1KGyesBK9P5dphlWRsty5ptWVZN17YsWwww17KstZZlPZPO97Py2nmSzmT8ZuItrwlAKWPMwZTPDwGl0nmOt702AE8gIyrpudr/RU/xUspw+vcZDHV60+vSBDhsjPk7g+977Gty2fuv2/5evCEA+xzLsgoAvwE9jDGnL/v2OmQItA7wOTDVxc3LjtuMMTcBdwEvWpZ1u7sblFOWZQUC7YGJ6Xzbm16TSxgZP/P6rQ6WZb0NJALjMniKN/xfHAlUBuoCB5HhW2/Whcx7vx75mmT2/uvqvxdvCMBRQLk0X5dNeSzd51iW5Q8UBo67pHXZZFlWAPLijzPGTL78+8aY08aYsymf/w4EWJZVwsXNzBJjTFTKv0eAKcjwWVpZee08xV3AOmPM4cu/4U2vSYrDqUP9Kf8eSec5XvPaWJb1ONAWeCjlDfIKWfi/6HbGmMPGmCRjTDLwLem30Stel5T32U7AhIye44mvSQbvv277e/GGALwauN6yrIopvZTOwPTLnjMdSF2Vdh+wMKM/VHdKmTMZBWwzxnycwXOuTZ2/tizrZuQ18ribCcuy8luWVTD1c2SxzObLnjYdeNQSjYBTaYZ6PE2Gd/Pe8pqkkfbv4TFgWjrPmQO0tCyraMpQaMuUxzyKZVmtgTeB9saYmAyek5X/i2532fqHe0i/jVl5v/MEdwLbjTH70/umJ74mmbz/uu/vxd0r07Lygaym3YmsDnw75bH3kD9KgGBk6HAXsAqo5O42Z3AdtyHDG38BG1I+7gaeA55Lec5LwBZk9eMK4FZ3tzuDa6mU0saNKe1NfV3SXosFfJHyum0CGri73RlcS34koBZO85hXvCbITcNBIAGZl3oSWf+wAPgbmA8US3luA+C7ND/7RMrfzC6gm4deyy5k7i317yV1t0MZ4PfM/i964LX8lPJ38Bfypl/68mtJ+fqK9ztPuo6Ux39I/ftI81xPf00yev9129+LZsJSSiml3MAbhqCVUkopn6MBWCmllHIDDcBKKaWUG2gAVkoppdxAA7BSSinlBhqAlVJKKTfQAKyUUkq5gQZgpZRSyg3+D500RJBjqRw7AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 576x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "prstd, iv_l, iv_u = wls_prediction_std(res)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(8,6))\n",
    "\n",
    "ax.plot(x, y, 'o', label=\"data\")\n",
    "ax.plot(x, y_true, 'b-', label=\"True\")\n",
    "ax.plot(x, res.fittedvalues, 'r--.', label=\"OLS\")\n",
    "ax.plot(x, iv_u, 'r--')\n",
    "ax.plot(x, iv_l, 'r--')\n",
    "ax.legend(loc='best')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Time-Series Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "from statsmodels.tsa.arima_process import arma_generate_sample"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Gerando dados\n",
    "np.random.seed(12345)\n",
    "arparams = np.array([.75, -.25])\n",
    "maparams = np.array([.65, .35])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Parâmetros\n",
    "arparams = np.r_[1, -arparams]\n",
    "maparam = np.r_[1, maparams]\n",
    "nobs = 250\n",
    "y = arma_generate_sample(arparams, maparams, nobs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/dmpm/opt/anaconda3/lib/python3.8/site-packages/statsmodels/tsa/base/tsa_model.py:524: ValueWarning: No frequency information was provided, so inferred frequency M will be used.\n",
      "  warnings.warn('No frequency information was'\n"
     ]
    }
   ],
   "source": [
    "dates = sm.tsa.datetools.dates_from_range('1980m1', length=nobs)\n",
    "y = pd.Series(y, index=dates)\n",
    "arma_mod = sm.tsa.ARMA(y, order=(2,2))\n",
    "arma_res = arma_mod.fit(trend='nc', disp=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                              ARMA Model Results                              \n",
      "==============================================================================\n",
      "Dep. Variable:                      y   No. Observations:                  250\n",
      "Model:                     ARMA(2, 2)   Log Likelihood                -245.887\n",
      "Method:                       css-mle   S.D. of innovations              0.645\n",
      "Date:                Fri, 09 Jul 2021   AIC                            501.773\n",
      "Time:                        20:24:41   BIC                            519.381\n",
      "Sample:                    01-31-1980   HQIC                           508.860\n",
      "                         - 10-31-2000                                         \n",
      "==============================================================================\n",
      "                 coef    std err          z      P>|z|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "ar.L1.y        0.8411      0.403      2.089      0.037       0.052       1.630\n",
      "ar.L2.y       -0.2693      0.247     -1.092      0.275      -0.753       0.214\n",
      "ma.L1.y        0.5352      0.412      1.299      0.194      -0.273       1.343\n",
      "ma.L2.y        0.0157      0.306      0.051      0.959      -0.585       0.616\n",
      "                                    Roots                                    \n",
      "=============================================================================\n",
      "                  Real          Imaginary           Modulus         Frequency\n",
      "-----------------------------------------------------------------------------\n",
      "AR.1            1.5618           -1.1289j            1.9271           -0.0996\n",
      "AR.2            1.5618           +1.1289j            1.9271            0.0996\n",
      "MA.1           -1.9835           +0.0000j            1.9835            0.5000\n",
      "MA.2          -32.1794           +0.0000j           32.1794            0.5000\n",
      "-----------------------------------------------------------------------------\n"
     ]
    }
   ],
   "source": [
    "print(arma_res.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Conheça a Formação Cientista de Dados, um programa completo, 100% online e 100% em português, com 400 horas, mais de 1.200 aulas em vídeos e dezenas de projetos, que vão ajudá-lo a se tornar um dos profissionais mais cobiçados do mercado de análise de dados. Clique no link abaixo, faça sua inscrição, comece hoje mesmo e aumente sua empregabilidade:\n",
    "\n",
    "https://www.datascienceacademy.com.br/bundle/formacao-cientista-de-dados"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fim"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "### Obrigado\n",
    "\n",
    "### Visite o Blog da Data Science Academy - <a href=\"http://blog.dsacademy.com.br\">Blog DSA</a>\n"
   ]
  }
 ],
 "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.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}