{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Prediction (out of sample)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import statsmodels.api as sm\n",
    "\n",
    "plt.rc(\"figure\", figsize=(16,8))\n",
    "plt.rc(\"font\", size=14)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Artificial data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "nsample = 50\n",
    "sig = 0.25\n",
    "x1 = np.linspace(0, 20, nsample)\n",
    "X = np.column_stack((x1, np.sin(x1), (x1-5)**2))\n",
    "X = sm.add_constant(X)\n",
    "beta = [5., 0.5, 0.5, -0.02]\n",
    "y_true = np.dot(X, beta)\n",
    "y = y_true + sig * np.random.normal(size=nsample)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Estimation "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                            OLS Regression Results                            \n",
      "==============================================================================\n",
      "Dep. Variable:                      y   R-squared:                       0.985\n",
      "Model:                            OLS   Adj. R-squared:                  0.984\n",
      "Method:                 Least Squares   F-statistic:                     1034.\n",
      "Date:                Mon, 26 Oct 2020   Prob (F-statistic):           3.38e-42\n",
      "Time:                        07:28:11   Log-Likelihood:                 2.1453\n",
      "No. Observations:                  50   AIC:                             3.709\n",
      "Df Residuals:                      46   BIC:                             11.36\n",
      "Df Model:                           3                                         \n",
      "Covariance Type:            nonrobust                                         \n",
      "==============================================================================\n",
      "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "const          4.8261      0.082     58.588      0.000       4.660       4.992\n",
      "x1             0.5217      0.013     41.064      0.000       0.496       0.547\n",
      "x2             0.4193      0.050      8.397      0.000       0.319       0.520\n",
      "x3            -0.0213      0.001    -19.125      0.000      -0.024      -0.019\n",
      "==============================================================================\n",
      "Omnibus:                        0.719   Durbin-Watson:                   2.493\n",
      "Prob(Omnibus):                  0.698   Jarque-Bera (JB):                0.760\n",
      "Skew:                          -0.087   Prob(JB):                        0.684\n",
      "Kurtosis:                       2.422   Cond. No.                         221.\n",
      "==============================================================================\n",
      "\n",
      "Notes:\n",
      "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n"
     ]
    }
   ],
   "source": [
    "olsmod = sm.OLS(y, X)\n",
    "olsres = olsmod.fit()\n",
    "print(olsres.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## In-sample prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 4.2927808   4.7556813   5.18412681  5.55526333  5.85448478  6.07783268\n",
      "  6.23264654  6.336358    6.41362686  6.49228952  6.5987856   6.75381391\n",
      "  6.96893174  7.24465603  7.57037841  7.92610816  8.2857566   8.6214234\n",
      "  8.90798052  9.1272007   9.27075209  9.34156721  9.35336112  9.32837807\n",
      "  9.29373684  9.2769752   9.30152552  9.38286492  9.52597242  9.72451088\n",
      "  9.9618682  10.2138866  10.45283175 10.65194935 10.78986086 10.8540759\n",
      " 10.84304438 10.76641119 10.64343156 10.4998075  10.36346501 10.25996592\n",
      " 10.20830818 10.21780451 10.28655237 10.40174654 10.541783   10.67980836\n",
      " 10.78813184 10.84277477]\n"
     ]
    }
   ],
   "source": [
    "ypred = olsres.predict(X)\n",
    "print(ypred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create a new sample of explanatory variables Xnew, predict and plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[10.81346852 10.67115085 10.43226692 10.13429316  9.82656174  9.55818234\n",
      "  9.36601836  9.26566152  9.24761441  9.27961584]\n"
     ]
    }
   ],
   "source": [
    "x1n = np.linspace(20.5,25, 10)\n",
    "Xnew = np.column_stack((x1n, np.sin(x1n), (x1n-5)**2))\n",
    "Xnew = sm.add_constant(Xnew)\n",
    "ynewpred =  olsres.predict(Xnew) # predict out of sample\n",
    "print(ynewpred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Plot comparison"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA6MAAAHWCAYAAACCFl9HAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAB1xElEQVR4nO3de3yO9R/H8de1GcYwx9gcw285VTJRFEJLClGhExEpHSiLUVIKGTnkUA6lEh2kkdMioRQaI+SUQ9icmeMw2/X742uy2Wa4d9/37r2fj8f9mF3XdV/70N229/39fj9fy7ZtRERERERERJzJy9UFiIiIiIiISM6jMCoiIiIiIiJOpzAqIiIiIiIiTqcwKiIiIiIiIk6nMCoiIiIiIiJOpzAqIiIiIiIiTpfL1QUUK1bMLl++vKvLEBERERERkSywevXqw7ZtF0993OVhtHz58kRFRbm6DBEREREREckClmX9m9ZxTdMVERERERERp1MYFREREREREadTGBURERERERGnUxgVERERERERp1MYFREREREREadzeTfdqzlx4gQHDx4kISHB1aWIB/Hx8aFEiRIULFjQ1aWIiIiIiORIbh1GT5w4wYEDBwgMDMTX1xfLslxdkngA27aJj48nJiYGQIFURERERMQF3Hqa7sGDBwkMDCRfvnwKouIwlmWRL18+AgMDOXjwoKvLERERERHJkdw6jCYkJODr6+vqMsRD+fr6avq3iIiIiIiLuHUYBTQiKllGry0REREREddx+zAqIiIiIiIinkdhVERERERERJxOYTQLdOzYEcuysCzr0hYijRo1YuzYsde0RnHJkiVYlsXhw4ezsFoRERERERHnUxjNIk2aNGHfvn3s2rWLn376iYcffpi3336be+65h9OnT7u6PBEREREREZfKEWE0IjqGekMWU6HPXOoNWUxEdEyWf808efJQsmRJAgMDuf3223nttddYsmQJa9asYejQoQBMnTqV2rVrU6BAAUqUKMFjjz12ae/LXbt20ahRIwCKFy+OZVl07NgRgAULFnDPPfdQuHBhihQpQkhICJs2bcryv5OIiIiIiIijeHwYjYiOIWzmemLi4rGBmLh4wmaud0ogTa169eo88MADfP/99wCcP3+ed955h3Xr1jFnzhwOHz5M+/btAShTpsyl6zZu3Mi+ffsYNWoUAKdPn6ZHjx6sWrWKJUuWUKhQIR5++GHOnz/v9L+TiIiII7jijWMREXGtXK4uIKuFR24hPiExxbH4hETCI7fQqmag0+upWrUqixYtAqBTp06Xjt98882MHz+eKlWqsHfvXkqXLk2RIkUAKFGiBMWKFbt0bZs2bVLc87PPPqNgwYKsWrWK+vXrO+FvISIi4jjJbxwn/7xOfuMYcMnPahERcQ6PHxmNjYu/puNZzbbtS/tbrlmzhpYtW1KuXDkKFChAcHAwALt3787wHtu3b+eJJ56gYsWKFCxYkJtuuomkpKSrPk9ERMQdZfTGsYiIeC6PD6MB/r7XdDyr/f3339x8882cPn2akJAQ8uXLx5dffsmff/7JggULAK463fbhhx/m0KFDfPLJJ6xcuZLo6Ghy5cqlaboiIpItudsbxyIi4hweH0ZDQ4Lw9fFOcczXx5vQkCCn17JhwwYWLFjAo48+yubNmzl8+DCDBg3i3nvv5ZZbbuHgwYMprs+dOzcAiYn/vVt85MgRNm3aRN++fWnSpAlVqlTh5MmTXLhwwal/FxEREUdxtzeORUTEOTw+jLaqGcjg1jUI9PfFAgL9fRncukaWr0E5d+4c+/fvJzY2lnXr1vHhhx/SsGFDatWqRa9evShbtix58uRhzJgx7Nixg7lz5/LWW2+luEe5cuWwLIu5c+dy6NAhTp06ReHChSlWrBgTJ07kn3/+YenSpXTr1o1cuTx++a+IiHgod3rjWEREnCdHJJhWNQOd3gBh0aJFlCpVCm9vb/z9/alevTpvv/02zz//PLlz5yZ//vx8/vnn9O3bl7Fjx3Lrrbfy4Ycf8sADD1y6R2BgIO+88w79+vXjueee45lnnmHKlCl88803vPLKK1SvXp1KlSoxfPjwK5oaiYiIZBfJP6PDI7cQGxdPgL8voSFBal4kIuLhLNu2XVpAcHCwHRUVlea5TZs2UaVKFSdXJDmJXmMiIiIiIlnLsqzVtm0Hpz7u8dN0RURERERExP0ojIqIiIiIiIjTKYyKiIiIiIiI0ymMioiIiIiIiNNlKoxalnWvZVmzLcuKsSzLtiyrY6rzrS3LirQs69DF8w2zoFYREREREZeLiI6h3pDFVOgzl3pDFhMRHePqkkSypcyOjPoBG4BXgfg0zucHfgdec1BdIiIiIiJuJyI6hrCZ64mJi8cGYuLiCZu5XoFU5Dpkap9R27bnAfMALMuaksb5Ly+eK+bI4kRERMR1IqJjtPenSCrhkVuIT0hMcSw+IZHwyC36/0PkGmUqjIqIiEjOkjz6k/xLd/LoD6BfuCVHi41La5Jg+sdFJH0uaWBkWVZXy7KiLMuKOnTokCtKEBERkQxkNPoj4kmudf1ngL/vNR0XkfS5JIzatj3Btu1g27aDixcv7ooSREREJAMa/ZGc4HrWf4aGBOHr442dBEnnvQHw9fEmNCTISVWLeA5t7eJglmVl+OjYsaOrSxQREbkqjf5ITnCtMwBsGwIvBFJ5R332f9KE48srE+jvy+DWNTR9XeQ6aM2og+3bt+/Sn+fMmUOXLl1SHPP1TflDPCEhAR8fH6fVJyIikhmhIUEp1oyCRn/E82R2BsDff8P06eaxfTvkzu1Hs2bQtWtFHnywojNKFfFImd1n1M+yrNsty7r94nPKXvy87MXzRS6eq37xKZUuni+ZFUW7s5IlS156+Pv7pzh29uxZ/P39mT59Ovfddx++vr588sknTJkyBT8/vxT3WbJkCZZlcfjw4UvHfv/9dxo0aEC+fPkIDAzkhRde4MSJE87864mISA7RqmYgg1vXINDfFwsyPfoTHw/z58MHH8DEiTB7NqxcCbt2mXMi7iSjGQC7dsGQIXDbbVCtGgwaBBUqwOTJsH8/RETAgw86tVwRj5PZkdFg4JfLPn/n4uNzoCPQAvjssvMTL7tuwA1V6IHCwsIYNmwYkydPxsfHh0WLFl31OevXr+f+++/nnXfeYdKkSRw9epQePXrQqVMnZsyY4YSqRUQkp2lVMzBTUw937IB588zjl1/g7Nn0ry1QAG66Ke1HQADcdx+ken9WJMukngGQeDo357cFcmRfRSqEmWvq1oVRo+Dxx6FkjhtmEclamd1ndAlgZXB+CjDFIRVdRY8esHatM77Sf26/HUaOdNz9Xn75ZR599NFrek54eDht27bl9ddfv3Rs/Pjx1KxZk4MHD1KiRAnHFSgiIpKBc+dg2bL/AujWreb4//4Hzz9vRovq1oXjx+HAgfQfmzfD0qVw5Mh/9y5YEDp2hBdfhCDNCJYs1qpmIIkXoM/Q4+z9szhndxcD26LkrTB4MLRta0ZDRSRraM2oCwQHB1/zc1avXs0///zDN998c+mYbdsAbN++XWFURESy1L//mum38+bBzz/DmTOQNy80bAgvvQTNmkGlShcvvnAB9uyhoK8vZW7xgzvygVf6K4MSEuDwYdiyBSZNgvHjYfRoaNIEuneHhx6CXA76jSUiOobwyC3ExsUT4O9LaEiQGs/kYH/8AQNfCOSfdYHcfDM80Q/at4eqVV1dmUjOkO3CqCNHKF0lf/78KT738vK6FCyTJSQkpPg8KSmJ5557jp49e15xv8BA/RAVERHHS0yEjz4yAXHjRnOsfHl49lkz+tmwIeTLd/HiXbvgk0j46SeTVo8fT3mz/PnNHF0/vysePn5+lPLzo5S/Pw0738eHQ+9h0pRcfPwxPPIIlCkD3brBc8/Bjbz3mryNR/KUzORtPAAF0hzm6FEIC4MJE6B0aZgxA1q3BivdeYCXSUiA6Gj4/XfTXjcgwDwCA6FUKfBVx2mRzMp2YdQTFS9enDNnznDixAkKFiwIwNpUc5HvuOMONm7cSKVLbzuLiIhknb//hk6dTPOh+vVh+HATQIOCLv7CfuoULF4CkRcDaPJc3dKl4dFHzTzdhARzXerHyZPm45EjZsg1+fiJE/Dee5QoXpy+jzxC7wmPMudUQ8Z84kO/fvDOO/DYY2a0tG7dTAaHy2S0jYfCaM5g2/Dll9Crlwmkr70GAwaY90nSdeIErFgBv/1mHitWZNyNq3BhE0yTQ2pyUL38zzfd5LjhfpFsTP8XuIE6deqQP39+wsLC6NmzJ+vWrWPcuHEprunduzd169alW7duPP/88xQoUIDNmzfz448/8sknn7iochER8TQJCRAeboJfgQIwbRq0aweWnQTr1sEHkSaALl9uLvb1NUOkL74I998Pt9xy7Skx2enTsGCBGab66iu8J0ygZdGitGzVij0TH2XE2vuY/GVuvvoKatY0X/KJJy4bnb2KzG7jIZ5p0yZ44QWzTrluXVi40HTKvUJs7H/B87ffzOs+KclMNa9ZE7p2Ne/Q1Ktn5qrHxppHTEzKj7GxZkrB/v1mmsHlLMsE0tq1oW9fU5BIDqQw6gaKFCnCV199RWhoKJ9++in33nsvAwcO5Omnn750za233sqyZct48803adCgAYmJidx888088sgjLqxcREQ8ybp1ZjR0zRozAvnRR3DT+kXw9BTzm/vBg+bC224zHQVDQv77hdwR8ueHNm3MIz7ehN4ZM+DbbykzeTIf+vvzwUMtWVT4Md5c0oQuXfIQGgpvvGFGuq62bXeAvy8xaQTP9Lb3EM9w5gy89x4MG2Zmhk+YAJ07X7aM+fRp+OYbWLLEhM+dO83xfPlMSHzzTRM+69ZNewi1cGGz90t6EhPh0KGUITUmBvbuNXsf3XUXPPAAvP22QqnkOFbqtYrOFhwcbEdFRaV5btOmTVSpUsXJFUlOoteYiAicPw/vv2/2USxSBMaNgzZ37DRzGCMioHhxM+oZEgJNmzp/f4uzZ00YnjEDZs2C48exCxbkYN0WTI57lHdWPUCV2/IweTLUqpX+bVKvGQXw9fHO1P6p4v7Sak7lExvISy+ZJc0dOphR/+LFLz7h8GHzjsuYMWbObvHicM89JnjWr2+2U7jaOxw36uRJ8z9ceLiZth4SYkLpXXdl7dcVcTLLslbbtn1FF1eFUcnR9BoTkZwuKso0JNqwAZ56CkYOjqfopA/ggw/A2xveesuMgubJ4+pSjfPnTYOk774zQfnYMc4UL0fY+XcYd+IpXgv1ZsCA9HvIqJuuZ0r9RsOFE3k5vrgap7aUpEoV06G5QYOLF+/aZRZBT55sRuBbtDDD63ffff1TzG/UqVMwdqxCqXgshVGRNOg1JiI51dmzpnFLeLgZ6PzkY5uHLkRAz56mqVD79uakO3dsT0gwU3kHDIDVq4n1r0r3uPfYULEVEydZNGzo6gLFWeoNWUxMXDx2ksXJqPLE/fY/sC3KNN7FttkVyZ0b+OsvGDoUvv7ahM6nnoLQUPfax+XUqf9GSg8fNjMSBgxQKJVsL70wmv6mXyIiIuKRfv/dzED84AMzKrrph8089FGI2duiQAGzdm7aNPcOomCmUD70EPz5J8yYQUDJJH6gNd/trcu7jRbTrduVO8yIZ4qNi+fCyTwcmF6HY79UJW/ZI5TqvBSv2zaR+4+lphX0bbeZad6vvmrWhX72mXsFUTCLWt94w9Q3dKhZwH333Wak9I8/XF2diMMpjIqIiOQQZ86YZaD165vZiYt+OMmkIm9QsF4NWLUKRo82+ydems/o3iKiY6g3ZDEVwuZRb1thZn2xACZP5tbi+1hMYx79pClPVP6TH390daWS1fIdDmDflHs4f6AQRR9aS4nWq2h+cDE/Tn/DdHuOijJdjHbvNlN0S5d2dckZ8/Mzo7ZphdLff3d1dSIOozAqIiKSzVwKYX3mUm/IYiKiY676nD17THOfESPghW42m9/6isYvBpnpgB06mH1CX3452+x9mLxGMCYuHhuIiYunz+xNRNQMwWvbVhgxggb+a5l76E7Ot2jDaw9uutQMWDxHYqLZhmjTp7eTK995yjy1lGcTP2PR5BeZ8MP7VEg6bdZi/vsv9OtnOt9mJ8mhdNcu8/9qdLTpYH3//Qql4hEURkVERLKRtEJY2Mz1GQbSnTvh3nvNjhIrPl7L2PX34tvlKTM6tHIlTJoEJUo47y/hAOGRW1J0xQWIT0gkPHKL2WqmRw98du/gQv93eCjPQsLnV2dhmU7MHLkbF7fLEAc5eBCaNTNLKp96ymLBm3NZMrcT4fNHk5QnL38OHkv+XdvNhrTpdbTKLvLnN/sX7dxpQum6dSaUNmsGBw64ujqR66YwKiIiko1kGMLSsHWr2a3idFwCW5v3pM6LtWDzZhNAV6yAO+90RtkOF5vGfqFXHC9QgFzv9CfP3h3EdejBownTaN6zMrMq9GDvGg2TZme//go1a8KyZTB9yL98fuZRGvVoTwW/XDBrFv+L2UrtPi9mm5H+TEsOpTt2mFC6dKn5f/ivv1xdmch1URgVERHJRjIVwi7asMGMiOY9G8f2oGbcNH0kPP+8SaidO4NX9v01IMA/7ZGuNI8XK0bRKcPJtWMb2+56hof+HUOBWv8juut4M89Tso2kJLOEslEj8M97lp2d36PdO1Ww5s2DgQNh40azVYurtmhxluRQ+uuvcOGCGSXV4mjJhrLvTyEREZEcKLMhbM0a07elvL2Tjf53U2DNMvj8c7NtRHZbN5eG0JAgfH28Uxzz9fEmNCQo3ed4ly9D9d8nsu+nDWwuWJuaE19kV+DdnFsRndXligMcPQqtWkHv3jDwzh/5K6kapca9Bc2bm9H+N980U7Rzklq1TDfpoCBo2RKGDUPz0CU7URiVG/bSSy/R8LLN3Dp27MhDDz10Q/ccMGAA1atXv8HKREQ8T2ZC2IoVcN99cE+uP1h+oQ55ju6HhQvhmWecXW6WaVUzkMGtaxDo74sFBPr7Mrh1DVrVvPp2NGWa3kLNgz8x9cFp5D3wL7nuCiauYw84cSLL687urqd5liP8+SfccQf8M38bO6o2J+yPFnj75oFFi+C776BsWafU4ZYCAsx85TZtTLOj556D8+ddXZVIpiiMZpGYmBi6du1K6dKlyZ07N4GBgXTp0oW9e/emuO5qwW3dunW0bNmSkiVLkjdvXsqWLUubNm34999/s/qvcN1GjRrF1KlTM3Xtrl27sCyLqKioFMd79erF0qVLs6I8EZFs7WohbNkyaNoUnsn7LTPjGuFdpBCsWEFEwUouCRFZqVXNQJb3uY+dQ5qzvM99mQqiyXLnsXhqbnuip23mszzdKPj5aM6UrwIzZmhkKR3X0zzrRtk2jBkDTe8+Ta9j/dhgVafCnl/N9izr1kHjxln2tbOVfPngm2+gf3/49FPzTeDwYVdXJXJVCqNZYOfOnQQHB7NhwwY+//xz/vnnH6ZOncrGjRupXbs2u3btytR9Dh06ROPGjfHz82Pu3Lls3ryZL7/8kooVK3LCwe/ennfgO2iFChXC39//hu7h5+dH0aJFHVOQiIiHSS+ELVwID4TYvJt3EKMPtMWqXRv++IOI0/mdHiKyi2bt/bl/61i6VF/BlmMl4LHHSHyguWkQIylca/OsG3XiBLRra7P05e/YlusWXjoxCK/27WDLFrNhro9PlnzdbMvLy+xzM22a6ZJdpw78/berqxLJkMJoFujevTteXl4sWrSIxo0bU7ZsWRo1asSiRYvw8vKie/fumbrP8uXLOXbsGJ999hm1atWifPnyNGjQgKFDh1KjRo10n5c82vree+9x00034efnx7PPPkt8/H/NLRo2bMgLL7xAr169KF68OPXq1QPg77//pnnz5hQoUIASJUrQvn179u/ff+l5iYmJ9OrVi8KFC1O4cGF69OhBYqrmD6lHe23bZvjw4VSuXJk8efJQunRpwsLCAKhQoQIAtWvXxrKsS9N9U0/TTUpKYuDAgZQpU4Y8efJQo0YNZs2adel88gjr999/T9OmTcmXLx9Vq1Zl4cKFmfq3FhHJ7ubMgUean2e6byd6Hu4HTz5ppjAWK+b0EJHdlC0LH6+5k297/cmrjOTsol9JqloNBg3SdMfLXEvzrBv111/QrsZGun3XmO94nGK3FIPffjPrnkuVcvjX8yjt25suu6dPw113wYIFrq5IJF0Kow529OhRFixYQPfu3cmXL1+Kc/ny5ePFF19k/vz5HDt27Kr3KlmyJElJScyYMQP7GqcMLV26lHXr1vHzzz/z/fff89NPP9G7d+8U10ydOhXbtvn111/54osv2LdvH/feey/Vq1dn1apVLFq0iFOnTtGiRQuSkpIAGD58OBMnTuSTTz7hjz/+IDExka+++irDWvr27cvAgQMJCwtj48aNfPfdd5QpUwaAVatWAbBgwQL27dvHzJkz07zHqFGjCA8P54MPPmD9+vU88sgjtG7dmrVr16a4rl+/frzyyiusW7eO2rVr065dO06dOnVN/3YiItnNjBnwbKtj/JL3AVoem2I2XvzyS8iTB3BuiHDVmsIb5eMDg8Nz8cC8V6lbaDOzLjwE/frB7bebX+zl2joY34Avx51k8R29mLX7duoXWAvjxmFFRZmOsZI5derAqlVQoYJp8PTRR5p+Lm4p+22+1KMHpAogWe7222HkyExdum3bNmzbpkqVKmmer1q1KrZts23bNu68yt5udevWpW/fvnTo0IHu3btTu3ZtGjZsyJNPPkm5cuUyfK63tzefffYZfn5+VK9enQ8++IDOnTszePBg8ufPD5hRyeHDh196Tv/+/bntttv44IMPLh374osvKFKkCFFRUdx5552MHDmSN954g8cffxwwITEyMjLdOk6dOsWIESMYOXIknTp1AqBSpUrcddddABQvXhyAokWLUrJkyXTvM2zYMHr16sUTTzwBwLvvvsuyZcsYNmxYivWpPXv25OGHHwZg0KBBfPHFF6xdu5b69etn+O8lIpJdTZ0K7zyznai8zSl7bqc58OSTKa4J8PclJo3g6egQkbymMHkUNnk6MHBN6zldqVkzuHV9IO3bf8fEX+fx+Z6XKN6wIXToYPZ1vPhzKycKDQlK8d8Xrt7B+FqcOW3zRfNveHjp65RiH2efeo58IwZBsWIOuX+OU7asGU1+6il45RXYtAlGjdL0ZnErGhnNIlY6+1slj3Cmdz61999/n/379zNhwgRq1KjB5MmTqVq1Kj///HOGz7v11lvx8/O79Pldd93F+fPn2b59+6VjtWrVSvGc1atXs2zZMvz8/C49kkcwt2/fzvHjx9m3b9+lIAng5eVFnTp10q3j77//5ty5czS+gQYDJ06cIDY29tJU4mT169fn71RrIW699dZLfw4ICADg4EFtbC4inmnSJPj46eVE5apDWd/DWIsWXRFE4fq2QbkenjIdODAQFi+GWm8+SPlTG/ikWF/sr6aZ7TPGjjX7OuZAN9LB+Gp2zfubDTc1ptvS9lCqFPbyP8j35QQF0Rvl5wczZ5r9cMaPN++2ZGJ2noizZL+R0UyOULpK5cqVsSyLjRs30qpVqyvOb9q0CcuyqFixYqbvWbRoUR577DEee+wxBg8eTM2aNRk4cOANBTzg0ghpsqSkJJo3b86wYcOuuPamm266NFX3Wlzr9OKMpBXgUx/zuezdvuRz11O3iIi7GzMGfn95Gr9Yz+JdrhzW/HlQqVKa1yaHhfDILcTGxRPg70toSJDDRyudOR04q+XKBQMHwr335uOpp97nE+8nmVPyZQJeegk++QRGjzYbueYwrWoGOvZ1c/Ikm58aSMXZIyhkFWDjS+OpNrILeHtf/bmSOV5eMGQIVKkCXbpA3bpmkXnlyq6uTEQjo45WpEgRQkJCGDduHGfOnElx7syZM4wdO5ZmzZpRpEiR67p/7ty5qVix4lXXQa5fv57Tp09f+nzFihWXnpueO+64g40bN1KuXDkqVaqU4lGgQAEKFSpEqVKlWLFixaXn2LZ9ad1nWqpWrUqePHnSHcnNnTs3wBVNkC5XsGBBAgIC+O2331Ic/+2336hatWq6zxMR8VSffWpz4OWBTONJvOvfhdfKFekG0WQ3sg1KZjlrTaEzNW1qVgf5312VwE2LGFZ3BolxJ6BRI2jbFnbvdnWJ2ZNtc2Hat8SVqsIts8OZX7wDZ9ZsodpH3RREs0qHDmbI/+hRs6Z08WJXVySiMJoVxowZw4ULF2jSpAmLFy9mz549LFmyhKZNm2LbNmPGjElx/YkTJ1i7dm2Kx65du5gzZw5PPfUUc+bMYevWrWzZsoVhw4Yxb948HnnkkQxruHDhAp06dWLjxo0sXLiQPn360KVLlytGQy/XvXt3jh8/Ttu2bVm5ciU7duxg0aJFdO3alZMnTwLw6quvMnToUGbMmMGWLVvo0aMH+/btS/eeBQoU4NVXXyUsLIzPPvuM7du3s2rVKsaPHw9AiRIl8PX1JTIykgMHDnD8+PE07xMaGsqwYcOYPn06W7dupX///vz666+8/vrrGf47iIh4mp8X2Rx5rjcD6U/SU8/gtegnuM43OB3NWdOBna1UKbNtzpAhFn1Xt6Fywib+efodmD0bbrnFDKHGZ7/RX5fZvJn4e5qS68m27DhdglFtf+eBvZMIvD3nrsd1mvr1TWOjgAAICTGj/CIupDCaBSpWrEhUVBTVqlXj6aef5uabb+aJJ56gSpUq/Pnnn5e2M0n266+/UrNmzRSPXr16UbVqVfz8/OjVqxc1a9bkzjvvZOrUqQwbNoy+fftmWEODBg2oVq0ajRo14pFHHuG+++5j6NChGT4nICCA5cuX4+XlxQMPPEC1atXo3r07efLkIc/Fjoyvv/46zz77LM899xx16tQhKSmJJ9NYn3S5wYMH07t3bwYOHEiVKlVo06YNe/fuBSBXrlyMHj2aSZMmERAQQMuWLdO8xyuvvEJoaChvvPEG1atX54cffuD777/n9ttvz/Bri4h4kr//hjUP9aeXHc65517E64spcHGGiTvIyjWFrubtbZbdrVwJeQv7UvnL/gxot5kLzR6C/v2halX44Qd1LM3IqVPQpw9JNW7l3O+reT3vWHZ+8yevfn2XO72MPV+FCvD773D//dCtGwwe7OqKJAezHLmm73oEBwfbUVFRaZ7btGlTul1pJX0dO3bk8OHDzJkzx9WluD29xkQku9i/H766ZSCvH+/PqXbP4ffVJ2YtmDhdfDz06WOWjVatChGv/kLlj16BDRugSRPTsbRqVSKiY7J8nW62kJQE332H3asX1t69fMqzfFl1CBNnlbja7HLJSomJZuruV1+ZniyvvurqisSDWZa12rbt4NTH9VNMRETEzZ05A9/VHsrrx/tzuHkHBVEX8/U1eTMy0jQmrfZSI4a2jyZp9BhYvRpuvZV/nuzCoGl/EBMXj81/29xkxb6rbru3a2IifPMN3HortGvH9mNFuZvl/PHcp8yLUhB1OW9vmDIFHnnEbJ04aZKrK5IcSD/JRERE3FhiInxddyQv7+3N3nvaU2zWZAVRN3H//bB+PbRsCb375aLhd93ZvWgrdOnCzdMnM2/cc7RdF4l3kmnSlxXb3CTv7eqM0JtpiYkwbRrUqAHt2nH8WCLdCnzF7Ymr6fb53UycaAK9uIFcuWD6dHjgAeja1YySijiRfpp5oClTpmiKroiIk2XV6FTE/ePotL4n225rQ+nFX6jTqJspWhS+/RY+/9x03a3esBhf3DWeh58ZyY4igXyw4COWftKFzn9G4HfujMO3uXGrvV0vXICpU6FaNXjySc5f8OL9Gl9TJHYDa6s+wcoob555xvllyVXkyWP2Im3Y0Ezb/eEHV1ckOYjCqIiIyA3KqtGpX56cRJvF3VlfoQWV/5xuRjHE7VgWPPMMrFsHt91mfp9fubodbR75kC6t3ySmUAneWjyJP8Z1YNDyKbBrl8O+tlvs7XrhgknjVavC00+T5JObme2/w3/3Xwz9ty0fjfVm+XKTUcVN+fqa7tB33mm2LJo/39UVSQ7h9mHU1Q2WxHPptSUijpIVo1PrXvucBtO6ElW8GVU3fAs+PjdapmSxChVgyRLTnDRuc3H2f9aA2bkeou0TQ3j4mREsqVyHtisioGJFePxxuGzf7uvl0r1dExLgs8/M9jYdO0L+/Gx6fya3Ja2lzfRHaf6wF5s2wYsvXvuAvtuug/Vkfn4wb56ZXt26tXkxi2Qxtw6jPj4+xGvfLski8fHx+OiXOxFxAEePTu0aPJ3qIzqxqkBjbvl7Jt758txIeeJE3t6m0+6qlRYli3tz8Ns6HJpVk5hidUj8cipeO3dCaKjZuPSuu+Duu+G778zo4nVwyd6u58+bZjdBQdCpExQqxKmvZtHtzjVU7fcIJ0558eOP5q8VEHDtt3fLdbA5hb+/6cxVsSI89BD88YerKxIP59ZhtESJEsTExHDmzBmNYonD2LbNmTNniImJoUSJEq4uR0Q8gCNHp45M+J7SfZ9mVZ57KBc9C79ieW+0PHGBO+6ArRt96NcPrD0B/DXqLuaMDWR3UmkYMgT27IGPPoKDB80oaaVKMGIEc37dfE0jgk7b2/X8eVi0yHRdrVQJunSBokWxZ//IN6FRVHqtBRMnWbz+OmzcaHLM9XKrdbA5UbFi5s2SUqWgWTNYs8bVFYkHc+t9RgFOnDjBwYMHSUhIcGJV4ul8fHwoUaIEBQsWdHUpIuIBkkdyLv8F2tfH+5pDwZlvZuPTrg1R3nUouHwB1er4ZUW54mQHDpipu+PHm89feAH69oUSJTCdZ3/8EUaMgGXLOJXbl69vvZ/ZVRuwqUQFcuXNmzXhMrOFz5sHc+aYcHLyJOTNC/fdBy+9xM6gB3ixu8WCBRAcDBMmQM2aN/5lK/SZS1q/nVrAziHNb/wLSObs3g333AOnT8PSpVr0KzckvX1G3T6MioiIZAcR0TGER24hNi6eAH9fQkOCrilAJM6ZT1KLlkTbNTnx/UKatNabZZ7m33/h3XfN1o6+vtCzJ7z+upkZCdD5pfE8/Ms3NN/8Gz5JiZzNlZv1N1Vi283VeaJHWzOtt1SprCswKQmio2HuXBNA//zTHA8MNEOdDz0E991Hgk8+RoyAAQPMtOT334fu3R3X6LnekMXEpDHFPdDfl+V97nPMF5HM2b7dBFLbhmXLoHJlV1ck2ZTCqIiIiJuyFy7iQrOH+CuxGus/XETHnoVdXZJkoS1boH9/syVM4cLQuze8/DJUe9eMCBY/dYzaezdyR8wm7ojdTLUD28mTeHFNablyJpTWrWs+3n475M59/cWcOmWm386ZY0ZB9+0z7YHr1oXmzU0AvfVWsCxOnDA1jxoFGzZAq1YwejSUKZPxl7jWN2ocNdNAHOTvv6FBA/MOyq+/mtegyDVSGBUREXFHy5aR0OQB/k6oTMTLi3l7dFFXVyROEh0Nb75pMmDJkuB75xYS/7cdyzvl72bl/bxZ0riQ6b77xx/msXevOZk3L9SqZcJjhQpw5ox5nD595Z/TOnbkiOmKW7AgPPCACaDNmkHx4sB/A2KffgozZpinVKkCgwaZMHo11xssb3SmgTjY2rXQqJHZWHfZsuvrTCU5msKoiIiIu4mOJqF+Q/45E8Dwh5cyIaIEXm7dWlCywm+/mTWkv/4KPoXOULD+VvJXjcXystMPbnv3mlCaHFBXrzZNhpLlzQv58plH/vwpP17+52LF4P77oV69FNsH7dljtg6dMsXM1CxQANq3N81z77zTDJ5mhqbcepCVK6FJEzMUvnTppTcsRDIjvTCq3bNFRERcYetWzjRsyqEz/jQv8R0lav/N7HUJGv3JgerXN7/bR0ZC95652DH3do4urIZ/heM0etCHMkmFSExMtSazdGl47DHzADh3DuLiTMj09b2uBZxnz8KsWWYUdOFCMyraqJFZG9q6tcmu18rR2x6JC9WpY6ZzN2tm3sBYvNjMMxe5AXr/VURExNn27uXUvY05ddKLZvl+5Pyj+9l/+rT2UszBLMvMkt22MTezZkHXZ30obhdjyohCBAdDkSLw8MMwfLjZaSMxMdUN8uSBm24CP79rCqK2bQZVX3rJzLxs1w42b4a33oIdO0zeeOqp6wui4Nhtj8QNNGgAP/xg1pE2a2Y6LIvcAE3TFRERcaYjR0i6515Obd5NI+tnYp86T55Sxy+d1vRFuVxsrBk1XbLEPLZuNccLFYJ774WGDc3o5a23psygCQlmoPTo0Ywf69ebR548ZvSzUyezc4ujpourGZGHmjUL2rQxw/rz5l3/uxWSY2jNqIiIiKudOgVNmpAQtZamiQvY8FAh/KrFprhEeylKRmJiUobTbdvMcX9/KF8ejh0zITOjASvLMtcXKWJ2bWnXzjyyasalmhF5qOnT4cknzZB+RMSNdXUWj6cwKiIi4krnzsHDD5P082IeSfqe1Q1uJ1fdDVdcppFRuRZ795pw+ssvsH+/CZhXexQq5Lg9QSWHmzQJunQx3a2mTnXckLp4HDUwEhERcZXERHj6aVi4kOe8ppDYvCWjBsTw5izvK6YvhoYEubBQyW5KlzaDU08+6epKJEd67jkzFN+7txlaHzMm862WRVAYFRERyVq2DS++CN99x1v5hrOiXAdWTIOCBQPx9kbTF0Uke3vjDbNf7dChZh/Sd991dUWSjWQqjFqWdS/QC6gFBADP2rY95bLzFvA20BUoDKwEutu2vdHRBYuIiGQr/frBhAlMKN6XsRdeY9VsKFjQnGpVM1DhU0SyvyFDzAjpwIFmLniPHq6uSLKJzI6M+gEbgC8uPlJ7A3gd6AhsAfoDCy3LCrJtWz2fRUQkZxo+HAYPZkH553lx93ssiIRKlVxdlFwrNeARuQrLgo8/Nh20evY0gfSZZ1xdlWQDmQqjtm3PA+YBWJY15fJzF0dFewBDbNv+/uKxDsBB4AngE8eVKyIikk1MmQK9erG+ymM03zSWkaMtmjRxdVFyrVJvTRITF0/YzPUACqQil/P2hq++guPHzR5BhQpBy5aurkrcnCNaXlUASgI/JR+wbTseWAbc7YD7i4iIZC+zZsFzz7G/RlOCN31Jp+e8eeklVxcl1yM8ckuKJlMA8QmJhEducVFFIm4sTx744QcIDoa2bc3+QyIZcEQYLXnx44FUxw9cdi4Fy7K6WpYVZVlW1KFDhxxQgoiIiJtYsgTatuVUlWBqbJvJnfXzMHasGkxmV7Fx8dd0XCTH8/ODuXOhYkVo0QJWr3Z1ReLGHLkZUOoNS600jpkLbXuCbdvBtm0HFy9e3IEliIiIuNDq1dCiBQnlKnLXkbn4Fvfj+++1F3x2FuDve03Hs4uI6BjqDVlMhT5zqTdkMRHRMa4uSTxJ0aLw00/m4wMPwObNrq5I3JQjwuj+ix9Tj4KW4MrRUhEREc/0119w//3YhYvQxu8ndhwvyuzZUKKEqwuTGxEaEoSvj3eKY9l9P9jkdbAxcfHY/LcOVoFUHCow0ARSLy9o2hR273Z1ReKGHBFGd2ICadPkA5Zl5QXuAX53wP1FRETc28aN0Lgxdr589K69mB/XBDJlCtx+u6sLkxvVqmYgg1vXINDfFwsI9PdlcOsa2bp5kdbBitNUrgyRkXDyJNx/P2h5nqSS2X1G/YDkZvReQFnLsm4Hjtq2vduyrJFAP8uyNgNbgTeBU8A0h1csIiLiTjZvhsaNwceHz55aTPiQm+nfHx57zNWFiaN42n6wWgcrTnX77TBnjgmjDzwAv/zy32bLkuNldmQ0GIi++PAF3rn453cvnh8KfAiMBaKAUsD92mNUREQ82rZtcN99AET2+YXnPqhM69bw9tsurkskA566DlbcWP36MGOGWc7QogXE640PMTIVRm3bXmLbtpXGo+PF87Zt2wNs2y5l23Ze27Yb2La9IUsrFxERcaXt26FRI7hwgTXDFtPyjSDq1oWpU80SKRF35YnrYCUbePBB+PxzWLYM2rWDCxdcXZG4Af24FBERuVa7dpkR0bNn2TnpZ5q+WpWyZWH2bPDV4JK4OU9cByvZxBNPwEcfmW+WnTtDUpKrKxIXy9SaUREREblo924zInryJEe+W0zjLjXw9ob586FYMVcXJ5I5nrYOVrKR7t3h6FHo3x8KF4YRI7QRcw6mMCoiIpJZe/eaIHrsGGd+/JmQnrdz4AAsWWL2dxcRkUx48004cgRGjYJChWDAAAXSHEphVEREJDNiY83U3EOHuDB/IY8NqkV0NMyaBbVru7o4EZFsxLLgww/hxAl49104cwaGDlUgzYEURkVExKNFRMcQHrmF2Lh4Avx9CQ0JuvbpiQcOmO1b9u3DXhDJi5/XYd48+PhjeOihrKlbRMSjeXnBpElmof2wYSaYjhsH3t5Xf654DIVRERHxWBHRMYTNXE98QiIAMXHxhM1cD5D5QHrokBkR3b0bFixg0JK7mTgRwsLg+eezqnIRkRzAywvGjDFTdQcPNoH0iy/Ax8fVlYmTqJuuiIh4rPDILZeCaLL4hETCI7dk7gZHjkCTJrBzJ8yZwxc77+HNN+Gpp+D997OgYBGRnMayYNAgGDIEvv4aHnlE+5DmIAqjIiLisWLj0v6FJr3jKRw7Bk2bwtatMHs2ixIb0bmzGSSdPFlLm0REHKp3bzNNd948syfpyZOurkicQGFUREQ8VoB/2pt+pnf8kv37zYjoxo0QEcFfJZrQujXccgvMnAm5c2dBsSIiOd0LL8CXX8Kvv5p1+keOuLoiyWIKoyIi4rFCQ4Lw9UnZDMPXx5vQkKD0n7R+Pdx5J2zeDD/8wJ6qITz4IBQsaN6wL1Qoi4sWEcnJnnzSvOv311/QsCHs2+fqiiQLKYyKiIjHalUzkMGtaxDo74sFBPr7Mrh1jfSbF82fD/XqQWIi/PorcXc/eGm22Lx5UKaMU8sXEcmZWrSAuXPNev177oFdu1xdkWQRy7ZtlxYQHBxsR0VFubQGERERxo6FV16B226DH3/kfPFAHnjAzBZbsMDMGBMRESdasQKaNYP8+WHRIrNWQrIly7JW27YdnPq4RkZFRCRnS0yEHj3gpZegeXNYtowLNwXSsSP88gt8+qmCqIiIS9StC0uXwoULZoR0zRpXVyQOpjAqIiI518mT0KoVjBoFPXvCDz9wzsePtm1h+nSz7d3TT7u6SBGRHOzWW80UlXz5oFEj+O03V1ckDqQwKiIi2UJEdAz1hiymQp+51BuymIjomBu74d695p32+fPNdgIffsipeG+aNze9M0aMgD59HFO7iIjcgMqVTQgtWRLuvx8iI11dkTiIwqiIiLi9iOgYwmauJyYuHhuIiYsnbOb66w+kq1ebjrk7dpgmGS+8wJEjZjeXJUtgyhQzc1dERNxEmTJmhDQoCB5+GL7/3tUViQMojIqIiNsLj9xCfEJiimPxCYmER2659pvNmgX33gs+PvD77xASQmwsNGgAa9ea3286dHBM3SIi4kAlSpjF/LVrw+OPw+efu7oiuUEKoyIi4vZi4+Kv6XiabBuGD4dHHoHq1WHlSqhenX/+Mbu5/PuvmbHbsqWDihYREcfz94effjKd5Tp2hJEjzfd3yZYURkVExO0F+Pte0/ErJCTACy9Ar17w6KNmLm7Jkvz1F9Svb/oYLV5semOIiIiby58ffvzRvLnYs6d5F3HfPldXJddBYVRERNxeaEgQvj7eKY75+ngTGhJ09Sfv22e2bPnkEwgLg6+/Bl9ffv/dTM3NlcssQ6pdO4uKFxERx8uTB777Dj78EBYuhGrVYNo0jZJmMwqjIiLi9lrVDGRw6xoE+vtiAYH+vgxuXYNWNQPTf1JiIowZYzZJX7rUbBg6aBB4eREZCU2bQvHisHw5VKnitL+KiIg4ire3GRldu9Y0NnrySTP75eBBV1cmmWTZLn73IDg42I6KinJpDSIi4mGioqBbN9M1t2lTGDvWbA0AfPstPPWUeRN9wQK46SYX1yoiIjcuMdH0BXjrLShYEMaPN8FU3IJlWatt2w5OfVwjoyIi4jni4uCll8y2LTExZkpuZOSlIDpxIrRrB3XqmGWj6QVRh+9pKiIiWcvbG954A9asgXLl4LHHoH17OHLE1ZVJBhRGRUQk+7NtmD7dTMkdP94E0s2boW1bsCwAPvgAunaFZs1MPi1UKO1bOXxPUxERcZ5q1eCPP2DgQLNXV7VqMHu2q6uSdCiMioi4OY3SXcXWrWYq7hNPmE3RV62C0aMvpc2kJOjdG/r0MW+SR0RAvnzp386he5qKiIjz+fjAm2/Cn39CyZKm2+4zz8CxY66uTFJRGBURcWMapcvA2bPw9ttQo4ZZIzp2LKxYAbVqXbpk40azdcvQoWZnl6lTze8oGXHInqYiIuJ6t91m3qDs39902q1e3WwoLW4jl6sLEBHxdLYN//5rfh6uXGlyk22bTq7pPYoVMx8zGqXLsJOsp4uMhO7dYft2MyI6fLh59/uic+fg/fdhyBDTx+LLL02TxYszdjMU4O9LTBrBM9N7moqIiPvInRveeQdatIAOHeDBB6FzZ/NzI731GuI0CqMiIg527JiZGZQcPlet+q/LfJ48ULOm+bh5s9nf8sgRM5U0LVbue/DOdx6vfOfxKXyavGWPkLf8YWLJoaN069fDe++Zlrj/+x8sWgSNG6e45NdfoUsX2LLFdM398EMT7DMrNCSIsJnrU7wJkOk9TUVExD3VqmU6rL/zjmki8NNPZsuvJk1cXVmOpjAqInIDzp+Hdev+C50rV5oljMmqVDENc+6803RwrVHDvEl7ucREE2APHfrvcfiw+Thu/n6OH/Mi6Uxu4ncW5/TG0gDkLXqGbnEmhzVqZEZSPdbOnaY50fTpsGGDSfLvvmu6JubJc+my48fN2tBPPoHy5c22LSEh1/7lkkecwyO3EBsXT4C/L6EhQTl7JFpExBPkyWP2m27Z0oySNm0Kzz8PvXpBpUquri5H0j6jIiLX4cwZCA83j9OnzbGbbjKBM/kRHJzODKBz58xw6MmTV33s2rWfTVtjyXv2NBcsb2J8SrEvoRznvasRtftm/jkbSAylqVSzAI0bm3B6zz2QP79T/zkcb/9+M/o5fbpZBwpm8Wf79qZdf6qhzpkzTQPdAwfM/ufvvOMB/wYiIpJ14uPNnqQjRpjpSXfeaZZ9tG2bYtmHOEZ6+4wqjIqIXAPbNhkpNBT27DH7abdta36GlSmTwZrE2FiYM8e0l1+0yATSjFgW+PlBgQKczO3L3oRccP4cpU4fxf/MiSsuP+1dgN2JpdlLILFWaazSgRS9NZAKDcoS9FBlvCtVuHrnHleLizOpcvp0WLzY/HJw++0mgLZrB2XLXvGU2FgTQn/4wVw6caJ5E0BERCRT9u41e1JPmwbR0eDlBffdZ4Jp69ZaV+ogCqMiIjdozRp49VX47TcTfEaNgnvvTedi2zbrG2fPNo8//zTHK1SAhx+GqlWhQIH0H/nymR+IaYmPNyksJsb8EL348cLuGE5tNp/7ndxHLv5b85hoeXOm5M343hZErqr/M+stkx8BAZdSdER0jHOnp8bHm5A+bRrMm2fmPVeqZAJo+/ZmnnMakpJgwgQzLff8eTMS2rOn++dtERFxY5s2mTdEp00zDfLy5IHmzU0wbd4c8uZ1dYXZlsKoiMh1OnDAbFc2ebJZm/n++9CpE3h7p7rw/HlYuvS/ALp7twl5deqYLn4tWpgQmpmWrjcqMZFjmw+weua/bJ+/lVNrtlL23FaCrK38z9pG3qTLGiDlzw+VK7O3RFlmn8nPtkKlOJKvEHG+BYgv4E+Px+vQ/J4q11d3YqL5B4yJ+S9AX/5YsQJOnYJSpczoZ/v2Zmgzg6+1aRN07WreFGjc2KwRrVjxOv6NRERE0mLb5k3kadPMqOmBA6Y1e+vWJpjed18avwRIRhRGRUSu0fnzMHo0DBxo1oi++qpZXpJixs7Ro2bPstmzzceTJ8HX1zRFaNHCvJPqBmtPEhJg+XL48UeYHZHE2R0xBLGFxqW30ihwK1VybeXkpnWUOrYfbzuN1r7e3lC4MBQpAkWLpvyY/Dh58sqwuX+/CaSXy5XLhM/AQNPRqX17M8ScwQ/2EyfMP29EhJnJmz+/6ZLboYNzsr2IiORQFy7AkiUmmH7/vfmBdNNN5g3UJ56A2rWd+4Po3DkzK2rvXrNeaM8e8wvL2287r4broDAqIpJJtg1z58Jrr8G2bSZPDh8OQZfv7BEfb4ZIw8PND4GSJc3024cfNsN1+fK5rP6rsW2zrczs2Sac/v67OeZdIJ6CN+8lqOQ6ShbcQ5GkY/jHn6Jw/AnevLukCd5HjpiPl//55MlL9z6ZNz92QCAFK5YzYTOtR4kS6U9Bvsz+/abGiAj4+Wfzz1y8uFmnO2CAuY2IiIjTxMebJSXTppklJufPw803Q7Vq4O//36Nw4ZSfX36sYMH0fwYmJJg3cvfsSRk2L/88ea+4y5UpY2ZjuTGFURGRTNi0yaw9jIw04XPECLM1SwoLFkD37rBjh9nI8uWXzdTSTAQsd3TwoPnZ+vqwgxzbWgQ7wez6lavwKXKXOk6JCmf4ondl7rjjyg61EdEx9P9uDblPHueMT17O5PbF18ebwa1rXNda023bTPj84Qczg9e2zc/5Rx6BVq3grrs0M0pERNxAXJz5YfX992YZSlyc2aft+HHzwys9lmUC6eVB9cwZEzb377/yuYUKmbBZurT5mNafs0H7eIVREZEMnDgB/fvDmDGmie2AASZvpmiIExMDPXrAjBkQFMRvr71L76PFPGYvyojoGHp/u5G4HYU4v8+fc/sLkbC/EBdO+gIma1epYmYkBQebR4+flrAveW+bywT6+7K8z31X/Zq2DVFRJoBGRMDff5vjd9xhwmerVlC9uqbiiohINpGUZGYMJYfTuLiUj9THjh0zy3vSC5wFCrju7+JACqMiIuk4cgRCQkxH9y5dzBrRFNtYXrgAY8eaLkYXLsCbbzK7SXt6z9lKfMJ/6yEzMyLo9G611yit+uqUDGT1ahMa//zTPA4duvgEryRyFz9J7pLH8cp7HpK8sJMssC2evLM8Fy6Q4pGYmPLzv/4yGd/b2ywbbdXK7EVerpwr/xVERETEkRRGRUTSsH+/6TW0bZuZadO8eaoLVq2Cbt1MUn3gARNKb76ZekMWExMXf8X9MhoRjIiOIWzm+msOsO7Gts1soqgoeHn0Lo7sys/5/YVISvDG8rLBsvH2hiIFfMiVy/Qr8vbm0p8v/7xMGRNAmzc3/ZBERETE86QXRnO5ohgREXewZ4/pNRQTYxoWNW582cm4OOjXD8aPN51fv/sO2rS5NF80No0gmtFxgPDILSmCKEB8QiLhkVuyVRi1rP9mEVnlfQibuTrbB2wRERFxPoVREcmRduww24QdOwY//QT16l08YdtmT7GePc1c1FdegXffNc0GLhPg75vmyGiAv2+6X/N6Aqy7Sw6c7jz1WERERNyTwqiI5DibN5tR0LNnYfFiqFXr4omtW+HFF80+IrVrm40ta9ZM8x6hIUFpTrkNDQlK83q4vgCbHbSqGXjN4dPd186KiIhI1sue+xCIiFyndetMo5zERLOH9aUgOmkS1KhhFkKOGwd//JFuEAUTwAa3rkGgvy8WZq3o1aamhoYE4euTcl+SqwVYT5S8djYmLh4biImLJ2zmeiKiY1xdmoiIiDiRRkZFJMdYtcr0IMqf3wx+/u9/F08MGwahoaal7pQpULJkpu53rSOCmtJqeMraWREREbkxCqMikiP8+qvp2FqsmJmaW748Zn1o//7w3nvw2GMwdSrkzp2ldVzPlFZP44lrZ0VEROTaaZquiHi8hQvNoGdgoAml5ctjNqV+9VUTRDt3hunTszyIipHeGtnsvnZWREREro3DwqhlWQUsyxppWda/lmXFW5b1u2VZtR11fxGR6/Hjj/DQQ1C5MixdagIpFy6YAPrRR6Zr7sSJZuNLcQqtnRURERFw7MjoJCAE6ADUAH4CFlmWlbPno4mIy3z7LbRuDbfdBr/8AiVKAOfOQbt2Zm3ogAEwfPilvUPFOa6n+ZOIiIh4Hsu27Ru/iWX5AieBNrZtz7rs+Gpgvm3bb6b33ODgYDsqKuqGaxARudznn0OnTmb/0DlzLm4TeuaMSaeRkTBiBPTo4eoyRURERDyeZVmrbdsOTn3cUQ2McgHewNlUx+OB+g76GiIimTJ5Mjz3HDRtCj/8YLrncvy4ma/7++/mgk6dLl2vPS9FREREnM8hYdS27ZOWZf0BvGlZ1gZgP9AeuAv4J/X1lmV1BboClC1b1hEliIgAJmt262a2cPnhB8ibFzh0yBz46y/TqOjxxy9dn7znZfJWI8l7XgIKpCIiIiJZyJFrRp8GkoC9wDngFWA6kJj6Qtu2J9i2HWzbdnDx4sUdWIKI5GQHD8LDrRLJVfAMG4MiaTxyMQsio+Dee+Hvv2HWrBRBFDLe81JEREREso7Dwqht29tt224A+AFlbNu+E/ABdjrqa4iIpCcxEZo8dJZjx6Bwi9V45b2A984dVG/7EAl79sKCBfDgg1c8T3teioiIiLiGw/cZtW37tG3b+yzLKozprjvras8REblR/fvD+j/zUuT+DeS+6QT/O7SL76b1Jv+5M3TrOBQaNEjzeZ6652VEdAz1hiymQp+51BuymIjoGFeXJCIiIpKCI/cZDbEsq5llWRUsy2oK/AJsAT5z1NcQEUnLjz/CoEHgd9tu/Grspca+bXwzLQyAx58YwmK/9Neme+Kel8nrYGPi4rH5bx2sAqmIiIi4E0d10wUoBAwGSgNHge+BfrZtJzjwa4hIDpKZLrc7dsDTT8Mdd0DuVjuw9h3ksxkDOJUnH0+0e589/iUJzGCUM/l+ntRNN6N1sNn57yUiIiKexWFh1Lbtb4FvHXU/EcnZMtPlNj4eHn0ULAtmzIAN+0pzc5su5L6QQNsnhrDHv2SmRjlb1Qz0qJCmdbAiIiKSHTh8zaiIiCNkpsvtyy9DdDRMnQoVyts8PKY/VQ9sZ2C7vuwoWoZAf18Gt67hUUEzMzx1HayIiIh4FkdO0xURcZirje5Nnmweb74JzZsDQ8Nh+nSsQYMIDwsj3Im1upvQkKAUo8qQ/dfBioiIiOfRyKiIuKWMRveio6F7d2jSBAYMAObNgz59oG1b8zGHa1UzkMGtaxDo74sFOXaEWERERNybZdu2SwsIDg62o6KiXFqDiLif1GtGwYzu9WtyK/2fDSAhAdasgeJHt0CdOlChAvz2G+TP78KqRURERCQ1y7JW27YdnPq4pumKiFtKq8vt602DmPJOAHv3wrJlUDz3cWjZEnx8ICJCQVREREQkG1EYFRG3lbrL7eDBZk/Rjz6CurUToeWTsH07LFoE5cq5sFIRERERuVYKoyKSLfz8s2lW1K6dWS9Kv7dg7lwYNw4aNHB1eSIiIiJyjdTASETcXkwMtG8PQUEwcSJY335jhkm7doVu3VxdnoiIiIhcB4VREXFr58/DY49BfDzMnAl+26Lh2WehXj0zX9eyXF2iiIiIiFwHTdMVEbf2xhvwxx/wzTdwS5GDULsVFC0K338PuXO7ujwRERERuU4KoyLithYvhlGj4JVX4PFW56HJo3DwoNnC5aabXF2eiIiIiNwAhVERcUunT0OXLlCpEgwZAvToAb/+Cl99BbVqubo8EREREblBCqMi4pbeegt27IClS8H3i09g/HgzZ/eJJ1xdmoiIiIg4gBoYiYjbWbECRo6EF16Ae61f4aWX4IEHYNAgV5cmIiIiIg6ikVERcSvnzkHnzlC6NHzw+kGo9xhUqADTp4O3t6vLExEREREHURgVEbcyaBD8/TfMnWNT4PWucOwYLFwI/v6uLk1EREREHEhhVETcxl9/mTD61FPw4MEpMGsWDBsGNWq4ujQRERERcTCFURFxCxcumOm5hQvDqB47oeEr0LAh9Ozp6tJEREREJAsojIqIWxgxAqKi4NvpiRTp2QG8vGDKFPNRRERERDyOwqiIOEVEdAzhkVuIjYsnwN+X0JAgWtUMBGDbNujfH1q2hEf/HW72E/38cyhXzsVVi4iIiEhWURgVkSwXER1D2Mz1xCckAhATF0/YzPUAtLgtkOeegzx5YEL3dVjN34Q2beDpp11ZsoiIiIhkMYVREcly4ZFbLgXRZPEJiYRHbmH/ykCWLYNPx5+jxOtPQ5Ei8PHHYFkuqlZEREREnEFhVESyXGxcfJrHd++2eWMQNG4MHbe/BevXw9y5UKyYkysUEREREWdTZxARyXIB/r5XHLNtOPXzbSQmwpfPLcUaPgyefx4efNAFFYqIiIiIsymMikiWCw0JwtfHO8WxhC2lidtajGH9T1CqTweoWNHsKSoiIiIiOYKm6YpIlkvumpvcTbdYroJsWVqdunXh+U2vwp49sHw5+Pm5uFIRERERcRaFURFxilY1Ay+F0nbt4K8z8G37H/B6dQq8+SbUrevaAkVERETEqTRNV0ScatYs+OYbGNJjP2UGdoVatcwmoyIiIiKSoyiMiojTxMXBCy/ArTVsXt3wHJw6BV9+CT4+ri5NRERERJxMYVREnKZXLzhwAGa3mITXvLnwwQdQpYqryxIRERERF1AYFRGnWLIEJk+GwZ3/odzInmZz0ZdecnVZIiIiIuIiCqMikuUSEqB7d6hY7gKv//WMmZY7ZQp46VuQiIiISE6l3wRFJMuNHAl//w0/3jMU75V/wNixULq0q8sSERERERdSGBWRLLV3L7zzDrx87zqqfP02tG0L7du7uiwRERERcTGFURHJUj17AhcuEH60ExQpYkZFLcvVZYmIiIiIi+VydQEi4rl++glmzICF948gz09r4LvvoGhRV5clIiIiIm5AI6MikiXOnTPNchuX3UbjZf2hVSto08bVZYmIiIiIm1AYFZEsMWwY/LMtiW8LdcHKk0fTc0VEREQkBU3TFRGH27UL3n8fxtecSJHopTBpEgQEuLosEREREXEjGhkVEYd79VUozV66/PMG3HcfdOrk6pJERERExM0ojIqIQ82ZA7Nn28wr/wJeFxJg4kRNzxURERGRK2iarog4THw8vPIKvBbwDZU2zYHhw+Hmm11dloiIiIi4IY2MiojDDBkCJ3YeZvDpl+HOO818XRERERGRNGhkVESuSUR0DOGRW4iNiyfA35fQkCBa1Qzkn3/ggw9gYfme5N4bZ5oWeXu7ulwRERERcVMKoyKSaRHRMYTNXE98QiIAMXHxhM1cj23DhH6BPOQ1j3t2TYX+/aFGDRdXKyIiIiLuTGFURDItPHLLpSCaLD4hkbAPjxK7oAB7/btBharQt6+LKhQRERGR7MIha0Yty/K2LGugZVk7Lcs6e/Hje5ZlKeyKeJDYuPgrjiWd9+af2ZWYUDSMAsf3wuTJkCePC6oTERERkezEUWGxN9Ad6ACsB24FPgfOAQMd9DVExMUC/H2JSRVIj/9eibtORtGOcaZhUd26LqpORERERLITR3XTvRv40bbtH23b3mXb9mxgNlDHQfcXETcQGhKEr89/TYkSDvtxblUA0/J3hvLl4b33XFeciIiIiGQrjgqjvwGNLMu6BcCyrKrAfcA8B91fRNxAq5qBDG5dg0B/X7Dh9JIaDPQZSJnT22DCBPDzc3WJIiIiIpJNWLZt3/hNLMsC3gPCgETM9N/3bdt+M53ruwJdAcqWLVvr33//veEaRMS5vv4aPmgfzWqv2nh1eAY+/dTVJYmIiIiIG7Isa7Vt28FXHHdQGG0HhAOhwEbgdmAUEGrb9uSMnhscHGxHRUXdcA0i4jwnTkC1oAssPH4nQQVjsTZtgsKFXV2WiIiIiLih9MKooxoYhQPDbNv++uLn6y3LKocZKc0wjIpI9vP22/Dk/uHcQjR8OUNBVERERESumaPCaD7M9NzLJeK4Naki4iaio2H+qK2s934bWraGNm1cXZKIiIiIZEOOCqM/An0sy9qJmaZbE3gN+MJB9xcRN5CUBC92S2KKd2e88/vCmDGuLklEREREsilHhdGXMfuJjgNKAPuAicC7Drq/iLiBiROh1qpx1OU3GPkZlCrl6pJEREREJJtySAOjG6EGRiLZw4ED0LTyLlaerk7epvWx5s8Hy3J1WSIiIiLi5tJrYKQ1nSKSKaG9bD481YXceS2sTz5REBURERGRG+Koaboi4sF++QVyTf2MJiyC8LFQrpyrSxIRERGRbE5hVEQydO4cvN0lljler5F49714d+vm6pJERERExANomq6IZGhYuE2v7d3I73Me788mg5e+bYiIiIjIjdNvlSKSrh07YOu7X9OCH/F+fyBUquTqkkRERETEQyiMikiabBv6dTnIhwkvc/72O6FHD1eXJCIiIiIeRGtGRSRNM2dCy8Wv4O99Au+pn4K3t6tLEhEREREPopFREbnCyZMwr2sE7fgG+veHatVcXZKIiIiIeBiFURG5wpDex3jv6Auc/t/teIf1dnU5IiIiIuKBNE1XRFJYuxYqj3+NEtYhvL+eBz4+ri5JRERERDyQRkZF5JKkJPisXSQdmcL5nr2hZk1XlyQiIiIiHkphVEQumTL6BK9v6cLxUrfg+/5bri5HRERERDyYpumKCAAHDoDduw+l2Ys1YznkzevqkkRERETEg2lkVEQAmPTMUjqfH09chx5Yd9/l6nJERERExMMpjIoIS+efoe1PnTnqfzNFxr3n6nJEREREJAfQNF2RHO7cOdj+ZH8asJ1z0xdDvnyuLklEREREcgCNjIrkcF+9spIOx0awu9nz5HmgkavLEREREZEcQmFUJAdbv+I09SZ04JhvAGW/HurqckREREQkB1EYFcmhzp+HLc1fozJbyTX1cyhY0NUliYiIiEgOojAqkkN998QPPHp0Attbh+Lf+j5XlyMiIiIiOYzCqEgOtObHGB74/jl2Fq1F5ekDXV2OiIiIiORA6qYrksOcOZXE2XYd8LXOknvBV5A7t6tLEhEREZEcSCOjIjnMwgeGc/eZn/m35ygKBAe5uhwRERERyaEURkVykFUfr6HZ8n6svbk1VYZ1dnU5IiIiIpKDKYyK5BDHY09T9OX2HM1Vgv8tnQiW5eqSRERERCQHUxgVySHWNupJhQvbODryS/KVLuLqckREREQkh1MYFckB/uz7Aw22TuTXum9QtXsjV5cjIiIiIqJuuiKeIiI6hvDILcTGxRPg70toSBCtagZydH0MFYc8x0bfWtT96V1XlykiIiIiAiiMiniEiOgYwmauJz4hEYCYuHjCZq7HTkyiUvOO3GyfxfvraeQpoG1cRERERMQ9aJquiAcIj9xyKYgmi09IZFvnkdQ4uJglrUdzS4v/uag6EREREZErKYyKeIDYuPgrjlXZ/i89/xrN4iJtCPm6kwuqEhERERFJn8KoiAcI8PdN8Xnec2cZNetD9lOSsvMnkMtH27iIiIiIiHtRGBXxAKEhQfj6eF/6vO9306iUsIOv242n0p3axkVERERE3I8aGIl4gFY1AwGzdrTa4mU8EzOTSSV78fpXD7m4MhERERGRtCmMiniIVjUDaVHU5tRbj7LGK5j7f30fr6vMfUhvOxgRERERkaymMCriKS5cIOa+pyly4Ry7Bn/FHZUy3sYlve1gAAVSEREREclyWjMq4glsm/2Pv0KZ7UuYePs4Hul99W1c0tsOJjxyS1ZVKSIiIiJyicKoiAc4+s5HlPxhPBP9Q3nm5w5YmWiem9Z2MBkdFxERERFxJIVRkWzuzIx5FHqnJ3NzteSe3wZTJJPNc1NvB3O14yIiIiIijqQwKpKNJa1bD+3bsY7byPv9V9xSzfvqT7oo9XYwAL4+3oSGBDm6TBERERGRK6iBkUh2deAAx+99mPgLfqx9ZzadWuS/pqdfvh2MuumKiIiIiLMpjIpkoSzbOuXsWQ7Va0X+Ewf5uPUy+rxV+rpu06pmoMKniIiIiLiEwqhIFsmyrVNsm0MtOlF8+woGVJ9Bv6+DM9WwSERERETEnWjNqEgWyaqtU+Jee5fiC6czrMggXlnaBh+fG7qdiIiIiIhLKIyKZJGs2Dol/tPp+I8cwDSfDjz0W59Md84VEREREXE3CqMiWcTRW6ck/b4C7y7P8iv3UHzmJ9xSRXNzRURERCT7ckgYtSxrl2VZdhqPuY64v0h25NCtU/79l9NNW7InKZDNg2bS9KE8DqpSRERERMQ1HNXAqDZw+W/dpYDVwLcOur9ItuOwrVNOnOBY/Yewzpxjatsl9O9TLAuqFRERERFxLoeEUdu2D13+uWVZnYETwHeOuL9IdnXDW6ckJnLsgfYU2LuJvrfP5/0vq6hzroiIiIh4BIdv7WJZlgV0Bqbatn3G0fcXyUlOdO1F4T/m0b/4ePr83FSdc0VERETEY2RFA6OmQAVgUnoXWJbV1bKsKMuyog4dOpTeZSI52tmRH1Pw05GMz/0qTyzrps65IiIiIuJRsiKMdgH+tG17bXoX2LY9wbbtYNu2g4sXL54FJYhkb0nTvsbntZeYS3MqRgznlltcXZGIiIiIiGM5NIxallUCaAlMdOR9RXKS86PGw5NP8Jtdjz1Dp3N/M++rP0lEREREJJtx9MhoR+Ac8LWD7yvi+Wyb030GkrvHi/zIw/w9fAHdQgu4uioRERERkSzhsAZGFxsXPQd8bdv2SUfdVyRHSEoi7tme+H8xmqlez+D3zWReeNTh/cVERERERNyGI3/bbQhUBp5y4D1FPF9CAodbdKLYgql8nLcHty0azl31smI5t4iIiIiI+3BYGLVt+xdAOyCKXIv4ePY3eJySf85heOH3ePiPvvwvSP8biYiIiIjn0zxAEVeJiyO2dgtK/vMbQ8qNp9OqbpQo4eqiREREREScQ2FUxAWS9h1g/+0hFDv4N+F3TOflZW3Jn9/VVYmIiIiIOI/CqIiTnd+6i6PBTSl0MpaPH/yR12eFkEv/J4qIiIhIDqNfgSVLRETHEB65hdi4eAL8fQkNCaJVzUBXl+VyJ//YwLmGIeQ+H8/33Rbx8ri7sLREVERERERyIIVRcbiI6BjCZq4nPiERgJi4eMJmrgfI1oH0RgP2gVkryNvmQc4n5mX1oGU8E1Y9C6sVEREREXFvCqPicOGRWy4F0WTxCYmER25xmzB6rcHyRgP2zo8juenF1uynFLFfLOThpys45i8iIiIiIpJNaTNDcbjYuPhrOu5sycEyJi4em/+CZUR0TLrPyShgZ8ROTOLPTuMIfOFhduaqzNmfl1NfQVRERERERGFUHC/A3/eajjvb9QTL6wnY0TN3srpIE2p/1p3VBRtRaM0Sqja66fqKFhERERHxMAqj4nChIUH4+ninOObr401oSJCLKkrpeoLltQTs3buS+LT2eCq3qUHQySiWd5xAnaMLKF3d/7rqFRERERHxRAqj4nCtagYyuHUNAv19sYBAf18Gt67hNutFr2fkNjMB+9QpGP7STnZUbEKnqBfZX+FuvDZuoN5nXfDyVstcEREREZHLqYGRZIlWNQPdJnymFhoSlKIZEVx95Db575JW06PERPj8syS2vv4xb554A+9cXhwZNIFKbzyH9m0REREREUmbwqjkOBkFy6s9L/U1S5ZA+Is7eX1TZzrxC3F17sfv24n4li2bVeWLiIiIiHgEhVHJkW505Paff+CNXkmUnPUx31pvkNvXC3vURPyf66zRUBERERGRTFAYFbkGR47AoEHw4+idTLI7cS9LSGwcgvfkCaDRUBERERGRTFMDI5Gr+PdfGD0aGjeGkiWSOPfhWP6iBvXzrYFJk/D+ab6CqIiIiIjINdLIqEgqtg3r1sGsWRARAWvXgheJvBI4k89LDKX0/ihoHAITJ0KZMllSQ0R0zDWvaRURERERyU4URkWAhAT49VcTQGfNMqOhlgWN6pxhYavPuHf1h+TeswMqVYIpU+CZZ7JsbWhEdEyKbr8xcfGEzVwPoEAqIiIiIh5DYVRyrJMnITLShM+5c+HYMcibF5o2hfd7HKLlnjH4fT4WVhyBunVhZDi0bAne3le/+Q0Ij9ySYtsZgPiERMIjtyiMioiIiIjHUBgVj3f4MGzeDJs2pXz8+685X6QItGhhcmZIxX/IN344hE2Bs2fNidBQqFfPaV1yY+Pir+m4iIiIiEh2pDAqHsG2Yc+eKwPnpk0mjCbz9YVbboG774bOneGee6B+fci1eiUMHQo//AA+PmYa7uuvm4udLMDfl5g0gmeAv6/TaxERERERySoKo5Jt2DYcOgTbtqV8bN1q9v08ffq/a4sUgSpVoFUr8zH5UbYseCX3kE5KMvNz7ws3C0b9/SEsDF5+GUqWdMHf0AgNCUqxZhTA18eb0JAgl9UkIiIiIuJoCqPido4eNQEzdejctg1OnPjvuly5oEIFqFwZGjb8L3DecgsUL57BrNoNG+Dbb2H6dJNiy5aFESPMUGmBAs74K2YoeV2ouumKiIiIiCdTGBWXO3ECliyBj786xS+LLc4ezn/pnJcXlC9vAuddd5mPyY9y5cyM2kzZvBm++caE0L//Njdu0ADeeQcee+wabuQcrWoGKnyKiIiIiEdTGBWnS0yE1avhp5/M448/4MIFsHzykrfMUfyr78an6Cn8ip/lg2dv5rE7rzOUbdtmwue338Jff5mh0nvugbFjoXVrl07FFRERERHJ6RRGxSn+/fe/8Pnzz2YbFYBataBXL/jx8BpOFjqAlSvp0nOSgJGLt1xbGN2xA777zoyCRkebY3ffDaNGwaOPQkCA4/5SIiIiIiJy3RRGJcssXw5ff20C6Nat5lhgoGkqdP/90LixWdsJ8HWffaS1xDMz25n8NHcFW8dO4Z61i7lt3zZzsE4dGD7cTMEtU8Yhfx8REREREXEchVFxuOho6NsXFiyAfPlMc6EXXzQB9JZb0m4sdM3bmezcCTNmcOzzr7h/4zruB9aVrMyghs/yc/UGvNypidZcioiIiIi4MYVRcZitW+Gtt8wSzcKF4YMP4KWXTCC9mkxtZ7JjB8yYYabhRkUBcCDwf3zSoCPzguqxu3CpS5eGR25RGBURERERcWMKo3LD9uwxTWmnTIG8eeHNN8060EKFMn+PdLczKXgWhgwxAXTNGnNxcLBJuo8+SrMJm7DTuF9mpveKiIiIiIjrKIzKdTt0CAYPhnHjwLbNKGhYGNx00/Xd79J2Jv/8Y8Jn517/NSG6804IDzdNiMqXv/ScAP9d1za9V0RERERE3ILCqFyzEydMb6APP4QzZ6BDB3j7bbPv53VLSoI5c2DECLPpKJgmRMOGmQCazs0zNb1XRERERETcjsKoZFp8vBkFHTwYjhwxGfHdd6FKlRu46enT8MUXJoRu22Y63w4ZAu3bQ9myV316utN7tV5URERERMStKYzKVdm2yYv9+kFMjOmKO2iQ2SP0usXGwtix8PHHcPQo1K5t9oFp3Rp8fK7pVpem94qIiIiISLahMCoZSkyEHj1gzBioWxemTjVbtVy3tWvNKOj06XDhAjzyCLz2Gtx9d9p7voiIiIiIiEdSGJV0nT4NTzwBs2eb7rgffABeXtdxo6QkmD/fLDJdvBjy54cXXoBXXoGKFR1et4iIiIiIuD+FUUnT/v3w8MNmN5WxY+HFF6/jJvHx/60H3bIFSpeGoUOhSxfw93d0ySIiIiIiko0ojMoVNm2CZs3M1i2zZsFDD13HTRYvhs6dYdcus7h02jTT8ega14OKiIiIiIhnup5Jl+LBfvnFLN88dw6WLbuOIHrypBlGbdzYBM9Fi+DPP013XAVRERERERG5SGFULpk6FUJCICAAVqy4jm65P/8MNWqYDrmvvw7r1plQqsZEIiIiIiKSisKoYNswcCA8/TTUrw/Ll0O5ctdwgxMnoFs3aNIE8uSB336DYcPA1zfLahYRERERkexNYTSHS0gwSzv79zdhdMGCa+wttHChGQ2dMMG03F271szzFRERERERyYAaGOVgx4+bnkKLFsHbb5tHWjNqI6JjCI/cQmxcPAH+voSGBNGqYgETPidOhKAgM5x6113O/0uIiIiIiEi2pDCaQ+3ZAw8+CJs3w2efQceOaV8XER1D2Mz1xCckAhATF8+P4VO4/5dx5Du4H0JD4Z13NCVXRERERESuicJoDhQdDc2bw+nTZlpu48bpXxseueVSEC1w7jR9F0+m/V8/sat4WcovXw516zqpahERERER8SQKoznM6tXQsCEULmxm1lavnvH1sXHxANyzcw0fzB/NTaeO8nGdNoyo/yRbFERFREREROQ6KYzmIAcPwiOPQJEi8McfZguXqwnw9+WuX3/kg/mj2VEkkDZPhbM2IIhAf03LFRERERGR6+ewMGpZVilgCPAgUADYAbxg2/ZSR30NuX4JCfD443DokBkRzUwQBRh3fAW3zRvJsvI1ef6RfsTnzouvjzehIUFZW7CIiIiIiHg0h4RRy7L8geXAb0Bz4BBwM3DQEfeXG9erFyxdClOnwh13ZPJJw4dz25B+7Lu3KW83eY2zpxMJTO6mWzMwS+sVERERERHP5qiR0TeAfbZtP3PZsZ0OurfcoC++gNGjoWdPePLJTDzBtuG998zmo48/TqmpU/nFxyfL6xQRERERkZzDy0H3aQWstCzrG8uyDlqWtdayrJcsK61dK8WZVq+Grl2hUSMYOjQTT7BtCAszQbRDB5g2DRRERURERETEwRwVRm8GXsSsEw0BRmHWj3ZP62LLsrpalhVlWVbUoUOHHFSCpJbcsOimm+CbbyDX1cbBk5Lg1Vfhgw/ghRfg00/B29sptYqIiIiISM7iqGm6XkCUbdthFz+PtiyrMiaMjkl9sW3bE4AJAMHBwbaDapDLpG5YVLz4VZ6QmAjPPw+TJ8Prr0N4OGhgW0REREREsoijRkb3AX+nOrYJKOug+8s1Sm5YNHFiJhoWJSTAM8+YINq/v4KoiIiIiIhkOUeNjC4HUu/18T/gXwfdX65BcsOiHj3gqaeucvG5c9C+PfzwAwwZAr17O6NEERERERHJ4RwVRkcAv1uW1Q/4BqgJvAL0ddD9JZNWrzazbRs1MgOcGYqPh9atYcECk15fftkpNYqIiIiIiDgkjNq2/adlWa2AQcBbwO6LH8c54v6SOckNi0qUyETDolOnoEULWLIEJk2Czp2dVaaIiIiIiIjDRkaxbXsuMNdR95Nrc00Ni+Li4MEHYdUqmDoVnnjCWWWKiIiIiIgADgyj4lqhoaZh0ZdfXqVh0fHj0LgxrF8P335rpumKiIiIiIg4mcKoB/jySxg1KhMNi5KSoEMH+OsvmDXLjI6KiIiIiIi4gMJoNrd6NXTtCg0bZqJh0dChJoSOHKkgKiIiIiIiLuWofUbFBQ4d+q9h0bffXqVh0c8/Q79+0LYtvPKK02oUERERERFJi0ZGs7Hu3U0H3d9/v0rDor17zV6it9xiOudaltNqFBERERERSYvCaDY1fz589x28995VGhadOwePPmr2FP3+e/Dzc1qNIiIiIiIi6VEYzYbOnIEXX4QqVUwX3Qy99hqsXGmS6y23OKU+ERERERGRq1EYzYYGDoRdu8xWLrlzZ3Dh1Kkwbhz06mVGR0VERERERNyEwmg2s2EDDBsGzz4L996bwYV//WXa7DZoAIMHXzocER1DeOQWYuPiCfD3JTQkiFY1A7O+cBERERERkcsojGYjSUnQrRsUKmR2aUlXXBy0aQP+/vD115fa7EZExxA2cz3xCYkAxMTFEzZzPYACqYiIiIiIOJW2dslGPv0Uli83I6PFiqVzUVISdOxo5vF+9x2ULHnpVHjklktBNFl8QiLhkVuyrGYREREREZG0aGQ0mzh4EN54w0zN7dAhgwuHDoVZs2DkSKhXL8Wp2Lj4NJ+S3nEREREREZGsopHRbKJXLzh1Cj7+OINtQn/+Gfr1g7Zt4ZVXrjgd4O+b5tPSOy4iIiIiIpJVFEazgcWL4csvzcholSrpXLRnD7RrZ7ZvmTQpzcQaGhKEr493imO+Pt6EhgRlQdUiIiIiIiLp0zRdN3fuHLzwAlSsaAY9073oscfMx5kzwc8vzcuSmxSpm66IiIiIiLiawqibGzIEtm6FyEjwTW827WuvwcqVMGMGBGU8ytmqZqDCp4iIiIiIuJym6bqxrVth0CAz+/b++9O5aOpUGDfOLCpt08ap9YmIiIiIiFwvhVE3Zdtmeq6vL4wYkc5F69dD167QoAEMHuzU+kRERERERG6Epum6qa++Mo2Lxo1LsVXofxIT4dlnoUAB+PpryKX/lCIiIiIikn0owbihY8fMMtA6deD559O5aMwYWL0apk9PJ62KiIiIiIi4L4VRN9SnDxw9CgsXgldaE6n37IE334QHHjB7ioqIiIiIiGQzWjPqZn7/HSZMgB494Lbb0rno5ZfNNN1x49LcT1RERERERMTdaWTUjSQkmGm5ZcrAgAHpXBQRAbNmwQcfQIUKTqxORERERETEcRRG3ciIEbBhg8mafn5pXHDiBLz0Etx6K/Ts6fT6REREREREHEVh1E3s2mVGQ1u1ghYt0rnozTchNha+/x58fJxXnIiIiIiIiINpzaib6NnTNCsaPTqdC/7803TQffFF02ZXREREREQkG1MYdQN//GGWgoaFmfWiV7hwAbp2NVu4vP++s8sTERERERFxOE3TdTHbhr59oUQJePXVdC4aNQrWroUZM6BQIWeWJyIiIiIikiUURl1s4UJYssRMz02zadGuXdC/Pzz8MLRu7eTqREREREREsobCqAslj4qWK2dm4aZ5QffuZi/RMWNS7CkaER1DeOQWYuPiCfD3JTQkiFY1A51XvIiIiIiIyA1QGHWh77+H1athyhTIkyeNC2bMgHnz4MMPoWzZS4cjomMIm7me+IREAGLi4gmbuR5AgVRERERERLIFy7ZtlxYQHBxsR0VFubQGV7hwAapXB29v+Osv8zGF48fhllsgIABWroRc/71vUG/IYmLi4q+4Z6C/L8v73JfFlYuIiIiIiGSeZVmrbdsOTn1cI6Mu8sUXsGULzJyZRhAF01r34EGYMydFEAWITSOIZnRcRERERETE3WhrFxc4exYGDIA774RWrdK44I8/4OOP4ZVXoFatK04H+Pumed/0jouIiIiIiLgbhVEX+Phj2LMHBg1K0ZPISEgw3YwCA+Hdd9N8fmhIEL4+KYdTfX28CQ0JyqKKRUREREREHEvTdJ3s5El4/31o3Ng8rjB8OGzYALNmQYECad4juUmRuumKiIiIiEh2pTDqZCNGwOHDZlT0Cjt2wDvvmP1EW7TI8D6tagYqfIqIiIiISLalabpOdPgwDBsGjzxi1oumYNvwwgvg4wOjR7ukPhEREREREWfRyKgTDRkCp0/De++lcfLrr+Gnn+Cjj8x6UREREREREQ+mkVEn2bsXxoyBp5+GqlVTnTxzBnr1guBgMzoqIiIiIiLi4TQy6iTvvgtJSWZLlyuMGgWxsWZ0NM1NR0VERERERDyLRkadYNs2+PRT6NYNypdPdfLwYTN/9+GH4Z57XFGeiIiIiIiI0ymMOkH//pAnD/Trl8bJ99+HU6dMIBUREREREckhFEaz2Nq1ZvZtjx5w002pTu7cCWPHwrPPprGQVERERERExHNpzWgW69cPCheG0NB0TubKxYLHujFwyGJi4+IJ8PclNCRIe4iKiIiIiIhHUxjNQr/9BvPmmRm4/v6pTq5eDdOns6XTy/T87TDxCYkAxMTFEzZzPYACqYiIiIiIeCyHTNO1LGuAZVl2qsd+R9w7u7JtCAuDUqXg5ZfTONm7NxQtystlmlwKosniExIJj9zivGJFRERERESczJEjo1uAhpd9npjOdTnC/PlmZHTcOMiXL9XJn36Cn3+GkSPZti/trVxi4+KzvkgREREREREXcWQDowu2be+/7HHIgffOVpKSoG9fuPlm6Nw5jZO9e0OFCtCtGwH+vmneI73jIiIiIiIinsCRYfRmy7JiLMvaaVnW15Zl3ezAe2cr334L69bBu+9C7typTn71lTn5/vuQJw+hIUH4+qQcHfX18SY0JMh5BYuIiIiIiDiZZdv2jd/EspoBBYDNQAngTeAWoJpt20fSuL4r0BWgbNmytf79998brsFdXLhgdmnJm9ds6+J1edw/exaCgqB4cVi16tLJiOgYwiO3qJuuiIiIiIh4HMuyVtu2HZz6uEPWjNq2PT/VF1sB7AA6AB+mcf0EYAJAcHDwjadhN/LNN7BtG8ycmSqIAowZA7t3w6efpjjZqmagwqeIiIiIiOQojpyme4lt26eAjUDlrLi/u0pKMrNvq1eHli1TnTx2DAYNgpAQaNzYJfWJiIiIiIi4iyzZZ9SyrLyYabq/ZMX93dUPP8CmTTB9ehqjooMHQ1wcfPCBK0oTERERERFxK47aZ3SYZVkNLMuqYFlWHWAGkB/43BH3zw5sG957DypXhsceS3Vy924YPRqefhpuu80l9YmIiIiIiLgTR42MlgamA8WAQ8AKoK5t257Tmegq5s0zDYs++wy8U28d2r+/+ThwoLPLEhERERERcUuOamDUzhH3ya5s2+TMcuXgySdTnfzrL/jiC3j9dShb1iX1iYiIiIiIuJssWTOa0yxeDCtXwvjx4OOT6mSfPlCoEISFuaQ2ERERERERd6Qw6gDvvQcBAdCxY6oTixfD/PkwdCgUKeKK0kRERERERNySwugN+u03WLIERoyAvHkvO5GUBG+8AWXKwMsvu6o8ERERERERt6QweoPefx+KF4cuXVKd+PZbWL0apkxJlVJFRERERETEIVu75FRRUbBgAbz2GuTPf9mJ8+ehXz+oUQOeespl9YmIiIiIiLgrjYzegPffB39/ePHFVCc+/hh27DDrRa/Y50VEREREREQ0Mnqd1q+HiAh49VUoWPCyEydPmn1e7rsPQkJcVZ6IiIiIiIhbUxi9ToMGgZ8fvPJKqhNjx8Lhw2bY1LJcUpuIiIiIiIi7Uxi9Dlu3mv5E3bun2rHl5EkID4dmzaBuXZfVJyIiIiIi4u4URq/DkCGQJ49pXJTCRx/B0aMwYIAryhIREREREck2FEav0a5d8OWX0LUrlChx2YkTJ2DYMHjwQbjzTleVJyIiIiIiki0ojF6joUPBywt69Up14qOP4NgxjYqKiIiIiIhkgsLoNYiNhcmT4dlnoXTpy04cPw7Dh8NDD0Ht2i6rT0REREREJLtQGL0Gw4ZBYiL07p3qxOjRGhUVERERERG5BgqjmXToEHz8MTz1FFSocNmJuDj48ENo0QJq1XJVeSIiIiIiItmKwmgmjRgBZ89CWFiqE6NHm0D69tuuKEtERERERCRbyuXqArKDY8dgzBh4/HEICrrsRPKoaMuWcMcdAERExxAeuYXYuHgC/H0JDQmiVc1Al9QtIiIiIiLirhRGM+Gjj+DkSejbN9WJkSNN86KLa0UjomMIm7me+IREAGLi4gmbuR5AgVREREREROQymqZ7FSdPmszZsiXceutlJ44dM3N3H3kEbr8dgPDILZeCaLL4hETCI7c4rV4REREREZHsQGH0KsaPN7mzX79UJ0aOhBMnUqwVjY2LT/Me6R0XERERERHJqRRGM3DmjNk+NCQk1fahx46ZMNq6Ndx226XDAf6+ad4nveMiIiIiIiI5lcJoBv780wx+XjEq+uGHV4yKAoSGBOHr453imK+PN6EhQYiIiIiIiMh/1MAoAw0awJ49UKzYZQePHoVRo+DRR1MtIv2vSZG66YqIiIiIiGRMYfQqUgRRMKOip06lu69oq5qBCp8iIiIiIiJXoWm61+LIETMq+thjUL26q6sRERERERHJthRGr8Xw4XD6NPTv7+pKREREREREsjWF0cw6fBg++ggefxyqVXN1NSIiIiIiItma1oxmICI65lIzondXfMVTp09jaVRURERERETkhimMpiMiOoawmeuJT0ikyJnjtP79B+ZWbcCFc4Vo5eriREREREREsjlN001HeOQW4hMSAei6aia+CecYUbct4ZFbXFyZiIiIiIhI9qeR0XTExsUDUPR0HM+smcPsqveyvVgZrIvHRURERERE5PppZDQdAf6+AHRZNZM8FxIYfXf7FMdFRERERETk+imMpiM0JIjA8yd5Jnous6o2YEfR0vj6eBMaEuTq0kRERERERLI9TdNNR6uagRT6uxAHCpdkzN3tCPT3JTQkiFY1A11dmoiIiIiISLanMJqBRk82gyd2stiyXF2KiIiIiIiIR9E03atREBUREREREXE4hVERERERERFxOoVRERERERERcTqFUREREREREXE6hVERERERERFxOoVRERERERERcTqFUREREREREXE6hVERERERERFxOoVRERERERERcTqFUREREREREXG6LAmjlmX1tSzLtixrTFbcX0RERERERLI3h4dRy7LqAl2Avxx9bxEREREREfEMDg2jlmUVAr4COgPHHHlvERERERER8RyOHhmdAMywbXuxg+8rIiIiIiIiHiSXo25kWVYXoBLwtKPuKSIiIiIiIp7JIWHUsqwgYBBwj23b5zNxfVegK0DZsmUdUYKIiIiIiIhkI46apnsXUAzYYFnWBcuyLgANgBcvfp7n8ott255g23awbdvBxYsXd1AJIiIiIiIikl04appuBBCV6thnwDbMiOlVR0tFREREREQk57Bs286aG1vWEmCDbdsvXeW6Q8C/WVKE4xQDDru6CMnx9DoUd6DXobgLvRbFHeh1KO4gO7wOy9m2fcWUWIc1MLpeaRXlbizLirJtO9jVdUjOptehuAO9DsVd6LUo7kCvQ3EH2fl1mGVh1Lbthll1bxEREREREcneHL3PqIiIiIiIiMhVKYxmzgRXFyCCXofiHvQ6FHeh16K4A70OxR1k29dhljUwEhEREREREUmPRkZFRERERETE6RRGRURERERExOkURjNgWdaLlmXttCzrrGVZqy3LusfVNUnOYVnWAMuy7FSP/a6uSzyfZVn3WpY127KsmIuvu46pzlsXX5+xlmXFW5a1xLKsai4qVzxUJl6HU9L4HrnCReWKh7IsK8yyrD8tyzphWdYhy7J+tCyreqpr9D1RslQmX4fZ8nuiwmg6LMtqC4wCBgE1gd+B+ZZllXVpYZLTbAFKXfao4dpyJIfwAzYArwLxaZx/A3gdeBmoDRwEFlqWVcBpFUpOcLXXIcAiUn6PfNA5pUkO0hAYB9wN3AdcABZZllXksmv0PVGyWkOu/jqEbPg9UQ2M0mFZ1krgL9u2u1x2bBsww7btMNdVJjmFZVkDgEdt265+tWtFsoplWaeAl2zbnnLxcwuIBcbYtv3+xWO+mF++etm2/YmrahXPlfp1ePHYFKCYbdsPuaouyXksy/IDjgOtbNv+Ud8TxRVSvw4vHptCNvyeqJHRNFiWlRuoBfyU6tRPmHckRJzl5otT1HZalvW1ZVk3u7ogyfEqACW57PujbdvxwDL0/VGcr75lWQcty9pqWdZEy7JKuLog8XgFML8/H7v4ub4niiukfh0my3bfExVG01YM8AYOpDp+APMNR8QZVgIdgWZAF8xr73fLsoq6sijJ8ZK/B+r7o7jaAuAZoDFmiuSdwGLLsvK4tCrxdKOAtcAfFz/X90RxhdSvQ8im3xNzuboAN5d6DrOVxjGRLGHb9vzLP7+4CH0H0AH40CVFifxH3x/FpWzb/vqyT9dblrUa+BdoDsx0TVXiySzL+hCoD9S3bTsx1Wl9TxSnSO91mF2/J2pkNG2HgUSufEerBFe+8yXiFLZtnwI2ApVdXYvkaMkdnfX9UdyKbduxwF70PVKygGVZI4D2wH22be+47JS+J4rTZPA6vEJ2+Z6oMJoG27bPA6uBpqlONcV01RVxOsuy8gK3APtcXYvkaDsxv3xd+v548bV5D/r+KC5kWVYxIBB9jxQHsyxrFPAEJgBsTnVa3xPFKa7yOkzr+mzxPVHTdNP3IfClZVmrgOVANyAA+NilVUmOYVnWMOBHYDfmHda3gPzA566sSzzfxS59lS5+6gWUtSzrduCobdu7LcsaCfSzLGszsBV4EzgFTHNBueKhMnodXnwMAL7H/KJVHhiM6WD6g5NLFQ9mWdZY4GmgFXDMsqzkEdBTtm2fsm3b1vdEyWpXex1e/H45gGz4PVFbu2TAsqwXMXtHlcLsddbTtu1lrq1KcgrLsr4G7sU01DoErADesm37b5cWJh7PsqyGwC9pnPrctu2OF7cyeBt4HiiMabbV3bbtDU4rUjxeRq9D4AUgArMPuD/ml69fMN8j9zilQMkRLMtK7xfld2zbHnDxGn1PlCx1tdfhxe2EIsiG3xMVRkVERERERMTptGZUREREREREnE5hVERERERERJxOYVREREREREScTmFUREREREREnE5hVERERERERJxOYVREREREREScTmFUREREREREnE5hVERERERERJxOYVRERERERESc7v/3zvp02FKOYQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1152x576 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(x1, y, 'o', label=\"Data\")\n",
    "ax.plot(x1, y_true, 'b-', label=\"True\")\n",
    "ax.plot(np.hstack((x1, x1n)), np.hstack((ypred, ynewpred)), 'r', label=\"OLS prediction\")\n",
    "ax.legend(loc=\"best\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Predicting with Formulas"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using formulas can make both estimation and prediction a lot easier"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from statsmodels.formula.api import ols\n",
    "\n",
    "data = {\"x1\" : x1, \"y\" : y}\n",
    "\n",
    "res = ols(\"y ~ x1 + np.sin(x1) + I((x1-5)**2)\", data=data).fit()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We use the `I` to indicate use of the Identity transform. Ie., we do not want any expansion magic from using `**2`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Intercept           4.826102\n",
       "x1                  0.521686\n",
       "np.sin(x1)          0.419346\n",
       "I((x1 - 5) ** 2)   -0.021333\n",
       "dtype: float64"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "res.params"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we only have to pass the single variable and we get the transformed right-hand side variables automatically"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0    10.813469\n",
       "1    10.671151\n",
       "2    10.432267\n",
       "3    10.134293\n",
       "4     9.826562\n",
       "5     9.558182\n",
       "6     9.366018\n",
       "7     9.265662\n",
       "8     9.247614\n",
       "9     9.279616\n",
       "dtype: float64"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "res.predict(exog=dict(x1=x1n))"
   ]
  }
 ],
 "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
