statsmodel.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "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:                 3.828e+08\n",
      "Date:                Fri, 25 Aug 2023   Prob (F-statistic):               0.00\n",
      "Time:                        14:57:13   Log-Likelihood:                -144.23\n",
      "No. Observations:                 100   AIC:                             294.5\n",
      "Df Residuals:                      97   BIC:                             302.3\n",
      "Df Model:                           2                                         \n",
      "Covariance Type:            nonrobust                                         \n",
      "==============================================================================\n",
      "                 coef    std err          t      P>|t|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------\n",
      "const          1.6107      0.273      5.900      0.000       1.069       2.153\n",
      "x1            -0.0602      0.089     -0.679      0.499      -0.236       0.116\n",
      "x2            10.0009      0.001   1.11e+04      0.000       9.999      10.003\n",
      "==============================================================================\n",
      "Omnibus:                        0.243   Durbin-Watson:                   2.171\n",
      "Prob(Omnibus):                  0.886   Jarque-Bera (JB):                0.411\n",
      "Skew:                          -0.077   Prob(JB):                        0.814\n",
      "Kurtosis:                       2.726   Cond. No.                     1.04e+03\n",
      "==============================================================================\n",
      "\n",
      "Notes:\n",
      "[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.\n",
      "[2] The condition number is large, 1.04e+03. This might indicate that there are\n",
      "strong multicollinearity or other numerical problems.\n"
     ]
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt \n",
    "import numpy as np\n",
    "import pandas as pd  \n",
    "import statsmodels.api as sm \n",
    "plt.style.use('ggplot')\n",
    "\n",
    "nsample = 100\n",
    "x = np.linspace(0, 10, 100)\n",
    "X = np.column_stack((x, x ** 3))\n",
    "beta = np.array([1, 0.1, 10])\n",
    "e = np.random.normal(size=nsample)\n",
    "\n",
    "X = sm.add_constant(X)\n",
    "y = np.dot(X, beta) + e\n",
    "\n",
    "model = sm.OLS(y, X)\n",
    "results = model.fit()\n",
    "print(results.summary())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "Multi-dimensional indexing (e.g. `obj[:, None]`) is no longer supported. Convert to a numpy array before indexing instead.",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "Cell \u001b[0;32mIn[4], line 11\u001b[0m\n\u001b[1;32m      9\u001b[0m ax\u001b[39m.\u001b[39mplot(x, y, \u001b[39m\"\u001b[39m\u001b[39mo\u001b[39m\u001b[39m\"\u001b[39m, label\u001b[39m=\u001b[39m\u001b[39m\"\u001b[39m\u001b[39mdata\u001b[39m\u001b[39m\"\u001b[39m, c\u001b[39m=\u001b[39m\u001b[39m'\u001b[39m\u001b[39mm\u001b[39m\u001b[39m'\u001b[39m)\n\u001b[1;32m     10\u001b[0m ax\u001b[39m.\u001b[39mplot(x, results\u001b[39m.\u001b[39mfittedvalues, \u001b[39m\"\u001b[39m\u001b[39mr--.\u001b[39m\u001b[39m\"\u001b[39m, label\u001b[39m=\u001b[39m\u001b[39m\"\u001b[39m\u001b[39mOLS\u001b[39m\u001b[39m\"\u001b[39m)\n\u001b[0;32m---> 11\u001b[0m ax\u001b[39m.\u001b[39;49mplot(x, iv_u, \u001b[39m\"\u001b[39;49m\u001b[39mr--\u001b[39;49m\u001b[39m\"\u001b[39;49m)\n\u001b[1;32m     12\u001b[0m ax\u001b[39m.\u001b[39mplot(x, iv_l, \u001b[39m\"\u001b[39m\u001b[39mr--\u001b[39m\u001b[39m\"\u001b[39m)\n\u001b[1;32m     13\u001b[0m ax\u001b[39m.\u001b[39mlegend(loc\u001b[39m=\u001b[39m\u001b[39m\"\u001b[39m\u001b[39mbest\u001b[39m\u001b[39m\"\u001b[39m)\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/matplotlib/axes/_axes.py:1632\u001b[0m, in \u001b[0;36mAxes.plot\u001b[0;34m(self, scalex, scaley, data, *args, **kwargs)\u001b[0m\n\u001b[1;32m   1390\u001b[0m \u001b[39m\"\"\"\u001b[39;00m\n\u001b[1;32m   1391\u001b[0m \u001b[39mPlot y versus x as lines and/or markers.\u001b[39;00m\n\u001b[1;32m   1392\u001b[0m \n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m   1629\u001b[0m \u001b[39m(``'green'``) or hex strings (``'#008000'``).\u001b[39;00m\n\u001b[1;32m   1630\u001b[0m \u001b[39m\"\"\"\u001b[39;00m\n\u001b[1;32m   1631\u001b[0m kwargs \u001b[39m=\u001b[39m cbook\u001b[39m.\u001b[39mnormalize_kwargs(kwargs, mlines\u001b[39m.\u001b[39mLine2D)\n\u001b[0;32m-> 1632\u001b[0m lines \u001b[39m=\u001b[39m [\u001b[39m*\u001b[39m\u001b[39mself\u001b[39m\u001b[39m.\u001b[39m_get_lines(\u001b[39m*\u001b[39margs, data\u001b[39m=\u001b[39mdata, \u001b[39m*\u001b[39m\u001b[39m*\u001b[39mkwargs)]\n\u001b[1;32m   1633\u001b[0m \u001b[39mfor\u001b[39;00m line \u001b[39min\u001b[39;00m lines:\n\u001b[1;32m   1634\u001b[0m     \u001b[39mself\u001b[39m\u001b[39m.\u001b[39madd_line(line)\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/matplotlib/axes/_base.py:312\u001b[0m, in \u001b[0;36m_process_plot_var_args.__call__\u001b[0;34m(self, data, *args, **kwargs)\u001b[0m\n\u001b[1;32m    310\u001b[0m     this \u001b[39m+\u001b[39m\u001b[39m=\u001b[39m args[\u001b[39m0\u001b[39m],\n\u001b[1;32m    311\u001b[0m     args \u001b[39m=\u001b[39m args[\u001b[39m1\u001b[39m:]\n\u001b[0;32m--> 312\u001b[0m \u001b[39myield from\u001b[39;00m \u001b[39mself\u001b[39;49m\u001b[39m.\u001b[39;49m_plot_args(this, kwargs)\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/matplotlib/axes/_base.py:488\u001b[0m, in \u001b[0;36m_process_plot_var_args._plot_args\u001b[0;34m(self, tup, kwargs, return_kwargs)\u001b[0m\n\u001b[1;32m    486\u001b[0m \u001b[39mif\u001b[39;00m \u001b[39mlen\u001b[39m(xy) \u001b[39m==\u001b[39m \u001b[39m2\u001b[39m:\n\u001b[1;32m    487\u001b[0m     x \u001b[39m=\u001b[39m _check_1d(xy[\u001b[39m0\u001b[39m])\n\u001b[0;32m--> 488\u001b[0m     y \u001b[39m=\u001b[39m _check_1d(xy[\u001b[39m1\u001b[39;49m])\n\u001b[1;32m    489\u001b[0m \u001b[39melse\u001b[39;00m:\n\u001b[1;32m    490\u001b[0m     x, y \u001b[39m=\u001b[39m index_of(xy[\u001b[39m-\u001b[39m\u001b[39m1\u001b[39m])\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/matplotlib/cbook/__init__.py:1327\u001b[0m, in \u001b[0;36m_check_1d\u001b[0;34m(x)\u001b[0m\n\u001b[1;32m   1321\u001b[0m \u001b[39mwith\u001b[39;00m warnings\u001b[39m.\u001b[39mcatch_warnings(record\u001b[39m=\u001b[39m\u001b[39mTrue\u001b[39;00m) \u001b[39mas\u001b[39;00m w:\n\u001b[1;32m   1322\u001b[0m     warnings\u001b[39m.\u001b[39mfilterwarnings(\n\u001b[1;32m   1323\u001b[0m         \u001b[39m\"\u001b[39m\u001b[39malways\u001b[39m\u001b[39m\"\u001b[39m,\n\u001b[1;32m   1324\u001b[0m         category\u001b[39m=\u001b[39m\u001b[39mWarning\u001b[39;00m,\n\u001b[1;32m   1325\u001b[0m         message\u001b[39m=\u001b[39m\u001b[39m'\u001b[39m\u001b[39mSupport for multi-dimensional indexing\u001b[39m\u001b[39m'\u001b[39m)\n\u001b[0;32m-> 1327\u001b[0m     ndim \u001b[39m=\u001b[39m x[:, \u001b[39mNone\u001b[39;49;00m]\u001b[39m.\u001b[39mndim\n\u001b[1;32m   1328\u001b[0m     \u001b[39m# we have definitely hit a pandas index or series object\u001b[39;00m\n\u001b[1;32m   1329\u001b[0m     \u001b[39m# cast to a numpy array.\u001b[39;00m\n\u001b[1;32m   1330\u001b[0m     \u001b[39mif\u001b[39;00m \u001b[39mlen\u001b[39m(w) \u001b[39m>\u001b[39m \u001b[39m0\u001b[39m:\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/pandas/core/series.py:1033\u001b[0m, in \u001b[0;36mSeries.__getitem__\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m   1030\u001b[0m     key \u001b[39m=\u001b[39m np\u001b[39m.\u001b[39masarray(key, dtype\u001b[39m=\u001b[39m\u001b[39mbool\u001b[39m)\n\u001b[1;32m   1031\u001b[0m     \u001b[39mreturn\u001b[39;00m \u001b[39mself\u001b[39m\u001b[39m.\u001b[39m_get_values(key)\n\u001b[0;32m-> 1033\u001b[0m \u001b[39mreturn\u001b[39;00m \u001b[39mself\u001b[39;49m\u001b[39m.\u001b[39;49m_get_with(key)\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/pandas/core/series.py:1048\u001b[0m, in \u001b[0;36mSeries._get_with\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m   1043\u001b[0m     \u001b[39mraise\u001b[39;00m \u001b[39mTypeError\u001b[39;00m(\n\u001b[1;32m   1044\u001b[0m         \u001b[39m\"\u001b[39m\u001b[39mIndexing a Series with DataFrame is not \u001b[39m\u001b[39m\"\u001b[39m\n\u001b[1;32m   1045\u001b[0m         \u001b[39m\"\u001b[39m\u001b[39msupported, use the appropriate DataFrame column\u001b[39m\u001b[39m\"\u001b[39m\n\u001b[1;32m   1046\u001b[0m     )\n\u001b[1;32m   1047\u001b[0m \u001b[39melif\u001b[39;00m \u001b[39misinstance\u001b[39m(key, \u001b[39mtuple\u001b[39m):\n\u001b[0;32m-> 1048\u001b[0m     \u001b[39mreturn\u001b[39;00m \u001b[39mself\u001b[39;49m\u001b[39m.\u001b[39;49m_get_values_tuple(key)\n\u001b[1;32m   1050\u001b[0m \u001b[39melif\u001b[39;00m \u001b[39mnot\u001b[39;00m is_list_like(key):\n\u001b[1;32m   1051\u001b[0m     \u001b[39m# e.g. scalars that aren't recognized by lib.is_scalar, GH#32684\u001b[39;00m\n\u001b[1;32m   1052\u001b[0m     \u001b[39mreturn\u001b[39;00m \u001b[39mself\u001b[39m\u001b[39m.\u001b[39mloc[key]\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/pandas/core/series.py:1082\u001b[0m, in \u001b[0;36mSeries._get_values_tuple\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m   1077\u001b[0m \u001b[39mif\u001b[39;00m com\u001b[39m.\u001b[39many_none(\u001b[39m*\u001b[39mkey):\n\u001b[1;32m   1078\u001b[0m     \u001b[39m# mpl compat if we look up e.g. ser[:, np.newaxis];\u001b[39;00m\n\u001b[1;32m   1079\u001b[0m     \u001b[39m#  see tests.series.timeseries.test_mpl_compat_hack\u001b[39;00m\n\u001b[1;32m   1080\u001b[0m     \u001b[39m# the asarray is needed to avoid returning a 2D DatetimeArray\u001b[39;00m\n\u001b[1;32m   1081\u001b[0m     result \u001b[39m=\u001b[39m np\u001b[39m.\u001b[39masarray(\u001b[39mself\u001b[39m\u001b[39m.\u001b[39m_values[key])\n\u001b[0;32m-> 1082\u001b[0m     disallow_ndim_indexing(result)\n\u001b[1;32m   1083\u001b[0m     \u001b[39mreturn\u001b[39;00m result\n\u001b[1;32m   1085\u001b[0m \u001b[39mif\u001b[39;00m \u001b[39mnot\u001b[39;00m \u001b[39misinstance\u001b[39m(\u001b[39mself\u001b[39m\u001b[39m.\u001b[39mindex, MultiIndex):\n",
      "File \u001b[0;32m/Library/Frameworks/Python.framework/Versions/3.10/lib/python3.10/site-packages/pandas/core/indexers/utils.py:343\u001b[0m, in \u001b[0;36mdisallow_ndim_indexing\u001b[0;34m(result)\u001b[0m\n\u001b[1;32m    335\u001b[0m \u001b[39m\"\"\"\u001b[39;00m\n\u001b[1;32m    336\u001b[0m \u001b[39mHelper function to disallow multi-dimensional indexing on 1D Series/Index.\u001b[39;00m\n\u001b[1;32m    337\u001b[0m \n\u001b[0;32m   (...)\u001b[0m\n\u001b[1;32m    340\u001b[0m \u001b[39min GH#30588.\u001b[39;00m\n\u001b[1;32m    341\u001b[0m \u001b[39m\"\"\"\u001b[39;00m\n\u001b[1;32m    342\u001b[0m \u001b[39mif\u001b[39;00m np\u001b[39m.\u001b[39mndim(result) \u001b[39m>\u001b[39m \u001b[39m1\u001b[39m:\n\u001b[0;32m--> 343\u001b[0m     \u001b[39mraise\u001b[39;00m \u001b[39mValueError\u001b[39;00m(\n\u001b[1;32m    344\u001b[0m         \u001b[39m\"\u001b[39m\u001b[39mMulti-dimensional indexing (e.g. `obj[:, None]`) is no longer \u001b[39m\u001b[39m\"\u001b[39m\n\u001b[1;32m    345\u001b[0m         \u001b[39m\"\u001b[39m\u001b[39msupported. Convert to a numpy array before indexing instead.\u001b[39m\u001b[39m\"\u001b[39m\n\u001b[1;32m    346\u001b[0m     )\n",
      "\u001b[0;31mValueError\u001b[0m: Multi-dimensional indexing (e.g. `obj[:, None]`) is no longer supported. Convert to a numpy array before indexing instead."
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAArUAAAH5CAYAAACf0sbLAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAA9hAAAPYQGoP6dpAABXKklEQVR4nO39e3yV5Z3v/7/uZCVhhQTCYcUEI4Jg0LYIOhUdtT+1u1Md8LHV2q1+tYfxgO3Qg1PbaZ1WrTpadU/BTmunasVaH9bxtPFQRVu1M3vXemytVdSCCCiHxCyUhEBCYCX374+wFlk5YIAka93J6/l4zGju+2atK1yNvNfFdX0+QRiGIZIkSVKEFeR6AJIkSdK+MtRKkiQp8gy1kiRJijxDrSRJkiLPUCtJkqTIM9RKkiQp8gy1kiRJijxDrSRJkiLPUCtJkqTIi+V6APlg06ZNpFKpIXmvRCJBMpkckvfS4HAOo835iz7nMPqcw2gb6vmLxWKMGzfuw58bgrHkvVQqxY4dOwb9fYIgyLyf3YmjyTmMNucv+pzD6HMOoy2f58/tB5IkSYo8Q60kSZIiz1ArSZKkyDPUSpIkKfIMtZIkSYo8Q60kSZIiz1ArSZKkyDPUSpIkKfIMtZIkSYo8Q60kSZIiz1ArSZKkyDPUSpIkKfIMtZIkSYq8WK4HIEmSpPwXtoe0vtjKe9vfo6W4hVFzRhEUBrkeVoahVpIkSbvVvLSZ5BVJUnWpzLVYdYzE1QnK55bncGS77HGofeONN3jkkUdYvXo1mzZt4lvf+hZz5szJ3A/DkPvuu4+nn36arVu3csghh3DhhRdSXV2deWbLli3cfvvt/OlPfyIIAo466ijOO+88Ro0alXnmnXfeYfHixbz99tuMGTOGk08+mVNPPTVrLM899xz33nsvyWSSqqoqzj33XI444oi9+X2QJElSL5qXNlN3UR2EUEKSOOtopYa2+kTn9VvJi2C7x3tq29ramDJlChdccEGv9x9++GEef/xx5s+fzw9+8ANKSkq49tpr2b59e+aZH//4x6xdu5bLLruMSy+9lDfffJNbbrklc7+lpYVrrrmGiRMncv311/O5z32O+++/n6eeeirzzPLly/n3f/93PvnJT3LDDTdw5JFH8m//9m+8++67e/otSZIkqRdhe0jyiiSEUMVjHM3ZzOYSjuZsqsLHAEh+P0nYHuZ4pHuxUnv44Ydz+OGH93ovDEOWLl3KZz7zGY488kgAvvrVrzJ//nxeeukljj32WNatW8crr7zCddddx7Rp0wA4//zzue666/j85z/P+PHjeeaZZ0ilUixYsIBYLMYBBxzAmjVrePTRR/nUpz4FwNKlS5k9ezb/83/+TwDOPvtsXnvtNZ544gkuuuiiXse3Y8cOduzYkfk6CALi8Xjm3wdb+j2G4r00OJzDaHP+os85jD7nMFpaX2wlVZeihCQzWERABwABHcxgEZvCObRtSLDtxW2UHlOa07EO6J7ahoYGGhsbOeywwzLXSktLmT59OitWrODYY49lxYoVjB49OhNoAWbOnEkQBKxcuZI5c+awYsUKDj30UGKxXcObNWsWDz/8MFu2bKGsrIwVK1ZwyimnZL3/rFmzeOmll/oc34MPPsgDDzyQ+Xrq1KnccMMNJBKJgfj2+62qqmpI308DzzmMNucv+pzD6HMOo+G97e8BEGddJtCmBXQQZz1tJCjbXsZ+1fvlYogZAxpqGxsbARg7dmzW9bFjx2buNTY2MmbMmKz7hYWFlJWVZT1TWVmZ9UxFRUXmXvrZ3b1Pb04//fSsIJz+lJhMJkmlUn39sgETBAFVVVXU19cThrlfpteecw6jzfmLPucw+pzDaGkpbgGglRpCCrKCbUgBrewPwJbiLXTUdfT6GvsqFov1awFyRFU/KCoqoqioqNd7Q/mDFYahP8gR5xxGm/MXfc5h9DmH0TBqzihi1THa6hM0hjOp4FUCQkIKWM4ltAUJYtUxRs0ZlfP5HNDmC+nV1KampqzrTU1NmXsVFRVs3rw56357eztbtmzJeqb7imv6667P7O59JEmStG+CwoDE1QlGhysZx18AeIPLeZ57qA/mAZC4KpEX9WoHNNRWVlZSUVHBa6+9lrnW0tLCypUrqa2tBaC2tpatW7eyatWqzDPLli0jDEOmT5+eeebNN9/M2hLw6quvMmnSJMrKyjLPdH2f9DMHH3zwQH5LkiRJI1r53HJqj34UgCQn0MAnaaNzhbb61uq8KOcFexFqt23bxpo1a1izZg3QeThszZo1bNy4kSAImDt3LkuWLOGPf/wj7777LjfddBPjxo3LVEOoqalh9uzZ3HLLLaxcuZK//vWv3H777RxzzDGMHz8egOOOO45YLMbNN9/M2rVrefbZZ3n88cez9sPOnTuXv/zlL/z6179m/fr13Hfffbz99tucfPLJA/DbIkmSJICCDRsY88fO8l3br/86h959KAc8cABTn5+aN4EWIAj3cAPE66+/zlVXXdXj+vHHH89XvvKVTPOFp556ipaWFg455BAuuOACJk2alHl2y5YtLF68OKv5wvnnn99n84Xy8nJOPvlkTjvttKz3fO6557jnnntIJpNUV1fvdfOFZDKZVeprsARBQHV1NXV1dTnfd6K94xxGm/MXfc5h9DmH0VN+7bWU/8d/0Hb00XywZMmQz19RUVG/Dortcagdjgy16i/nMNqcv+hzDqPPOYyWYMsW9jvySAo2b+b9X/yC7SedlLehdkRVP5AkSdKHC9tDWl9opfyeX1CweTM7pk2j7VOfIvfHwfpmqJUkSVJG89JmklckSdWlGMskYszhg43H0/HEVsbMG/PhL5AjA1r9QJIkSdHVvLSZuovqMq1xA1Ks4Fusb/p76i6qo3lpc66H2CdDrSRJkgjbQ5JXJCGEKh7jaM5mNpdwNGdTxVIAGq5oIGzPz73QhlpJkiTR+kJrZoV2BgszLXEDOpjBIkrCJKkNKRp/35jbgfbBUCtJkiRSDZ1Nr+KsIyB7NTaggzjrAdhet33Ix9YfhlpJkiQRq0zXD0j1uBdSQCv7A1BcXTyEo+o/Q60kSZKIHxUnVh2jiicBMmu1IQUs5xLaggSxSTEqPlGRszHujiW9JEmSRFAYUH3xDiov/R0Ar/EDOojTyv60BZ3NDyqvriQozM9qtYZaSZIkAbDf2/9JAe00Fh/BB9v/NnM9Vh0jcVWC8rnlORzd7hlqJUmSRNDYSOmvfgVA6rZvUBOvIdWQIlYZI35UPG9XaNMMtZIkSSJob6fl7LMpfuUVtn/yREqD/A6x3RlqJUmSRqiwPeysT9uQIlYZp/3KqwmCECIWaMFQK0mSNCI1L20meUWSVN2uEl6x6hiJq/N772xfLOklSZI0wjQvbabuorqdHcTq+SjfYwJ/IFWfou6iOpqXNud6iHvMUCtJkjSChO0hySuSEEIVj3E055DgWT7GZUwKHwYg+f0kYXv4Ia+UXwy1kiRJI0jrC607V2iTzGBRpiVuABzMjykJk6Q2pGh9oTW3A91DhlpJkqQRJNXQuYc2zjoCOrLuBXQQZ33Wc1FhqJUkSRpBYpWddQJaqaH7BoOQAlrZP+u5qDDUSpIkjSDxo+LEqmMU8QEBZIJtSAHLuYS2IEFsUmfDhSiJVgSXJEnSPgkKAxJXJxg3/z8B2MhxrOcMWtmftiABQOKqRN53EOvOlVpJkqQRpnxuOeE5x7OtsJo1nEcjs2kjQaw6RvWt1ZGsU+tKrSRJ0gjU/m9f4v1rL2DcH7dT3pAiVtm55SBqK7RphlpJkqQRKiiOUXrM8IiDw+O7kCRJUp/C9rCzPm1Digkv3UVs5kRazzgdiopyPbQBY6iVJEkaxpqXNpO8IkmqLkWMzRzMQgrZxpa15RT+89/nengDxoNikiRJw1Tz0mbqLqrLdBCbxs86Ay3TePtHB9G8tDnXQxwwhlpJkqRhKGwPSV6RhBCqeIyjOYtqngBgM4dCEJD8fpKwvXsLhmgy1EqSJA1DrS+0ZlZoZ7CIoEv/sGqWUhImSW1I0fpCaw5HOXAMtZIkScNQqiEFQJx1BHRk3QvoIM76rOeizlArSZI0DMUqO+sBtFJDSHbt2ZACWtk/67moM9RKkiQNQ/Gj4sSqY7QFCd7l7Mzmg5AClnMJbUGC2KTOhgvDwfCI5pIkScoSFAYkrk5Qd1Edq4OL2BCeTpy1tHIAbUECgMRVich2EOvOlVpJkqRhqnxuOdW3VhOritFGgkaOoI0EseoY1bdWUz63PNdDHDCu1EqSJA1TJU8/TflfXmTM0vlsXVlKqiFFrLJzy8FwWaFNM9RKkiQNR2FI+cKFFP/lL1BQQPid7+R6RIPKUCtJkjSMhO0hrS+0Uvx//y/Ff/kLHfE4Wy+8MNfDGnSGWkmSpGGieWkzySuSpOpSzObHANQVnMKWF4opn5vjwQ0yD4pJkiQNA81Lm6m7qI5UXYqJ/D8qeJUOYryz9X9Rd1EdzUubcz3EQWWolSRJiriwPSR5RRJCqOIxPsr3AQhIMZ4XAUh+P0nYHu7uZSLNUCtJkhRxrS+0kqpLUUKSGSzM9A8LgBksoiRMktqQovWF1lwOc1AZaiVJkiIu1ZACIM46ArJXYwM6iLM+67nhyFArSZIUcbHKzrP/rdQQdot3IQW0sn/Wc8ORoVaSJCni4kfFiVXHaAsSLOeSTLANKWA5l9AWJIhN6my6MFwN37guSZI0QgSFAft/qZGyK/837/B5nuce4qynlf1pCxIAJK5KDLsuYl0ZaiVJkoaB/V5ZTCm/JxhVwLJtV9JGZ5iNVcdIXJWgfG55jkc4uAy1kiRJEVe4ciXxhx8GoOPB71KzpYZUQ4pYZeeWg+G8QptmqJUkSYq48p/8hCAMaf30p0kd9jFKcz2gHDDUSpIkRVDYHtL6QivBsrepXvIgAFv+6Z9yO6gcMtRKkiRFTPPSZpJXJEnVpZjBvxPQzgclR/PB+oMon5Xr0eWGJb0kSZIipHlpM3UX1ZGqSzGG19mP3wCwpu1z1F1UR/PS5hyPMDcMtZIkSRERtockr0hCCFU8xuF8jQI6CIFS1gCQ/H6SsD3c7esMR4ZaSZKkiGh9oZVUXYoSksxgUaYlbgDMYBElYZLUhhStL7TmdqA5YKiVJEmKiFRDCoA46wjoyLoX0EGc9VnPjSSGWkmSpIiIVXae8W+nhO4bDEIKaGX/rOdGEkOtJElSRMSPihOrjlHNEwSQCbYhBSznEtqCBLFJnQ0XRpqRF+MlSZIiKigMqL54O1WXPg7A61xJirG0sj9tQWdb3MRViRHRQaw7Q60kSVKEVL12OwWkaCw+go3bj89cj1XHSFyVoHxueQ5HlzuGWkmSpIgoXLeO0nvvBWDHr75HDTWkGlLEKju3HIzEFdo0Q60kSVIeS7fDTTWkqF6yiCCVou2449hxzNGU5npwecRQK0mSlKe6tsMtIUktDwDw3jH/aIjrxt8PSZKkPJRuh0sIJSSJs5blfJsy3mLVv02i+uDmEbt/tjeGWkmSpDzTvR1uZ/ewjkzpLuhsh1t2UtmI3kfblXVqJUmS8kzPdrid3cMCOkZ8O9y+GGolSZLyjO1w95yhVpIkKc+k29y2UmM73H4y1EqSJOWZdDvcQlqBXXtmbYfbN+O9JElSngkKAxJXJxg3/woCQj7g47zLubbD3Q1DrSRJUh4aN3UdlfwXAG/zZbYyDbAdbl8MtZIkSXmofOFCAFrmncK4fziectvh7pahVpIkKU+kW+KyfC3Vv32SMAjY8q1vUlprQ9wPY6iVJEnKA11b4kLA+/yCCWNeoX1lNeW1uR5d/rP6gSRJUo6lW+KmGy5U8Gc6KGFd8ynUXVRH89LmXA8x7w34Sm1HRwf33Xcfv//972lsbGT8+PEcf/zxnHHGGQRB5/6PMAy57777ePrpp9m6dSuHHHIIF154IdXV1ZnX2bJlC7fffjt/+tOfCIKAo446ivPOO49Ro0ZlnnnnnXdYvHgxb7/9NmPGjOHkk0/m1FNPHehvSZIkadD0bIm7kICws3xXeAn1wTxb4vbDgK/UPvTQQzz55JNccMEF3HjjjZx77rk88sgjPP7445lnHn74YR5//HHmz5/PD37wA0pKSrj22mvZvn175pkf//jHrF27lssuu4xLL72UN998k1tuuSVzv6WlhWuuuYaJEydy/fXX87nPfY7777+fp556aqC/JUmSpEGT3RK3M9CCLXH31ICv1K5YsYKPf/zjHHHEEQBUVlbyzDPPsHLlSqBzlXbp0qV85jOf4cgjjwTgq1/9KvPnz+ell17i2GOPZd26dbzyyitcd911TJvWWb7i/PPP57rrruPzn/8848eP55lnniGVSrFgwQJisRgHHHAAa9as4dFHH+VTn/pUr2PbsWMHO3bsyHwdBAHxeDzz74Mt/R5D8V4aHM5htDl/0eccRp9z2FN7QzsAcdZmAm1auiVuGwnaG9pz/vuWz/M34KG2traWp59+mg0bNjBp0iTWrFnD8uXL+cIXvgBAQ0MDjY2NHHbYYZlfU1payvTp01mxYgXHHnssK1asYPTo0ZlACzBz5kyCIGDlypXMmTOHFStWcOihhxKL7foWZs2axcMPP8yWLVsoKyvrMbYHH3yQBx54IPP11KlTueGGG0gkEgP927BbVVVVQ/p+GnjOYbQ5f9HnHEafc7jLqI+Moo46StjY417XlriVH6lkXPW4oR5er/Jx/gY81J522mm0trbyjW98g4KCAjo6Ojj77LP5xCc+AUBjYyMAY8eOzfp1Y8eOzdxrbGxkzJgxWfcLCwspKyvLeqaysjLrmYqKisy93kLt6aefzimnnJL5Ov0pI5lMkkql9ur73RNBEFBVVUV9fT1h2L2Ts6LAOYw25y/6nMPocw57CqeHxKoKqKnvXHgLCXbtqU23xK2O0Tq9lW1123I61lzMXywW69cC5ICH2ueee45nnnmGr3/965ktAXfccQfjxo3jhBNOGOi32yNFRUUUFRX1em8of7DCMPQHOeKcw2hz/qLPOYw+57CLAph66p8pv+UtUpTyZ35MEc09WuJSMLR5ZXfycf4GPNTeddddnHrqqRx77LEATJ48mWQyyUMPPcQJJ5yQWU1tampi3LhdS+hNTU1MmTIF6Fxx3bx5c9brtre3s2XLlsyvr6ioyKzapqW/Tj8jSZKU91IpKp/6CQDry85m65Zd2y9tidt/A179oK2tjYKC7JctKCjIpPnKykoqKip47bXXMvdbWlpYuXIltbWdlYVra2vZunUrq1atyjyzbNkywjBk+vTpmWfefPPNrG0Dr776KpMmTep164EkSVI+ir39NgUffED7uHEUvfDP1NxfQ9VPq6i5v4apz0810PbTgK/U/s3f/A1Llixh4sSJ1NTUZCoSnHjiiUDnXoy5c+eyZMkSqqurqays5J577mHcuHGZagg1NTXMnj2bW265hfnz55NKpbj99ts55phjGD9+PADHHXcc999/PzfffDOnnnoqa9eu5fHHH+eLX/ziQH9LkiRJAyrdDjfVkCJWeQCpPzxH0VvLoWIMpcfkenTRFIQDvCGitbWVe++9lxdffJGmpibGjx/Psccey2c/+9lMpYJ084WnnnqKlpYWDjnkEC644AImTZqUeZ0tW7awePHirOYL559/fp/NF8rLyzn55JM57bTT9njMyWQyq9TXYAmCgOrqaurq6vJuH4r6xzmMNucv+pzD6HMOu7fD7RSrjpG4Ov+3GeRi/oqKivp1UGzAQ20UGWrVX85htDl/0eccRt9In8N0O1zCzrq0lfwX9ZxMW9BZ0an61uq8Drb5HGoHfE+tJEmSeureDncOX2Qqv+BozqIqfAyA5PeThO0jL+wPBEOtJEnSEOi7HS62wx0AhlpJkqQhkGro3EMbZ12f7XC7Pqc9Y6iVJEkaArHKzgPz7cTpvsGgazvc9HPaM4ZaSZKkIRA/Kk6sOkY1jxJAJthmtcOdFCN+VDyXw4wsPwpIkiQNgaAwYP9/bKb6iscBWMbVtFPeox1uUBjkcpiRZaiVJEkaImUfh7aaWrYmJ/J+2ycy122Hu+8MtZIkSUNkx6xZfPDck7CpiZrlJTs7inVuOXCFdt8YaiVJkgZJdjvcdHgtgAnjbIc7wAy1kiRJg6BrO9xx/Ili3uSd/c5iwjWT3WYwCAy1kiRJA6xrO9wS3qOWRcTZQMF721lz0flwKwbbAWZJL0mSpAHUvR3u0fx/xNlACOxgDGA73MFgqJUkSRpA2e1wF2W1w53Oz2yHO0gMtZIkSQMoux1uR9Y92+EOHkOtJEnSAEq3ud3OONvhDiFDrSRJ0gBKt8Ot4re2wx1CfkSQJEkaQEFhQOLqBHXz5zKKDbzP0bRRZTvcQWaolSRJGmDlc8vh50ey4ooppOp27Z21He7gMdRKkiTto6zOYRMgfkw55XPLKTuprJeOYq7QDgZDrSRJ0j7o2jkMYCbfIRUfT/PV/0LpOdMpPaY0xyMcGQy1kiRJeym7c1iSBE8xgRfpaI2x+p8/z7iK/dxqMEQMtZIkSXuhe+ewzkYLnXVpG5nFtmASye8nKTupzC0HQ8CSXpIkSXuhZ+ewXY0WxvFnO4cNMUOtJEnSXrBzWH4x1EqSJO2FdEewVmoIyd5eYOewoWeolSRJ2gvpzmHbGc8OxmSu2zksN/zoIEmStBcyncMuquOF8FfU8ABNHEYrNXYOywFXaiVJkvZS+dxyqm+tJqgeyzt8kUYOp40EseoY1bdWW85rCLlSK0mStAe6dg8rb/gjZeefYOewPGColSRJ6qeu3cPG8yL78x0arzuc937yK8pPGZvr4Y1obj+QJEnqh3T3sM7atPUczKLO69sPpu7L79G8tDnHIxzZDLWSJEkfonv3sKM5hzjvEQLbqAQg+f0kYXuY24GOYIZaSZKkD5HdPWwhAZ3hNQCm8x92D8sDhlpJkqQPkd09LHs11u5h+cFQK0mS9CHSXcE6KKT7BgO7h+UHQ60kSdKHSHcP66CUNiZmgq3dw/KHHyckSZI+xK7uYSmeD+9hNGsooplW9rd7WJ5wpVaSJKkf0t3DYtUlbGUajcy2e1gecaVWkiSpF107h41789fERyXhy1+i7IWpdg/LQ4ZaSZKkbrp2DitkCwfxA4pppKV+NB03zKf0mNJcD1HduP1AkiSpi+zOYUlqWdgZaDmA5Xcda+ewPGWolSRJ2qln57Cz2I//BuB9jiQMYnYOy1OGWkmSpJ2yO4ctymq0UMNDdg7LY4ZaSZKknbI7h3Vk3bNzWH4z1EqSJO2U7gjWyv52DosYQ60kSdJO6c5hbUElq7jQzmER4scMSZKknXZ1DqtjbXAuDeGnibPezmER4EqtJElSFxWHNXV2DquK0UbCzmER4UqtJEka0bp2Dhvd9DrVl51J2ZlnUvbsDbT+cbudwyLCUCtJkkasrp3DoIMjuJyADlLvtBIUxyg9xqgUFW4/kCRJI1L3zmFT+TljeJMUcV5/7vN2DosYP35IkqQRp3vnsBkszDRa+IA5bA8mkPx+krKTytxyEBGu1EqSpBFnd53DEvzezmERZKiVJEkjjp3Dhh9DrSRJGnF2dQ6rISR7e4Gdw6LJUCtJkkacXZ3DEiznm4Q7I5Gdw6LLjx+SJGnE6do5rD6Yx6Zwjp3DIs5QK0mSRoyujRYmPnc75RdOYtWv/4a2+gRtdIbZWHWMxFUJO4dFjKFWkiSNCF0bLZSyhloWUUA7O/75V+yYM8fOYRFnqJUkScNeutECIZTQwKFcSwHtbORYVv9wEtW3tjPmtDG5Hqb2gQfFJEnSsNa90cLRnE05KwmBJj4KQPL7ScL2cPcvpLxmqJUkScNadqOFXZ3DAuAgbrPRwjBhqJUkScNadqOF7NVYGy0MH4ZaSZI0rKUbKLQziu4bDGy0MHwYaiVJ0rCWbrTQHBzKWs7MdBCz0cLw4kcSSZI0rHVttLAq+EfWh5+10cIw5EqtJEka1go++ICxszZTfWs1saoYbSRoZDZtJIhVx6i+tdpGC8OAK7WSJGnY6do5bNL936fsxScoueEGyl44PXPdRgvDi6FWkiQNK107h43lL4xhCSEBm9cnKCoMKD2mNNdD1CBw+4EkSRo20p3DOuvSbuAQrgOgjnmsvmE/mpc253iEGiyGWkmSNCz07Bz2OeK8RwhsZQpg57DhzFArSZKGhd11DpvOf9g5bJgz1EqSpGFhV+ewtXYOG4EG5aDYBx98wF133cUrr7xCW1sbVVVVLFiwgGnTpgEQhiH33XcfTz/9NFu3buWQQw7hwgsvpLq6OvMaW7Zs4fbbb+dPf/oTQRBw1FFHcd555zFq1KjMM++88w6LFy/m7bffZsyYMZx88smceuqpg/EtSZKkPJfuCNbKJEKga00DO4cNfwO+UrtlyxYuv/xyYrEY3/3ud7nxxhv5whe+wOjRozPPPPzwwzz++OPMnz+fH/zgB5SUlHDttdeyffv2zDM//vGPWbt2LZdddhmXXnopb775JrfcckvmfktLC9dccw0TJ07k+uuv53Of+xz3338/Tz311EB/S5IkKQLSncPagiqW8y3CnTHHzmEjw4B/VHn44YeZMGECCxYsyFyrrKzM/HsYhixdupTPfOYzHHnkkQB89atfZf78+bz00ksce+yxrFu3jldeeYXrrrsus7p7/vnnc9111/H5z3+e8ePH88wzz5BKpViwYAGxWIwDDjiANWvW8Oijj/KpT32q17Ht2LGDHTt2ZL4OgoB4PJ7598GWfo+heC8NDucw2py/6HMOo28w5zCIBVRenWDDRXXUB/PYFM7p0Tms8upKCmLuvtxb+fwzOOCh9o9//COzZs1i0aJFvPHGG4wfP55Pf/rTmaDZ0NBAY2Mjhx12WObXlJaWMn36dFasWMGxxx7LihUrGD16dCbQAsycOZMgCFi5ciVz5sxhxYoVHHroocRiu76FWbNm8fDDD7NlyxbKysp6jO3BBx/kgQceyHw9depUbrjhBhKJxED/NuxWVVXVkL6fBp5zGG3OX/Q5h9E3UHMYtoc0/r6R7XXbKU3+iapfXc3Ef/sBy38Up21dgjY6/4wvqSlh+o+mk/jM0P6ZP1zl48/ggIfahoYGnnzySebNm8fpp5/O22+/zS9+8QtisRgnnHACjY2NAIwdOzbr140dOzZzr7GxkTFjxmTdLywspKysLOuZrivAABUVFZl7vYXa008/nVNOOSXzdfpTRjKZJJUa/E3jQRBQVVVFfX09YWg5kShyDqPN+Ys+5zD6BnIOm5c203B5A6m6FAE7+DgLCHiHwpKbOfC5H3VWQ3gvRWy/zi0HqcIUdXV1A/SdjEy5+BmMxWL9WoAc8FDb0dHBtGnTOOecc4DO1dB3332XJ598khNOOGGg326PFBUVUVRU1Ou9ofyPYxiG/sc44pzDaHP+os85jL59ncN0kwVCKCHJFH7BaN5hO+N47fefI/FEM+Vzy3u8pwZGPv4MDvimknHjxlFTU5N1raamho0bNwK7VlObmpqynmlqasrcq6ioYPPmzVn329vb2bJlS9Yz6VXbtPTX6WckSdLw07PJwllU8zgAGzmGVDDGJgsj0ICH2hkzZrBhw4asaxs2bMgsG1dWVlJRUcFrr72Wud/S0sLKlSupra0FoLa2lq1bt7Jq1arMM8uWLSMMQ6ZPn5555s0338zaNvDqq68yadKkXrceSJKk4SG7ycKirJq01Txuk4URasBD7bx583jrrbdYsmQJ9fX1PPPMMzz99NOcdNJJQOdejLlz57JkyRL++Mc/8u6773LTTTcxbty4TDWEmpoaZs+ezS233MLKlSv561//yu23384xxxzD+PHjATjuuOOIxWLcfPPNrF27lmeffZbHH388a8+sJEkafnY1WVhHQEfWPZssjFwDvqd2+vTpfOtb3+Luu+/m//yf/0NlZSVf/OIX+cQnPpF55tRTT6WtrY1bbrmFlpYWDjnkEL773e9SXFyceebrX/86ixcv5uqrr840Xzj//PMz90tLS7nssstYvHgxl156KeXl5Zxxxhl9lvOSJEnDw64mC/vbZEEZQZhvu3xzIJlMZtWvHSxBEFBdXU1dXV3eba5W/ziH0eb8RZ9zGH0DMYdhe8jqo1aTqk9RFT6a2YKQbrJQH8wjVh1j6vNTCQrzr55qlOXiZ7CoqCg31Q8kSZIGU1AYkLg6Qd1FddQHp7ApPKpHk4XEVQkD7QhjSw1JkhQJYXtIy7MtbLl7LYnf/pRJP64gVhWjjQSNzKaNBLHqGNW3Vvco56Xhz5VaSZKU95qXNpO8IkmqLsUM/jfjeRwe+RMdP15MbHyMVEOKWGVnkwVXaEcmV2olSVJeSzdaSNWlSPA01TxOCKxpO5v6L9fT3tjOmNPGUHpMqYF2BDPUSpKkvNW10cIkHuYjXJO5V8q7ADZaEGColSRJeaxro4WD+VGmfFcAzGCRjRaUYaiVJEl5K91AYQK/p/vGAhstqCtDrSRJyludDRRCqlna456NFtSVoVaSJOWt+FFxYtVFLOM6mvgI4c7okm600BYkiE3qrHqgkc2PNZIkKa+E7WHnXtqdZboSVyao+3KKP/NTSsKkjRbUK0OtJEnKG13r0QZsZyyvs6X6SCq+XMGWh7bQVpegjc4wG6uOkbgqYaMFAYZaSZKUJ9L1aAmhhCRTuZUqnuLdurNYdfOXqbq5ykYL6pOhVpIk5VzXerRVPMYMFhLQWXu2kM5yXRuv2sjU56caZNUrD4pJkqSc61qPdgaLMoEWYBKPWo9WH8pQK0mSci5dZzbOOgI6su5Zj1b9YaiVJEk5l64zW8A2uje8tR6t+sNQK0mSci5+VJyiKpjGLQSQCbbWo1V/+XFHkiTlXFAYMPFfq1k9/wKmcAdv8F2KabYerfrNUCtJknKie5OFspPK2PLzM3nl8hNI1bfTsvM569GqPwy1kiRpyDUvbabh8gZSdSkKaCVGKx3VlSSuTjD1xYOywq71aNUfhlpJkjSkkkuSbJi/IdNkYTo/poJX+GvdpdRddCzVt1a7Kqs95kExSZI0ZML2kJUXr8w0WTias0jwDEVsoYJXAEh+P0nY3r0GgrR7hlpJkjRkWl9opW1d284mCwuzmizUsMQmC9prhlpJkjRkUu91bbKQvRprkwXtC0OtJEkaMrH9Oo/zjGJtj3s2WdC+MNRKkqQhEz8qTln1FqaxGICQYOc/bbKgfePHIEmSNOjSNWnbG9qZeP5UGq49gTG8yTKuJE6DTRa0zwy1kiRpUDUvbSZ5RZJU3a59sgUV36SAbaQai2ljEmCTBe0bQ60kSRo0zUubqbuoDkIoZRXFbKKVybQ1JegIixn/rfEUTy22yYL2maFWkiQNirA9JHlFEkKo5lFqWUhA5z7a5eE3qQ/msfnuzUx9fqphVvvMg2KSJGlQtL7QSqouRQlJallEOrYGhMxgkTVpNaAMtZIkaVCka81O4A/WpNWgM9RKkqRBEauMUUgLB3BPj3vWpNVAM9RKkqRBET8qTm38JuK8xw7GEO6MHdak1WDwo5EkSRow6Xq0qYYUFSt/y36tjxNSwDL+lW1UE2e9NWk1KAy1kiRpQHSvR7uFUsbEDmTzYX/H1rojSNWlaKMzzFqTVgPNUCtJkvZZ13q0JSSJs45Wavhj6mbCl2Psd8tEYuNjtDe0U/mRSlqnt7oJUgPKUCtJkvZJ13q0VTzGDBYSEGb2ztYH89h41UamPj+VglgB46rHsa1uG2EYfviLS/3kZyRJkrRPutajTQda6CzbZT1aDRVDrSRJ2ifpOrOjWWk9WuWMoVaSJO2TzjqzIdU81uOe9Wg1VAy1kiRpn8SPijNp7FMk+AMhAeHOhrjWo9VQ8iOTJEnaY13r0cZTa5m+7UcArOZ83uMk69FqyBlqJUnSHulajzZgO0fwVQpooXn/OWxo/zyp+tB6tBpyhlpJktRvPevRrmUjx1DMJl5b/89MvGU/YuNjpBpSxCo7txy4QquhYKiVJEn90rMe7SICOggp4C2+zvZgYqYerUFWQ82DYpIkqV961qPtADrLdh3Mj61Hq5wy1EqSpH5J15kt5R3r0SrvGGolSVK/pOvMTuSZHvesR6tc8391kiSpT11LdxVOKCQx7iX23/Rw5z0CAsLserTV1qNVbhhqJUlSr7qW7gIooYGPB9cCsI7TWcv/Zz1a5Q1DrSRJ6qF76a5S3mEqP6cobKKZg1ldsYD2xpj1aJU3DLWSJClLX6W7ANop5nW+D/ES9r+livaN7dajVV4w1EqSpCzZpbt2BVqAAnYQUkx7XTtBQcCY08bkcKTSLlY/kCRJWdIlueKsywq0AAGhpbuUlwy1kiQpS7okVyuVhGRvKbB0l/KVoVaSJGWJHxUnVh1jMvfvLNnVGWyzSndNsnSX8osfsSRJEpBdk/bAw59h/7qHCQn4K5fSRqWlu5TXDLWSJCmrJm0pqzmYqwB4d9TneW/bpzPPWbpL+cpQK0nSCNe1Jm2cdzmM71DINj7g46ze9gXGf2s8xVOLLd2lvGaolSRpBOtZk/aHBEAIvM9REBSy+e7NTH1+qmFWec2DYpIkjWDZNWkXZmodBMB0fkZJmCS1IUXrC625HKb0oQy1kiSNYNk1acOsewEd1qRVZBhqJUkawXbVpK0h7BYLrEmrKDHUSpI0woTtIS3PtrD5oc2wrY1D4guBdpZzSSbYWpNWUePHLkmSRpCupbsAallIFY8yhld4kTvYxBzirLcmrSLHUCtJ0gjRtXRXCUn25wEm8SghAW/xNQrGFdG2KUEbnWHWmrSKEkOtJEkjQM/SXQszB8M2chybgjkUjgrY/979ad/Ybk1aRY6hVpKkEaBn6a5dlQ4m8gdKwiRtdQmCgoAxp43J4UilveNBMUmSRoB0Sa5S1li6S8OSoVaSpBEgXZJrOxXdIq2luzQ8GGolSRqmupbuCttDYtUxtgYH8xb/ZOkuDTt+HJMkaRjqXrorRjMdFWMhhA3BqbwfHmPpLg0rgx5qH3roIe6++27mzp3LP/zDPwCwfft27rzzTp599ll27NjBrFmzuPDCC6moqMj8uo0bN/Lzn/+c119/nVGjRnH88cdzzjnnUFhYmHnm9ddf584772Tt2rVMmDCBM844gxNOOGGwvyVJkvJa99JdY3iVg/kJ6xtP5x0+T0FFzNJdGnYGNdSuXLmSJ598kgMPPDDr+i9/+UtefvllLrnkEkpLS1m8eDELFy7kX//1XwHo6Ojguuuuo6KigmuuuYZNmzZx0003UVhYyDnnnANAQ0MD119/PX/3d3/H1772NZYtW8bNN99MRUUFs2fPHsxvS5KkvLW70l3VLGUtZxOMKrJ0l4adQQu127Zt4yc/+Qlf+tKXWLJkSeZ6S0sLv/vd77j44ov52Mc+BsCCBQv4xje+wYoVK6itreUvf/kL69at4/LLL6eiooIpU6Zw1lln8atf/YozzzyTWCzGb3/7WyorK/nCF74AQE1NDX/961957LHH+gy1O3bsYMeOHZmvgyAgHo9n/n2wpd9jKN5Lg8M5jDbnL/qcww/X+mLX0l2LsiodlLCRIjbTVldCQWEBZaeXDfn4nMNoy+f5G7RQe9ttt3H44Ydz2GGHZYXaVatW0d7ezsyZMzPX9t9/fyZOnJgJtStWrGDy5MlZ2xFmz57Nbbfdxtq1a5k6dSpvvfVW1msAzJo1izvuuKPPMT344IM88MADma+nTp3KDTfcQCKR2PdveA9UVVUN6ftp4DmH0eb8RZ9z2Lf3tr8HQJx1BHRk3UuX7mojQdn2Mvar3i8XQwScw6jLx/kblFD7hz/8gdWrV3Pdddf1uNfY2EgsFmP06NFZ18eOHUtjY2Pmma6BNn0/fS/9z/S1rs+0trayfft2iouLe7z36aefzimnnJL5Ov0pI5lMkkoNfl2+IAioqqqivr6eMOxeUEVR4BxGm/MXfc7hh2spbgEgzjs97nUt3bWleAsddR09nhlszmG05WL+YrFYvxYgBzzUbty4kTvuuIPLLrus12CZS0VFRRQVFfV6byh/sMIw9Ac54pzDaHP+os857ClsD2l9oZUd9TuIV3zA9MafdV4nICDMLt1VHWPUnFE5/T10DqMtH+dvwEPtqlWraGpq4jvf+U7mWkdHB2+++SZPPPEE3/ve90ilUmzdujVrtbapqSmzOltRUcHKlSuzXrepqSlzL/3P9LWuz8Tj8bwL05IkDabu5btgPG/xDSbwLCtZQJx6S3dp2BvwUDtz5kx++MMfZl372c9+xqRJkzj11FOZOHEihYWFvPbaaxx99NEAbNiwgY0bN1JbWwtAbW0tS5YsoampKbPF4NVXXyUej1NTUwPAwQcfzJ///Oes93n11VczryFJ0kjQvXxXnHW0UkM9J1PPSUBAG537Hy3dpeFswENtPB5n8uTJWddKSkooLy/PXP/kJz/JnXfeSVlZGaWlpdx+++3U1tZmAumsWbOoqanhpptu4txzz6WxsZF77rmHk046KbN94NOf/jS/+c1vuOuuuzjxxBNZtmwZzz33HJdeeulAf0uSJOWl7PJdj2aqHaS3GtQzj4IJBVReWUmsytJdGt5y0lHsi1/8IkEQsHDhQlKpVKb5QlpBQQGXXnopt912G5dddhklJSUcf/zxnHXWWZlnKisrufTSS/nlL3/J0qVLmTBhAl/+8petUStJGjFaX+havmsh6bga0MEMFrGJObS9nyBWFaP0mNKcjlUabEGYb7t8cyCZTGbVrx0sQRBQXV1NXV1d3m2uVv84h9Hm/EWfc5ht80Obqf9KPQdxM5O5t8f9V7iRRmZT9dMqxpw2Jgcj7Mk5jLZczF9RUVFuqh9IkqTBk65ykGpIkUqmqOBlari/53NdynfFKv3jXsOf/yuXJCkiulc5KOUdDuf7FNBBE4cwhhUEdPQo3xU/Kp7jkUuDz1ArSVIEdK9yUM6bTOOnFLGFJj7KX1hEEU3EWW/5Lo1IhlpJkvJcdpWDx3ZWOeggBLYzlmVcQ0dBMW0dCdroDLOW79JIY6iVJCnPZVc56Ay0AAFQRDMF7IAOmHjlRGKJGLFKy3dp5DHUSpKU51INnXto46zLBNq0gA7irKeNBLFELG+qHEhDrSDXA5AkSbuXrl4wltfoXkTJKgdSJ//XL0lSHupauqtwQiH7jXuWKZvuIABCgqzOYVY5kAy1kiTlne6lu8awjFlcRUDIev4n73IucTZY5UDqwlArSVIe6V66q4KXmc5NFLKdjRzD2xUX09FYQBuVgFUOpDRDrSRJeaJn6a6FBDt30bZSxRtcRhAvYv9bqmjf2G6VA6kLQ60kSXmiZ+muXcfCRtFAEVtoq4sTFARWOZC6sfqBJEl5oj+lu7o+J2kXV2olScqxdKWDthVtQAet7EdIQVawtXSXtHv+VEiSlEPZlQ5CpnMTpazlLb7GwfxkZztcS3dJH8ZQK0lSjnSvdHAgv2ASjwOwjjN4nnuIs97SXVI/GGolScqB3VU6eI9P8gFHA9BGZ5i1dJe0e4ZaSZJyILvSwcKsSgeV/Der+DJtJBh/8XhKjyu1dJf0Iax+IElSDqQrGCT4XVaghexKB8W1xZQeU2qglT6EK7WSJA2RdJWDVEOKVDJFAds5gHt7PmelA2mP+ZMiSdIQyK5ysFNBMa913MDB3MgYllvpQNoHhlpJkgZZ9yoHcdbSygG0dSTYwsH8mf/Yed1KB9LeMtRKkjSI+qpyEBKwnG9SzzwogLaOhJUOpH1gqJUkaRD1VeUgIGQGi9jEHNo6Eky8ciKxRIxYZcxKB9JeMNRKkjSI0lUOyljRZ5WDNhLEEjHGnDYmF0OUhgVDrSRJA6x7lYNCWpjCHT2fs8qBNGD8CZIkaQB1r3JQwDYO418oZyUpRlHIdqscSIPAUCtJ0gDpWeVgHRN4hgpeJcVoXmEhOxhvlQNpEBhqJUkaAD2rHCzKrMg2cShv8xW2FMyADqxyIA0CQ60kSQMgu8pBZ6CFzsNgY1hOG5XQgVUOpEFiqJUkaQCkqxyU8m4m0KZZ5UAafIZaSZL2QbrSQduKNgLa2Z//0/MZqxxIg86fLEmS9lLXSgcB7RzCD5jIc3RQQEC4s3OYVQ6koWColSRpL2RXOqhnBv+b8fyZDmK8zpVsodYqB9IQMtRKkrSHsisdPLqz/S2EQB2n8D7HAlY5kIaSoVaSpD3Us9JBpwCYxCO8yzm0kWD8xeMpPa7UKgfSEDDUSpLUD11b325fsR2AOOsICLOe61rpoLi2mNJjSnMxXGnEMdRKkvQhure+DUhRRDOt1BBSkFXCy0oHUm4U5HoAkiTls/SBsPR2g3G8yEz+hcO5mJBClnMJ4c4/TrMqHUyy0oE0lPwIKUlSH/pqfQvQQSGlvEM989jEHCsdSDlmqJUkqQ/ZB8IWZu2fDeiglRqgs8qBlQ6k3DLUSpLUh3Tr29Gs7OVAWJg5EDbu4nGU1JYQq4xZ6UDKEUOtJElddK1ykEqmiNHEVG7r+VyXA2GjjxttlQMpxwy1kiTt1L3KAUBx0E5h2EaKURSynYAOW99KechQK0kS3dveJomzjlZqaAsT/IWFFNJKO6M9ECblKUOtJGnE61nloPNQWHpFtp55nUUwO2x9K+UrQ60kacTrq8pBQAczWMQm5tDWkWDilROJJWIeCJPykKFWkjQi9db2diL/vdu2t7FEjDGnjcnFcCV9CEOtJGnE6e1A2AT+wDRu7fGsbW+laPCnU5I0ovR2IKyMt5jGLQR00Mw0ylhtlQMpYgy1kqQRY3dtbwHqOJkVfItiPrDKgRQxhlpJ0oiRfSAsO9CGBKzhPEIKbXsrRZChVpI07KUPhTUvbQaglDVZgRbSbW830EalbW+lCDLUSpKGte6HwgrZyoHcSQh0jaq2vZWizVArSRq2uh8KK+d1pnIHo3mHdooooN0DYdIwYaiVJA1LfXUJA0hRyivcyA7GeSBMGiYMtZKkYamvLmEAhWxjB+M8ECYNI4ZaSdKw0VuXsCoe222XsLHnjaV8brkHwqSIM9RKkoaF3rqEAWzmI7s9FFY+t9xDYdIwUJDrAUiStK/SB8LS2w0q+DMlJAHYxBxWcx7hzj/ysg6FTfJQmDRcuFIrSYq0vg6EhQQs55vUM493+QLv8fceCpOGMVdqJUmR1teBsICQGSzKrNi2kaCR2bTRWbar+tZqD4VJw4grtZKkyAnbQ1qeb8k6EJbgv3d7IMwuYdLwZqiVJEVKckmSVV9dlXUgrJpHOYibezxrlzBp5DDUSpIio3lpMxvmb8h0CIvzLpX8jkksBWAzMyjnLbuESSOQoVaSFAlhe0jD5Q1dDoQt2hleO63mH3iHL1DCRg+ESSOQoVaSlNfSDRW2PrO1y4GwzkALnfVnQwLqmQsEdgmTRihDrSQpb/XWUKGS32YCbVpA6IEwaYQz1EqS8lK6ocKu/bNrGccfmcx/9njWA2GSDLWSpLzTV0OFtEYOYyzLPBAmKcNQK0nKO301VIDO/bNvchmAB8IkZRhqJUl5IX0gLLuhwlO9NFTo3D+b7g4GHgiTZKiVJOWB3g6EAbRyACGdFQ7Suu6fHX/xeEqPK/VAmKSBD7UPPvggL774IuvXr6e4uJja2lo+97nPMWnSpMwz27dv58477+TZZ59lx44dzJo1iwsvvJCKiorMMxs3buTnP/85r7/+OqNGjeL444/nnHPOobCwMPPM66+/zp133snatWuZMGECZ5xxBieccMJAf0uSpEHU80DYu7QymTYSvM9xvMtZTOb+XvfPTvjmBMOsJGAQQu0bb7zBSSedxLRp02hvb+c///M/ueaaa1i0aBGjRo0C4Je//CUvv/wyl1xyCaWlpSxevJiFCxfyr//6rwB0dHRw3XXXUVFRwTXXXMOmTZu46aabKCws5JxzzgGgoaGB66+/nr/7u7/ja1/7GsuWLePmm2+moqKC2bNnD/S3JUkaBH0dCAsJWM43qWceq/kyGzjD/bOSdmvAQ+33vve9rK+/8pWvcOGFF7Jq1So+8pGP0NLSwu9+9zsuvvhiPvaxjwGwYMECvvGNb7BixQpqa2v5y1/+wrp167j88supqKhgypQpnHXWWfzqV7/izDPPJBaL8dvf/pbKykq+8IUvAFBTU8Nf//pXHnvssT5D7Y4dO9ixY0fm6yAIiMfjmX8fbOn3GIr30uBwDqPN+csPmb2z76VIJVM7D4Q1ZB0ICwiZwSI2MSfTTCG9f7akpoSJ359I2dyyXH4b2kv+HEZbPs/foO+pbWlpAaCsrPM/PqtWraK9vZ2ZM2dmntl///2ZOHFiJtSuWLGCyZMnZ21HmD17Nrfddhtr165l6tSpvPXWW1mvATBr1izuuOOOPsfy4IMP8sADD2S+njp1KjfccAOJRGIAvtP+q6qqGtL308BzDqPN+cud5JIkKy9eSdu6tsy1QlqYwfW9HAjryDRUmHzZZEZ/ZDTF1cVUfKLCFdphwJ/DaMvH+RvUUNvR0cEdd9zBjBkzmDx5MgCNjY3EYjFGjx6d9ezYsWNpbGzMPNM10Kbvp++l/5m+1vWZ1tZWtm/fTnFxcY/xnH766ZxyyimZr9OfMpLJJKlUqsfzAy0IAqqqqqivrycMww//Bco7zmG0OX+51by0mQ3zN3TZO7sOSFHLTyhl7W4PhIWzQzqO6aAtaCMoDJzDCPPnMNpyMX+xWKxfC5CDGmoXL17M2rVrufrqqwfzbfqtqKiIoqKiXu8N5Q9WGIb+IEeccxhtzt/QC9tDGi5v6LJ3dtHOg1+dQbaNidTzd0zm3l4PhI2aMyprzpzD6HMOoy0f52/QQu3ixYt5+eWXueqqq5gwYULmekVFBalUiq1bt2at1jY1NWVWZysqKli5cmXW6zU1NWXupf+Zvtb1mXg83usqrSRp6KX3z259ZmuXZgqdgRY6A20IvMr1bGUaGzjdA2GS9krBQL9gGIYsXryYF198kSuuuILKysqs+wcddBCFhYW89tprmWsbNmxg48aN1NbWAlBbW8u7776bFVpfffVV4vE4NTU1ABx88MFZr5F+Jv0akqTcal7azOqjVrPuf61j079vAmA8z2cCbVoAFNEMQBuJTFOFWHWM6lurbaggqV8GfKV28eLFPPPMM3z7298mHo9n9sCWlpZSXFxMaWkpn/zkJ7nzzjspKyujtLSU22+/ndra2kwgnTVrFjU1Ndx0002ce+65NDY2cs8993DSSSdltg98+tOf5je/+Q133XUXJ554IsuWLeO5557j0ksvHehvSZK0h3qrPVvBK0zm7t3unZ145URiiRixypgNFSTtkSAc4A0RZ555Zq/XFyxYkGmMkG6+8Ic//IFUKtVr84VkMsltt93G66+/TklJCccffzznnntuj+YLv/zlL1m3bt0+NV9IJpNZpb4GSxAEVFdXU1dXl3f7UNQ/zmG0OX9DI2wPWX3UalJ1qazas2mbqaWclVl7Z+uDecSqY0x9fupug6xzGH3OYbTlYv6Kior6dVBswENtFBlq1V/OYbQ5f4MnU3u2obP27MYrN1JCkqM5KyvQhgQ8zz1A0GPvbH+2GjiH0eccRls+h9pBr1MrSRrempc2k7wiSapuV2nEArYxnR/1Uns2JM6GzL5ZgFh1jMRVCffOStonhlpJ0l7ruXd2Ha3UEBJQwas9nu+6f3b8xeMpPa7UvbOSBoShVpK0V8L2kOQVyV5qz3buk/0r36WCl6lhSa+1Zyd8c4JhVtKAMdRKkvqt+97ZXbVndx0GC+hgBot4nnt4n79lHWdae1bSoDPUSpL6pbe9sxBSw7297J3tIM562khk/g/cPytp8BhqJUkfqre9s+2UMJXbGc+fejxv7VlJQ81QK0narb73znY2UWinhI0cSyX/3eve2XHnjzPIShp0hlpJUq/S+2e3PrO1y97ZRZk2twEQAq9yA03MYhVfdu+spJwx1EqSeuht/2ycdZlAmxZAZj+te2cl5ZKhVpKUpfv+2bG8QhW/YTUXEFKQFWzdOyspXxhqJWmE61qmq3BCYZf9s4/u3G7QuRJbxGaWc0mPerTunZWUDwy1kjSC9V6mC8bx4s7as7uU8TabmMPz3OPeWUl5x1ArSSNUb2W6tjGRSTzGAdxH94iarj3byGz3zkrKO4ZaSRqB+i7TFWS2G6RLdmV+TZf9s+MvHk/pcaXunZWUNwy1kjRC9N3itmuZrpAQ+CvfIaC9z/2zE745wTArKa8YaiVpBOht72xAisnc1WuZrjaqaGQ2m5jj/llJkWColaRhrre9szGamcIvKWNVj+e7bjOw9qykqDDUStIw03eJrsd2VjQIM8/uYAzv87fsx5M9thkUVhVS9aMq2je2W3tWUt4z1ErSMNJXia4S6nsE2pCAP7OIFqaxmgt6bDOovLqS0ceNHtLxS9LeMtRK0jDR2zaDVmpoI0GcDVmBFjoPhRXTTAtuM5AUfYZaSYqw9FaDHfU72Hjlxl5LdC3nm2xiTla5LrDFraThxVArSRHV21aDzhJdu7YZBITMYBHPcw/L+aYtbiUNW4ZaSYqg3rcaVO3sBNZ9m0FnJ7B65lmiS9KwZaiVpAjYfUWD9OorPVrbgiW6JI0MhlpJynN9VzTo3g2ss7Vtkk+Q4A+W6JI0ohhqJSmP9V3RYAJx1vXaDWwDn+FtvmaJLkkjiqFWkvJIf7cZNFPL61xDSEFWsE1vNXCbgaSRxlArSXli99sMulY0gDGsYBRrWc4lPSsakKBgQgGVV1YSq3KbgaSRwVArSXmgr20GKeJM5bYeFQ0AAgr6rGiw3/X7uTIraUQx1EpSjvSncUI7xcRo6/lrrWggSVkMtZKUA303TuhazSCkkDZa2J9NHMEkHrOigST1wVArSUOst60GbUyghvt7rWawgm/QyN/wLp+3ooEk9cFQK0mDbPcVDToPgO2+ccJkwG0GkrQ7hlpJGkR9VTQYzdvM4IeZIJtunNDACVTy/9xmIEl7yFArSQOk64psrDJG+wft1H25Z0WDgHaO4B97rMwGQB2nsooFbjOQpD1kqJWkAdDrimwBvWwz6Fx93cpUylmRFWxtnCBJe89QK0n7qM9Wth0Tmcj/7bbNoIMZLOIlfs4Y3rRxgiQNEEOtJO2F/tSY3UYlcd7r8WsDOihms40TJGkAGWolaQ/1XWO2ayvbkDjv0U6MAlK9bjMAKxpI0kAx1ErSbvT38Fecdb22sn2D71NEU89tBlY0kKQBZaiVpD7s/vDXozuDaufhr7eZT0hBVvOEkAK2MIM2Er1uM7CigSQNHEOtJNH/FdmOjmL2534O5FdZh7+m8XNWcREHcWuPg18UQFuH2wwkaTAZaiWNeB9ejiu9dQBCCimgvcdrBHTQzAye554eK7JVP6siNj6WCcxuM5CkgWeolTSi9V2OK7Hz8NeizJaCAAhoZysHUNptD601ZiUptwy1kkacsD2k5fmW3ZbjWsfpvM9xWXtk097iEkax3hqzkpRHDLWSRpTkkiSrvrrqQ8tx1bCEBj7V6+GvVvankdnWmJWkPGKolTRs9Tj8tamdui/t2mpQympGsYEqftOjHFcAFNHEci7pdUXWw1+SlF8MtZIir3t4jR8VZ8tvtvTr8Fefr0kBW5nGBxzd64qsh78kKb8YaiVFWm+VCwoqCuho7Ays6cNfISHjO/7MJmb2CLQhAe9yFpO5zxVZSYooQ62kyOhvLdm2xs4AegD37Kwbu2trQZz/X48V2oCQTRzFBj7jiqwkRZShVlIk9L+WbAENHE+c9YxhRdZrhEATHyXBM70e/rIclyRFl6FWUl7quiq7Y/UO3l/4/m5qyXatXNBBJf9Fb2upAbCV2j4Pf1mOS5Kiy1ArKe/0uipLz1qyST7BaNb0WrmgjpOo4re9NkiwHJckDT+GWkk51d99stDRo5ZsJf+v8zUga2U2pIA1XEBTl0NhHv6SpOHNUCtpSOxL2a2QApo4tMeKLMBaPksbCaZxS4/wWs+8Xldkq2+upnBcoYe/JGkYMdRKGnT9LbvVSg1BR4pKnmIqt2dWXwM6GMsbva7IruNM2kiQ5MTs8BpCwbgC2jZlr8jO+MkMUn+bIgx7BmRJUnQZaiUNuH4d8tpZdmsSD3Mw/05A2CO0dhUQ8i6f5QCW9Hs7QdlJZVmrw6VHl5KoSVBXVzc0vxGSpCFjqJW01/q9pYCe2wmWcwnN1HIwP+qyItu5P7aZaZSzqschr/WcyXrO3KNasqXHlGZeIwjcYiBJw5WhVtKH6m947X1Lwf6Mor5H2a0ZLOINLuuxMhsAq/gqo1jf+yEv8ICXJKkHQ62k3dqT/bC7OnndzUHcttstBZ3NDwJCgv6X3Qph/LfGUzy12ANekqQshlpJWfZkP2xvWwo2MYeD+HmPLQXpf8+8DwVs5qMs55uW3ZIk7TNDrTRC7f1+2IC1nElIwGTuyapQsLstBe9yFgdwf7/Lbu1un6wkSd0ZaqVhbt/2w9YAZAItdG4YmMy9vb7Xri0FBZnnIX3I6wzWc0a/ym65IitJ2lOGWmmY2Nfw2tY4kRLeYzL/ySQe2bltIGAt/ysroKY1cihj+WuP/bCdWwou6fchr+5lt1yRlSTtDUOtlMd6C6pBYdBra9nklXt2mKtzS8HCzGGudkYRY1vW+weEHMD9vR7mepOrGMeLvYbXXrcU7OaQV9eyW5Ik7Q1DrZQH+rvKGquOUXZaGVse2tJj3yv01dzgQQ7mJ5nwuo7PsI6zd4bRdIktiLGNDgIKurWi7Wx6sAf7Yd1SIEnKAUOtNMT2aZW1LkHjzxp7XifR7TAXbGUqMZoZxcbM6wZADQ+ymY/1uqXgr/wLh3L9Pu+HdUuBJGmoGWqlfdTfLQK7qy4A/S+Z9QFHUcMDHMB9mdXX9zmKCbzU5TAXlLG61/EGOyvH9naYq4nZA7Yf1i0FkqShZKiVerG7oLrpvzex+Y3NFFYW9rrK2tcWgT5XX3usshbwNhexmY92qzrQwQx+2Gu5rAm80GuDg9WcxxTu2KPDXO6HlSRFkaFWw86erJwCA7+Xta8tAr2uvgY7V1mfz6r3Op2baWJGjy0C6UYGvQXb3g5z1fP3tDFhj8Kr+2ElSVFkqFVO7Gvw3JO/3t/dyikBdGzqyLq2b3tZC1jJAho5nCqWUsOSzBaBrRxIM4dQxZNZNV8n8nyvv0cB7b1uEXiZH3EE/9Tj+iou4iBu3efw6n5YSVIUGWqHSNge0vpiK+9tf4+W4hZGzRkF9D+s7en1gXjtwXqNfQ2ee/zX+x+yctqfvawruJgP+Fuq+TUHctfOoBrQyEeo4I0uVQQ6mM5Nva6klvEOhWzr9YBW99XXkAJe5we9lsxqZmavWwfqmUcDnxyQ8OqWAklS1ARhGIYf/tjwlkwm2bFjx6C9fvPS5l4L4Pc3rO3p9YF47cF6jd3tK03r7fqePNvb4ap65vW4vpovUsi2na1eO0PqBxxJitFU8t/d/iq/51/57047JRTS1uP6Ks5jKr/s1yprPfO6fI87g2qP773n9bTYJFdeuwuCgOrqaurq6vA/fdHkHEafcxhtuZi/oqIiEomef851F/lQ+8QTT/DrX/+axsZGDjzwQM4//3ymT5++R68xmKG2eWkzdRfVQbjvYW1Pr+fza/Q3eC7nEoCsPairmM/7/C2VPMWB3L0zkHY+O4o6DuRX3VY9YRtVjKK+x/U9iXa91XDt7XU6twj8lCP4So/w+jz39Lr6Ws+87JDaZZU164PFpBhlp/bywWJSjInfn0hsfMzwuhv+YRp9zmH0OYfRZqgdJM8++yw33XQT8+fP5+CDD+axxx7j+eef50c/+hFjx47t9+sMVqgN20NWH7WaVF2qH2Ft98FuT68PxGt3vsbCzCrmSr5MkhPZj6c4iJ9nhcl6/p6p/JzJ3Jt5fgOnENBONY9nPbuJIzmas3oEzCY+ylhe73Y9yPz/D/39poB2int0xdpT9ZzIfj1Wagt4mZ9wBF/r9yprX7/XsO+rrH1tAdHu+Ydp9DmH0eccRpuhdpB897vfZdq0aVxwwQUAdHR08I//+I/8/d//PaeddlqP53fs2JEVXoMgIB6Pk0wmSaV6nmjfVy3PtrD2s2spIcnRnN0tDO0cQ5fnQwJe5j96WeGDbVQyioY+VgQX9Pir8p6v3fezHRRRwI5ewiRZz/YlfXjpb/h6v559g8v4KFd/6LO709df7yf5WybyXI/vZQVfp5Yfd/vee36PH7aa2ldQ7dcWgd2sviauTHSusr6XIrafIXWwBEFAVVUV9fX1/mEaUc5h9DmH0ZaL+YvFYv0KtZE9KJZKpVi1alVWeC0oKGDmzJmsWLGi11/z4IMP8sADD2S+njp1KjfccEO/fqP2xnvb3wMgzrpeSzN1FxAyltd6fTZOQy/Pd+x8vntb095eu+9nC+m5St1XmO29jWoHY3mz1+d7G8eu8v+7hASs4Rym7NxO0PV69/GEFPAKP+r1r/dX8g3e59heg2dIUb9WtHutGLAzqPZ1vY1dB7FKDiih8uxKGv6zgbZ1Xa7XlDD9R9OZeOpEGn/fyPa67RRXF1PxiQoD7BCrqqrK9RC0j5zD6HMOoy0f5y+yoXbz5s10dHRQUVGRdb2iooINGzb0+mtOP/10TjnllMzXQdAZJAZtpba4BYBWanopzQT0Ule0iZm9PBuwnH9iBj/q4/nur9N7EOz9tQtYxlV8jCv69Rp/7mOvaBOH9Tru3l6js/D/t3oNnm1U9yt4NnNI/5sHfEgg7Ve5q0kxyk8tp/mhZtrqsq/3tcp64NcP7Nwi0OV6qjBFfUM9zABmwDa2dX6tIeEKUfQ5h9HnHEabK7V5oqioiKKiol7vDcbEjJozilh1jLb6BMvDngEM+h/WOvdiFvbx/Df36bXf57h9fo2+rvf2GvscPPt6djeBNH24qmsgLRhXALBH5a4S300QXxmn4Y0GCisL+9wmEIYhFED8b+M9ryvnwjB0LiLOOYw+5zDa8nH+Ihtqx4wZQ0FBAY2NjVnXGxsbe6ze5kpQGJC4OkHdRXXUB/PYFO5lWNvL63n5GgMRPPt49kMD6b8k+l2LF3qv1RoUBow7YRzbZmzLux9mSZJGssiG2lgsxkEHHcSyZcuYM2cO0HlQbNmyZZx88sk5Ht0u5XPL4VZIXpHc+7C2h9cH4rUH7TUGMHj29uyHBdLerttoQJKk6It09YNnn32Wn/70p8yfP5/p06ezdOlSnnvuOW688cY9Wq0d7OYL0Fnea9uL2yjbXsaW4i0juqNYlA9FWYom2py/6HMOo885jDZLeg2iJ554gkceeYTGxkamTJnCeeedx8EHH7xHrzEUoRb8QR4OnMNoc/6izzmMPucw2vI51EZ2+0HaySefnFfbDSRJkjT0CnI9AEmSJGlfGWolSZIUeYZaSZIkRZ6hVpIkSZFnqJUkSVLkGWolSZIUeYZaSZIkRZ6hVpIkSZFnqJUkSVLkGWolSZIUeYZaSZIkRZ6hVpIkSZEXy/UA8kEsNrS/DUP9fhp4zmG0OX/R5xxGn3MYbUM5f/19ryAMw3CQxyJJkiQNKrcfDKHW1la+853v0NramuuhaC85h9Hm/EWfcxh9zmG05fP8GWqHUBiGrF69GhfHo8s5jDbnL/qcw+hzDqMtn+fPUCtJkqTIM9RKkiQp8gy1Q6ioqIjPfvazFBUV5Xoo2kvOYbQ5f9HnHEafcxht+Tx/Vj+QJElS5LlSK0mSpMgz1EqSJCnyDLWSJEmKPEOtJEmSIs9QK0mSpMiL5XoAI8UTTzzBr3/9axobGznwwAM5//zzmT59eq6HpX548MEHefHFF1m/fj3FxcXU1tbyuc99jkmTJuV6aNpLDz30EHfffTdz587lH/7hH3I9HPXTBx98wF133cUrr7xCW1sbVVVVLFiwgGnTpuV6aPoQHR0d3Hffffz+97+nsbGR8ePHc/zxx3PGGWcQBEGuh6devPHGGzzyyCOsXr2aTZs28a1vfYs5c+Zk7odhyH333cfTTz/N1q1bOeSQQ7jwwguprq7O2ZhdqR0Czz77LHfeeSef/exnueGGGzjwwAO59tpraWpqyvXQ1A9vvPEGJ510Etdeey2XXXYZ7e3tXHPNNWzbti3XQ9NeWLlyJU8++SQHHnhgroeiPbBlyxYuv/xyYrEY3/3ud7nxxhv5whe+wOjRo3M9NPXDQw89xJNPPskFF1zAjTfeyLnnnssjjzzC448/nuuhqQ9tbW1MmTKFCy64oNf7Dz/8MI8//jjz58/nBz/4ASUlJVx77bVs3759iEe6i6F2CDz66KP8j//xPzjxxBOpqalh/vz5FBcX81//9V+5Hpr64Xvf+x4nnHACBxxwAFOmTOErX/kKGzduZNWqVbkemvbQtm3b+MlPfsKXvvQlw1DEPPzww0yYMIEFCxYwffp0KisrmTVrFlVVVbkemvphxYoVfPzjH+eII46gsrKSo48+msMOO4yVK1fmemjqw+GHH87ZZ5+dtTqbFoYhS5cu5TOf+QxHHnkkBx54IF/96lfZtGkTL730Ug5G28lQO8hSqRSrVq1i5syZmWsFBQXMnDmTFStW5HBk2lstLS0AlJWV5Xgk2lO33XYbhx9+OIcddliuh6I99Mc//pGDDjqIRYsWceGFF/Ltb3+bp556KtfDUj/V1taybNkyNmzYAMCaNWtYvnw5hx9+eI5Hpr3R0NBAY2Nj1n9LS0tLmT59ek6zjXtqB9nmzZvp6OigoqIi63pFRUXmh1vR0dHRwR133MGMGTOYPHlyroejPfCHP/yB1atXc9111+V6KNoLDQ0NPPnkk8ybN4/TTz+dt99+m1/84hfEYjFOOOGEXA9PH+K0006jtbWVb3zjGxQUFNDR0cHZZ5/NJz7xiVwPTXuhsbERgLFjx2ZdHzt2bOZeLhhqpT2wePFi1q5dy9VXX53roWgPbNy4kTvuuIPLLruM4uLiXA9He6Gjo4Np06ZxzjnnADB16lTeffddnnzySUNtBDz33HM888wzfP3rX+eAAw5gzZo13HHHHYwbN87504Ax1A6yMWPGUFBQ0OOTS2NjY4/VW+W3xYsX8/LLL3PVVVcxYcKEXA9He2DVqlU0NTXxne98J3Oto6ODN998kyeeeIK7776bggJ3Y+WzcePGUVNTk3WtpqaGF154IUcj0p646667OPXUUzn22GMBmDx5MslkkoceeshQG0Hp/NLU1MS4ceMy15uampgyZUpuBoWhdtDFYjEOOuggli1bltls3dHRwbJlyzj55JNzPDr1RxiG3H777bz44otceeWVVFZW5npI2kMzZ87khz/8Yda1n/3sZ0yaNIlTTz3VQBsBM2bM6LFla8OGDSQSiRyNSHuira2tx89ZQUEBYRjmaETaF5WVlVRUVPDaa69lQmxLSwsrV67k05/+dM7GZagdAqeccgo//elPOeigg5g+fTpLly6lra3NT6cRsXjxYp555hm+/e1vE4/HM6vupaWl/lV2RMTj8R57oEtKSigvL3dvdETMmzePyy+/nCVLlnDMMcewcuVKnn76aS666KJcD0398Dd/8zcsWbKEiRMnUlNTw5o1a3j00Uc58cQTcz009WHbtm3U19dnvm5oaGDNmjWUlZUxceJE5s6dy5IlS6iurqayspJ77rmHcePGceSRR+ZszEHox6Qh8cQTT/DII4/Q2NjIlClTOO+88zj44INzPSz1w5lnntnr9QULFvjBJMKuvPJKpkyZYvOFCPnTn/7E3XffTX19PZWVlcybN49PfepTuR6W+qG1tZV7772XF198kaamJsaPH8+xxx7LZz/7WWIx19fy0euvv85VV13V4/rxxx/PV77ylUzzhaeeeoqWlhYOOeQQLrjggpw2JjLUSpIkKfLcSCZJkqTIM9RKkiQp8gy1kiRJijxDrSRJkiLPUCtJkqTIM9RKkiQp8gy1kiRJijxDrSRJkiLPUCtJkqTIM9RKkiQp8gy1kiRJirz/PwecWChGQEP9AAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 800x600 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "pred_ols = results.get_prediction()\n",
    "iv_l = pred_ols.summary_frame()[\"obs_ci_lower\"]\n",
    "iv_u = pred_ols.summary_frame()[\"obs_ci_upper\"]\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(8, 6))\n",
    "\n",
    "ax.plot(x, y, \"o\", label=\"data\", c='m')\n",
    "ax.plot(x, results.fittedvalues, \"r--.\", label=\"OLS\")\n",
    "ax.plot(x, iv_u, \"r--\")\n",
    "ax.plot(x, iv_l, \"r--\")\n",
    "ax.legend(loc=\"best\")\n",
    "plt.show()"
   ]
  }
 ],
 "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.10.0rc1"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}