{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<!--NOTEBOOK_HEADER-->\n",
    "*This notebook contains material from [cbe67701-uncertainty-quantification](https://ndcbe.github.io/cbe67701-uncertainty-quantification);\n",
    "content is available [on Github](https://github.com/ndcbe/cbe67701-uncertainty-quantification.git).*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<!--NAVIGATION-->\n",
    "< [5.1 Ridge Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.01-Contributed-Example.html) | [Contents](toc.html) | [5.3 Elastic Net Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.03-Contributed-Example.html)<p><a href=\"https://colab.research.google.com/github/ndcbe/cbe67701-uncertainty-quantification/blob/master/docs/05.02-Contributed-Example.ipynb\"> <img align=\"left\" src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open in Colab\" title=\"Open in Google Colaboratory\"></a><p><a href=\"https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.ipynb\"> <img align=\"left\" src=\"https://img.shields.io/badge/Github-Download-blue.svg\" alt=\"Download\" title=\"Download Notebook\"></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpages": {
     "level": 1,
     "link": "[5.2 Lasso Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2-Lasso-Regression)",
     "section": "5.2 Lasso Regression"
    }
   },
   "source": [
    "# 5.2 Lasso Regression\n",
    "\n",
    "Created by Haimeng Wang (hwang22@nd.edu)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "nbpages": {
     "level": 1,
     "link": "[5.2 Lasso Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2-Lasso-Regression)",
     "section": "5.2 Lasso Regression"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.linear_model import Lasso\n",
    "from sklearn.linear_model import LassoCV\n",
    "from matplotlib.ticker import MultipleLocator, FormatStrFormatter"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpages": {
     "level": 1,
     "link": "[5.2 Lasso Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2-Lasso-Regression)",
     "section": "5.2 Lasso Regression"
    }
   },
   "source": [
    "This example was adapted from:\n",
    "\n",
    "McClarren, Ryan G (2018). Uncertainty Quantification and Predictive Computational Science: A Foundation for Physical Scientists and Engineers, Chapter 4: Local Sensitivity Analysis Based on Derivative Approximations, Springer, https://doi.org/10.1007/978-3-319-99525-0_4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpages": {
     "level": 2,
     "link": "[5.2.1 Lasso Regression Basics](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2.1-Lasso-Regression-Basics)",
     "section": "5.2.1 Lasso Regression Basics"
    }
   },
   "source": [
    "## 5.2.1 Lasso Regression Basics\n",
    "\n",
    "Lasso stands for the Least Absolute Shrinkage and Selection Operator. Unlike the ridge regression, the lasso regression make the penealty to be *1-norm* ($L_1$) of the coefficient\n",
    "\n",
    "$$\\hat{\\beta}_{lasso} = \\min_{\\beta} \\sum^{I}_{i=1}(y_i-\\boldsymbol{\\beta} \\cdot \\boldsymbol{X}_i)^2 + \\lambda \\Vert \\boldsymbol{\\beta} \\Vert_1 $$\n",
    "\n",
    "Using the $L_1$ penalty tends to make some of the coefficients to be zero, which is also called a sparse model. In a sparse model, many of the coefficients are close to zero, and there are a few large non-zero coeffcients. In other words, this model does not include variables that are not important (variables with zero coeffcients)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpages": {
     "level": 2,
     "link": "[5.2.2 An example of the lasso regression ](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2.2-An-example-of-the-lasso-regression)",
     "section": "5.2.2 An example of the lasso regression "
    }
   },
   "source": [
    "## 5.2.2 An example of the lasso regression ##\n",
    "This is an example adapted from the example in Section 5.3 of the book. Let us assume that we have a simulation that has 200 input variables and only 120 simulations can be afforded. The *sklearn.linbear_model* is used to do the fit.\n",
    "\n",
    "Firstly, we used a model with 5 large non-zero ceofficients (between 5~30) and the rest coefficients are 0.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "nbpages": {
     "level": 2,
     "link": "[5.2.2 An example of the lasso regression ](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2.2-An-example-of-the-lasso-regression)",
     "section": "5.2.2 An example of the lasso regression "
    },
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Residue of the fiting\n",
      "0.9999789441481239\n",
      "Weight\n",
      "0.01873762187420794\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAq8AAAIlCAYAAAD7ShEjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3Xl8VNX9//H3ISSEBAJJoEFUgpbFBQQF6lYkICB1B9G6/OpC1dYNEeFLW1pBRSuFVkFt1fZbsbhvRQVXICwKylYUVKB8BRWFQEjYAiQk+fz+mMVMlslkMpPMwOv5eMxjknvPvedzh1HfHu49x5mZAAAAgHjQpLELAAAAAEJFeAUAAEDcILwCAAAgbhBeAQAAEDcIrwAAAIgbTRu7gHjinGNqBgAAgAZiZq7yNkZeAQAAEDcYeQ1DQ8+Nm5OTI0lasGBBg/Yby/hMAvF5BOLzqIrPJBCfRyA+j6r4TAI19OfhXJUBVz9GXgEAABA3CK8AAACIG4RXAAAAxA3CKwAAAOJGTIZX51ymc+5G59y/nXMbnXMHnHO7nXMfOud+6ZxrUql9R+ecBXm92FjXAgAAgMiJ1dkGLpf0N0lbJeVK+kZSlqRhkv4h6WfOucut6mP/n0qaVc351kaxVgAAADSQWA2vGyRdLGmOmZX7NjrnfidpmaTL5Amyr1U6brWZTWyoIgEAANCwYvK2ATObb2ZvVQyu3u3bJD3h/TWnwQsDAABAo4rVkddgDnnfS6vZ19459ytJmZJ2SlpqZp81WGUAAACIqrgKr865ppKu9f76bjVNBnlfFY9ZIOk6M/smUnX4VpmoK1bpAAAAR5JwM1MwMXnbQBAPSeom6W0ze6/C9v2S7pfUS1K699VPnoe9ciTNc86lNmypAAAAiLS4GXl1zo2UdLekdZJ+UXGfmW2XdE+lQxY55wZL+lDS6ZJulDQtErUwggoAAFC7cDOTc67GfXEx8uqcu02e4PmFpP5mVhDKcWZWKs/UWpJ0TpTKAwAAQANxVadKjS3OuVGSHpZnrtZzvaOsdTn+Ennmfn3PzIbUsxaTpFj/zAAAAOKZb+TVzKoMwcb0yKtzbpw8wXW1PCOudQquXmd437+KWGEAAABoFDEbXp1zf5DnAa2V8oy45gdpe7pzLqma7QMk3eX99dmoFAoAAIAGE5O3DTjnrpM0Q1KZpEcl7a6m2WYzm+Ftv0DSyZIWSNri3X+KpAHen/9gZpMiUBe3DQAAAERZsNsGYnW2geO87wmSRtXQZqE8AVeSZkoaKqmPpJ9JSpSUJ+llSY+Z2eKoVdpYinZK62ZL+/KkFlnSCRdKqZmNXRUAAEBUxeTIa6yKiZFXM2nxVGnhFKms+IftCc2kfmOlvmOkINNLAAAAxLp4HHlFTRZPleZXcwdEWfEP288Z27A1AQAANJCYfWAL1Sja6RlxDWbhFE87AAAOU2amxx57TD179lRKSoqcc3LOafPmzf42y5cv10UXXaQ2bdqoSZMmcs5p4sSJkjxLljrnNGPGjIjU07FjRznnWMSogRBe48m62YG3ClSnrFhaP6dh6gGAGHb99dfLOReVtdXRuB588EHdcccd+vTTT2VmysrKUlZWlhISEiRJ//3vf5WTk6PZs2ersLBQbdq0UVZWllq0aNHIlTeM1atXa+LEiREL57GG2wbiyb680Nrt3RbdOgAAaETTpnlWe//LX/6iUaNGVVlK9KmnntL+/fvVt29fvfnmm2rdunXA/g4dOqhr165q1apVROr58Y9/rOTkZKWkpETkfPW1evVq3XvvverXr5+uv/76xi4n4giv8aRFVmjtWraLbh0AADSS7du3a8eOHZKkm266qUpwlaTPP/9cknTFFVdUCa6S9K9//SuiNc2bNy+i50Nw3DYQRwo7DFaxJQZtU2yJKjx2UANVBABAwzpw4ID/55puA/C1OVJuEzjSEF7jyHubDml66dCgbaaXDtX7mw81UEUAcHgpKytTbm6u7rzzTvXq1UtZWVlKSkpS+/btNXToUM2fP7/GY8vLyzVjxgz1799fmZmZSkxMVNu2bXXyySdrxIgRevfdd6scs2nTJt1yyy3q0qWLmjdvrpSUFGVnZysnJ0d//OMflZ9f/eKSubm5GjZsmNq1a6ekpCS1a9eu1vpCcejQIT311FM699xz1bZtWzVr1kzZ2dkaPHiwnnrqKRUVFVU5pri4WH/5y190+umnq1WrVmrevLm6du2q0aNHa9u24LexlZSU6LHHHlPfvn2VkZHh72/EiBH68ssvA9ouWLBAzjl17NjRv833oJbvYazKD07dcMMN/v0Vj6vtgS0z00svvaQLLrhA7dq1U7NmzXT00UfrnHPO0cMPP6ydOwMfjK7tga26XKeP757tiRMnqqysTI888oh69OihlJQUZWRk6MILL9SKFSuqHOec0w033CBJWrhwYcBnVLnGvXv36v7771evXr3UsmVL/3e9d+/eGjt2rNauXVttbY3OzHiF+JJkno+scUybu8Gyx71lf/rdjXbwnkyzCWn+18F7Mu1Pv7vRsse9ZdPmbmi0GgEgVlx33XUmyfr16xfyMWvWrDHfv+slWbNmzSw1NTVg2wMPPFDtsVdffXVAu1atWllSUpL/99NPPz2g/cqVK61ly5b+/YmJida6deuAc7zzzjtV+hk/frx/v3POWrdubd55yE2S/eY3v6nT5+SzZcsW69mzp/88TZo0qXLu3NzcgGO2b99up556asDnVfGa0tPTbenSpdX29/3331uPHj0C+qt4bHJysr322mv+9h999JFlZWVZmzZt/G2ysrL8rylTpljv3r0tKyvLEhMTTZKlpaX59/fu3dt/rn79+pkke/rpp6vUtWvXLhs4cGCVz7hJkyb+bZWPy87OrvbzCec6fXzf3/Hjx9uQIUP835EWLVoEHLtkyZKA47KysiwtLc3fvuJnlJWVZR999JH/Ok866aSAutLT0wOuc9y4cdX+2TWECpmrah6rbiOv2AyvL3zytWWPm23Z42Zbz3HP29jf3W1TfzfCxv7ubus57nn/vheXfd1oNQI4vOzcV2wvfPK1TZu7wV745Gvbua+4sUsKWTjhdf369Xb55ZfbW2+9Zdu2bbPy8nIzM8vLy7P777/fEhISzDlnH3/8ccBxCxcu9AeAhx9+2Pbs2WNmZuXl5fb999/bjBkz7O677w44pn///v5Qu2rVKv/2oqIiW758uY0aNapKMHnhhRf8/1G//fbbbceOHWZmlp+fb3fccYd/38yZM0O+ZjOzgwcP2mmnnWaSrE2bNvbMM8/Yvn37zMxs//79/noqX7cvVKWnp9vLL79spaWlZma2fPly6969uz9g+ur0KSkpsT59+pgkO+ecc2zRokVWXOz5bm3bts3uvvtuk2QpKSm2cePGgGM3bdpktf33OFg4rW3/BRdcYJKsefPmNm3aNCssLDQzs+LiYluzZo3dc889NmvWrIBjagqv9blO3/e3devWlpGRYS+99JL/2E8//dS6detmkqxPnz5VruHpp5+u9bt/7733miRr27atzZ492w4dOuSvecOGDfbQQw/ZU089VePx0UZ4PUzC6859xdZ5/Nv+kFrdq/P4t+PqPy4AYlN5ebk9Om9DlX/ndB7/tj06b4M/1MWycMJrbe677z6TZNdff33A9smTJ5skGzJkSMjnat68uUmqEghrUl5ebp06dTJJduWVV1bb5qqrrjJJlp2dbWVlZSHX8vjjj/tHTj/99NOQjlm0aFHQEeJt27ZZenq6SbI//OEPAfv+/ve/+4PXwYMHqz3/LbfcYpLstttuC9gezfA6Z84c/2hrdddUk5rCa32u0/f9lWSLFy+uctyKFSv8+zdv3hywL5Tw+rOf/cwk2UMPPRTaRTawYOGVe17jSEZqkkYO6BS0zcgBnZSRmtRAFQE4XD2eu1FT39+gktLygO0lpeWa+v4GPZ67sZEqa1wXXXSRJOmjjz4K2J6WlibJ8yR8eXl5leOq4ztm69atIbVfvXq1Nm70fO6///3vq20zYcIESdLXX3+tZcuWhXRe6Yen72+44QadcsopIR3z6quvSpJ69+6tIUOGVNmflZWlX//615Kkl19+OWDfM888I0m67bbb1KxZs2rPf/XVV0uSPvjgg5DqiQTf53DeeedVe011FYnr7Nu3r376059W2d6rVy8dc8wxkn6YXaEu6vr9iyWE1zhzW/9OGjO4i5KaBv7RJTVtojGDu+i2/sHDLQDUpqCoRNPnBw+n0+dvVEFRSQNV1LAOHDighx9+WDk5OfrRj36kxMRE/8Mup556qiTp+++/Dzhm4MCBSkpK0qpVq5STk6Nnn322SpvKzj//fEnStddeq9/85jf6+OOPdehQzQ/crlq1SpL8D4FVp2vXrjr66KMD2tfm0KFDWrlyZUBNofCdv3///jW2GTBggCRpw4YN/oe9SktL/cF69OjRateuXbWvoUM9Dyh/++23IddUXx9//LGkun0ONYnUdfbp06fGPnx/1oWFhXWuz3eN06dP1y9+8Qu988472rt3b53P0xiY5zXOOOd0+4DOuvr0bH3wxTbl7SlWVlozDTqpHSOuACLi/c+3VRlxrayktFwffLFNP+/ToYGqahhbt25VTk6ONmzY4N+Wmpqq9PR0NWnSRGVlZcrPz6/y1H2nTp30t7/9TbfffrsWL16sxYsXS/I8hT5kyBDdfPPN/uDrM2XKFK1fv15LlizR5MmTNXnyZCUnJ+vMM8/U5Zdfruuvv17Nmzf3t/fNbeoLLDU55phj9N133/nb16agoEClpaWSPJP3hyqUenwjg2am/Px8paamqqCgQCUlJf6+a1Nxaqxoy8vzLAZUl8+hJpG6zpYtW9Z4THJysiQF/Z+emlx77bX66KOP9NRTT+nZZ5/Vs88+qyZNmuiUU07RRRddpFtuuUVHHXVUnc/bEBh5jVMZqUn6eZ8OGnluZ/28TweCK4CI2b63lmWovfL2hNYunowaNUobNmzQ8ccfr9dee00FBQXat2+ftm/frm3btvlH5qozYsQIbdq0SY888oguueQSZWZmavPmzXriiSfUq1cvPfjggwHtMzMz9eGHH+qDDz7QyJEjdeqpp6qkpES5ubm69dZb1a1bN23ZsqVKP8XFkf3cPbcXhq+u9VS8rcK3vGttr3gUD9f55JNPau3atbrnnnuUk5OjZs2aafXq1br//vvVuXPnBr1loy4IrwCAAD9qWf29eZVlpYXWLl6UlJTojTfekCQ999xzGjZsmNLT0wPa+EbmapKVlaU777xTs2bN0o4dO7Rs2TINHTpUZqY//OEP+uyzzwLaO+c0cOBATZs2TatWrVJ+fr6efPJJZWRk6KuvvtJdd93lb9u2bVtJ0jfffBO0Bl/g9bWvTWZmppo29fxF7Ndffx3SMRXPH+wYXy3OObVp08bfX0JCgiTpiy++CLm/hpCV5VnJsi6fQ01i+TorOvnkk3XvvfcqNzdXu3bt0ltvvaXu3burqKhI1113XVijutFGeAUABBh8crsq99VXltS0iQaddHgtRZ2fn+8fRaz8V/w+c+fODfl8zjn16dNHr7zyio455hiVl5frww8/DHpMenq6br75Zv8o7cKFC/37TjvtNElSUVFRjQ9jbdiwQd99911A+9okJiaqV69ekqS33347pGMqnn/hwoU1jhr6Fk3o0qWLUlNT/f317t1bkvT666+H3F9DOOOMMyTV7XOoSWNeZ5Mmnn9+6zqam5SUpAsvvFCvvPKKJM9tNP/9738jXl99EV4BAAGO1JlN0tLS5JyTJK1Zs6bK/q1bt+rRRx+t9ljfvY3VSUhIUGKiZ2lvXzguLy/332daHd+9rhX/Sr5nz57q1Mnz51L5FgSfiRMnSvLca/uTn/ykxvNXdu2110qSZsyYUWV0uCbDhw+X5HnS3TdiXVFeXp6eeOIJSdIVV1wRsO/666+XJL322mvKzc0N2k84DyOFy/c5vP/++9WuiFZXjXWdvpkEdu3aVWObYN/ZivdaR/o2lUggvAIAqjicZjY5dOiQ8vPzg74OHTqkFi1a+EfeRowYodWrV0vyBM158+apX79+NY5k/e53v9Pw4cM1a9asgIdz8vLyNHLkSG3atEnOOQ0aNEiStGfPHnXq1EkPPPCA1qxZo7KysoC+xo8fL8kzZZOPc06TJk2SJL3xxhu64447/MuU7ty5UyNHjtQLL7wgSZo0aZJ/9C0Uv/zlL9WzZ08VFxfr3HPP1cyZM7V//35JngeJli1bpptuukmffPKJ/5i+ffv6p5MaMWKEXn31Vf91rFy5UoMHD1ZhYaH/VorK/Z1xxhkqLy/XhRdeqGnTpgV8btu3b9cLL7ygnJwcTZs2LeTrqK+f/exn+tnPfiYz02WXXaZHH33UHwBLSkq0Zs0a3X333Zo1a1ZI52us6/TNRvHFF18E/JlVNHDgQI0cOVKLFi0KeFjs888/94fuo446St27d49YXRETyg3EvBpxkYJ9+WYrZpgtmOx535ffsP0DOKLt3FdsLy7zrLD14rL4XGErlJdvcvmPP/7Yv3iAJEtNTfX/npGRYbNmzap2gvw777wz4HxpaWkBS4BKgcvKFhYWBuxLTEy0jIwMS0hI8G87/vjj7dtvv61yXRWXh61uSc9wl4f95ptv/Ks2SbKEhARLT0+vdXnYikvKJicnV1ketvIqYT55eXl29tln+9s65yw9PT1g+VNJNnHixIDjor3CVmFhoX9/TZ9xXZaHDfc6fd/fCRMmhHWd55xzjv/cGRkZlp2dbdnZ2f7leisvWZuenm7Jycn+bSkpKTZ37twa+462Cn/GVfIYU2XFKjNp8VRp4RSprMKQ/dtjpX5jpb5jJO9fbwFAtPhmNjlSnH766Vq6dKkmTpyohQsXqqioSEcddZSGDBmi8ePH+0cWK7vrrrv04x//WPPmzdOXX36prVu3qri4WMcee6zOOuss3Xbbberbt6+/fVpammbPnq25c+dqyZIl2rJli3bs2KHU1FR17dpVl156qe64445qp0maNGmSBgwYoOnTp2vp0qUqLCxUZmamzjzzTI0cOVLnnntuWNd+7LHHasWKFXrqqaf08ssva+3atdq/f786dOigE044QcOHD69yK0Lbtm21dOlS/fWvf9ULL7yg9evXq6SkRJ07d9YFF1yg//mf/6lxuqUf/ehHWrhwoV566SU999xzWrlypQoKCpSUlKQTTjhBZ599ti677DINHDgwrOsJV+vWrTV//nw9++yzmjlzplavXq3du3frqKOOUqdOnTR06FBdfPHFIZ+vsa7z9ddf1z333KN33nlH3333nX/E9+DBg5Kkf/zjH3r77be1YMECbdq0Sdu2bZMknXDCCRo4cKBGjx6t4447LqI1RYqzOJ2CojF4/++zYaazWDRFmj+p5v0Dfi+dMzb6dQAAADQw3/3nZlZlpI7wWgcNFl6Ldsr+cqJcWc03SVtCM7nRX0qpmdGtBQAAoIEFC688sBWL1s0OGlwlefavn9NABQEAAMQGwmsM2l8QfD1sf7ud30W5EgAAgNhCeI1Ba/ckh9iuee2NAAAADiOE1xi0OvWnKrbEoG2KLVGrU89uoIoAAABiA+E1BqVlZGl66dCgbaaXDlWrzKwGqggAACA2MM9rDBp8cjud8eZQ6ZA0sum/1cwd8u8rtkRNLx2qv7uh+vgwW1ccAACgNkyVVQcNOc/rY/P/q6nvb1C69mhQwkplqVB5StcHZb1UqDSNGdxFtw/oHPU6AAAAGhrzvEZIQ4ZXM9PjuRs1ff5GlZSW+7cnNW2ikQM66bb+nfx/sAAAAIcTwmuENOgKW14FRSX64IttyttTrKy0Zhp0UjtlpCY1WP8AAAANjfAaIY0RXgEAAI40rLAFAACAwwLhFQAAAHGD8AoAAIC4QXgFAABA3CC8AgAAIG4QXgEAOELk5OTIOacZM2Y0dilRtXfvXo0ePVo//vGPlZSUJOecOnbsGNDmhRde0JlnnqmWLVvKOSfnnBYsWCBJ/t83b95c71o2b97sPx8ig+VhAQBHhFmzZmno0KGSpEGDBun999+PeB+7du3SI488IkmaOHFixM+P0AwbNkxz586VJKWlpSkjI0Nt27b173/++ed1zTXXSJISExOVlZUlSUpKOjLmUZ81a5ZWr16tnJwc5eTkNHY5dUZ4BQAcEZ555hn/z/PmzdOWLVt0zDHHRLSPXbt26d5775VEeG0sn3/+uebOnavExEQtWrRIZ5xxRpU2Dz/8sCTprrvu0p/+9Cc1bRoYh7p27SrJE2zrKzEx0X++WDFr1iz/Pw/xGF65bQAAcNjbuXOn5syZo5SUFF199dUqLy/Xs88+29hlIQo+//xzSdIpp5xSbXCt2GbEiBFVgqskrVu3TuvWrdPRRx9d73qOPvpo//kQGYRXAMBh7/nnn9ehQ4d0ySWX6Fe/+pWkwJFYHD4OHDggSWrRokW92iB2EV4BAIc9X1C95ppr1LdvX3Xo0EHr1q3TsmXLaj22qKhIU6dO1VlnnaWMjAwlJyfr+OOP18UXX6znnntOhw4dkuT569fjjjvOf5zvIR3fq+JtBB07dgx4QKg6NT00VFJSojlz5uimm25Sjx491KZNGyUnJys7O1vXXHONVq5cGfoHU0c7d+7UhAkT1KtXL7Vu3VopKSnq0qWLrrzySr3xxhvVHpOXl6e7775bJ5xwglJSUtSqVSv95Cc/0Z///GcVFxcH7W/Hjh367W9/q+7du6tFixZKTU1Vt27dNH78eBUUFAS0nThxopxzuv766yVJCxcuDPj8Z8yYUeXBqeOOO86/zXecVPsDW6F+J6TQHtiqy3X6VPwOFRQUaPTo0TruuOPUrFkzHX300brpppu0devWgGMWLFgg55z/n4d77723yve0ok2bNumWW25Rly5d1Lx5c6WkpCg7O1s5OTn64x//qPz8/BqvKarMjFeIL0nm+cgAAPFi7dq1JskyMzOtpKTEzMzGjRtnkuzWW28Neuznn39uHTt2NN+//5s2bWqtW7f2/y7JNm3aZGZmQ4cOtTZt2vi3Z2VlBbymTJniP292drZJstzc3Br7rnx+n7feeiug/5SUFEtOTg6o8V//+le15+zXr59JsqeffrrWz62yRYsWWWZmpr+fpKQkS0tLC6ilsk8++cQyMjL8+1u2bBlQa48ePSwvL6/a/hYvXhxwbFJSkjVv3tz/+7HHHmvr1q3zt58yZYplZWX5a0pMTAz4/BcsWOD/2XeONm3a+LeNHDmy1s/erG7fCTOzTZs21fj5hHOdPr7v0MyZM/0/p6SkWLNmzfzHduzY0QoKCvzHfPTRR5aVleX/M0hNTa3yPfVZuXKltWzZ0n+uxMTEKtf5zjvvVHtNkVDhM6uax6rbyIvwCgBmZrYv32zFDLMFkz3v+/Ibu6I6Gzt2rEmyW265xb/ts88+M0mWkZFhxcXF1R63c+dOO/bYY02SHXfccTZr1ix/2z179tjixYvthhtusG+//dZ/TG1Bxac+4TU3N9duuOEGmzdvnuXn//Dn8fXXX9uoUaNMkiUnJ9vXX39d5ZzhhteNGzf6Q2HPnj1t/vz5VlpaamZmBQUF9t5779mwYcMCjikoKLCjjjrKJFn37t1t2bJlZmZWWlpqr7zyiqWnp5skGzhwYJX+Nm/e7A9KN954o61bt87KysqsvLzc1q5da0OGDDFJdtJJJ/nr8Hn66adNkvXr16/G6wkWToPtj/R3oj7X6fsOtW7d2nr27GlLliwxM7NDhw7ZG2+84T/v2LFjq/R73XXXmSSbMGFCjZ9R//79TZKdfvrptmrVKv/2oqIiW758uY0aNcrfZzQQXgmvAFA35eVmC/9kdl9bswlpP7zua+vZXl7e2BWGpLS01B+gFi9eHLCve/fuJsleffXVao/1hd42bdrYli1bQuqvIcJrbUaMGGGSbOLEiVX2hRteL7/8cpNkXbp0sT179oR0zH333ecPV1u3bq2y/7333vNf47x58wL2XXPNNSYpYDS0ouLiYuvRo4dJsldeeSVgXzTDa6S/E/W5Tt93KCsrK+B/YnymTp3qD9mVhRJefaO/H3/8cS1XGB3Bwiv3vAIAqlo8VZo/SSqrdE9iWbFn++KpjVNXHb3//vvaunWrsrOzdfbZZwfs883zWdODWzNnzpQkjRkzJiJPnTeUiy66SJL00UcfReR8+/bt07///W9J0n333aeWLVuGdNyrr74qSbrxxhvVrl27KvsHDx6sM888U5L08ssv+7cfOHBAr7zyiiRp9OjR1Z47KSlJw4cPlyR98MEHIV5J/UXyOxGp67z55puVmZlZZfull14qyXPfalFRUZ3rS0tLk6Qq983GAuZ5BQAEKtopLZwSvM3CKVKvEVJq1f9oxhJfML3qqquqPIxy1VVX6be//a3eeecd7dixI2AS+82bN2vbtm2SpPPPP7/hCg5RQUGBHn/8cb3zzjtav369du/erbKysoA233//fUT6WrFihUpLS+Wc05AhQ0I6pqSkRGvXrpUk9e/fv8Z2AwYM0NKlS7Vq1aqA/kpKSiRJp59+eo3H+mYM+Pbbb0Oqqb4i/Z2I1HX26dOn2u0Vw/WuXbuUmppap/rOP/98Pf3007r22mt166236tJLL1WvXr0iMvdtfRFeAQCB1s2uOuJaWVmxtH6OdNq1DVNTGHbv3u1/Av7qq6+usr9Dhw7q27evFi1apOeff1533nmnf19eXl5Au1jyxRdfaMCAAQE1tmzZUs2bN5dzTiUlJSosLAxrtK06vn5atWqlVq1ahXRMQUGBysvLJSnoCKVvkYgdO3b4t1Uc6at4jTXZv39/SDXVV6S/E5G6zppGwpOTk/0/V5z9IFRTpkzR+vXrtWTJEk2ePFmTJ09WcnKyzjzzTF1++eW6/vrr1bx58zqfNxK4bQAAEGhf7f8hlSTt3RbdOurppZde0sGDByV5JqyvPCWQc06LFi2SVPXWAc8td7HphhtuUF5enk477TS9++672rt3r/bs2aO8vDxt27bN/1fRkbqG+p6ntumwKvOF3vT09JCeRwk23VgkRfo7EavX6ZOZmakPP/xQH3zwgUaOHKlSfEu0AAAgAElEQVRTTz1VJSUlys3N1a233qpu3bppy5YtDVqTD+EVABCoRVZo7VpWvY8xltRlEYL//Oc/WrNmjf/3ivdofv311xGtS5J/VSdfuK5s9+7d1W7/5ptvtGzZMiUkJOjNN9/UeeedV2Wi/VBG8erC91ns3r27xroqy8jIUJMmnogR7PPzhZ+Kt2xkZXm+f4WFhf6/po8Fkf5OxOp1VuSc08CBAzVt2jStWrVK+fn5evLJJ5WRkaGvvvpKd911V6PURXgFAAQ64UIpoVnwNgnNpK4XNEw9Ydi4caOWLFkiSVq9erUKCwtrfPkecKoYdjt27OgPK2+//XbI/foCmxR8pK5169aSVOPI1fLly6vdXjHs1fTX8XPnzg2p1lD17t1bTZs2lZnpnXfeCemYpKQkdevWTZKUm5tbY7v58+dLkk477bQq/UnS66+/Hm7ZERfud6ImjXmdvu9pXUeT09PTdfPNN+vBBx+U5FkIojEQXgEAgVIzpX5jg7fpNzamH9byBdEePXqoR48eat26dY2vyy+/XJL03HPPBTz09Itf/EKS9Oc//1nfffddSP36ntCWPA/J1KR79+6SVO2qVGamyZMnV3uc757TvLw8bd++vcr+NWvW6Pnnnw+p1lC1aNFCQ4cOlSRNmDBBe/fuDek431PyM2bMqPaJ9ffff19Lly6VJF1xxRX+7S1bttRll10mSZo0aVLQkeTS0lLt27cvtAuJgHC+EzVpzOv0fU9r+o6Wl5ertLS0xuN997rW9ZaQSCG8AgCq6jtGGvD7qiOwCc082/uOaZy6QmBm/imNhg0bVmv7iy66SImJidq2bZvee+89//Zx48bp6KOPVn5+vvr27as333zT/3T4vn37tGDBAl155ZUBo6etW7dW+/btJUlPP/10jX36wtqcOXM0efJk/8NVmzdv1lVXXaUVK1ZUe9yJJ56oY445Rmamn//859q4caMkzwM5r7/+ugYNGlTlNoJIePDBB9WyZUtt2LBB55xzjnJzc/33bO7atUtz5szRBRcEjsTffvvtOuqoo3TgwAENGTLEf01lZWV67bXXdOWVV0qSBg4cqAEDBgQc+9BDDykjI0Nbt27VWWedpX//+98BQWnjxo165JFHdOKJJ9b4WUVDON+JYBrrOk8++WRJ0rvvvlvt/1js2bNHnTp10gMPPKA1a9b4/6euvLxc8+bN0/jx4yVJ5513XsRqqpNQbhLmxSIFAI5Q+/LNVj7jWWFr5TNxscLW/Pnz/ROcr127NqRjzjvvPJNkV1xxRcD2zz77zI455pigS2RWnsj+nnvu8e9LTU217Oxsy87Otocffjig3bBhw/ztmjRp4j9vcnJywAT+lc//+uuvW5MmTfz7W7ZsaUlJSSbJOnToYDNnzjRJlp2dXeU667M87Pz58wOuvVmzZtaqVauAz6KyTz75xL+Slq/WisvDnnLKKTUuD7ts2TJr3769v23Tpk0tMzMzYPlTSbZgwYKA46K5SIFZ3b8TtS1cEe511mehix07dviXpG3SpIm1a9fO/z01MyssLAzoOzEx0TIyMiwhIcG/7fjjjw9YSSzSKnxmVfNYdRt5EV4BIF75Vg/q0qVLyMc89dRT/kBWWFgYsG/37t02adIk6927t6WlpVlycrIdf/zxdumll9oLL7xghw4dCmhfWlpqkydPtlNOOcVSUlL8/xGuvJpRcXGxPfDAA9a1a1dLSkqytm3b2mWXXWaffvqpmQUPUAsWLLBBgwb5w2Dnzp1tzJgxtnPnTsvNzY1KeDUzy8vLs3Hjxlm3bt0sNTXVUlNTrUuXLnbVVVfZm2++We0xW7dutbvuusu6dOliycnJ1qJFC+vdu7dNmTLFDhw4ELS/PXv22OTJk+2ss86y9PR0S0hIsNatW1vv3r1t3Lhxtnz58irHRDu8mtXtOxHKqmvhXGd9V2n79NNPbdiwYZaVlRUQSs3MysrKbPbs2TZq1Cj7yU9+Yu3bt7fExERLS0uzPn362AMPPBDySmvhChZenWc/QuGc8yRYPjMAAICo8S0qYmau8j7ueQUAAEDcILwCAAAgbhBeAQAAEDcIrwAAAIgbhFcAAADEDcIrAAAA4gbhFQAAAHGD8AoAAIC4QXgFAABA3CC8AgAAIG4QXgEAABA3CK8AAACIG4RXAAAAxA3CKwAAAOIG4RUAAABxIybDq3Mu0zl3o3Pu3865jc65A8653c65D51zv3TOVVu3c+4s59zbzrkC59x+59xnzrlRzrmEhr4GAAAARJ4zs8auoQrn3K8l/U3SVkm5kr6RlCVpmKRWkl6TdLlVKN45d4l3+0FJL0kqkHSRpK6SXjWzyyNQl0lSLH5mAAAAhwvnnCTJzFyVfbEYxJxzAySlSppjZuUVtreTtEzSsZKGm9lr3u1pkjbKE2zPNrMV3u3JkuZLOlPSVWb2Yj3rIrwCAABEWbDwGpO3DZjZfDN7q2Jw9W7fJukJ7685FXYNl9RW0ou+4Optf1DS772/3hK9igEAANAQYjK81uKQ9720wrYB3vd3q2m/SNJ+SWc555pFszAAAABEV9PGLqAunHNNJV3r/bViUO3qfd9Q+RgzK3XObZJ0sqTjJX1Z3zpycnLCOm7BggX17RoAACBuhJuZgom3kdeHJHWT9LaZvVdheyvv++4ajvNtbx2twgAAABB9MfnAVnWccyMlTZO0Tp6Hsgoq7NsgqbOkzma2sZpjl8jz0NaZZvZxPWrggS0AAIAoi7sHtipzzt0mT3D9QlL/isHVyzey2krVS6vUDgAAAHEo5sOrc26UpMckrZUnuG6rptl673uXao5vKuk4eR7w+ipadQIAACD6Yjq8OufGSXpY0mp5guv2GprO974PqWbfOZJSJC0xs+LIVwkAAICGErPh1Tn3B3ke0Fop6Vwzyw/S/FVJ+ZKudM71rnCOZEmTvL/+LVq1AgAAoGHE5ANbzrnrJM2QVCbpUVV/r+pmM5tR4ZhL5QmxByW9KM/ysBfLuzyspCusnhfLA1sAAADRF4/Lw06UNKGWZgvNLKfScWdLGi/PzALJ8iwZ+09J082sLAJ1EV4BAACiLO7Ca6wivAIAAERf3E+VBQAAAEiEVwAAAMQRwisAAADiBuEVAAAAcYPwCgAAgLhBeAUAAEDcILwCAAAgbhBeAQAAEDcIrwAAAIgbhFcAAADEDcIrAAAA4gbhFQAAAHGD8AoAAIC4QXgFAABA3CC8AgAAIG4QXgEAABA3CK8AAACIG4RXAAAAxA3CKwAAAOIG4RUAAABxg/AKAACAuEF4BQAAQNwgvAIAACBuEF4BAAAQNwivAAAAiBuEVwAAAMQNwisAAADiBuEVAAAAcYPwCgAAgLhBeAUAAEDcILwCAAAgbhBeAQAAEDcIrwAAAIgbhFcAAADEDcIrAAAA4gbhFQAAAHGD8AoAAIC4QXgFAABA3CC8AgAAIG4QXgEAABA3CK8AAACIG4RXAAAAxA3CKwAAAOIG4RUAAABxg/AKAACAuEF4BQAAQNwgvAIAACBuEF4BAAAQNwivAAAAiBuEVwAAAMQNwisAAADiRtPGLgARVLRTWjdb2pcntciSTrhQSs1s7KoAAAAixplZY9cQN5xzJkkx95mZSYunyhZOkSsr/mFzQjO5fmOlvmMk5xqxQAAAgNA5b24xsyoBhpHXw4AtmiqXO0mV/3RdWbE0f5LM5AmxAAAAcY6R1zqIyZHXop0qnXqCmlpJjU1KXZKajlnHLQQAACAuBBt55YGtOFf02aygwVWSmlqJ9n32RgNVBAAAED2E1zj31Vf/F1K7TZs2RrkSAACA6CO8xrntSg+pXZ6F1g4AACCWEV7j3O7swSq2xKBtii1RezsObqCKAAAAoofwGudyTj1Rj5cPC9rm8fJh6tfzxAaqCAAAIHoIr3EuIzVJSTljNOXQFVVGYIstUVMOXaGknDHKSE1qpAoBAAAiJ6ypspxzrSWdImmvmf2n0r6jJD0qaZCkMklzJN1tZtvrX27jismpsuSp5/HcjZo5f5X62XJlqVB5StdC10e/GHCabuvfyT/lBAAAQKwLNlVWuOH1bkl/kvRXM7ujwvamkv4j6STJP2e+SfpCUi+zWuZ0inGxGl59CopK9MEX25S3p1hZac006KR2jLgCAIC4E43w+p6kgZL6mtmSCtuvkTRT0gFJf/G+j5WUJulOM3ssjPpjRqyHVwAAgMNBNJaH7eR9X1Np+xXyjLROMLOp3s43SnpR0nBJcR1eAQAA0LjCHXndLUlm1qrS9kJ5RlmPMbOt3m1J8ozAFphZ23pX3IgYeQUAAIi+aCwPm1z5WOdcV0mtJP3XF1y9nZZI8oVaAAAAIGzhhtftklKcc+0qbBvofV9STfvmknaH2RcAAAAgKfzwutz7PlqSnHMpkn4tz/2u8yo2dM4dLU943SoAAACgHsINr0/KMxXW3c65LyVtkHSypB2SXq/Utr/3vfLDXQAAAECdhBVezew9SRPlGWntKqm9pHxJ15jZgUrNr/a+54ZZIwAAACApzNkG/Ac710HS6ZJ2SVpmZrsr7U+SNE6ekPykmW2rR62NjtkGAAAAoi/iixQcqQivAAAA0RfxqbKcc/c450bXof1I59w9dexjuHPuUefcYufcHuecOeeeraFtR+/+ml4v1qVvAAAAxKZwFykol7TNzNqH2H6TpA5mllCHPlZL6iFpn6Qtkk6Q9JyZ/b9q2naUtEnSp5JmVXO6tWb2aqh9B6mJkVcAAIAoi8bysA3hLnlC60ZJ/RTaA1+rzWxiNIsCAABA42mo8Joh6WBdDjAzf1j1pW8AAAAc2aIeXp1zl0tqKWl9tPuS1N459ytJmZJ2SlpqZp81QL8AAABoACGFV+fcnZLurLS5rXPuq2CHSWotKU2e+WDnhFVh3Qzyvn4owrkFkq4zs28i1UlOTk5Yxy1YsCBSJQAAAMS8cDNTMKGOvLaW1LHStoRqttVknqT7Qmwbjv2S7pfnYS1foD5FnoUU+kua55zraWZFUawBAAAAURbSbAPOuR6Sevp+lfRPSbsljQpyWLmkPfI86f9/9SrSuRx5HtiqdraBIMc1lfShPAspjDKzafWsg9kGAAAAoqzesw2Y2afyTEPlO+E/JR0ws2ciVGNUmFmpc+4f8oTXcyTVK7wCAACgcYX1wJaZhbW4QSPZ4X1PbdQqAAAAUG/xFELDdYb3PdjDZQAAAIgD9Z4qyznXRFJneeZyTQzW1swW1be/Gmo4XdJ/zKyk0vYB8ix2IEnVLi0LAACA+BF2eHXOHSXpj5KGS2oewiFWl/6cc5dKutT7azvv+5nOuRnen/PNbIz358mSTvZOi7XFu+0USQO8P//BzJaE2jcAAABiU0izDVQ5yLn2kj6R1F6e2QdCUpd7ZZ1zEyVNCNLkazPr6G37S0lDJXWT1EaeEeA8SUslPWZmi0Ptt5aamG0AAAAgyoLNNhBueH1K0o2S9koaL+kNSd+bWVm9Ko1xhFcAAIDoi0Z4/VaeUdefm9mr9S0wXhBeAQAAoi8a4fWgPLcLpJpZaX0LjBeEVwAAgOgLFl7DnSpruzyLFBwxwRUAAACNL9zwOldSS+dc50gWAwAAAAQTbnh9UFKRPFNUAQAAAA0irPBqZhslXSypn3PuA+dcf+ccy68CAAAgqsJ9YCucKbHMzOq9oldj4oEtAACA6Av2wFa4YTLkhQkAAACASAk3vPaPaBUAAABACMK6beBIxW0DAAAA0ReNeV4BAACABheR8Oo82jjnOkTifAAAAEB16hVenXOnOedel7RbUp6kryrtT3fOPemce8I5l1SfvgAAAICww6tz7heSlkq6VFILeWYgCLgvwcwKJR0n6SZJg8IvEwAAAAgzvDrnTpT0d0mJkqZL6i0pv4bm/5In1F4STl8AAACAT7hTZY2WlCTpcTMbJQVduGC+9/3MMPsCAAAAJIW/wtb/SeooKdvMtni3bZX0IzNLqKb9PkllZtaqfuU2LqbKAgAAiL5oTJXVXlKRL7iG4ICk5mH2BQAAAEgKP7wWS0pyvlgchHOuuaTW8sxIAAAAAIQt3PC6WZ6HtTqH0PZ8SQmSvgizLwAAAEBS+OH1XXlmELgzWCPnXKakP0kySXPC7AsAAACQFH54fVjSPkm/ds5NcM61rLjTOdfcOXe1pBXyzPO6U9IT9aoUAAAAR7ywZhuQJOfchZJelef2gUPyBOEESeskHS/PVFpOnvtjLzSzeZEouDEx2wAAAED0RWO2AZnZbEnnSFopT1BtKk9YPVFSM+/P/5F0zuEQXAEAAND4wh55DTiJc6dI+qk8U2glSNom6SMzW1Hvk8cQRl4BAACiL9jIa0TC65GC8AoAABB9UbltAAAAAGhohFcAAADEjaa1NXDO/dP741YzG19pW12Ymf0yjOMAAAAASSHc8+qcK5dnkYH1ZnZSpW21Lg9boZ2ZWUL9ym1c3PMKAAAQfcHuea115FXSv+QJoFur2QYAAAA0GGYbqANGXgEAAKKP2QYAAABwWCC8AgAAIG6Ecs9rFc65H0m6UtIOM3uhlrbXSMqU9LyZ5YfTH7yKdkrrZkv78qQWWdIJF0qpmY1dFQAAQIMJK7xK+n+SpkiaGELbHpLu9v48Pcz+jmxm0uKp0sIpUlnxD9vfHiv1Gyv1HSO5UCZ+AAAAiG/h3jZwsff9tRDazpRnqqxLwuwLi6dK8ycFBlfJ8/v8SZ79AAAAR4CwZhtwzn0r6ShJzc3sUC1tkyQdkPStmXUMp8hY0SizDRTtlP3lRLnKwbUCS2gmN/pLbiEAAACHhWjMNtBW0q7agqu30xJJuyRlhdnXkW3d7KDBVZJn//o5DVQQAABA4wk3vO6V1Mo5l1xbQ2+bNEn7w+zriLa/4PvQ2u38LsqVAAAANL5ww+vn3mMvDKHtRZISJK0Ls68j2to9tf7/gbdd8yhXAgAA0PjCDa9vyvMQ1lTnXPuaGjnnjpY0VZ6lZGeF2dcRbXXqT1VsiUHbFFuiVqee3UAVAQAANJ5ww+sTkrZIOlbSaufcXc65zs65JO+rs3NutKT/eNt8J+mvkSn5yJKWkaXppUODtpleOlStMrmlGAAAHP7CmufVzPY75y6V9K6kNvKMrlY3X5OTlC/pYjMrCrvKI9jgk9vpjDeHSoekkU3/rWbuh2fkii1R00uH6u9uqD4+qV0jVgkAANAwwpoqy3+wc8dI+qOkKyRV/rvtEkkvShpvZofF00SNMlWWpMfm/1dT39+gdO3RoISVylKh8pSuD8p6qVBpGjO4i24f0LlBawIAAIiWYFNl1Su8VuggRVJvSb7hv62SVpjZgXqfPIY0Vng1Mz2eu1HT529USWm5f3tS0yYaOaCTbuvfyf+HDAAAEO+iHl6PFI0VXn0Kikr0wRfblLenWFlpzTTopHbKSE1qlFoAAACihfAaIY0dXgEAAI4E0VhhCwAAAGhwtc424Jz7yvvjRjMbXGlbXZiZ/TiM4wAAAABJoU2V1dH7frCabXXB37UDAACgXkIJrzd433dXsw0AAABoMDywVQc8sAUAABB99Xpgyzn3unPuH5W2dXDOHR2xCgEAAIAQ1Dry6pwrl7TNzNpX2rbVzI6oAMvIKwAAQPTVd6qsckkJ1Z23fmUBAAAAdRNKeC2QlOmcaxXtYgAAAIBgQpltYLmkIZLecs69KGmfd3tz59y1denMzP5Vx/oAAAAAv1Duee0raZ48QdfX2Knu87aamYUSlmMW97wCAABEX7B7XkOaKss5d4akOyV1l5QizyIFZZK21KUQMzuuLu1jDeEVAAAg+uodXqs5YZUZCI4EhFcAAIDoq+9sAwAAAEBMqPUeVOfcKkk7zOy8Cpv7SyqOWlUAAABANVikoA64bQAAACD66nvbQJmkxOrOW7+yAAAAgLoJJbzukJThnDsq2sUAAAAAwYQy7+qHkoZLWuCce0M/LFLQwjl3T106M7P76lgfAAAA4BfKPa/dJH0kqaXqt0iBzCyhrsfEEu55BQAAiL5g97zWOvJqZmudcz0k/Uo/LFKQI+mQpKWRLBQAAAAIhkUK6oCRVwAAgOir18hrDb6RlFePmgAAAIA6C2vk9UjFyCsAAED0RX15WOfRxjnXIULnG+6ce9Q5t9g5t8c5Z865Z2s55izn3NvOuQLn3H7n3GfOuVHOubh+SAwAAAA/qFd4dc6d5px7XdJueW4j+KrS/nTn3JPOuSecc0l1OPXvJd0uqaek70Ko4xJJiySdI+nfkh6XlCTpYUkv1qFfAAAAxLCww6tz7hfyzDZwqaQW8kyfFTC0a2aFko6TdJOkQXU4/V2SukhKk3RLLXWkSfq7PCuB5ZjZL81srDzBd6mk4c65K+vQNwAAAGJUWOHVOXeiPIExUdJ0Sb0l5dfQ/F/yhNpLQj2/meWa2X8ttJtLh0tqK+lFM1tR4RwH5RnBlWoJwAAAAIgP4c42MFqev5Z/3MxGSZJzrqyGtvO972eG2VdtBnjf361m3yJJ+yWd5ZxrZmbFUaoBAAAADSDc8DpAnhW2JtfW0My+d87tlxSRh7mq0dX7vqGavkudc5sknSzpeElfRqLDnJycsI5bsGBBJLoHAACIC+FmpmDCvee1vaQiM9sSYvsDkpqH2VdtWnnfd9ew37e9dZT6BwAAQAMJd+S1WFKyc87Vdl+qc665PMFxV5h91ZfvIbKITc7KCCoAAEDtws1MvnleqxPuyOtmeR7W6hxC2/MlJUj6Isy+auMbWW1Vw/60Su0AAAAQp8INr+/KM6J5Z7BGzrlMSX+SZ9RzTph91Wa9971LNf03lWeqrlJVmoMWAAAA8Sfc8PqwpH2Sfu2cm+Cca1lxp3OuuXPuakkr5AmPOyU9Ua9Ka+abzWBINfvOkZQiaQkzDQAAAMS/sMKrmeVJulrSIUn3SNohKVOSnHOfSyqQNFNStjz3x15lZnsiUXA1XpVnjtkrnXO9fRudc8mSJnl//VuU+gYAAEADcqGtA1DDwc79RNJj8ixSUJ3/SPq1mS2v43kvlWflLklqJ+k8ef7af7F3W76ZjanU/lVJB+VZDrZA0sXyTKP1qqQrQlzwoLa6TJIicCoAAADUwPfAlplVeXKrXuG1QgenSPqpPFNoJUjaJumjiite1fF8EyVNCNLkazPrWOmYsyWNl2cxhGRJGyX9U9J0M6tpAYW61kV4BQAAiLKoh9cjBeEVAAAg+oKF13Af2AIAAAAaXLiLFPg555IkDZLnvtcfyTMt1g5JyyXNNbOS+vYBAAAASPUMr865myXdL6lNDU3ynXO/N7O/16cfAAAAQKrHPa/OucmSxuiH5Ve/k7TF+/Mxko72/mySppjZb+pRZ0zgnlcAAIDoi/gDW865fpJyvb++JukPZrauUpuu8ozKDpcnwOaY2WLFMcIrAABA9EXjga3bvO//a2aXVw6u3s7Wm9kVkv5XntHZ28PsCwAAAJAU/sjrFnkWD2hvZttraZsl6XtJW83smLCqjBGMvAIAAERfNG4bOCipyMwyQ2y/U1KqmSXXubMYQngFAACIvmjcNrBXUkvnXK1h1DnXXFJLSfvC7AsAAACQFH54/UyeZWBHhNB2hDxTcn0aZl8AAACApPDD63PyPIT1Z+fcL2tq5Jy7UdKf5ZltYGaYfQEAAACSwr/ntYmkeZL6yRNMt8gzddZ33t+PldRfnrlenaQFks61OL9ZlHteAQAAoi/iD2x5T5om6Z+Shnk3VT6Rr7PXJP3SzPaE1VEMIbwCAABEX1TCa4WT95F0paTekn7k3bxd0gpJL5rZ8np1EEMIrwAAANEX1fB6JCG8AgAARF+w8No0zBMmSTpBUkl1q2tVanuCpCRJX5rZoXD6AwAAAKTwZxv4uaT/SBoVQtvx3rbDw+wLAAAAkBR+eL3M+x7K9Ff/K8/DW4RXAAAA1Eu44bWb9z2UhQdWet+7h9kXAAAAICn8eV73SzpoZhkhti+QlGhmLevcWQzhgS0AAIDoC/bAVrgjryWSmofYufO2JfEBAACgXsINr5skJTnnzgyh7VmSmkn6Osy+AAAAAEnhh9cP5HkI6yHnXI3TbXn3/VGeUdf3w+wLAAAAkBR+eJ0u6aCkn0qa65w7tXID59xpkuZ52xRLmhZukQAAAIBUjxW2nHO/kDSjwqZt8twaYJKOk5Qlz+isSbrOzJ6tV6UxgAe2AAAAoi9qy8M6586X9JikjjU0+UrS7Wb2btidxBDCKwAAQPRFLbx6T54gqb88D2a1827eKmmJpFwzK69XBzGE8AoAABB9UQ2vRxLCKwAAQPRFY55XAAAAoMERXgEAABA3CK8AAACIG4RXAAAAxA3CKwAAAOIG4RUAAABxg/AKAACAuEF4BQAAQNwgvAIAACBuEF4BAAAQNwivAAAAiBuEVwAAAMQNwisAAADiBuEVAAAAcYPwCgAAgLhBeAUAAEDcILwCAAAgbhBeAQAAEDcIrwAAAIgbhFcAAADEDcIrAAAA4gbhFQAAAHGD8AoAAIC4QXgFAABA3CC8AgAAIG4QXgEAABA3CK8AAACIG00buwDUU9FOad1saV+e1CJLOuFCKTWzsasCAACICmdmjV1D3HDOmSTFxGdmJi2eKls4Ra6s+IfNCc3k+o2V+o6RnGvEAgEAAMLjvBnGzKqEGUZe45QtmiqXO0mV/0RdWbE0f5LM5AmxAAAAhxFGXusgZkZei3aqdOoJamolNTYpdUlqOmYdtxAAAIC4E2zklQe24lDRZ7OCBldJamol2vfZGw1UEQAAQMMgvMahr776v5Dabdq0McqVAAAANCzCaxzarvSQ2uVZaO0AAADiBeE1Du3OHqxiSwzaptgStbfj4AaqCAAAoGEQXuNQzqkn6vHyYUHbPF4+TP16nthAFQEAADQMwmscykhNUlLOGE05dEWVEdhiS9SUQ1coKWeMMlKTGqlCAACA6GCqrDqImamyvDU8nqK5VykAACAASURBVLtRM+evUj9briwVKk/pWuj+f3v3HmZXVR98/PubSTIkkyu5EIRggCQgggqi3ASSKGBfKNXWC+0DVaut9KFaWy+9PPV9tbbvYx9t66u0Sq2VV2yV1r6FilBAIAG8C1JECCSQcM09k8tMSIbMrPePvU/mzJlzzpyZzG2ffD/PM8+Zs69rr7P23r+zztprvY6rVp7BNSuWHOxmQpIkqUjqdZVl8DoEEyl4LdnR1c2dj25i8+79HDWzjYtOWWiNqyRJKjSD1xEyEYNXSZKkZuMgBZIkSWoKBq+SJEkqDINXSZIkFYbBqyRJkgqjqYLXiNgQEanG36bxTp8kSZIOzaTxTsAo2AV8rsr0zrFOiCRJkkZWU3WVFREbAFJKi0dp+3aVJUmSNMrsKkuSJElNoRmbDbRFxJXAcUAX8DBwb0qpZ3yTJUmSpEPVjMHrQuCGimnrI+I9KaXVI7GD5cuXD2u9VatWjcTuJUmSCmG4MVM9zdZs4KvAG8kC2HbgNOA6YDFwW0S8evySJkmSpEPVVA9s1RIRnwU+DNyUUnrrIWzHB7YkSZJGWb0Htg6X4HUJsBbYkVKaewjbMXiVJEkaZfY2AFvy1/ZxTYUkSZIOyeESvJ6Tvz41rqmQJEnSIWma4DUiXhkRR1aZ/nLg2vzt18c2VZIkSRpJzdRV1tuBP46Ie4D1wB7gROBS4AjgVuCz45c8SZIkHapmCl7vAU4CTidrJtAO7ATuJ+v39Ybkk1aSJEmFdlj0NjBS7G1AkiRp9NnbgCRJkpqCwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwJo13AqSRtqOrmzt+sYkte/azYEYbF79yIUe2TxnvZEmSpBEQKaXxTkNhREQCMM/GV63gNKXE392zjs/fvY7uA70Hl58yqYUPrlzCNSuWEBHjmPLhMRjvz/wononwmU2ENEgTSfk50T6lFQK69vdMmPOjdL9OKQ24cTdV8BoRxwJ/DrwZmAtsBG4CPplS6hiB7Rc7eO3aDmtugc7NMP0oOPkyaJ/bf/rkdoiAzi2wZyPMOJquKXN4fFMn+/fuhsntREDq7qJ15kKWXnAFc+YfXXcfHXu7WXvvN+nZvYm2aTM4aeFM2tnXt6/uzr70QNU0dmzdyBOrv8EzT6/n8Y5Eb0q0s48tzGFVvI7fPv8Ejt+2ioceXUMnUwmy+aX/57GLFcf08oqly6B9/sH9dtGWHdvOTUzq2syB9qNonXfCwOMCOrZuPHgcVY+9kbyuNr9Gnq/ZuIennt988HjnsYsFsZPtMYdTli7hDUvnEd1d1fOx1r6qfe471h/cb3neNJzu7s6hpSFfdm/HRjY/v55tMZd90xfx9IIV7Egza144S19Obrj7QZann7CADjqZyqQWeOMJ7Zx10nHZxa5GGkqfX3SsZ17awcJjjmfqUUuqnwdVyuPBsrJ3d//Pfzif9WBlvsZ5MNj5VZ6u8hvTsW17uaT1p7R3b2+8jJSnodpnXVZea5Wb1LWNe276Kk+seZh5qYMtaTZbmd33mb3qFKIyLSOhrGynPRt5oOMIbt/QS09vb7/rxlUrz+Cas+YQa77TWP4P8pl2TZnL7T1n8tz+aQfzPDo2HCznac7i+tfMWudU6Zzp3MUju4/gofY3MG/6lKqf6ZCuU43mYyP5MJz5jey3Wnk78vihn2cjka7hprHatbba/KGWtyr3r7bZR/XdY2tsK625hR89/Ch3PbX34DlRft8sv6++d96jxKHk0yE4LILXiDgR+D6wALgZWAO8HlgBPA6cl1Lafoj7KE7wWl7YJ02D9atJT95FpJ6DiyRaiTmLYNdz0HtgWLt5KbWyru0VtM1dzKLW7Ux+4Sf9ttVLQIKWaCTPWiASlOVvopVdRyyk/cVNTI6eqmv15OW6taF9NOal1MrGOWew6PVvAeDZH9/E0R0P9kvDgdTCs0eezeLXXUq8tDe7qACsvZ309PeIsnxI0UqcuBJOWA5r74CnvzfsPB9UTIJFr4dZi2DXs/Dcjyv21QpzFpF2PtuvPNTc1uLzYMnFsH41PHUX9A6yDkC0whCP96XUyo96T+au3tOZ1BL9A9LOLTz2xOPs3vQkZ7Ssq1kWqqUhHX8hz/7kZo7ueIDJ0TtgsdrnwcDy2D+9LeyYtID5vVtpKc/HimOvLAsViSRFEKkvXYkWYsB+s8+sMo3Vzq+XUiuPH3EaN+89jZ7eXla2PMRZLWv651m0sn/xcn4x5XTmbV7Nsbse7H8Mgxz74LL09nQ8QysD87x/FuRlbOkl9YPiQQLorilzeHzjHuZtWsWxOx+gZZD9vpRaeD7NY1Hrdlr7HfvAz6QnWnl+5hlsW3ghp3b/jClPr6p5HpS2e0xsq1reDl5bXvcrxPp7Gz+nylS75tW6Vh5ILWyYfTa7jj6P7hf39KuAOFiZ0L2jf54DrL194Hk7SNlOMYmod+2p81n3+/wGlMcqKq5Nlfe3ast2HbeSXY9+l4Vbv08LZffCmESUp6vel7RaeVNV/Wtt6dqTdj03tHwcgt5o5bmysjt5w6rBr/vUuK+2tsGFH4XzP5LlzSg7XILX24GLgQ+mlL5QNv1vgD8ArkspXX2I+5j4wWtKcN9nYfVnoGf/eKdGkiQ1k5V/Bhd8dNR30/TBa0ScADwJbABOTKnvK3NEzCBrPhDAgpRS1yHsZ+IHr/d+Bu7+i/FOhSRJakatbfCHj416E4J6wWuzdJW1Mn+9ozxwBUgp7QG+B0wDzh7rhI2pru1ZjaskSdJo6NkPj39nXJPQLF1lnZS/PlFj/lqyJgXLgLsOdWfLly8f1nqrVq061F3Xt+YWmwpIkqTRtWdTw4sON2aqp1lqXmflr7tqzC9Nnz0GaRk/nZvHOwWSJKnZzVg4rrtvlprXwZTaS4xIY9VRr0EdrulHjXcKJElSM2ttg5MubXjx4cZM9fplb5aa11LN6qwa82dWLNecTr4sK1SSJEmj4cKPjml/r9U0S/D6eP66rMb8pflrrTaxzaF9blaoJEmSRlJrW9ZN1vkfGe+UNE1XWScC66jfVVYLML/pu8pqoJ/Xqh1cxySenXU6X9v+Cnp6E/PZyfzYyY6Yw4knnMicaZPp3rsbpkyHlJj+zF2cvP/hIXX63kMrm+afw+xXXkJ77Du4reodPkfWCXJZ5xHlnS3P2vx9Xt7xQyaV7b/Uyf13e89gUsDJRwZrOhhwPK9YuoTzl84nujvp5Aie2LyH/R2baO3azMLezRzb9QgtqXqH0L0EAUQDLVB6UrZsrQEaSuldxWv5gzcty0bKyUcy6ZxyZJaurt20TZtBx94DrH3qSeamDram2WxhNpNagjedMC3ryB+yfNxwP1TtgLp6J/cvpVYe6F3Cs2lBX94smUcM0gn3gdTC6p7TuD+9ium8yOuWLeKC1l/U7XC9ehnr4LjYxmtb1zKJxjtqr9zWdF5kD0dwQcsjXNDycL9yUa43Wnlh+ml0z1jEopZtAwbVqJVPWaf+K3ik7TVMeeq7nLz/5/06gu9JQUvULhd9ZfN0Lmx9hAtbft6vk/S+ju23D9gu9D9Xh7JsvzyrOvhBeRpb6Jp6NLO6N/cfVKHUMf2JK7PO26dMz6Z3bukbEWj6Au59Yist624fMCBCqYzNPHoJpyxdBtMX1DnvR87QyvbAQRlq5XM/BzvtX1FjcJJJvPSyM3kuzWdK53O8bM/Pa15bemjl3t7TuLfntLw8Ty27jiXO5ecDynata+1QrlODG3ywjnqDMtQ8pwbRN1DMWw8OVMKeTbDz6bqd9peuTT+I0/jNuY+zaPfPag4ScsgD3FQddKF6GvvK43yOi60DBlsZLB+rX/Om0hpw6qx9HDNpDz3TFtA2eyHLFkxnx39/e8DAOtXz6VW86YR2zlq2iKg2EE3pGE/9NUfYGg0OUlCha3vWlcWeTQdvNtnQglN5qP085s+YwkUtDzC9e1vW8PqkS6F9Lju6urnz0U1s3r2fo2a2cdEptcc3zoYhvJHoeIq5qYOjj1nM1IVLD26rXxrK9jFoekvLQt31S/vv2b2RA9MW8Mz8FWxPM/qleyjHUy/v6O6snq7y+WX/d06Zx529r2Xrnm6OePI2tr7wNHuYSgDTeZHNzOHOntfSwUw+cvEyfm/l0gHJqNTQsZTSvv2pvsBi7glVP5NSGp/dN6369sqWTZPbuW/ddh5c+wwv9M4+mPYpk1r44MolXLNiSXahGSzvapUx9tTP0/JAqex4yre1YMYUntnxIv9+/39zYfoJR9GRBQClIL/aUKS1yuggZbe87LXOPJqlF7yTOdOmVD2Gqvlcdrzl8xcdsZeLWh6gpXMzf3HvDm576QwALmp9gKPo6FdujprUyV2XdjG9extp+lF8ZdspfPm+p+ofe1nZPfjFrWt33zGUhrtt9LwtUz58bykNm5nD6tIwrKUyUq28ludZRVA8oCzUKBflX/hqXQ/qfu7leVP2mcyN3Ry3dTWT9m6mbdrMfOjNF6vnzWB5N8j5sYMZVc/xUjnftX0zr+n6HqfOfJFpc4+pX15rfNZMmU4Aqbuz73jKvjj3y/NBrnnVyu70zmcGvfYM9vn1K4/VDHJ/mzX3qL7PvPyal+dD1+6dvNAzm81Hv5FFc6f23QsHu/ZUy5sG7mmlfHpidxtb9uxjwYwjWDZz/8H9DiUfh3OP7tm98eBnHd2d1fOp0TI8hg6X4LVyeNjHgLPIhod9Ajj3sBoeVhNG6ab++bvX0X2g7xv1gMCvAIb1ZWCMFSGNjbj27rV89o7aLZ2qfemZCMc+EdIgqfgOi+AVICIWAX8OvBmYS9Zc4CbgkymlHSOwfYNXDZs3dQ1FM33pkaShOmyC19Fm8CpprPmlR9LhyOB1hBi8SpIkjb56wWuzdJXV1JYvXz4qw6sVmXnSn/nRn/kxkHnSn/nRn/kxkHnS30TKD4NXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhREppvNNQGBFhZkmSJI2RlFJUTrPmVZIkSYVhzaskSZIKw5pXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJ7CIODYi/ikiXoiI/RGxISI+FxFzxjttoyEi5kbE+yLiPyJiXUS8GBG7IuL+iHhvRLRULL84IlKdv2+O17GMpPxzr3WMm2qsc25E3BoROyJib0Q8HBEfiojWsU7/SIqIdw/ymaeI6ClbvmnKSES8LSK+EBH3RcTuPP1fH2SdIZeDiLgsIlbl515nRPwoIt418kd0aIaSHxGxNCL+KCLujohnI6I7IjZHxM0RsaLGOoOVtatH9wiHZoj5MezzIiLeFRE/zsvGrrysXDZ6RzZ8Q8yT6xu4ttxVsU5hykgM8f5att6EvIZMGsmNaeRExInA94EFwM3AGuD1wO8Db46I81JK28cxiaPh7cAXgY3APcAzwFHArwL/CPxSRLw9pZQq1vtv4KYq23tkFNM61nYBn6syvbNyQkT8CvDvwD7gRmAH8MvA3wLnkeVzUT0EfLLGvPOBlcBtVeY1Qxn5M+DVZJ/5c8DJ9RYeTjmIiN8DvgBsB74OdANvA66PiNNSSh8ZqYMZAUPJj08B7wQeBW4ly4uTgMuByyPi91NKn6+x7s1k5a7ST4eZ7tEypPKRG9J5ERGfBT6cb//LwBTgCuDbEfGBlNK1w0j3aBpKntwEbKgx7yrgBKpfW6AYZWTI99cJfQ1JKfk3Af+A24EEfKBi+t/k07803mkchWNemZ8YLRXTF+YnWgJ+rWz64nza9eOd9lHOlw3AhgaXnQlsAfYDZ5ZNP4Lsy1ACrhjvYxqlfPpBfnyXN2MZAVYAS4EAlufH9fWRKgd5Xu0ju+ksLps+B1iXr3POeOfDMPPj3cDpVaZfSHZz3Q8cXWWdBLx7vI91FPJjyOcFcG6+zjpgTsW2tudlZ/GhHMN45kmdbcwG9uZlZF5RywhDv79O6GuIzQYmoIg4AbiYLGj5u4rZ/wvoAq6KiPYxTtqoSindnVL6dkqpt2L6JuBL+dvlY56wYnkbMB/4Zkrp4Lf+lNI+sloIgN8dj4SNpog4FTgbeB74zjgnZ1SklO5JKa1N+d1gEMMpB78FtAHXppQ2lK3TAfzv/O2E+Rl0KPmRUro+pfSzKtNXA6vIahDPHflUjp0hlo/hKH32f5mXidJ+N5Ddp9qA94zSvodlhPLkKmAq8P9SSttGKGljbhj31wl9DbHZwMS0Mn+9o0pB2xMR3yMLbs8G7qpcuUm9lL8eqDLvZRHxfmAu2Te+H6SUHh6zlI2Ntoi4EjiO7MvLw8C9KaWeiuVKZee/qmzjXrIahHMjoi2ltH/UUjv23p+/fqVKnsDhUUbKDacc1Fvntoplmkm9awvAayLiQ2Q1Ts8D96SUnhuTlI2+oZwXg5WPj+fL/K8RT+X4+u389R/qLFP0MlLtHJjQ1xCD14nppPz1iRrz15IFr8s4DILXiJgE/Gb+ttpJcVH+V77OKuBdKaVnRjd1Y2YhcEPFtPUR8Z689qikZtlJKR2IiPXAK8nabz02KikdYxExFbgS6CVru1XN4VBGyg2nHNRbZ2NEdAHHRsS0lNLeUUjzmIuIlwNvJLsR31tjsd+veN8TEf8IfCivhSqyhs6L/Fe+Y4DOlNLGKttZm78uG6V0jouIOAc4DXgipXRPnUULW0bq3F8n9DXEZgMT06z8dVeN+aXps8cgLRPBp4FTgVtTSreXTd9L9iDGa8na1Mwha8N2D9nPH3c1SdOKr5LdYBcC7WQX0+vI2hfdFhGvLlv2cCw77yA7nttSSs9WzDtcykil4ZSDRteZVWN+oUREG/DPZD9zfqL8p/DceuADZDfkduBlZGVtA1lN/z+NWWJH3lDPi8PxugLwO/nrl2vMb4YyUuv+OqGvIQavxRT562i1bZowIuKDZE+3riFre3RQSmlLSul/ppQeTCntzP/uJauV/hGwBHjfmCd6hKWUPpm3V9qcUtqbUnokpXQ12cN7U4FPDGFzzVh2SjeY6ypnHC5lZBiGUw6apuzk3fzcQPbE9I3AZyuXSSmtTildm1J6Ij/vNqaU/o3sIaAO4NcrvjgWxiieF4UvGyURMYssEO0Grq+2TNHLSL37ayOr56/jcg0xeJ2YBvt2MrNiuaYUEdcA/4ese5sVKaUdjayXUjpA38/HF4xS8iaCUiP78mM8rMpORJxC9qDNc2RdIDXkMCgjwykHja6z+xDSNe7ywPXrZN38/Ctw5VAe6Mlr90tlranKTp3zYrCyMViNWxFdCUxjGA9qFaGMNHB/ndDXEIPXienx/LVW+6Gl+WutNrGFlzd+v5asv8EV+RORQ7E1f23Gn4RLtuSv5cdYs+zkbZuOJ2uU/9ToJm3MDPagVj3NXEaGUw7qrXM0WT49V+T2rvmxf4Osb9J/AX4jD9iGqpnLzoBjSyl1kT2IND0vC5Wa8Z5UelBrwC86DZqwZaTB++uEvoYYvE5MpYbhF1eOehERM8h+6noR+OFYJ2wsRMQfkXWC/BDZibVlkFWqOTt/bZYgrZpz8tfyY7w7f31zleUvIKtJ+H4z9DQQEUeQ/dTVC3xlGJto5jIynHJQb51fqlimcCJiCvAtshrXrwFXDeMLT8lZ+Wszlp1a50VTl49yEXEW2eAGT6SUVg1zMxOyjAzh/jqxryFpAnSe61/VDoUPu0EK8uP7eH58PwWOHGTZs4ApVaavJOsoOQHnjvcxHWJ+vLJaPgAvJ3vCNwF/WjZ9Jtk3/qYfpIAscE3Atw+3MkJjgxQMqRyQ1aQUZpCCIeZHG1n/v4nsZ/GWBrZ5fpVpAfxJvp2twMzxPvZh5seQzwsKOEjBUPKkYtmv5Mt+uJnKyBDvrxP6GhL5hjXBVBke9jGyC84Ksp9mzk1NNjxsPvbx9UAP2fBy1dpPbUgpXZ8vv4osuFtF1uYR4FX09SP38ZTSX4xagsdARHwC+GOy2vj1wB7gROBSsovIrcBbU0rdZeu8hayGaR/wTbIh/S4neyL2W8A7UhOc+BFxH/AGshG1vl1jmVU0SRnJP9e35G8XApeQ1ercl0/blsqGXhxOOYiIDwCfJ7v53Ejf0I7HAn+dJtDwsEPJj4j4KtloSNuAv6f6AyOrUlktW0QksmvtT8h+Mp9F9qvXqWRP6781pXTHiB7UIRhifqxiGOdFRPw18If5Ot8iG9zhnWT9xE644WGHes7k68wEXgAmA8ekOu1di1RGhnp/zdeZuNeQ8f4m4F/dbz6LyLpJ2pgXgKfJGljX/cZU1D+yp+bTIH+rypZ/L3ALWbcknWTfEJ/JT5gB34iL+EfWfc03yJ4G3UnWmfRW4E6yvvmixnrnkQW2HWRNTH4O/AHQOt7HNEL58oq8PDxb75iaqYw0cH5sGIlyQDaE5GqyL0pdZDfmd4338R9KfpAFaYNdWz5Rsf3P5PnwAtnNe29+Hl4LnDDex3+I+THs8wJ4V14muvIyshq4bLyP/1DzpGyd383nfaOB7RemjDSQF/3ur2XrTchriDWvkiRJKgwf2JIkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJBRARyyMi5aP6jOR2F5e2GxGLx3p9SRoqg1dJkiQVxqTxToAkqSF7gcfHOxGSNN4MXiWpAFJKPwZOHu90SNJ4s9mAJEmSCsPgVdJhLyIWRMRL+UNHlw+y7Kfy5daVTTsuIq6JiO9ExBMR0RURnRHxaER8LiKOq7O9Vfn2PhERkyPiwxHx04jYmU9fni9X84GtiGiJiPMi4tMR8cOIeC4iuiNie0SsjoirI2Jyg3mxNCKuz7exPyKeiYgvRcQxjaxfZ7tviYibIuKFPG0dEXHvYGmLiHdExG0RsTn/jHZGxNqI+M88z484lHRJKp5IaUQfXJWkQoqIW4BLgW+llN5eY5kAngSOBz6RUvpkPn0VcGHZoruAGfRVEOwCLksp3V9lm6V1/wo4HzgXOADsAeYAK1JKq/Ig9h6AlFJUbGMxsL5s0gGyNrIzy6bdB1ySUnqxzrpXAF/O094JtAJT83k7gItSSg/WWf/4lNKGivnTgW8Al5VN3p3vo3QcPwAuTSl1VKz7FeC3yiZ1kuXptLJpA/YpqblZ8ypJma/lr78cEbNrLHMeWeAKcEPZ9EeAPwZOAaallGYDbcBZwH8Bs4AbI2IqtV0DvAp4DzAzpXQkMA94uIG0HwBuBt4JHAO0pZRmkQWI7wFeIAuM/3KQ7VxHFoielVKaAbQDlwDPAEcC/xERMxpIT7kbyALXdcBvkB3bLLIA9FeAp4BzgH8qXyki3kAWuPYCfwTMTSnNSCm1k+XLJcD/BbqHmB5JBWfNqyQB+c/Pm8gCzfenlP6hyjLXAb8D3J9SOr/B7bYCD5IFplellL5eMX8VfbW2l6eUvl1jO8upUfPaQBrOBH4CdAHzUkr7yuYtpq/mdDtwSkppS8X6rwAeAqYAH0spfabG+v1qQSPiUuAWsnw9M6X0fJW0HQusIQuUT08pPZRP/xhZbfQdKaVLhnK8kpqbNa+SBOQB3b/lb6+qnB8RbcA78rc3VM6vs90estpXgDfUWfQXtQLXQ5VS+imwhSxAfE2dRb9UGbjm6z8GfCt/e8UQdv2+/PWGaoFrvu3nyINystrUkp356/z8C4AkAQavklSu1HTgvIg4vmLeZcBsYD/wr5UrRsT5+YNOa/KHtUqjTiXgY/lix9bZ9/cOJeERMSV/+OmO/KGofRVpWNBAGu5uYN6rGn34i75g/XciYlOtP+BN+XIvL1v3u8A+4HTgvoh4b5XPRNJhyH5eJanP/WQ/gR8PXAl8qmxeqTb2P1NKO8tXioi/oi9ABegBOuhrjzmdrNazvc6+B9R4NioiFpAFe6eVTd4HbMvTAjCfrMKiXhqq1o5WzJtE1v518yBpmkzWNhWyphiz6i2fO/ggVkrpqYh4H/Alsjax5+Tb3UpWU/svZJ+Fbd+kw4w1r5KUywOhUpvUg00HImIu8D/yt18rXyciLqIvcP17sgCyLaV0ZEppYUppIfC3pcXr7L6nzrzB/G2+3+1kDzkdnVKamlKaX5aGFxpIw0gGguU/9V+RUooG/t7dLzEp/TNZbezVwI3As2RB+DuAm4DVEVHeo4Kkw4DBqyT1VwpOl0bE2fn/7wQmA1vpa79aUmoDentK6ZqU0iN5wDQo6AAAAwNJREFUO9dyC0cnqQdrOH81f/t7KaWvppQ2VSzTSl8taD31mhSU+nk9QNZtVl15G+Jd+dvT6i07yHZ2pJSuSyldkVI6DlgCfJos0D4f+MRwty2pmAxeJalMSmkdWb+j0Ff7Wnr9RkrpQMUqi/LXn1XbXt437MoRTWR/84FSR/1V00DW9rSRzvxXNDDv4ZTSSw2mrdSO9+0RMSL3m5TSkymlPyFrNgBw0UhsV1JxGLxK0kCl2td3RsQpwNkV08uVahdfXWNbVwMnjGDaKu2m7+f+AWmIiEkM3r9rydURMaCGNiJOAt6Wv71xCGkrdTe2DPhovQUjoj0ippS9bxtk26XBFg6luYWkAjJ4laSBbiR72GouWUf4AI+llB6osmypGcEvRcTHI6IdICJmR8SfAl8ga4s6KlJKnfTVcP5NRKws1XJGxKnArcCZZH28DmYycGdEvC5fPyLiTcDtZIMuPEv2AFWjabsZ+I/87acj4osRsaw0P+8h4az8gben6esRAeDaiPjXiPi1/IG00jrTI+Jq4DfzSbc2mh5JzcHgVZIq5MOU3pK/PTN/rVbrWpp+X/7/nwN7ImIHWcD6l2TB7RdHKaklHyILTo8B7gL2RsRu4OdkP/f/NlnPA4N5P3Ai8OOI2EM2HOudZA9N7QR+NaW0e4hpuxL4Zv7/1cDjeVdiO8hqT39I9sDbXPo/MDYZeDtZ/7KbI2JPRHSQDZv7RbIBE+6n8VplSU3C4FWSqisPVnvp64Wgn7z958XAJ4EngJfInuj/MfC7wOWM8k/beY3w68n6n91Gdm3fk78/N6XU6KAKPyIL1r9G1hxiElkXWV8GTssHOxhq2vamlH6dLIi+gWw42Bay7sO2kPUf+zFgacVABp8CPkhWc7uG7EGx0jp3kvWqsDyl1EiNsqQm4vCwkiRJKgxrXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBXG/weKp2icsQcBlgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 864x648 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Generate some sparse data\n",
    "np.random.seed(711)\n",
    "\n",
    "# number of samples is 120 and the dimension of varaibles os 200\n",
    "n_samples, n_features = 100, 200\n",
    "X = np.random.randn(n_samples, n_features)\n",
    "\n",
    "# generate coefficients\n",
    "LargeCoef = 5\n",
    "SenCoef = np.zeros(n_features)\n",
    "for i in range(LargeCoef):\n",
    "    SenCoef[i] = 5 + 20*np.random.rand()\n",
    "for i in range(LargeCoef, n_features):\n",
    "    SenCoef[i] = 0.1*np.random.rand()\n",
    "    \n",
    "# generate y\n",
    "y = np.dot(X, SenCoef) + np.random.normal(loc=0, scale=0.001, size=n_samples)\n",
    "\n",
    "# cross-validation and fit\n",
    "reg = LassoCV(cv=50).fit(X, y)\n",
    "print(\"Residue of the fiting\")\n",
    "print(reg.score(X, y))\n",
    "print(\"Weight\")\n",
    "print(reg.alpha_)\n",
    "\n",
    "# plot the results\n",
    "fig, ax = plt.subplots(figsize=(12, 9))\n",
    "fig.subplots_adjust(bottom=0.15, left=0.2)\n",
    "\n",
    "ax.scatter(np.where(reg.coef_)[0], reg.coef_[reg.coef_ != 0], label='Lasso coefficients', linewidth=4)\n",
    "ax.scatter(np.where(SenCoef)[0], SenCoef[SenCoef != 0], label='Actual coefficients', linewidth=4)\n",
    "\n",
    "# ax.set_xlim([0,n_features])\n",
    "# ax.set_ylim([0, 0])\n",
    "ax.set_xlabel('variables', fontsize=25)\n",
    "ax.set_ylabel('coefficients', fontsize=25)\n",
    "ax.legend(fontsize=24,loc='upper right',frameon=False)\n",
    "\n",
    "ax.tick_params(which = 'both',direction='in',colors='black',\n",
    "               bottom = True,top=True,left=True, right=True,pad=15)\n",
    "ax.tick_params(which = 'major',direction='in',length=15,labelsize=20,width=1.5)\n",
    "ax.tick_params(which = 'minor',direction='in',length=6,width = 1.5)\n",
    "\n",
    "for axis in ['top','bottom','left','right']:\n",
    "    ax.spines[axis].set_linewidth(2.0)\n",
    "    \n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "nbpages": {
     "level": 2,
     "link": "[5.2.2 An example of the lasso regression ](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2.2-An-example-of-the-lasso-regression)",
     "section": "5.2.2 An example of the lasso regression "
    }
   },
   "source": [
    "Then, we used a model with 20 large non-zero ceofficients (between 5~30) and the rest coefficients are 0.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "nbpages": {
     "level": 2,
     "link": "[5.2.2 An example of the lasso regression ](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2.2-An-example-of-the-lasso-regression)",
     "section": "5.2.2 An example of the lasso regression "
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Residue of the fiting\n",
      "0.9999855603124613\n",
      "Weight\n",
      "0.033012951782999095\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAq8AAAIlCAYAAAD7ShEjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzs3Xl8lNXZ//HPyQ4hgSRgwIWgZXFDrYCKPkCIgFi1igt1+VWRqk+tFhGh1NoqWrQi1AWlrbbPI+67xX1BdgVEoSiogDwCikJYErYAWa/fH5MZM1kmM5OZZAa+79drXpPc9zn3ue6bAa6cnMWZGSIiIiIi8SChpQMQEREREQmWklcRERERiRtKXkVEREQkbih5FREREZG4oeRVREREROJGUksHEE+cc1qaQURERKSZmJmrfUw9ryIiIiISN9TzGobmXhs3Pz8fgLlz5zZru7FMz8Sfnoc/PY+69Ez86Xn40/OoS8/EX3M/D+fqdLj6qOdVREREROKGklcRERERiRtKXkVEREQkbih5FREREZG4oeRVREREROKGklcRERERiRtKXkVEREQkbih5FREREZG4oeRVREREROJGTCavzrkc59w1zrl/O+fWOuf2Oed2Ouc+dM79yjmXUKt8F+ecBXg931L3IiIiIiKRE6vbw14C/B3YBMwBvgVygQuBfwFnO+cusbr7tH4GzKjneiujGKuIiIiINJNYTV7XAD8H3jKzKu9B59wfgCXARXgS2Vdq1VtuZhOaK0gRERERaV4xOWzAzGab2Rs1E9fq45uBf1R/m9/sgYmIiIhIi3J1f/Me25xz44D7gAfN7ObqY12AdcBMPL2xOcB2YJGZfR7Btg0g3p6ZiIiISDxxzgFgZq72uVgdNlAv51wScGX1t+/WU2Rw9atmnbnAVWb2baTiyM/PD6ve3LlzIxWCiIiISMwLN2cKJCaHDQRwL3A88LaZvVfj+F7gz0AvIKv6NQDPZK98YJZzLr15QxURERGRSIubYQPOuVHAQ8Aq4AwzKwqiThLwIXAqMNrMHmpiDBo2ICIiIhJlgYYNxEXPq3PuBjyJ65fAwGASVwAzq8CztBZA/yiF13xKtsPSJ2DefZ73ku0tHZGIiIhIs4r5Ma/OudHAA3jWaj3TzLaEeImt1e/xO2zADBZMgXmTobL0x+Nvj4MB46DfWHB1fjAREREROeDEdM+rc248nsR1OZ4e11ATV4DTqt+/iVhgzW3BFJg90T9xBc/3syd6zouIiIgcBGI2eXXO/QnPBK2leHpctwUoe6pzLqWe4wXAzdXfPh2VQKOtZLunxzWQeZM1hEBERA4aZsYjjzzCSSedROvWrXHO4Zxj/fr1vjKffPIJ5513Hu3btychIQHnHBMmTAA8M+Cdc0yfPj0i8XTp0gXnnFYVaiYxOWzAOXcVcBdQCSwARrm6vxZfb2bTq7+eBBxXvSzWxupjJwAF1V//ycwWRjPmqFn1Zt0e19oqS2H1W3DylYHLiYgcREaMGMETTzzBgAEDlFQcYO655x7++Mc/ApCWlkZubi4AiYmJAHz99dfk5+ezd+9eEhISfAlsmzZtWizm5rR8+XJmzJhBly5dGDFiREuHE3ExmbwCR1a/JwKjGygzD5he/fVTwDCgD3A2kAwUAi8Cj5jZgqhFGm17CoMrt3tzdOMQERGJEQ895Fk86P7772f06NHU7uB67LHH2Lt3L/369eP111+nXbt2fuc7d+5Mjx49aNu2bUTi+clPfkJaWhqtW7eOyPWaavny5dx5550MGDBAyWtzMbMJwIQQyv8P8D/RiqcllaTkBDXTbE9Kew6OnydFRORgtmXLFrZu9czFvvbaa+skrgBffPEFAMOHD6+TuAI8+eSTEY1p1qxZEb2eBBazY17F473K3pRacsAypZbMzKpezRSRiIhIy9m3b5/v64aGAXjLHCzDBA42Sl5j3MbS1kytGBawzNSKYXy3PzZ+VSEiEs8qKyuZM2cON910E7169SI3N5eUlBQOPfRQhg0bxuzZsxusW1VVxfTp0xk4cCA5OTkkJyfToUMHjjvuOEaOHMm779bd1XzdunVcf/31dO/enVatWtG6dWvy8vLIz8/nL3/5C9u21T9Xec6cOVx44YV07NiRlJQUOnbs2Gh8wSgvL+exxx7jzDPPpEOHDqSmppKXl8eQIUN47LHHKCkpqVOntLSU+++/n1NPPZW2bdvSqlUrevTowZgxY9i8OfCQtrKyMh555BH69etHdna2r72RI0fy1Vdf+ZWdO3cuzjm6dOniO+adqOWdjFV74tTVV1/tO1+zXmMTtsyMF154gXPOOYeOHTuSmprKYYcdRv/+/XnggQfYvt1/knRjE7ZCuU+vESNG+O6rsrKSBx98kBNPPJHWrVuTnZ3Nueeey6efflqnnnOOq6++GoB58+b5PaPaMe7evZs///nP9OrVi4yMDN9nvXfv3owbN46VK1fWG1uLMzO9gnwB5nlkzee5jzdY3vg37L4/XGP7b88xuyPT99p/e47d94drLG/8G/b8kg0Raa9oyw/28cv328L//Z19/PL9VrTlB7M928w+nW42d5Lnfc+2iLQlIhJNV111lQE2YMCAoOusWLHCvP/WA5aammrp6el+x+6+++56615++eV+5dq2bWspKSm+70899VS/8kuXLrWMjAzf+eTkZGvXrp3fNd5555067dx2222+8845a9eunVXvAGmA/f73vw/pOXlt3LjRTjrpJN91EhIS6lx7zpw5fnW2bNliP/3pT/2eV817ysrKskWLFtXb3g8//GAnnniiX3s166alpdkrr7ziK//RRx9Zbm6utW/f3lcmNzfX95o8ebL17t3bcnNzLTk52QDLzMz0ne/du7fvWgMGDDDAHn/88Tpx7dixwwYNGlTnGSckJPiO1a6Xl5dX7/MJ5z69vJ/f2267zYYOHer7jLRp08av7sKFC/3q5ebmWmZmpq98zWeUm5trH330ke8+jz32WL+4srKy/O5z/Pjx9f7ZNYcaOVfdfKy+g3rFTvK6fU+pdbvtbcsb/6adNP5ZG/eHW2zKH0bauD/cYieNf9byxr9p3W5727bvKW1SO1WVlbbw8fF1EuTy29tZ5R1Zfsfsrg5m8+4zq6qK0F2KSKzavqfUnvt4gz30wRp77uMNTf63pjmFk7yuXr3aLrnkEnvjjTds8+bNVlX971xhYaH9+c9/tsTERHPO2eLFi/3qzZs3z5cAPPDAA7Zr1y4zM6uqqrIffvjBpk+fbrfccotfnYEDB/qS2mXLlvmOl5SU2CeffGKjR4+uk5g899xzvv/Ub7zxRtu6dauZmW3bts1++9vf+s499dRTQd+zmdn+/fvt5JNPNsDat29vTzzxhO3Zs8fMzPbu3euLp/Z9e5OqrKwse/HFF62iosLMzD755BPr2bOnL8H0xulVVlZmffr0McD69+9v8+fPt9JSz2dr8+bNdssttxhgrVu3trVr1/rVXbdunTX2/3Gg5LSx8+ecc44B1qpVK3vooYesuLjYzMxKS0ttxYoVdvvtt9uMGTP86jSUvDblPr2f33bt2ll2dra98MILvrqfffaZHX/88QZYnz596tzD448/3uhn/8477zTAOnToYG+++aaVl5f7Yl6zZo3de++99thjjzVYP9qUvMZx8mpm9vCsNZY3/s0GXw/PWtPkNhY+Pt4/QQ3mNe++CNydiMSiqqoqe3jWGt8Pz95Xt9vetodnrfEldbEsnOS1MXfddZcBNmLECL/jkyZNMsCGDh0a9LVatWplQJ2EsCFVVVXWtWtXA+zSSy+tt8xll11mgOXl5VllZWXQsUybNs3Xc/rZZ58FVWf+/PkBe4g3b95sWVlZBtif/vQnv3P//Oc/fYnX/v37673+9ddfb4DdcMMNfsejmby+9dZbvt7W+u6pIQ0lr025T+/nF7AFCxbUqffpp5/6zq9fv97vXDDJ69lnn22A3XvvvcHdZDMLlLxqzGscuGFgV8YO6U5Kkv8fV0pSAmOHdOeGgV2bdP3irZs4ed2/Qq5n2hxB5IA1bc5apry/hrKKKr/jZRVVTHl/DdPmrG2hyFrWeeedB8BHH33kdzwzMxPwzISvqqqqU68+3jqbNm0Kqvzy5ctZu9bz3L1rnNZ2xx13ALBhwwaWLFkS1HXhx9n3V199NSeccEJQdV5++WUAevfuzdChQ+ucz83N5de//jUAL774ot+5J554AoAbbriB1NTUeq9/+eWXAzBz5syg4okE73M466yz6r2nUEXiPvv168d//dd/1Tneq1cvDj/8cODH1RVCEernL5YoeY0DzjluLOjG4lvPZNJFPRkzuDuTLurJ4lvP5MaCbvUuExKKr+c/T6orDz0u7+YIInJAKSopY+rswMnp1NlrKSopa6aImte+fft44IEHyM/P55BDDiE5Odk32eWnP/0pAD/88INfnUGDBpGSksKyZcvIz8/n6aefrlOmtp/97GcAXHnllfz+979n8eLFlJc3/G/xsmXLAHyTwOrTo0cPDjvsML/yjSkvL2fp0qV+MQXDe/2BAwc2WKagwLNX0Jo1a3yTvSoqKnyJ9ZgxY+jYsWO9r2HDPJOVv/vuu6BjaqrFixcDoT2HhkTqPvv06dNgG94/6+Li4pDj897j1KlT+eUvf8k777zD7t27Q75OS4jJdV6lftnpKfyiT+eIX7dyV/gbHOzd/j1a50DkwPL+F5vr9LjWVlZRxcwvN0fl36SWtGnTJvLz81mzZo3vWHp6OllZWSQkJFBZWcm2bdvqzLrv2rUrf//737nxxhtZsGABCxZ49sbp0qULQ4cO5brrrvMlvl6TJ09m9erVLFy4kEmTJjFp0iTS0tLo27cvl1xyCSNGjKBVq1a+8t61Tb0JS0MOP/xwvv/+e1/5xhQVFVFRUQF4Fu8PVjDxeHsGzYxt27aRnp5OUVERZWVlvrYbU3NprGgrLPRsDBTKc2hIpO4zIyOjwTppaWkAAX/oaciVV17JRx99xGOPPcbTTz/N008/TUJCAieccALnnXce119/PZ06dQr5us1BPa9CYmbHsOuu3NWq8UIiEle27G5kS+pqhbuCKxdPRo8ezZo1azjqqKN45ZVXKCoqYs+ePWzZsoXNmzf7eubqM3LkSNatW8eDDz7I+eefT05ODuvXr+cf//gHvXr14p577vErn5OTw4cffsjMmTMZNWoUP/3pTykrK2POnDn85je/4fjjj2fjxo112iktjexz9wwvDF+o8dQcVvHZZ58FO+ck7sTDfT766KOsXLmS22+/nfz8fFJTU1m+fDl//vOf6datW7MO2QiFklehW/9LG90IoT6llszy9DOiEJGItKRDMuofm1dbbmZw5eJFWVkZr732GgDPPPMMF154IVlZWX5lvD1zDcnNzeWmm25ixowZbN26lSVLljBs2DDMjD/96U98/vnnfuWdcwwaNIiHHnqIZcuWsW3bNh599FGys7P55ptvuPnmm31lO3ToAMC3334bMAZvwust35icnBySkjy/iN2wYUNQdWpeP1AdbyzOOdq3b+9rLzExEYAvv/wy6PaaQ25uLhDac2hILN9nTccddxx33nknc+bMYceOHbzxxhv07NmTkpISrrrqqrB6daNNyauQ1aETy468JuR6UyuG0TYnNwoRiUhLGnJcxzoTRGtLSUpg8LHh/9YmFm3bts3Xi1j7V/xeH3zwQdDXc87Rp08fXnrpJQ4//HCqqqr48MMPA9bJysriuuuu8/XSzps3z3fu5JNPBqCkpKTByVhr1qzh+++/9yvfmOTkZHr18uzS+PbbbwdVp+b1582b12CvoXfThO7du5Oenu5rr3fv3gC8+uqrQbfXHE477TQgtOfQkJa8z4QEz9/fUHtzU1JSOPfcc3nppZcAzzCar7/+OuLxNZWSVwHgtCvvYVGX6+v0wFZYAhXm/zEptWQmlw/nn27YAfefl4h4xtePKgi8ismogq5kp6c0U0TNIzMz0zcBdsWKFXXOb9q0iYcffrjeut6xjfVJTEwkOdnzb6s3Oa6qqvKNM62Pd6xrzV/Jn3TSSXTt6vlzqT0EwWvChAmAZ6ztKaec0uD1a7vyyisBmD59ep3e4YZcfPHFgGemu7fHuqbCwkL+8Y9/ADB8+HC/cyNGjADglVdeYc6cOQHbCWcyUri8z+H999+vd0e0ULXUfXpXEtixY0eDZQJ9ZmuOtY70MJVIUPIqALiEBPqOuJe9N65gSc87WZT3a57scAt9Sv9Gn9K/8bvya/lr+cX8rvxaTit9mGmVFzCqoNsB95+XiHhEe4m+5lReXs62bdsCvsrLy2nTpo2v523kyJEsX74c8CSas2bNYsCAAQ32ZP3hD3/g4osvZsaMGX6TcwoLCxk1ahTr1q3DOcfgwYMB2LVrF127duXuu+9mxYoVVFZW+rV12223AZ4lm7ycc0ycOBGA1157jd/+9re+bUq3b9/OqFGjeO655wCYOHGir/ctGL/61a846aSTKC0t5cwzz+Spp55i7969gGci0ZIlS7j22mv5+OOPfXX69evnW05q5MiRvPzyy777WLp0KUOGDKG4uNg3lKJ2e6eddhpVVVWce+65PPTQQ37PbcuWLTz33HPk5+fz0EMPBX0fTXX22Wdz9tlnY2ZcdNFFPPzww74EsKysjBUrVnDLLbcwY8aMoK7XUvfpXY3iyy+/9Pszq2nQoEGMGjWK+fPn+00W++KLL3xJd6dOnejZs2fE4oqYYAYQ69VymxTUu11rMzkQFikXkabZvqfUnl/i2WHr+SXxucNWMC/v4vKLFy/2bR4AWHp6uu/77OxsmzFjRr0L5N90001+18vMzPTbAhT8t5UtLi72O5ecnGzZ2dmWmJjoO3bUUUfZd999V+e+am4PW9+WnuFuD/vtt9/6dm0CLDEx0bKyshrdHrbmlrJpaWl1toetvUuYV2FhoZ1xxhm+ss45y8rK8tv+FLAJEyb41Yv2DlvFxcW+8w0941C2hw33Pr2f3zvuuCOs++zfv7/v2tnZ2ZaXl2d5eXm+7Xprb1mblZVlaWlpvmOtW7e2Dz74oMG2o63Gn3GdfExLZcUoq6pi8ZN/4OR1/+KUGmuwln5+N4uOvIbTrrwHF8JP1eHwri97+al5zPxyM4W7SsnNTGXwsR3V4ypykIjWEn2x6tRTT2XRokVMmDCBefPmUVJSQqdOnRg6dCi33Xabr2extptvvpmf/OQnzJo1i6+++opNmzZRWlrKEUccwemnn84NN9xAv379fOUzMzN58803+eCDD1i4cCEbN25k69atpKen06NHDy644AJ++9vf1rtM0sSJEykoKGDq1KksWrSI4uJicnJy6Nu3L6NGjeLMM88M696POOIIPv30Ux577DFefPFFVq5cyd69e+ncuTNHH300F198cZ2hCB06dGDRokX87W9/47nnnmP16tWUlZXRrVs3zjnnHH73u981uNzSIYccwrx583jhhRd45plnWLp0KUVFRaSkpHD00UdzxhlncNFFFzFo0KCw7idc7dq1Y/bs2Tz99NM89dRTLF++nJ07d9KpUye6du3KsGHD+PnPfx709VrqPl999VVuv/123nnnHb7//ntfj+/+/fsB+Ne//sXbb7/N3LlzWbduHZs3e5bNPProoxk0aBBjxozhyCOPjGhMkeIsTpegaAnVP302y3IWi6b/nr7r/97w+S7X03fEvVGPQ0RERKS5ecefm1mdnZg05jUGBbNd68nr/kXx1vjb0k1ERESkKZS8xqBgtmtNdeV8Pf+FZopIREREJDYoeY1BwW7XWrlLPa8iIiJycFHyGoOC3a41MTM29xwWERERiRYlrzEomO1aSy2Zbv1/0UwRiYiIiMQGJa8xKJjtWpcdeQ1ZHdTzKiIiIgcXJa8xqqHtWkstmUVdrue0K+vfGlBERETkQKZ1XkPQnOu8ehVv3cTX81+gctcmEjM70a3/L9TjKiIiIge0QOu8KnkNQUskryIiIiIHG21SICIiIiIHBCWvIiIiIhI3lLyKiIiISNxQ8ioiIiIicUPJq4iIiIjEDSWvIiIiIhI3lLyKiIiISNxQ8ioiInKQyM/PxznH9OnTWzqUqNq9ezdjxozhJz/5CSkpKTjn6NKli1+Z5557jr59+5KRkYFzDuccc+fOBfB9v379+ibHsn79et/1JDKSWjoAERGR5jBjxgyGDRsGwODBg3n//fcj3saOHTt48MEHAZgwYULEry/BufDCC/nggw8AyMzMJDs7mw4dOvjOP/vss1xxxRUAJCcnk5ubC0BKSkrzB9sCZsyYwfLly8nPzyc/P7+lwwmZklcRETkoPPHEE76vZ82axcaNGzn88MMj2saOHTu48847ASWvLeWLL77ggw8+IDk5mfnz53PaaafVKfPAAw8AcPPNN3PfffeRlOSfDvXo0QPwJLZNlZyc7LterJgxY4bv70M8Jq8aNiAiIge87du389Zbb9G6dWsuv/xyqqqqePrpp1s6LImCL774AoATTjih3sS1ZpmRI0fWSVwBVq1axapVqzjssMOaHM9hhx3mu55EhpJXERE54D377LOUl5dz/vnn89///d+Af0+sHDj27dsHQJs2bZpURmKXklcRETngeRPVK664gn79+tG5c2dWrVrFkiVLGq1bUlLClClTOP3008nOziYtLY2jjjqKn//85zzzzDOUl5cDnl+/Hnnkkb563kk63lfNYQRdunTxmyBUn4YmDZWVlfHWW29x7bXXcuKJJ9K+fXvS0tLIy8vjiiuuYOnSpcE/mBBt376dO+64g169etGuXTtat25N9+7dufTSS3nttdfqrVNYWMgtt9zC0UcfTevWrWnbti2nnHIKf/3rXyktLQ3Y3tatW7n11lvp2bMnbdq0IT09neOPP57bbruNoqIiv7ITJkzAOceIESMAmDdvnt/znz59ep2JU0ceeaTvmLceND5hK9jPBAQ3YSuU+/Sq+RkqKipizJgxHHnkkaSmpnLYYYdx7bXXsmnTJr86c+fOxTnn+/tw55131vmc1rRu3Tquv/56unfvTqtWrWjdujV5eXnk5+fzl7/8hW3btjV4T1FlZnoF+QLM88hERCRerFy50gDLycmxsrIyMzMbP368Afab3/wmYN0vvvjCunTpYt5//5OSkqxdu3a+7wFbt26dmZkNGzbM2rdv7zuem5vr95o8ebLvunl5eQbYnDlzGmy79vW93njjDb/2W7dubWlpaX4xPvnkk/Vec8CAAQbY448/3uhzq23+/PmWk5PjayclJcUyMzP9Yqnt448/tuzsbN/5jIwMv1hPPPFEKywsrLe9BQsW+NVNSUmxVq1a+b4/4ogjbNWqVb7ykydPttzcXF9MycnJfs9/7ty5vq+912jfvr3v2KhRoxp99mahfSbMzNatW9fg8wnnPr28n6GnnnrK93Xr1q0tNTXVV7dLly5WVFTkq/PRRx9Zbm6u788gPT29zufUa+nSpZaRkeG7VnJycp37fOedd+q9p0io8czq5mP1HdRLyauIiJmZ7dlm9ul0s7mTPO97trV0RCEbN26cAXb99df7jn3++ecGWHZ2tpWWltZbb/v27XbEEUcYYEceeaTNmDHDV3bXrl22YMECu/rqq+27777z1WksUfFqSvI6Z84cu/rqq23WrFm2bduPfx4bNmyw0aNHG2BpaWm2YcOGOtcMN3ldu3atLyk86aSTbPbs2VZRUWFmZkVFRfbee+/ZhRde6FenqKjIOnXqZID17NnTlixZYmZmFRUV9tJLL1lWVpYBNmjQoDrtrV+/3pcoXXPNNbZq1SqrrKy0qqoqW7lypQ0dOtQAO/bYY31xeD3++OMG2IABAxq8n0DJaaDzkf5MNOU+vZ+hdu3a2UknnWQLFy40M7Py8nJ77bXXfNcdN25cnXavuuoqA+yOO+5o8BkNHDjQADv11FNt2bJlvuMlJSX2ySef2OjRo31tRoOSVyWvIiKhqaoym3ef2V0dzO7I/PF1VwfP8aqqlo4wKBUVFb4EasGCBX7nevbsaYC9/PLL9db1Jr3t27e3jRs3BtVecySvjRk5cqQBNmHChDrnwk1eL7nkEgOse/futmvXrqDq3HXXXb7katOmTXXOv/fee757nDVrlt+5K664wgC/3tCaSktL7cQTTzTAXnrpJb9z0UxeI/2ZaMp9ej9Dubm5fj/EeE2ZMsWXZNcWTPLq7f1dvHhxI3cYHYGSV415FRGRuhZMgdkTobLWmMTKUs/xBVNaJq4Qvf/++2zatIm8vDzOOOMMv3PedT4bmrj11FNPATB27NiIzDpvLueddx4AH330UUSut2fPHv79738DcNddd5GRkRFUvZdffhmAa665ho4dO9Y5P2TIEPr27QvAiy++6Du+b98+XnrpJQDGjBlT77VTUlK4+OKLAZg5c2aQd9J0kfxMROo+r7vuOnJycuocv+CCCwDPuNWSkpKQ48vMzASoM242FmidVxER8VeyHeZNDlxm3mToNRLS6/6nGUu8ielll11WZzLKZZddxq233so777zD1q1b/RaxX79+PZs3bwbgZz/7WfMFHKSioiKmTZvGO++8w+rVq9m5cyeVlZV+ZX744YeItPXpp59SUVGBc46hQ4cGVaesrIyVK1cCMHDgwAbLFRQUsGjRIpYtW+bXXllZGQCnnnpqg3W9KwZ89913QcXUVJH+TETqPvv06VPv8ZrJ9Y4dO0hPTw8pvp/97Gc8/vjjXHnllfzmN7/hggsuoFevXhFZ+7aplLyKiIi/VW/W7XGtrbIUVr8FJ1/ZPDGFYefOnb4Z8Jdffnmd8507d6Zfv37Mnz+fZ599lptuusl3rrCw0K9cLPnyyy8pKCjwizEjI4NWrVrhnKOsrIzi4uKwetvq422nbdu2tG3bNqg6RUVFVFVVAQTsofRuErF161bfsZo9fTXvsSF79+4NKqamivRnIlL32VBPeFpamu/rmqsfBGvy5MmsXr2ahQsXMmnSJCZNmkRaWhp9+/blkksuYcSIEbRq1Srk60aChg2IiIi/PY3/RwrA7s3RjaOJXnjhBfbv3w94FqyvvSSQc4758+cDdYcOeIbcxaarr76awsJCTj75LR0dAAAgAElEQVT5ZN599112797Nrl27KCwsZPPmzb5fRUfqHpp6ncaWw6rNm/RmZWUFNR8l0HJjkRTpz0Ss3qdXTk4OH374ITNnzmTUqFH89Kc/paysjDlz5vCb3/yG448/no0bNzZrTF5KXuNdyXZY+gTMu8/zXrK9pSMSkXjXJje4chl1xzHGklA2IfjPf/7DihUrfN/XHKO5YcOGiMYF+HZ18ibXte3cubPe499++y1LliwhMTGR119/nbPOOqvOQvvB9OKFwvssdu7c2WBctWVnZ5OQ4EkxAj0/b/JTc8hGbq7n81dcXOz7NX0siPRnIlbvsybnHIMGDeKhhx5i2bJlbNu2jUcffZTs7Gy++eYbbr755haJS8lrvDKD+ZPh/mPgjVEw527P+/3HeI7HcK+BiMS4o8+FxNTAZRJTocc5zRNPGNauXcvChQsBWL58OcXFxQ2+vBOcaia7Xbp08SUrb7/9dtDtehM2CNxT165dO4AGe64++eSTeo/XTPYa+nX8Bx98EFSswerduzdJSUmYGe+8805QdVJSUjj++OMBmDNnToPlZs+eDcDJJ59cpz2AV199NdywIy7cz0RDWvI+vZ/TUHuTs7KyuO6667jnnnsAz0YQLUHJa7w6QGYCi0gMSs+BAeMClxkwLqYna3kT0RNPPJETTzyRdu3aNfi65JJLAHjmmWf8Jj398pe/BOCvf/0r33//fVDtemdog2eSTEN69uwJUO+uVGbGpEmT6q3nHXNaWFjIli1b6pxfsWIFzz77bFCxBqtNmzYMGzYMgDvuuIPdu3cHVc87S3769On1zlh///33WbRoEQDDhw/3Hc/IyOCiiy4CYOLEiQF7kisqKtizZ09wNxIB4XwmGtKS9+n9nDb0Ga2qqqKioqLB+t6xrqEOCYkUJa/xKNiZwBpCICLh6jcWCv5Ytwc2MdVzvN/YlokrCGbmW9LowgsvbLT8eeedR3JyMps3b+a9997zHR8/fjyHHXYY27Zto1+/frz++uu+2eF79uxh7ty5XHrppX69p+3atePQQw8F4PHHH2+wTW+y9tZbbzFp0iTf5Kr169dz2WWX8emnn9Zb75hjjuHwww/HzPjFL37B2rVrAc+EnFdffZXBgwfXGUYQCffccw8ZGRmsWbOG/v37M2fOHN+YzR07dvDWW29xzjn+PfE33ngjnTp1Yt++fQwdOtR3T5WVlbzyyitceumlAAwaNIiCggK/uvfeey/Z2dls2rSJ008/nX//+99+idLatWt58MEHOeaYYxp8VtEQzmcikJa6z+OOOw6Ad999t94fLHbt2kXXrl25++67WbFihe+HuqqqKmbNmsVtt90GwFlnnRWxmEISzCBhvWJsk4JPp/svGt7Qa+kTLR2piMS7Pds8/5bMneR5j4MdtmbPnu1b4HzlypVB1TnrrLMMsOHDh/sd//zzz+3www8PuEVm7YXsb7/9dt+59PR0y8vLs7y8PHvggQf8yl144YW+cgkJCb7rpqWl+S3gX/v6r776qiUkJPjOZ2RkWEpKigHWuXNne+qppwywvLy8OvfZlO1hZ8+e7Xfvqamp1rZtW79nUdvHH3/s20nLG2vN7WFPOOGEBreHXbJkiR166KG+sklJSZaTk+O3/Slgc+fO9asXzU0KzEL/TDS2cUW499mUjS62bt3q25I2ISHBOnbs6PucmpkVFxf7tZ2cnGzZ2dmWmJjoO3bUUUf57SQWaTWemTYpOCAcIDOBRSQOpOd4lsMa8DvPewwPFfDyDhno3r27r4epMd5f37722mt+v0rt2bMnX3zxBRMnTqR37960atWK/fv3c9RRR3HBBRfw3HPP+ZZ78rr99tuZNGkSJ5xwAmbGhg0b2LBhQ51f0T733HPcfffd9OjRg6SkJJKTk7nooov4+OOPGTJkSIOxDhs2jNmzZzN48GAyMjIoLy8nLy+PsWPH8p///KdOPJEycOBAVq9ezfjx4zn++ONJSkqioqKC7t27c9lll/H666/XqXPKKafw5ZdfcvPNN9O9e3fKy8tJSkqid+/eTJ48mY8//phDDjmk3vb69OnDqlWrmDRpEqeffjoZGRns2LGDVq1a0bt3b8aPH88nn3zCgAEDonK/DQnnMxFIS9xn+/btmTNnDhdeeCEdOnRg69atvs8peIYVvPnmm4wePZpTTjmFDh06sHv3btLT0+nTpw933303y5cvj9pnrTHONLEnaM45T/drSz+zpU94Jmc15ucPx/QajCIiIiL18W4qYmau9jn1vMah4s5DKLXAO1yUWjLFRwxupohEREREmoeS1zj03rpyplYMC1hmasUw3l8f+o4aIiIiIrFMyWsc2rK7lGmV5zO5fHidHthSS2Zy+XCmVZ5P4a6WWcJCREREJFqSWjoACd0hGamAY1rlBTxbWcDgxKXkUkwhWcys7EUxnvXbcjMbWWRcREREJM5owlYIYmXCVlFJGaf9ZRZlFVUNlklJSmDxrWeSnZ7SjJGJiIiINJ0mbB1gstNTGFXQNWCZUQVdlbiKiIjIAUfDBuLUDQM9yevU2Wv9emBTkhIYVdDVd15ERETkQKJhAyGIlWEDNRWVlDHzy80U7iolNzOVwcd2VI+riIiIxLVAwwaUvIYgFpNXERERkQONxryKiIiIyAFByauIiIiIxA0lryIiIiISN5S8ioiIiEjcUPIqIiIiInFDyauIiIiIxA0lryIiIiISN5S8ioiIiEjcUPIqIiIiInFDyauIiIiIxI2YTF6dcznOuWucc/92zq11zu1zzu10zn3onPuVc67euJ1zpzvn3nbOFTnn9jrnPnfOjXbOJTb3PYiIiIhI5Dkza+kY6nDO/Rr4O7AJmAN8C+QCFwJtgVeAS6xG8M6586uP7wdeAIqA84AewMtmdkkE4jKAWHxmIiIiIgcK5xwAZubqnIvFRMw5VwCkA2+ZWVWN4x2BJcARwMVm9kr18UxgLZ7E9gwz+7T6eBowG+gLXGZmzzcxLiWvIiIiIlEWKHmNyWEDZjbbzN6ombhWH98M/KP62/wapy4GOgDPexPX6vL7gT9Wf3t99CIWERERkeYQk8lrI8qr3ytqHCuofn+3nvLzgb3A6c651GgGJiIiIiLRldTSAYTCOZcEXFn9bc1EtUf1+5radcyswjm3DjgOOAr4qqlx5Ofnh1Vv7ty5TW1aREREJG6EmzMFEm89r/cCxwNvm9l7NY63rX7f2UA97/F20QpMRERERKIvJids1cc5Nwp4CFiFZ1JWUY1za4BuQDczW1tP3YV4Jm31NbPFTYhBE7ZEREREoizuJmzV5py7AU/i+iUwsGbiWs3bs9qW+mXWKiciIiIicSjmk1fn3GjgEWAlnsR1cz3FVle/d6+nfhJwJJ4JXt9EK04RERERib6YTl6dc+OBB4DleBLXLQ0UnV39PrSec/2B1sBCMyuNfJQiIiIi0lxiNnl1zv0JzwStpcCZZrYtQPGXgW3Apc653jWukQZMrP7279GKVURERESaR0xO2HLOXQVMByqBh6l/rOp6M5teo84FeJLY/cDzeLaH/TnV28MCw62JN6sJWyIiIiLRF4/bw04A7mik2Dwzy69V7wzgNjwrC6Th2TL2f4GpZlYZgbiUvIqIiIhEWdwlr7FKyauIiIhI9MX9UlkiIiIiIqDkVURERETiiJJXEREREYkbSl5FREREJG4oeRURERGRuKHkVURERETihpJXEREREYkbSl5FREREJG4oeRURERGRuKHkVURERETihpJXEREREYkbSl5FREREJG4oeRURERGRuKHkVURERETihpJXEREREYkbSl5FREREJG4oeRURERGRuKHkVURERETihpJXEREREYkbSS0dgLSQku2w6k3YUwhtcuHocyE9p6WjEhEREQnImVlLxxA3nHMGENfPzAwWTMHmTcZVlv54ODEVN2Ac9BsLzrVggCIiInKwc9W5iJnVSUrU83qQsflTcHMmUvuT4CpLYfZEzPAksSIiIiIxSD2vIYj7nteS7VRMOZokK2uwSIVLIWnsKg0hEBERkRYTqOdVE7YOIiWfzwiYuAIkWRl7Pn+tmSISERERCY2S14PIN9/8X1Dl1q1bG+VIRERERMKj5PUgsoWsoMoVWnDlRERERJqbJmwdRHbmDaF0zSRSXXmDZUotmd1dhkQvCC3RJSIiIk2gCVshiPcJW0UlZUy/90bGJL7QYJn7K3/BiN8/QnZ6SmQbr16ii3mTocYSXSSmgpboEhERkRo0YUsAyE5PISV/LJPLh1NqyX7nSi2ZyeXDSckfG/nEFTyJ6+yJ/okreL6fPdFzXkRERKQR6nkNQbz3vIIn9mlz1vLU7GUMsE/IpZhCspjn+vDLgpO5YWBX3087EVOyHe4/pm7iWlNiKoz5SkMIREREJGDPq5LXEBwIyatXUUkZM7/cTOGuUnIzUxl8bMfo9LgCLH0C3hjVeLmfPwwnXxmdGERERCRuaIctqSM7PYVf9OncPI3tKQyu3O7N0Y1DRERE4p7GvEr0tckNrlxGx+jGISIiInFPyatEXXHnIXUmiNVWaskUHzG4mSISERGReKXkVaLuvXXlTK0YFrDM1IphvL++4fVnRUREREDJqzSDLbtLmVZ5fsAluqZVnk/hrgCrEYiIiIigCVvSDA7JSAUc0yov4NnKAgYnLvUt0TWzshfFZAKQm5nasoGKiIhIzNNSWSE4kJbKak5FJWWc9pdZlFVUNVgmJSmBxbeeGb3lukRERCRuaIct8SjZ7llzdd59nveS7c3SbHZ6CqMKugYsM6qgqxJXERERaZR6XkMQtz2vZp7tV+dN9t/lKjEVBoyDfmMh0rtq1QnBs7PX1Nlr/XpgU5ISGFXQNTo7e4mIiEhc0g5bERK3yev8yTB7YsPnC/4I/cc1SyjNurOXiIiIxCUlrxESl8lryXa4/xj/HtfaElNhzFeQntN8cYmIiIg0QGNeD2ar3gycuILn/Oq3miceERERkSbQUlkHkOKtm/h6/vNU7tpMYmZHuvW/lKw9hcFV3r05usGJiIiIRICS1wOAVVWx+Mk/cPK6f3GK+3GXqtLP7+b/2p/OT4K5SEbHqMUnIiIiEika8xqCWB3zumj67+m7/u8Nnq8ikQQqGzxviak4jXkVERGRGKExrwew4q2bOHndvwKWqWok2V58xEglriIiIhIXlLzGua/nP09qjaEC9UlyVbxb2ZtSS/Y7XmrJTC4fzlVr+1NUUhaZgFpoIwQRERE5OGjMa5yr3BXcRKsvqrpwa/k1DE5cSi7FFJLFzMpeFJMJGDO/3Mwv+nQOP5CGNkJ4e1yzbYQgIiIiBz4lr3EuMTO4iVaFZFFMJi9WDqz//K5GltNqzIIp9W+EUFn64/Fm2ghBREREDlwaNhDnuvW/tM5wgNpKLZmZlb0ClsnNTA0/iJLt2LzJAYvYvMkaQiAiIiJNpuQ1zmV16MSyI68JWOaRymHVwwPql5KUwOBjm7BU1qo3cY1shOC0EYKIiIhEgJLXA8BpV97Doi7X1zsha1GX60kZGPjX9aMKupKdnhJ2+3uLfgiu3Pbvw25DREREBDTm9YDgEhLoO+JeirfexGfzX6By1yYSMzvRrf8v6NuhE6eZ4Zxj6uy1lFVU+eqlJCUwqqArNwzs2qT2V+5K45SgyrUKqpyIiIhIQ7RJQQhidZOCYBWVlDHzy80U7iolNzOVwcd2bFKPq9dj737CVYvODrhkV6kl80Tfd7huaJ8mtyciIiIHtkCbFKjn9SCSnZ7StOWwGpCZncvUimGMS36xwTJTK4bROSc34m2LiIjIwUXJqzTZkOM6ctrrw6AcRiX9268HttSSmVoxjH+6YSxuyqQwEREREcIcNuCcawecAOw2s//UOtcJeBgYDFQCbwG3mNmWpofbsuJ92EA0PTL7a6a8v4YsdtW7EcLYId25saBbS4cpIiIicSDQsIFwk9dbgPuAv5nZb2scTwL+AxwLeBsz4Eugl5lFaA/SlqHktWFmxrQ5awNOCnPaYUtERESCEI3k9T1gENDPzBbWOH4F8BSwD7i/+n0ckAncZGaPhBF/zFDy2rhoTQoTERGRg0c0ktf/A7oA7cxsd43jrwHnAuPNbEr1seHA88B8M8sPPfzYoeRVREREJPqikbzurL5g21rHi/H0sh5uZpuqj6Xg6YEtMrMOITcWQ5S8ioiIiERfoOQ13B220mrXdc71ANoCX3sT1+pGywBvUisiIiIiErZwk9ctQGvnXM21jwZVvy+sp3wrYGeYbYmIiIiIAOEnr59Uv48BcM61Bn6NZ2WBWTULOucOw5O8bkJEREREpAnCTV4fxbMU1i3Oua+ANcBxwFbg1VplB1a/rwizLRERERERIMzk1czeAybg6WntARwKbAOuMLN9tYpfXv0+J8wYRURERESAMFcb8FV2rjNwKrADWGJmO2udTwHG40mSHzWzzU2ItcVptQERERGR6Iv4UlkHKyWvIiIiItEX8aWynHO3O+fGhFB+lHPu9nDaEhERERHxCneTgipgs5kdGmT5dUBnM0sMubEYop5XERERkeiLxiYFIiIiIiLNrrmS12xgfygVnHMXO+ceds4tcM7tcs6Zc+7pBsp2qT7f0Ov5iNyFiIiIiLSopGg34Jy7BMgAVodY9Y/AicAeYCNwdBB1PgNm1HN8ZYhti4iIiEgMCip5dc7dBNxU63AH59w3gaoB7YBMPOvBvhVibDfjSVrXAgMIbp3Y5WY2IcR2RERERCROBNvz2g7oUutYYj3HGjILuCvIsgCYmS9Z9Q7aFREREZGDW7DJ6wxgffXXDvhfYCcwOkCdKmAXsNLM/i/cAEN0qHPuv4EcYDuwyMw+b6a2RURERCTKgkpezewzPONJAXDO/S+wz8yeiFZgYRpc/fJxzs0FrjKzbyPVSH5+flj15s6dG6kQRERERGJeuDlTIGGtNmBmCcGu8dpM9gJ/BnoBWdUv7zjZfGCWcy69xaITERERkYiIi+1hnXP5eBLRZ8zs/4VQLwn4EDgVGG1mDzUxDm1SICIiIhJlgTYpaPJSWc65BKAbnrVckwOVNbP5TW0vFGZW4Zz7F57ktT/QpORVRERERFpW2Mmrc64T8BfgYqBVEFWsKe01wdbqdw0bEBEREYlzYSWTzrlDgY+BQ/GsPhBUtXDaioDTqt8DrUkrIiIiInEg3O1hJwCH4dn9ahSQByRXT+Rq8BWhmOtwzp3qnEup53gBns0OAOrdWlZERERE4ke4v8Y/G88wgF+Z2csRjMfHOXcBcEH1tx2r3/s656ZXf73NzMZWfz0JOK56WayN1cdOAAqqv/6TmS2MRpwiIiIi0nzCWm3AObcfzzCAdDOriHhUnjYmAHcEKLLBzLpUl/0VMAw4HmiPZ+JYIbAIeMTMFkQoJq02ICIiIhJlgVYbCDd5/RbINLN2TY4ujih5FREREYm+QMlruONQPwAynHPdmhCXiIiIiEhIwk1e7wFK8Iw1FRERERFpFuFuD7sW+DkwwDk30zk3UNuvioiIiEi0hTvmtTKMtszMWmKTgojRmFcRERGR6IvG9rAtteGAiIiIiBzEwk1eB0Y0Cok9Jdth1ZuwpxDa5MLR50J6TtPLioiIiDRBWMMGDlYHxbABM1gwBeZNhsrSH48npsKAcdBvLFR35YdUVkRERCRI0Rg2IAeqBVNg9sS6xytLfzzef1zoZUVEREQiICI9r86THucArc3s2yZfMEYd8D2vJdvh/mP8e1FrS0yFMV8BYPcfgwtQ1hJTcWO+0hACERERCUk0NinwXvhk59yrwE4827F+U+t8lnPuUefcP5xzKU1pS5rBqjcDJ67gOb/6LVj1ZsDEFfCcX/1WBAMUERGRg13Ywwacc78E/gUkN1TGzIqdc0cCZwJvAMpkYtmewuDK7d7M3rJKWgdRdO/274MqJyIiIhKMsHpenXPHAP/Ek7hOBXoD2xoo/iSepbXOD6ctaUZtcoMrl9GRlbvSgiq6clerJgQkIiIi4i/cYQNjgBRgmpmNNrNlQEMbF8yufu8bZlvSTIo7D6HUGuxIB6DUkik+YjDL0/8rqLLL08+IZIgiIiJykAs3eS0ADJjUWEEz+wHYC3QOsy1pJu+tK2dqxbCAZaZWDOP99eVkZucGVbZtTpC9uSIiIiJBCDd5PRQoMbONQZbfB+j3xzFuy+5SplWez+Ty4XV6VUstmcnlw5lWeT6Fu0oZclxH/umGBSz7TzeMwcd2bM5bEBERkQNcuBO2SoE055yzRtaNcs61AtoBO8JsS5rJIRmpgGNa5QU8W1nA4MSl5FJMIVnMrOxFMZkA5Gamkp2ewqiCbkx5v+GyY4d0Iztdi0yIiIhI5IS1zqtzbjnQEzjGzNZUH9sEHGJmibXKXgS8BCwwswFND7nlHOjrvBaVlHHaX2ZRVlHVYJmUpAQW33om2ekpmBnT5qxl6uy1fnVSkhIYVdCVGwZ29a3TJiIiIhKsaOyw9S5wAnATcEOAhnOA+/CMj9UyWTHO05valSnvr2mwzKiCrr7eVOccNxZ04/JT85j55WYKd5WSm5nK4GM7qsdVREREoiLcntdc4GsgHbgLuB9YQ3XPa/VQgWHA3UAenmW0uprZrkgF3hIO9J5XQL2pIiIi0uIC9byGvT2sc+5c4GU8a72W45n8lQisAo7Cs5SWwzM+9lwzmxVWQzHkYEhevYpKytSbKiIiIi0iKslr9YVPAR7Bs0lBff4D/NrMPgm7kRhyMCWvIiIiIi0laslrjQZOAP4LzxJaicBm4CMz+7TJF48hSl5FREREoi/qyevBQsmriIiISPQFSl7D3aRARERERKTZhbtUlgiUbIdVb8KeQmiTC0efC+k5LR2ViIiIHMAaHTbgnPvf6i83mdlttY6FwszsV2HUixkaNlDNDBZMgXmTobL0x+OJqTBgHPQbC1pOS0RERMLUpDGvzrkqPJsMrDazY2sdCyZD8Zaz2rtvxRslr9XmT4bZExs+X/BH6D+u+eIRERGRA0pTd9h6Ek8CuqmeY3KwKdnu6XENZN5k6DVSQwhEREQk4hpNXs1sRDDH5CCx6k3/oQL1qSyF1W/ByVc2T0wiIiJy0NBqAxKaPYXBldu9ObpxiIiIyEFJyauEpk1ucOUyOkY3DhERETkohbVJgXPuEOBSYKuZPddI2SuAHOBZM9sWVpQxIm4nbEVySauS7dj9x+ACDB2wxFTcmK805lVERETCEo1NCv4f8ADQNYiyJ1aXvTzMtiRcZp6VAe4/Bt4YBXPu9rzff4zneDhJeHoOi48YGbDI4iM0WUtERESiI9zk9efV768EUfYpPEtlnR9mWxKuBVM8S1rV7iWtLPUcXzAl5EsWlZRx1dr+TC4fTqkl+50rtWQmlw/nqrX9KSopa0rkIiIiIvUKd9jAd0AnoJWZlTdSNgXYB3xnZl3CCTJWxNWwgZLtnh7WQCsDJKZCiL/ef37Jt/z+1RUAZLGLwYlLyaWYQrKYWdmLYjIBmHRRT37Rp3OTbkFEREQOTk1d57U+HYAdjSWu1Y2WOed2AEHO9JGIiNKSVlt2/3jNYjJ5sXJgveUKdzXStoiIiEgYwh02sBto65xLa6xgdZlMYG+YbUk4orSk1SEZqUGVy80MrpyIiIhIKMJNXr+orntuEGXPAxKBVWG2JeGI0pJWQ47rSEpS4I9NSlICg4/VUlkiIiISeeEmr6/jmYQ1xTl3aEOFnHOHAVPwbCU7I8y2JBxHn+sZ0xpIYir0OCeky2anpzCqIPAiE6MKupKdnhLSdUVERESCEW7y+g9gI3AEsNw5d7NzrptzLqX61c05Nwb4T3WZ74G/RSZkCUp6DtZ/XMAi1n9cWEta3TCwK2OHdK/TA5uSlMDYId25YWAwK6iJiIiIhC6s1QYAnHMnA+8C7fH0rNZbDNgGDDGz5WE1FEPiarUB4JFZa9g3ezKjkv5Nqvtxbl2pJTO1YhitCsZx45ndw75+UUkZM7/cTOGuUnIzUxl8bEf1uIqIiEiTBVptIOzktfrChwN/AYYDybVOlwHPA7eZ2fdhNxJD4il5LSop47S/zKKsoqrBJa1SkhJYfOuZ0U84I7nDl4iIiBzwopa81migNdAb8M7S2QR8amb7mnzxGBJPyWvN9VgDiep6rGaejRDmTfZftisxFQaMg35jwdX5TIqIiMhBLhrrvPoxs73A/EhcSyKj5nqsgUR1PVbvDl+1eXf4AmhkXK6IiIhITeFO2JIY1+LrsZZs9/S4BjJvsqeciIiISJCUvB6gWnw91lB2+BIREREJUqPDBpxz31R/udbMhtQ6Fgozs5+EUU/C4F2Pdcr7axosE9X1WKO0w5eIiIgc3IIZ89ql+n1/PcdCEfuznA4w3vVWp85eS1lFle94SlICowq6Rnc91ijt8CUiIiIHt0ZXG3DOXVX95U4zm1HrWEjM7Ilw6sWKeFptoKYWWY+1ZDt2/zG4AEMHLDEVN+YrLZslIiIifqK+VNbBIl6T15ayaPrv6bv+7w2f73I9fUfc24wRiYiISDwIlLw2OmHLOfeqc+5ftY51ds4dFrEI5YBTVFLGVWv7M7l8OKXmv39FqSUzuXw4V63tT1FJWQtFKCIiIvEomDGvFwC1Z9Wsx7MRgRJYqdf7X2ymrMKYxgU8W1lQ7w5fYMz8cnP0NkkQERGRA04wyWsVkFjPcW2NJA2quUlCMZm8WDmw3nJR3SRBREREDjjBrPNaBOQ459pGOxg5cLT4JgkiIiJyQAqm5/UTYCjwhnPueWBP9fFWzrkrQ2nMzJ4MMT6JU0OO68jtr3/ht0RXbVHdJEFEREQOSMEsldUPmIUn0fUWdoS+bquZWTDJcszSagOheWT21wE3SRg7pDs3FnRrxohEREQkHjR5qSzn3GnATUBPoDWeTQoqgY2hBGJmR4ZSPtYoeQ2NmTFtztqAmyR4P5wiIiIiXhFf59U5V3VunikAACAASURBVAVsNrNDmxxdHFHyGp4W2SRBRERE4paS1whR8ioiIiISfYGS10bHoDrnlgFbzeysGocHAlrjSERERESaVTATtur0slYf22RmB9UmBep5FREREYm+Jm0Pi2diVnI9xzXTRkRERESaVTDJ61Yg2znXKdrBiIiIiIgEEsy6qx8CFwNznXOv8eMmBW2cc7eH0piZ3RVifCIiIiIiPsGMeT0e+AjIoGmbFGBmiaHWiSUa8yoiIiISfU1abcDMVjrnTgT+mx83KcgHyoFFkQxURERERCQQrfMaAvW8ioiIiERfk3peG/AtUNiEmEREREREQhZWz+vBSj2vIiIiItHX1HVeg2nAOefaO+c6R+J6IiIiIiL1aVLy6pw72Tn3KrATzzCCb2qdz3LOPeqc+4dzLqUpbYmIiIiIhJ28Oud+iWe1gQuANniWz/Lr2jWzYuBI4FpgcAjXvtg597BzboFzbpdzzpxzTzdS53Tn3NvOuSLn3F7n3OfOudHOubhenktEREREfhRW8uqcOwb4J55tY6cCvYFtDRR/Ek9Se34ITfwRuBE4Cfg+iHjOB+YD/YF/A9OAFOAB4PkQ2hURERGRGBbuagNj8CSH08xsNIBzrrKBsrOr3/uGcP2bgY3AWmAAMKehgs65TDyJdCWQb2afVh//U3XbFzvnLjUzJbEiIiIicS7cYQMFeHbYmtRYQTP7AdgLBD2Zy8zmmNnXFty0/ouBDsDz3sS1+hr78fTgAlwfbNsiIiIiErvCTV4PBUrMbGOQ5fcBrcJsqzEF1e/v1nNuPp7E+XTnXGqU2hcRERGRZhLusIFSIM055xrrHXXOtQLaATvCbKsxParf19Q+YWYVzrl1wHHAUcBXkWgwPz8/rHpz586NRPMiIiIicSHcnCmQcHte1+OZrNUtiLI/AxKBL8NsqzFtq993NnDee7xdlNoXERERkWYSbs/ru8AJwE3ADQ0Vcs7l/P/27jxOsqo++P/n2z3TDdOzMMyKLA6yiqIRUfZlhgAm8KgYF5Kfuz4JvjTERDAmjyYYjY8RjP4Qoz5uREwUH40bagAZZgDFBQQRcIYZ2WEWmH2Rnpme8/xxb01Xd1dVV1V3dfet+bznNa/quvfcc8+9de693zp17znAx8juj/1Bk+saqVL3XaM2LJYtqJIkScNrNmYqjbBVSbMtr58AtgIXRcQ/RsS0QSvcNyL+DLiDrJ/XdcBnm1zXcEotqzOqzJ8+KJ0kSZIKqqngNaW0BvgzYCfwD8BTwCyAiLgPWA9cAzyb7P7YP00pbR6NAlewPH89cvCMiJhEFjzvYtDoX5IkSSqepkfYSildRzYowJ1kfb5OIvuJ/rlAd/73XcDpKaWbRl7Uqkr9yL6swrzTgSnAT1NKvS0sgyRJksZAs/e8ApBS+gXw0oh4AXAqWRdancBq4Cfl/a620DfJ+pu9MCI+VTZIwT7Ah/M0nxmDckiSJKnFor5xAMZWRLwSeGX+dj5wLtnP/rfm055OKV0yKP03gWfIhoNdD7ycrButbwKvrXPAg+HKlQAm4j6TJElqF6UHtlJKQ57cmqjB62XAP9ZI8khKacGgZU4B/hfZMLT7kA0t+yXgypRStaFrGy2XwaskSVKLtTR4jYgu4GzgeGAuWZdUTwG/BH6cUtoxohVMIAavkiRJrVcreB3RPa8R8efAh4DZVZI8HRHvTyl9fiTrkSRJkmAELa8R8S/AJfQPAvAE8Hj+90HAgfnfCbg8pfS+EZRzQrDlVZIkqfVG/baBiDgDuDl/+y3gAymlZYPSHEXWKvtqsgD2zJTSrRSYwaskSVLr1Qpem+3ntTQk7BdTSq8ZHLjmK1ueUnot8EWy1tl3NbkuSZIkCWi+5fVxsi6snpVSWjtM2nnAk8CqlNJBTZVygrDlVZIkqfVacdvAM8C2lNKsOtOvA3pSSvs0vLIJxOC1QdvWwbLrYOsamDoPjj4feuqqMpIkaS/WiuD1KWAGMD2l9MwwafcFNgGbU0rVeiUoBIPXOqUEt14BSy+HvrJReTu74YxL4bRLIIbURUmSJKA197zeQzYM7FvrSPtWsi65ft3kulQ0t14Biz88MHCF7P3iD2fzJUmSmtBs8PofZA9hfTwi3lYtUUS8Hfg4WW8D1zS5LrXatnVw57/D0o9lr9vWjSyvpZfXTrP08pGtQ5Ik7bWavW2gA7gJOIMsMH2crOusJ/L3BwMLyfp6DWAJcFYq+O/tbXfbQCt+3r/z3+H7Fw+f7uWfguPe2FjekiRprzDqI2yllHZHxCuALwGvIgtW3zB4vfnrt4C3FT1wbUuln/cHK/28D3D6pY3luXVNfem2rG4sX0mSJEYwPGxKaTPw6oh4CXAhcDwwN5+9FrgD+HpK6ZcjLqVG37Z1pKWXU6tdNS29nHjxWxvqIWBb1yx66ki3tWs2U+vOVZIkKdN08FqSB6cGqEWz7Dpi8ANVg0RfLyz/QUM/71/fdzznpcl0x86qaXrTZG7c/WIuqDtXSZKkTFMPbEVEV0S8ICKOriPt0Xnayc2sS62xff2T9aVb90RD+T7eO4Urd9UOS6/cdQGPPTOloXwlSZKg+d4GXgfcBby7jrT/K0/76ibXpRa4d3N940Xcu3nfhvKdO62bT/e9gst3vpbeNPD7Sm+azOU7X8un+17BvOndDeUrSZIEzd828Cf5az3dX30R+P/IgtevNbk+jbK7e07lhXX8vH93zym8tIF8z3nefP7he/fx6V2v5D/7FnF2553MYwNrmMmNfS9mA9PpmtTB2cfMH/lGSJKkvU6zLa/Pz1/rGXjgzvz12CbXpRaYvv+8un7enzFrXkP57t/TxcWLDgdgA9P5Rt9CPtX3Kr7Rt5ANTAfg4kWHs39PV3MFlyRJe7VmW16fBWxKKW0dLmFKaUtEbAQOaHJdaoFznjefE793AeyEiyd9e0ALbG+azJW7LuDzcQE/a6KF9J0Ls+D1ysUr2bFr957pXZM6uHjR4XvmS5IkNarZQQo2At0ppWFviIysl9ntwM6U0vTGizhxtNsgBVctXsEVNzzATDZX/Hn/knOO5F2Ljmg6//XbdnDj/atZs7mXedO7OfuY+ba4SpKkYdUapKDZ4PUu4AXAqSml24dJewpwK3BfSqnQtw60W/CaUuLTN6+s2UIajY6wJUmSNEKtCF4/BlxCFpSelVLaVSXdJGAxcArwyZTSexpe2QTSbsFriS2kkiRpImlF8HoQ8ADQTRbA/nVK6a5BaY4DPgGcBjwDHJ1SerThlU0g7Rq8SpIkTSSjHrzmmb4BuLps0mrgESABhwLzgMjfvyml9NWmVjSBGLxKkiS1XkuC1zzjPwauAhZUSfIg8K6U0n83vZIJxOBVkiSp9VoWvOaZdwILgZOBUr9Kq4CfAjenlHZXW7ZoDF4lSZJar6XB696kUMHrtnWw7DrYugamzoOjz4eeWeNdKkmSpGEZvI6SQgSvKcGtV8DSy6Gvt396ZzeccSmcdgmMpPsrg2JJktRiBq+jpBDB6y2Xw+IPV5+/6P1w+qWN59vqoFiSJCln8DpKJnzwum0d/OtzBwaXg3V2w9/8tvHW0lYFxZIkSYPUCl47xrw0ap1l19UOXCGbv/wHjeW7bV3W4lrL0suzdJIkSS1k8NpOtq6pL92W1Y3l26qgWJIkqUEGr+1k6rz60k2bP3yacq0KiiVJkhpk8NpOjj4/u6e1ls5uOOq8xvJtVVAsSZLUIIPXdtIzK3vyv5YzLm38Ya1WBcWSJEkNMnhtM+nU93D7gnfQmyYPmN6bJnP7gneQTn1P45n2zCIN05NAOr2JoFiSJKlBk8a7ABpdn17yO65YdhozeSFnd97JPDawhpnc2PdiNiybziVLfse7Fh3ReL67XsHvdz7AxZO+TXfs3DO9N03myl0XsO+uV/Cu0dwQSZKkCuzntQETvZ/X9dt2cOL/vokdu3ZXTdM1qYOf/d1Z7N/T1VS+M9k8NChmelP5SpIkVVKrn1dbXtvIDfetrhm4AuzYtZsb71/N615ySFP5bmA63+hbOCr5SpIkNcp7XtvI2i3D9MWaW7O5vnStzleSJKlRBq9tZO60YXoEyM2bXl+6VucrSZLUKIPXNnLO8+bTNan2R9o1qYOzj2msP9ZW5StJktQog9c2sn9PFxcvOrxmmosXHd7wQ1WtyleSJKlRPrDVZt65MAsyr1y8csDDW12TOrh40eF75k+UfCVJkhphV1kNmOhdZZVbv20HN96/mjWbe5k3vZuzj5k/Ki2jrcpXkiSppFZXWQavDShS8CpJklRUtYJX73mVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYUxabwLoOZseGoVK275On2bV9M5fT5HnH4hM+ccMN7FkiRJaqlIKY13GQojIhLAeO6ztHs3P/vK33PcQ1+gO3bumd6bJvOrQ9/OiW/8CNFhg7okSSquiAAgpRRD5hm81m8iBK+3X/0+Tnr4M9XnL3gHJ735o2NYIkmSpNFVK3htqya6iHg4IlKV/6vHu3wjteGpVRz30BdqpjnuoS+w4alVY1QiSZKksdWO97xuAj5ZYfrWsS7IaFtxy9d5admtApV0x05+fcu1vPRP3j1GpZIkSRo77Ri8bkwpXTbehWiFvs31NR73bbblVZIktae2um2g3XVOn19nOnsdkCRJ7amtHtiKiIeBbuBS4BBgG3APcEtKqW8U8h/XB7Y2PLWKKVcdO6CXgcF602S2v+s3dpslSZIKq9YDW+1428B84JpB0x6KiLeklJaOxgrOPPPMppZbsmTJiNY7c84B3H7o22v2NvCrQ9/OSQaukiRpAmg2Zqql3W4b+DJwFlkA2wMcC3wOWAD8KCJeOH5FGx0nvvEj3L7gHfSmyQOm96bJ3L7gHZz4xo+MU8kkSZJar61uG6gmIq4A3gN8J6V0wQjyGfd+XkuyEbaupW/zKjqnH8ARp7/OWwUkSVJb2OsHKYiIw4EVwPqU0qwR5DNhgldJkqR2tdcMUlDD2vy1Z1xLIUmSpBHZW4LXk/LXB8e1FJIkSRqRtgleI+J5EbF/henPBq7K3351bEslSZKk0dROXWW9BnhfRNwMPARsAQ4DzgP2AX4IXDF+xZMkSdJItVPwejNwFPAistsEeoCNwG1k/b5ek3zSSpIkqdD2it4GRou9DUiSJLWevQ1IkiSpLRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBiTxrsAklQ067ft4Ib7VrN2Sy9zp3VzzvPms39P13gXS5L2CpFSGu8yFEZEJAD32egwANB4GEm9Synx6ZtXcuXilezYtXvP9K5JHVy86HDeufBwIqJVRR+Wx5SkdlE6l6aUhpxUDV4bYPA6OkYrAPBCnXE/DFRtf4xGvbtq8QquuOGBqvMvOedI3rXoiFHblsFauW0qLs8BreO+HT8Gr6PE4HV0jDQAaNcLdaMnyXbdD80abn8k4OMjqHfrt+3gxP9904C8B+ua1MHP/u6sUb+4tXrbJoLRDBL2loDDc0DrjNe+3Vvqbj32muA1Ig4C/gl4GTALWAV8B/hgSmnDKORv8DpC5QHATDZzTuedzGUDa5nJDX0vBuBlk+7k0pNnsP9++0ME7NgKk3v2/L348eA9vzmIDUyvuI4BF+pt62DZdWxf/yT3bt6Hu3tOpbNnFgRs6+3joO7tnNt5Bz071sHUeXD0+dlyy66DrWsGrHfA38OlrXd+zyzStqe5+Ttf5oFl9zA7bWBt2o9HmM+SeAlvWHRc/0ky3xa2roGp8/jC08fw4ZvXDtiPW9mXINHDM7zomKNZ9Mq3Qs+syh9GKb/1D8GWVTDtANj/0D3lKl/ftq5ZXN93PI/3Tuk/obJl6LZtXdufV8+cyvuuUtoq6y3fT0PKPWj+x/57Gf+25HcV69UGpjO7YwtnxR1DppfMm7SVxedtzepChTLevraDG+9fQw/PDNjPpb9ns4m5sZFjjjqSY455YeXtqVIHqn42+TK3rnyaXyx/dMB6y7ehI2B3ouq2l4JqYM+FcUDdr7dc9dSbBmx4ahUPLP0ajz7yEMs3JHan/m2rWf+rHUt5wHHN4l9xZvrlnv0wJK9mVdr2avU8L9eG7TtYccvX6du8ms7p8zni9AuZOeeAoXnWOt9UOWZqnQNms4mFB+7muUccOWwZgf5jnW6Wr95K7/bN/eWd0lV/Gaud/8rrVem8vGYla554iKdjFmnmgv59M9w5oFWfX9m6vnD9L1lxy7U8m9XMjY17zsulY6ridaZsezqnzeGo+dPp4Zk9+2b7hlVVtzctu46f33M/339wNz/ceRwA53TeyQEdG3nB0Uex8JVvIXpmD7891c6lw52jGzjXjpW9IniNiMOAnwJzge8Cy4CXAguB5cApKaV1I1yHwesIff0Xj/K+/7qHd3Z+l4snfZvu2LlnXl8KAuiI4ffvztTJnbsP5/E0h7VpJk+x356L+jMdU3j3WYfT88hNpEd+QuzetWe5XamDpX3H8pP0fBZ13M0JHcuYHH1lOQcpgkjVW9fK0xIBVdImOoBEUG17OmHmwfRteJROhuaxM3Xy891Hkw4/m9M674cHb4LdfWXzO3gizebAWDdoG8rKEJ3EYYvgsIWwY1t2AgNYcT3p4duIVGm5rFxsehzK9l3/Pp/NIR3reHHnCjorLt+c3XSwufsApu9YTUd5vjEJDn4pzDgYNj1GevwXAz7TFJN4ZNqL+Pd1R3Jax72c3vEbJkX//qy2n4avCyNVeT+WG/L5TJoCDy0d8llXU/6ZHBxPc1zHygHbUKpDN+1+EUfPn8aGtU8wO63n4HhqSNoBSvt8v0P6L3IAK66Hh2+DSp97TIIFp7Dt2WexfE0W/HRPmZZdwHes33PB3NY1k2VPbqbroR9zdO9vqpah9Pl0HbGQ07gLHvlJ1f1IdMJhi7i173nEyhuHfJalvFbPOonn7Ed/ucoCi5pfvIbb9ip2E5AGntN2pk5WzTyOg1/yCuKhW0i/u6nKcTi8es4Bw+uASFDlutaXgo6gxnmsmkrnx+yYSBsfq7jNO1MH2/Y9gBk71gw5xnc+63ieSHN4umNgYFgKtvu2rGF2Ws/8Aw9l3/3mj+Dz62TXjINIGx9jclQ/Ly/lON591hHZdabqubS2attb6Vq4i06enHEcTx9wJkfNn0bPIzfV2J7hzz0VDTrX8vgvBi7f2Q1nXAqnXZLt3xbbW4LX64FzgItTSp8qm/6vwF8Dn0spXTTCdRi8jtCVN62gd/HHuHTyN8a7KJIkqVGL3g+nX9ry1bR98BoRzwF+BzwMHJZS/1e9iJhGdvtAAHNTSttGsB6D1xH6r9t+zXk3njWgxVWSJBVD6uwm/ua3Lb+FoFbw2i6DFCzKX28oD1wBUkpbgJ8AU4ATx7pgGujczjsMXCVJKqjo64XlPxjXMrTLIAVH5a/VHrddQXZLwZHATSNd2ZlnntnUckuWLBnpqguvZ8eIbjuWJEnjbPu6J5hSZ9pmY6Za2qXldUb+uqnK/NL0/cagLKpl6rzxLoEkSRqBezfvO67rb5eW1+GU7pcYlZtVbUEdgaPPhx9eCn29410SSZLUoN40mbt7TuGldaZvNmaq1a1du7S8llpWZ1SZP31QOo2XnllZVxuSJKlwrtx1ATNmje+vqO3S8ro8fz2yyvzS0DLVh6DR2Dntkux16eW2wEqSVAC9aTJX7rqAz8cF/OyY+eNalnbpKuswYCW1u8rqAObYVdYEUhq9495vVR1M4PZ4AX/4nB5OOOqQ7N6PJjoKhw5SjcECEh3EoE66Sx1R/3j3cUzqCP7wOVM44ciDiQodyJenDRJnddxVtZP0n3Isb5q9nIM33zVge0udzT+W5nDIcB3IV+r8O+8gnsPPaaqT+8fSHJ7d8XSFgQdqdypeWm864lxuXbmO+1esYFbawFNpP9ayHwFM5fdsYd8Bf0/qCF727OC4/XtZ9chy5m68Z1Dn+nnn3b2rq3ZmPuwADUT2s1PZfhp2IIzSfjzi3P6O6/PRZbZ27c8Da7bQu20z3VOm552k/x66pmbLbl0LW1bDxkeGdO6dYhKPzXgRX1l3NCfzG07vuGfAQApV5Z3vc9iirMP1rqlZPa1xHCQ6iTo7KE8xiShtL9R/fJUPYrDx0aGdmTdgZ+rknknP4/vPvIDTO+6tuW92RycdC06tq55XGiBgNJSOmSdiHn984h+w78x8tKzBn8+gQRV2E/m/yttW8ZyXd7SfJveUHV/rK58jatRdps6tq4xEJ70LFnJf94uYvXoJB226i440cD6HLSIdtpCfL3uUHz+4nb7dKT+u9+GMzns5o+M3dDD0HDKkPsYkdjzreB5Pc+ja+jjP2vKbAevqP8afrjhYQLP6z3lz95ynFh05h9O5c8j+aOS8vGfAEOZy4EGH0r3ffI6cN417H3piz36aw0bmxIaKeVU7p+0mSAk6a52zah2Lw5yjJ3XAG2ct4+BNv6pyrs2ub9/bfTI3VhpdrIXavp9XcJCCwtu2Dpb/gO3rnuDezftyd88pzJg1j7OPqTCuc56WdQ9WPzmX/p42H446L3u//AdZcDHM/K1ds7lx94t57JkpzJvePbQMZUP5fe/B3fxo53F7hhntCDho5r78fuNaFsYdzGMDa5jJ4nQ8Lz76cD7yqmOZNbW7Zh5dkzp472mzedvs+4mta4bfntK0wcP6Vdresr8rbidbKudbaZ/Pes6Q9a7ftoMb71/Nms29TO3uBIKtvbsG/F1pn254ahUrbrmWvs2r6Jx+AEec/rr+YSIHrXfr1EP2lPvgfbZzdsedTN3xdN2f61NbdvAH237C86f/nilTZwxdZrSGoaywH0v7Z9O6NZXLUGkbag3TWuszqVQHygOaCp9f1XzLj69h6ttW9hka5O9YV/FLQOmzTlNmDRg2+uzOO5nHhgFfetZ17M+lF//N0KFVaxzX2dCsWb2iayoBpB1bq3/5KNveW1Y8zS+XPzKgDGuYWf8FvFIdgCHl3b510/DnvFz58TWg7jdbd6vU03rnl5en7nNIHesqPzeVtnPn2hUsvuNe1uyeUfGL8Rw2MrdjE2e95FhmzjkoyzM/phLw82WPVDzXXrzo8IFDEFcow9xpXTy6/vd867Zfc0b65ZC6OVy9KN9P9z25ievvWzOgnpcv//6Fc3n7nN8O2E/l9XhA3R2tz6yna8Bxn7au5s713fzXw12191eL7S3B6+DhYX8LnEA2POwDwMkOD6vRVu1EUPUE0UAe0t7kqsUruOKG6nd2jVVrT0lKiU/fvJIrF69kx67+lr+xvoCr30jryEjPteu37eCG+1bx/V+v4ucPrWfX7v5YoN56UaR6Nd7Xpr0ieAWIiIOBfwJeBswiu13gO8AHU0rrRyF/g1dJaoGJelEf7wu4+k2kOjIagbD1qra9JnhtNYNXSWotL+oajnVk72DwOkoMXiVJklqvVvDaLv28trUzzzyzJcOrFZn7ZCD3x0Duj6HcJwO5PwZyfwzlPhloIu0Pg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmFESmm8y1AYEeHOkiRJGiMppRg8zZZXSZIkFYYtr5IkSSoMW14lSZJUGAavkiRJKgyDV0mSJBWGwaskSZIKw+BVkiRJhWHwKkmSpMIweJUkSVJhGLxKkiSpMAxeJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHidwCLioIj4UkQ8GRG9EfFwRHwyImaOd9laISJmRcTbI+LbEbEyIn4fEZsi4raIeFtEdAxKvyAiUo3/Xx+vbRlN+edebRtXV1nm5Ij4YUSsj4jtEXFPRLw7IjrHuvyjKSLePMxnniKiryx929SRiHh1RHwqIm6NiM15+b86zDIN14OIOD8iluTH3taI+HlEvGn0t2hkGtkfEXFERPxtRCyOiMciYkdErImI70bEwirLDFfXLmrtFjamwf3R9HEREW+KiF/kdWNTXlfOb92WNa/BfXJ1HeeWmwYtU5g6Eg1eX8uWm5DnkEmjmZlGT0QcBvwUmAt8F1gGvBT4K+BlEXFKSmndOBaxFV4DfAZYBdwMPArMA14FfAH4o4h4TUopDVru18B3KuR3bwvLOtY2AZ+sMH3r4AkRs9ZyfAAAEM1JREFU8QrgW8AzwLXAeuB/AJ8ATiHbz0V1N/DBKvNOAxYBP6owrx3qyPuBF5J95o8DR9dK3Ew9iIh3AZ8C1gFfBXYArwaujohjU0qXjNbGjIJG9seHgNcB9wM/JNsXRwEvB14eEX+VUrqyyrLfJat3g93RZLlbpaH6kWvouIiIK4D35Pl/HugCLgS+HxF/mVK6qolyt1Ij++Q7wMNV5r0BeA6Vzy1QjDrS8PV1Qp9DUkr+n4D/geuBBPzloOn/mk//7HiXsQXbvCg/MDoGTZ+fH2gJ+JOy6QvyaVePd9lbvF8eBh6uM+10YC3QCxxfNn0fsi9DCbhwvLepRfvp9nz7Xt6OdQRYCBwBBHBmvl1fHa16kO+rZ8guOgvKps8EVubLnDTe+6HJ/fFm4EUVpp9BdnHtBQ6osEwC3jze29qC/dHwcQGcnC+zEpg5KK91ed1ZMJJtGM99UiOP/YDteR2ZXdQ6QuPX1wl9DvG2gQkoIp4DnEMWtHx60Ox/BLYBb4iInjEuWkullBanlL6fUto9aPpq4LP52zPHvGDF8mpgDvD1lNKeb/0ppWfIWiEA3jEeBWuliHg+cCLwBPCDcS5OS6SUbk4prUj51WAYzdSDtwLdwFUppYfLltkAfCR/O2F+Bm1kf6SUrk4p3VVh+lJgCVkL4smjX8qx02D9aEbps//nvE6U1vsw2XWqG3hLi9bdlFHaJ28A9gX+K6X09CgVbcw1cX2d0OcQbxuYmBblrzdUqGhbIuInZMHticBNgxduUzvz110V5j0rIv4CmEX2je/2lNI9Y1aysdEdEa8HDiH78nIPcEtKqW9QulLd+e8KedxC1oJwckR0p5R6W1basfcX+esXK+wT2DvqSLlm6kGtZX40KE07qXVuAfiDiHg3WYvTE8DNKaXHx6RkrdfIcTFc/fhAnuYfR72U4+t/5q//p0aaoteRSsfAhD6HGLxOTEflrw9Umb+CLHg9kr0geI2IScAb87eVDoqz8//lyywB3pRSerS1pRsz84FrBk17KCLekrcelVStOymlXRHxEPA8svu3ftuSko6xiNgXeD2wm+zerUr2hjpSrpl6UGuZVRGxDTgoIqaklLa3oMxjLiKeDZxFdiG+pUqyvxr0vi8ivgC8O2+FKrK6jov8V74Dga0ppVUV8lmRvx7ZonKOi4g4CTgWeCCldHONpIWtIzWurxP6HOJtAxPTjPx1U5X5pen7jUFZJoKPAs8HfphSur5s+nayBzFeTHZPzUyye9huJvv546Y2ubXiy2QX2PlAD9nJ9HNk9xf9KCJeWJZ2b6w7ryXbnh+llB4bNG9vqSODNVMP6l1mRpX5hRIR3cB/kP3MeVn5T+G5h4C/JLsg9wDPIqtrD5O19H9pzAo7+ho9LvbG8wrAn+evn68yvx3qSLXr64Q+hxi8FlPkr626t2nCiIiLyZ5uXUZ279EeKaW1KaV/SCn9KqW0Mf9/C1mr9M+Bw4G3j3mhR1lK6YP5/UprUkrbU0r3ppQuInt4b1/gsgaya8e6U7rAfG7wjL2ljjShmXrQNnUn7+bnGrInpq8FrhicJqW0NKV0VUrpgfy4W5VS+r9kDwFtAP500BfHwmjhcVH4ulESETPIAtEdwNWV0hS9jtS6vtazeP46LucQg9eJabhvJ9MHpWtLEfFO4P8n695mYUppfT3LpZR20f/z8ektKt5EULrJvnwb96q6ExHHkD1o8zhZF0h12QvqSDP1oN5lNo+gXOMuD1y/StbNzzeA1zfyQE/eul+qa21Vd2ocF8PVjeFa3Iro9cAUmnhQqwh1pI7r64Q+hxi8TkzL89dq9w8dkb9Wuye28PKb368i629wYf5EZCOeyl/b8SfhkrX5a/k2Vq07+b1Nh5LdlP9ga4s2ZoZ7UKuWdq4jzdSDWsscQLafHi/y/a75tn+NrG/S/wT+LA/YGtXOdWfItqWUtpE9iDQ1rwuDteM1qfSg1pBfdOo0YetIndfXCX0OMXidmEo3hp8zeNSLiJhG9lPX74GfjXXBxkJE/C1ZJ8h3kx1Ya4dZpJIT89d2CdIqOSl/Ld/GxfnryyqkP52sJeGn7dDTQETsQ/ZT127gi01k0c51pJl6UGuZPxqUpnAiogv4JlmL61eANzTxhafkhPy1HetOteOiretHuYg4gWxwgwdSSkuazGZC1pEGrq8T+xySJkDnuf6v2KHwXjdIQb59H8i37w5g/2HSngB0VZi+iKyj5AScPN7bNML98bxK+wF4NtkTvgn4+7Lp08m+8bf9IAVkgWsCvr+31RHqG6SgoXpA1pJSmEEKGtwf3WT9/yayn8U76sjztArTAvi7PJ+ngOnjve1N7o+GjwsKOEhBI/tkUNov5mnf0051pMHr64Q+h0SesSaYCsPD/pbshLOQ7KeZk1ObDQ+bj318NdBHNrxcpfunHk4pXZ2nX0IW3C0hu+cR4AX09yP3gZTSh1tW4DEQEZcB7yNrjX8I2AIcBpxHdhL5IXBBSmlH2TKvJGthegb4OtmQfi8neyL2m8BrUxsc+BFxK3Aq2Yha36+SZgltUkfyz/WV+dv5wLlkrTq35tOeTmVDLzZTDyLiL4EryS4+19I/tONBwMfTBBoetpH9ERFfJhsN6Wng36j8wMiSVNbKFhGJ7Fz7S7KfzGeQ/er1fLKn9S9IKd0wqhs1Ag3ujyU0cVxExMeBv8mX+SbZ4A6vI+sndsIND9voMZMvMx14EpgMHJhq3O9apDrS6PU1X2binkPG+5uA/2t+8zmYrJukVXkFeITsBuua35iK+p/sqfk0zP8lZenfBlxH1i3JVrJviI/mB8yQb8RF/E/Wfc3XyJ4G3UjWmfRTwI1kffNFleVOIQtsN5DdYvIb4K+BzvHeplHaL8/N68NjtbapnepIHcfHw6NRD8iGkFxK9kVpG9mF+U3jvf0j2R9kQdpw55bLBuV/eb4fniS7eG/Pj8OrgOeM9/aPcH80fVwAb8rrxLa8jiwFzh/v7R/pPilb5h35vK/VkX9h6kgd+2LA9bVsuQl5DrHlVZIkSYXhA1uSJEkqDINXSZIkFYbBqyRJkgrD4FWSJEmFYfAqSZKkwjB4lSRJUmEYvEqSJKkwDF4lqQAi4syISPmoPqOZ74JSvhGxYKyXl6RGGbxKkiSpMCaNdwEkSXXZDiwf70JI0ngzeJWkAkgp/QI4erzLIUnjzdsGJEmSVBgGr5L2ehExNyJ25g8dvXyYtB/K060sm3ZIRLwzIn4QEQ9ExLaI2BoR90fEJyPikBr5LcnzuywiJkfEeyLijojYmE8/M09X9YGtiOiIiFMi4qMR8bOIeDwidkTEuohYGhEXRcTkOvfFERFxdZ5Hb0Q8GhGfjYgD61m+Rr6vjIjvRMSTedk2RMQtw5UtIl4bET+KiDX5Z7QxIlZExPfyfb7PSMolqXgipVF9cFWSCikirgPOA76ZUnpNlTQB/A44FLgspfTBfPoS4IyypJuAafQ3EGwCzk8p3VYhz9Ky/wKcBpwM7AK2ADOBhSmlJXkQezNASikG5bEAeKhs0i6ye2Snl027FTg3pfT7GsteCHw+L/tWoBPYN5+3Hjg7pfSrGssfmlJ6eND8qcDXgPPLJm/O11HajtuB81JKGwYt+0XgrWWTtpLt0yll04asU1J7s+VVkjJfyV//R0TsVyXNKWSBK8A1ZdPvBd4HHANMSSntB3QDJwD/DcwAro2IfanuncALgLcA01NK+wOzgXvqKPsu4LvA64ADge6U0gyyAPEtwJNkgfE/D5PP58gC0RNSStOAHuBc4FFgf+DbETGtjvKUu4YscF0J/BnZts0gC0BfATwInAR8qXyhiDiVLHDdDfwtMCulNC2l1EO2X84F/h3Y0WB5JBWcLa+SBOQ/P68mCzT/IqX0fyqk+Rzw58BtKaXT6sy3E/gVWWD6hpTSVwfNX0J/q+3LU0rfr5LPmVRpea2jDMcDvwS2AbNTSs+UzVtAf8vpOuCYlNLaQcs/F7gb6ALem1K6vMryA1pBI+I84Dqy/Xp8SumJCmU7CFhGFii/KKV0dz79vWSt0TeklM5tZHsltTdbXiUJyAO6/5u/fcPg+RHRDbw2f3vN4Pk18u0ja30FOLVG0vuqBa4jlVK6A1hLFiD+QY2knx0cuObL/xb4Zv72wgZW/fb89ZpKgWue9+PkQTlZa2rJxvx1Tv4FQJIAg1dJKle6deCUiDh00Lzzgf2AXuAbgxeMiNPyB52W5Q9rlUadSsB782QH1Vj3T0ZS8Ijoyh9+uiF/KOqZQWWYW0cZFtcx7wX1PvxFf7D+5xGxutp/4A/zdM8uW/bHwDPAi4BbI+JtFT4TSXsh+3mVpH63kf0EfijweuBDZfNKrbHfSyltLF8oIv6F/gAVoA/YQP/9mFPJWj17aqx7SItnvSJiLlmwd2zZ5GeAp/OyAMwha7CoVYaKraOD5k0iu/91zTBlmkx2bypkt2LMqJU+t+dBrJTSgxHxduCzZPfEnpTn+xRZS+1/kn0W3vsm7WVseZWkXB4Ile5J3XPrQETMAv44f/uV8mUi4mz6A9d/Iwsgu1NK+6eU5qeU5gOfKCWvsfq+GvOG84l8vevIHnI6IKW0b0ppTlkZnqyjDKMZCJb/1H9hSinq+P/mAYVJ6T/IWmMvAq4FHiMLwl8LfAdYGhHlPSpI2gsYvErSQKXg9IiIODH/+3XAZOAp+u9fLSndA3p9SumdKaV78/tcy81vTVH3tHC+Kn/7rpTSl1NKqwel6aS/FbSWWrcUlPp53UXWbVZN+T3Em/K3x9ZKO0w+61NKn0spXZhSOgQ4HPgoWaB9GnBZs3lLKiaDV0kqk1JaSdbvKPS3vpZev5ZS2jVokYPz17sq5Zf3DbtoVAs50Byg1FF/xTKQ3XtaT2f+C+uYd09KaWedZSvdx/uaiBiV601K6Xcppb8ju20A4OzRyFdScRi8StJQpdbX10XEMcCJg6aXK7UuvrBKXhcBzxnFsg22mf6f+4eUISImMXz/riUXRcSQFtqIOAp4df722gbKVupu7Ejg0loJI6InIrrK3ncPk3dpsIWR3G4hqYAMXiVpqGvJHraaRdYRPsBvU0p3Vkhbuo3gjyLiAxHRAxAR+0XE3wOfIrsXtSVSSlvpb+H814hYVGrljIjnAz8Ejifr43U4k4EbI+Il+fIREX8IXE826MJjZA9Q1Vu27wLfzt9+NCI+ExFHlubnPSSckD/w9gj9PSIAXBUR34iIP8kfSCstMzUiLgLemE/6Yb3lkdQeDF4laZB8mNLr8rfH56+VWl1L02/N//4nYEtErCcLWP+ZLLj9TIuKWvJusuD0QOAmYHtEbAZ+Q/Zz//8k63lgOH8BHAb8IiK2kA3HeiPZQ1MbgVellDY3WLbXA1/P/74IWJ53JbaerPX0Z2QPvM1i4ANjk4HXkPUvuyYitkTEBrJhcz9DNmDCbdTfqiypTRi8SlJl5cHqbvp7IRggv//zHOCDwAPATrIn+n8BvAN4OS3+aTtvEX4pWf+zT5Od27fk709OKdU7qMLPyYL1r5DdDjGJrIuszwPH5oMdNFq27SmlPyULoq8hGw62g6z7sLVk/ce+Fzhi0EAGHwIuJmu5XUb2oFhpmRvJelU4M6VUT4uypDbi8LCSJEkqDFteJUmSVBgGr5IkSSoMg1dJkiQVhsGrJEmSCsPgVZIkSYVh8CpJkqTCMHiVJElSYRi8SpIkqTAMXiVJklQYBq+SJEkqDINXSZIkFcb/A5uEfrp0XuIwAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 864x648 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Generate some sparse data\n",
    "np.random.seed(6250)\n",
    "\n",
    "# number of samples is 120 and the dimension of varaibles os 200\n",
    "n_samples, n_features = 100, 200\n",
    "X = np.random.randn(n_samples, n_features)\n",
    "\n",
    "# generate coefficients\n",
    "LargeCoef = 20\n",
    "SenCoef = np.zeros(n_features)\n",
    "for i in range(LargeCoef):\n",
    "    SenCoef[i] = 5 + 20*np.random.rand()\n",
    "for i in range(LargeCoef, n_features):\n",
    "    SenCoef[i] = 0.1*np.random.rand()\n",
    "\n",
    "# generate y\n",
    "y = np.dot(X, SenCoef) + np.random.normal(loc=0, scale=0.001, size=n_samples)\n",
    "\n",
    "# cross-validation and fit\n",
    "reg = LassoCV(cv=50).fit(X, y)\n",
    "print(\"Residue of the fiting\")\n",
    "print(reg.score(X, y))\n",
    "print(\"Weight\")\n",
    "print(reg.alpha_)\n",
    "\n",
    "# plot the results\n",
    "fig, ax = plt.subplots(figsize=(12, 9))\n",
    "fig.subplots_adjust(bottom=0.15, left=0.2)\n",
    "\n",
    "ax.scatter(np.where(reg.coef_)[0], reg.coef_[reg.coef_ != 0], label='Lasso coefficients', linewidth=4)\n",
    "ax.scatter(np.where(SenCoef)[0], SenCoef[SenCoef != 0], label='Actual coefficients', linewidth=4)\n",
    "\n",
    "# ax.set_xlim([0,n_features])\n",
    "# ax.set_ylim([0, 0])\n",
    "ax.set_xlabel('variables', fontsize=25)\n",
    "ax.set_ylabel('coefficients', fontsize=25)\n",
    "ax.legend(fontsize=24,loc='upper right',frameon=False)\n",
    "\n",
    "ax.tick_params(which = 'both',direction='in',colors='black',\n",
    "               bottom = True,top=True,left=True, right=True,pad=15)\n",
    "ax.tick_params(which = 'major',direction='in',length=15,labelsize=20,width=1.5)\n",
    "ax.tick_params(which = 'minor',direction='in',length=6,width = 1.5)\n",
    "\n",
    "for axis in ['top','bottom','left','right']:\n",
    "    ax.spines[axis].set_linewidth(2.0)\n",
    "    \n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "nbpages": {
     "level": 2,
     "link": "[5.2.2 An example of the lasso regression ](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.html#5.2.2-An-example-of-the-lasso-regression)",
     "section": "5.2.2 An example of the lasso regression "
    }
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<!--NAVIGATION-->\n",
    "< [5.1 Ridge Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.01-Contributed-Example.html) | [Contents](toc.html) | [5.3 Elastic Net Regression](https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.03-Contributed-Example.html)<p><a href=\"https://colab.research.google.com/github/ndcbe/cbe67701-uncertainty-quantification/blob/master/docs/05.02-Contributed-Example.ipynb\"> <img align=\"left\" src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open in Colab\" title=\"Open in Google Colaboratory\"></a><p><a href=\"https://ndcbe.github.io/cbe67701-uncertainty-quantification/05.02-Contributed-Example.ipynb\"> <img align=\"left\" src=\"https://img.shields.io/badge/Github-Download-blue.svg\" alt=\"Download\" title=\"Download Notebook\"></a>"
   ]
  }
 ],
 "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}