group79_other_algorithms2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "---\n",
    "**Model 3**: LassoRegressor with upsampling and  feature engineering.\n",
    "\n",
    "**To check the performance, simply run the following code cell, and result should be print out immediately.**\n",
    "\n",
    "the code for plotting and visualisation are also followed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "File loaded locally! \n",
      "Start training Lasso regressor...\n",
      "Training finished! Start predicting...\n",
      "The Performance of Our Lasso regressor:\n",
      "MSE:   105.419\n",
      "RMSE:  10.267\n",
      "MAD:   2.364\n",
      "MAE:   6.34\n",
      "NLL:   698.689\n",
      "RR:    -0.084\n",
      "AUC:   0.656\n",
      "Explained variance regression score: -0.077\n",
      "Time taken to run Logistic Regression = 1.8126661777496338 seconds\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAwUAAAGDCAYAAAB6E7yVAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdfZyNdf7H8dfHPYtmRSkRS1vYzVa23CWFMIu0Sm52sUm3NimxlZ9abUpR29Ym1OqOStG6HWUrRaR0s1SUyM1IKUlW7sZ8fn9cZ8w5c09z5pqZ834+Hucx1/W9ruuct5lj5nyu6/u9vubuiIiIiIhI4ioTdgAREREREQmXigIRERERkQSnokBEREREJMGpKBARERERSXAqCkREREREEpyKAhERERGRBKeiQEREROLGzG41s8fCzhEWMxtoZkvDzlEYzKyemf3PzMoexbGPmtn/HeXr3m1mNxzNsUfxWseb2Rozq1gUr1ecqCgQEREpxsxso5ntjXwY+8rMnjCzqmHnKih3H+vuVxT280Y+bB+KfF+iHycW9msVIEsnM3vTzHab2Tdm9oaZdS/qHHkxs3ZmlvpTnsPdN7t7VXc/lM9rZSuE3P1qd7/zSF/TzGoB/YFJkfV2ZpYe+VnvNrNPzexPWY5xM9uT5X0xImr7L83sBTP71sx2mdkqM7vRzMq6+9fA68CVR5q1pFNRICIiUvx1c/eqwG+AM4Bb4vEiR3MGOGTLIx9Sox9fZt3JzMoVpC0/OX1/zOwS4AXgKeAk4HhgNNDtSJ+/AK9/xJlLwWsPBBa4+96oti8j/x+qA8OAKWZ2apbjmmV5X9wLYGYNgRXAFuDX7n4McCnQHKgWOXYacFXc/kXFlIoCERGREsLdvwJeJigOADCzimY23sw2m9nXkW4alaO2jzCzbWb2pZldETmL2iiy7Qkzm2hmC8xsD3B+Xs9nZjXNbJ6ZfW9m35nZEjMrE9k20sy2Rp29bR9pv8PMnonK093MPo48x2Izaxy1baOZDY+cud1lZs+bWaWj+V5Fnmukma0C9phZuVzaGkdyfB/J1T3qObJ9f7K8hgH3A3e6+2Puvsvd0939DXcfnGXf8Wa208y+MLMuUe1/inRX2W1mG8zsqqht7cwsNZL5K2Cqmf088jP4JvJ888zspKhjapjZ1MjPe6eZ/dvMfgakACdGnTk/0czKmNlfzGy9me0wsxlmViPyPPUj75VBZrYZeC2qrVxkn4GRzLsj/65+kZ/no0DLyOt8H/W9/FtUzovM7EMz+yHy+p1z+VF2Ad7IaYMHFgDfAafncnxWfwWWufuN7r4t8jyfuntfd/8+ss8K4BdmdnIBn7NUUFEgIiJSQkQ+/HUBPo9qHgf8kqBQaATUIThTTeSD1o1Ah8i283J42r7AXQRnSZfm9XzATUAqUIvgjPitgEfO0g4Bfuvu1YBOwMYc8v8SeBa4IfIcC4C5ZlYhardeQGegAcEHvYH5f2dy1Qf4HZDk7mlZ2wAD5gKvAMcBfwamWexZ56zfn2inAnWBF/PJcQ7wKVATuBd4PFJQAGwHuhKc9f4T8ICZnRl1bG2gBnAyQZeWMsDUyHo9YC/wcNT+TwNVgKaRf9MD7r6H4H3zZZYrKtcDPQjeFycCO4F/Zsl+HtCY4Gd6WKTQ+AfQJfIzbwV86O5rgKvJvIqTlPWbYWZnE1xZuZng59CWHN4vEb+OfO+yiRQ13Qm+r5/ntE8OOpDPzyvyXvkcaFbA5ywVVBSIiIgUf/82s90EXR62A7fD4TPVg4Fh7v6du+8GxgK9I8f1Aqa6+8fu/iPBWdKsZrv7W+6eDuzP5/kOAicAJ7v7QXdf4u4OHAIqAk3MrLy7b3T39Tm81mXAfHdf5O4HgfFAZYIPlBn+4e5fuvt3BB/Yf5PD82RoETnDn/HI+pr/cPctWbqeRLe1AKoC97j7AXd/DZhHUDhk+/64+74sz39s5Ou2PDICbHL3KZG++E8SfA+PB3D3+e6+PnLW+w2CAuXcqGPTgdvdfb+773X3He4+091/jPx87iJS7JnZCQQf/q92952Rn1GOZ9kjrgJuc/dUd98P3AFcYrFdhe5w9z1ZvofR2X5lZpXdfZu7f5zP9yHDIOBfkfdBurtvdfe1ueybBOzO0nZi5ArEXuAl4EZ3/yDLPu9neW9kFDXHkv/Pi8hrZitoSjMVBSIiIsVfj8jZ2HbAaQRnRiE4214FeC/jww+wMNIOwdnfLVHPE72cU1t+z3cfwRnUVyLdRv4C4O6fE5z9vwPYbmbPWc4Dfk8ENmWsRAqRLQRXIzJ8FbX8I8GH9ty87e5JUY+Gefzbcmo7EdgSyZFhU5Y8OT1Hhh2RryfksQ9E/ZsixRlE/l1m1sXM3ragO9b3QDKZP1+Ab6KLETOrYmaTzGyTmf0AvAkkWTDeoS7wnbvvzCdPhpOBl6J+1msICrzjo/bJ8d8fufpwGcFVgW1mNt/MTivg69YFcioac7KTzL7+Gb6MXIGoTnC14oIcjjszy3vj5Uj7DvL/eRF5ze/z3asUUVEgIiJSQkTO+j5BcIYd4FuCs6VNoz78HBMZhAnBGdGTop6ibk5PG7Wc5/O5+253v8ndf0EwkPZGi4wdcPfp7t6G4IOmE3RDyurLyHbg8JWOusDWgn8Xjojn0/YlUNci4yIi6mXJk9NzZPiU4ENzz6MJZ8FtL2cS/DyPj3zQXUDQrSm317+JoNvSOe5enaDrDZFjtgA1zCynM9w5/Tu2EHT/if7wXMndC/Tvd/eX3b0jwYfstcCU/I6Jet2sBVxuVhF0Z8vp9fcDI4Ffm1mPAj7ff8jn5xW5UtII+G8Bn7NUUFEgIiJSsvwd6Ghmv4mc4Z5C0A/9OAAzqxPVVWIG8CcLBtNWIXNsQI7yez4z62pmjSIf5n8gOKt8yMxONbMLIh9y9xEUFjndtnIG8Dsza29m5Qk+4O4Hlv2E78dPsQLYA4wws/Jm1o6g2HmuIAdHuk7dCPyfBQOGq0f6ubcxs8kFeIoKBN2uvgHSLBiAfGE+x1Qj+P5+HxkUfHtUnm0EA4ofsWBAcnkzyygavgaONbNjop7rUeCujAG1ZlbLzC4qQO6M+/l3j4wt2A/8j8yf+dfASVnGikR7nOB92T7y/aqTx1WGBeQ8FgYAdz8ATCCf93aU24FWZnafmdWO/FsamdkzUcXU2cBGd9+U67OUQioKREREShB3/4ZgkGbGRFAjCbr0vB3pTvIfgjPJuHsKQfeK1yP7LI8csz+Pl8j1+YBTIuv/izzXI+6+mOCD7T0EVxq+IhjgemsO2T8F/gA8FNm3G8HtVg8cyfcgSsYdbqIfvy3owZHX7U7QD/9b4BGgfx7923N6jhcJutFcTnDl4Wvgb8DsAhy7m2Cw7wyCbjJ9gTn5HPZ3gnEY3wJvE3TvivZHgrEfawnGn9wQea21BIO8N0S6C50IPBh5vVciY1beJhgUXRBlCIq6Lwnu/nMecG1k22vAx8BXZvZtDv/ud4gMqgZ2EdxdKLc7/TwFJFvUHbVy8C+gnplF3wb2v1neF3+PvPZ6oCVQH/jYzHYRXK1ZSebYhX4EBVNCsaDIFRERkdLOgttFfgRU9My78YgUa2Y2Ftju7n8vgtc6jqBIOSOHgeWlmooCERGRUszMLgbmAz8juPNNursXtP+1iCQIdR8SEREp3a4i6LO+nqDP9zXhxhGR4khXCkREREREEpyuFIiIiIiIJDgVBSIiIiIiCa5c/ruEr2bNml6/fv2wY4iIFFvvvffet+5eK/89Sy/9rRARyVtefytKRFFQv359Vq5cGXYMEZFiy8wSapKdnOhvhYhI3vL6W6HuQyIiIiIiCU5FgYiIiIhIglNRICIiIiKS4FQUiIiIiIgkOBUFIiIiIiIJTkWBiIiIiEiCi+stSc1sI7AbOASkuXtzM6sBPA/UBzYCvdx9ZzxziIiIiIhI7oriSsH57v4bd28eWf8L8Kq7nwK8GlkXEREREZGQhNF96CLgycjyk0CPEDKIiIiIiEhEvIsCB14xs/fM7MpI2/Huvg0g8vW4nA40syvNbKWZrfzmm2/iHFNEREREJHHFdUwB0NrdvzSz44BFZra2oAe6+2RgMkDz5s09XgFFRERERI6GO+zbl/v2ChWgbNnYtr17C/78FStCmahT+O6Qnp79OQtDXK8UuPuXka/bgZeAs4GvzewEgMjX7fHMICIiIiJS2NauhVNOgSpVcn8sW5b9uGOPzfuY6McXX8Qeu2tXzs9ZGOJWFJjZz8ysWsYycCHwETAHGBDZbQAwO14ZREREREQKW2oqXHghrF8fdpLCE8/uQ8cDL5lZxutMd/eFZvYuMMPMBgGbgUvjmEFEpNjZt28f48ePZ/jw4VSqVCnsOCIicgS++w46dYItW4J1s6CbT07K5HD6vVKloBtQQQQfo/N/zsIQt6LA3TcAzXJo3wG0j9friogUZzt37qRHjx68+eabrF69mmeffZYy8foNLyIihW7GDPjkk2C5fHmYPx86diz48d99d/SvnZQErVsf/fF5ifdAYxERidi8eTOdO3dmzZo1AMyYMYO6desyfvz4kJOJiMj+/fD3v8O6ddm3/fnP0Cxyqvvqq2HPHhgxAp588sgKguJMRYGISBH48MMPSU5OZtu2bTHt8+bNY/To0VSvXj2kZCIicugQ/OEP8OKLOW/v0SOzKAC46Sbo2hVOPbVo8hUFFQUiInG2aNEievbsye7du2PaW7duzZw5c1QQiIiEwB1eew02bIA338y9IMhNaSoIQEWBiEhcPfXUUwwaNIi0tLSY9p49e/LMM89ooLGISEgmTYJrrsne3r8/nHtubNvppxdNpjCpKBARiQN35+677+a2227Ltm3o0KFMmDCBsvGYfUZEpBRITYUdO7K3n3ACHHdcbNvnnwd9/AuiXj34+c+D5b59g6sDr76aub1PH5g6NX53+CnOVBSIiBSytLQ0hgwZwqRJk7Jtu//++xk2bFgIqURESob77w/67Ofkvvtg+PDYtsGDYfHigj33tGlBMQBQvXpw56BevaBWraA70NChiVkQgIoCEZFCtWfPHnr37s28efNi2itUqMDTTz9Nr169QkomIlL8PfFE7gVBPFSsCLM1jS6gokBEpNBs376dbt268c4778S0JyUlMXv2bNq2bRtSMhGR+Fm4EP74R/jhh4LtP3w43HVXbFufPjBrFhw4kNl24olQs2bsfrVqZX++hg0Lfu//pKSC7ZeIVBSIiBSSbdu28UnGjDYR9erVIyUlhSZNmoSUSkQkvtq2hRYtIMsF0lwdOpS9LS0ttiBo1gzeeAOOOSb/53vssYK9ruQtQXtNiYgUvmbNmjFz5kzKlSt3eH358uUqCESkVKtSBV54AX71q8J5vjPPhJSUghUEUnh0pUBEpBBdeOGFPPbYY0ybNo0XX3xRcxCISEKoUAFWrizYvjndeG36dHjmmWC5YsXCyyUFp6JARKSQDRgwgD/+8Y+USdRbWIhIiZaaClOm5H6bz++/hzp14K9/zWwrU+anfZgvX/7oj5XCoaJAROQIpaen8/DDDzNw4MBcrwSoIBCRkmjrVmjdGjZvzn/fAwdg7Fgwi38uiT/91RIROQIHDhygf//+DB06lEsuuYSDBw+GHUlEpNAMHVqwggDggQdg7dr45pGioysFIiIFtGvXLnr27MmrkekvFy1axODBg5k6dSqmU2UiUgo8+ihs2QLvvw+33grVquW8X5ky0KEDNG5ctPkkflQUiIgUQGpqKsnJyaxevTqmfcGCBaSmplK3bt2QkomI5M0dtm+Hffti2w8eDMYHNG+e2VazJrz6Krz9dvChXxKHigIRkXx89NFHdOnShdTU1Jj2Ro0asXDhQhUEIlJspafDoEHBTME5qVQpuP1nu3aZbVWrqiBIRBpTICKSh8WLF9OmTZtsBcE555zDsmXLaNiwYUjJRETy5h6MEcitIIDg6kHnzvDRR0UWS4opXSkQEcnFs88+y8CBAzkQPc0m0L17d5599lmqVKkSUjIRkfyNHQsPP5y5Xrt2MJ9AtHLloHdvaNq0aLNJ8aOiQEQkC3dn/PjxjBgxItu2a665hoceeoiyOc2+IyJSTHz/Pfzzn5nrvXrBs88GA4RFcqKiQEQkyqFDhxg2bBgPPfRQtm133303I0eO1J2GRCRUL78cTC6W5SImzZvD6NHBclISvPUWXHghnHwyPPWUCgLJm4oCEZGIvXv30q9fP1566aWY9vLly/Ovf/2LP/zhDyElExEJLFwI3bpBWlr2benpsesNGgSFQeXKP222YUkMKgpERCK++uorlixZEtNWvXp1Zs2aRfv27UNKJSKJbuXK4OrAwYNw3305FwS5Oe64+OWS0kVFgYhIRIMGDZg3bx7nn38+e/fupU6dOixYsIDTTz897GgiksBWrIBRo2Lb6tULZhQuF/VJ7vjjizaXlC4qCkREopxzzjk899xzjB49mrlz52oOAhGJm//9D959N/PM/+bNkJoajAvIa+jSsccGVw5OO61ockpiUFEgIpJF9+7dSU5Oplw5/YoUkfj4+mto2xY++yz7tm3bgjsHZdzk7Kyz4NZbg+WKFaFPHzjllKLLKolBf/FEJCE9++yztG/fnuNy6XCrgkBECtuBA8FkYT/+CMnJORcEAFOnwjXXQLNmwXqLFsFDJJ70V09EEoq7c8cddzBmzBjOPvtsXn/9dU1CJiJFYsYM+OMfY9vKloV27TK7C1WqFMxCnFEQiBQVFQUikjAOHjzIVVddxdSpUwF455136N27N7NmzdKVAREJxeTJcPnlYacQAU1jISIJYffu3XTr1u1wQZBh0aJFrFq1KqRUIpJIypeHatWCx3HHwUMPqSCQ4kNFgYiUetu2beO8887j5Zdfjmk/9thjee211zjzzDNDSpa4zKyzmX1qZp+b2V9y2F7PzF43sw/MbJWZJYeRU6QwXXYZ/PBD8Pj6axgyJOxEIplUFIhIqbZmzRpatmzJBx98ENPeoEEDli1bRsuWLUNKlrjMrCzwT6AL0AToY2ZNsuw2Cpjh7mcAvYFHijalyE939dVwyy3w5ZdhJxHJn4oCESm1li5dSuvWrdm0aVNMe/PmzVm+fDm//OUvQ0qW8M4GPnf3De5+AHgOuCjLPg5UjywfA+hjlZQoX3wBU6bAPfdAgwawdWvYiUTypqJAREqlmTNn0qFDB3bu3BnTnpyczOLFizleU3+GqQ6wJWo9NdIW7Q7gD2aWCiwA/lw00UQKx9//DunpwfJ550GdrO9wkWJGRYGIlDoPPvggl156Kfv3749pHzx4MLNnz+ZnP/tZSMkkIqe5Wj3Leh/gCXc/CUgGnjazbH+zzOxKM1tpZiu/+eabOEQVOXLffQePPZa5Pnx4eFlECkr34BORUiM9PZ2bb76Z+++/P9u2O++8k9tuuw2znD6PShFLBepGrZ9E9u5Bg4DOAO6+3MwqATWB7dE7uftkYDJA8+bNsxYWIkVi82bYHvXOfO65YIIygNNPh44dw8klciRUFIhIqfHdd9/x4osvxrSVK1eOKVOmMHDgwHBCSU7eBU4xswbAVoKBxH2z7LMZaA88YWaNgUqALgVIsfTAA0F3oZwMH545MZlIcabuQyJSatSsWZOUlBSSkpIAqFq1KvPmzVNBUMy4exowBHgZWENwl6GPzWyMmXWP7HYTMNjM/gs8Cwx0d10JkBKlTp3gNqQiJYGuFIhIqdKkSRNmz57NgAEDmDVrFmeccUbYkSQH7r6AYABxdNvoqOVPgNZFnUvkaNStC2edFduWlAR//StUqBBOJpEjpaJAREqdtm3b8umnn1JBf41FpAjceGPwECnJ1H1IREqk1157jY0bN+a6XQWBiBS2hx+G1q2DMQK33grffht2IpHCo6JAREqcp59+mk6dOtGlSxe+++67sOOISAJYuxb+/GdYtixYv/tu+OUvM+8yJFLSqSgQkRLD3Rk7diz9+/cnLS2NtWvXctFFF7Fv376wo4lIKbdiRfa2a6+FKlWKPotIPKgoEJESIS0tjWuvvZbbbrstpn3p0qX85z//CSmViCSK//43c/nqq+HTT+HOO8PLI1LY4l4UmFlZM/vAzOZF1huY2QozW2dmz5uZOv6KSJ727NnD73//ex599NGY9goVKvD888/TtWvXkJKJSKKILgo6dQq6Dmn+ASlNiuJKwVCC+1BnGAc84O6nADsJZq0UEcnR9u3bueCCC5g7d25Me1JSEosWLaJXr14hJRORROEeWxScfnp4WUTiJa5FgZmdBPwOeCyybsAFQMaUo08CPeKZQURKrs8//5xWrVrxzjvvxLTXq1ePt956i7Zt24aUTEQSybZtsGNHsFytGtSvH2ockbiI95WCvwMjgPTI+rHA95HZLAFSgTpxziAiJdCKFSto2bIl69evj2lv1qwZy5cvp0mTJiElE5FEs2pV5vKvfw1lNCJTSqG4va3NrCuw3d3fi27OYdccp603syvNbKWZrfzmm2/iklFEiqc5c+Zw/vnn822Wm4B37NiRN998kxNPPDGkZCKSiNR1SBJBPGc0bg10N7NkoBJQneDKQZKZlYtcLTgJ+DKng919MjAZoHnz5jkWDiJS+kycOJEhQ4aQnp4e096/f38ee+wxypcvH1IyEUlU3bpB5cpBcdCxY9hpROIjbkWBu98C3AJgZu2A4e7ez8xeAC4BngMGALPjlUFESpZ9+/bx8MMPZysIRo0axZgxYzDd6kNEQtCkSfAQKc3C6BU3ErjRzD4nGGPweAgZRKQYqlSpEikpKYe7B5UpU4ZJkyZx5513qiAQERGJo3h2HzrM3RcDiyPLG4Czi+J1RaTkqVevHgsWLKBz585MmTJFcxCISKgOHIAKmlFJEkCRFAUiIkeiWbNmbNiwgcqVK4cdRUQSwOuvwxtvBPMRRFu6FLZvh0WLoHbtcLKJFBUVBSISio8++ggzo2nTpjluV0EgIkXh2Wehb9+892ndGl57DU4+uWgyiYRBd9oVkSK3ePFi2rRpQ5cuXdi6dWvYcUQkQb3yCgwYkP9+NWoED5HSTFcKRKRIPffccwwYMIADBw6wa9cukpOTWbJkCdWrVw87mogkkC1boGdPOHgwWG/cGC67LPt+xxwDl18ezGQsUpqpKBCRIuHuTJgwgZtvvjmmfdWqVTz99NNcd911ISUTkURUrRrccgv84x/BQOJXXoGTTgo7lUh41H1IROLu0KFDDB06NFtBAHDPPfdw7bXXhpBKRBJZUhLceits3AhLlqggENGVAhGJq71799KvXz9eeumlmPby5cszdepU+vXrF1IyERGoVEkDiEVARYGIxNGOHTvo3r07y5Yti2mvXr06s2bNon379iElE5FE5Q6aC1EkO3UfEpG4+OKLL2jVqlW2gqBOnTosWbJEBYGIFLm0tGAw8VNPhZ1EpPhRUSAihe69996jRYsWfPbZZzHtTZs2Zfny5Zx++ukhJRORROUO114LL7wQ3IZ0/PiwE4kULyoKRKRQpaSkcN5557F9+/aY9nbt2rF06VLq1q0bUjIRSWSjR8OUKZnr27aFl0WkONKYAhEpNO7O3XffzZ49e2La+/Tpw9SpU6lYsWJIyUQkUfz3v/D667FtGzfCgw9mrv/hD3DffUUaS6TYU1EgIoXGzJg1axatWrVi3bp1AIwYMYK7776bMmV0YVJE4u+tt2DYsNy3d+kC//oX6FeSSCz9lxCRQlWzZk1SUlKoXbs2Dz30EOPGjVNBICLFQosWwZiC8uXDTiJS/OhKgYgUuoYNG/LZZ59RrVq1sKOISII5/XQYOjR7e82a8Oc/w89+VvSZREoCFQUiclS++uortm7dyllnnZXjdhUEIhKGNm2Ch4gcGV3TF5EjtnbtWlq0aEGnTp2y3XZURKQobdwI06fDwYNhJxEp2VQUiMgRWbp0Ka1atWLTpk3s2LGDLl26ZLv9qIhIUfjuO/jd76BfP2jUCGbMCDuRSMmlokBECmzmzJl06NCBnTt3Hm7bsGED4zULkIgUsR9/hG7d4JNPgvVt26B27XAziZRkKgpEpEAefPBBLr30Uvbv3x/TfsUVVzB27NiQUolIIjp4EHr1gmXLMtuefBLatg0vk0hJp6JARPKUnp7O8OHDueGGG3D3mG1jxoxh8uTJlCunexaISNFwhyuvhPnzM9sefBD69Akvk0hpoL/kIpKrffv2MWDAAGZk6ahbrlw5pkyZwsCBA8MJJiIJa/58eOKJzPVbboHrrw8tjkipoaJARHK0c+dOevTowZtvvhnTXrVqVV588UU6deoUUjIRSWQvvJC53Lcv3HVXeFlEShMVBSKSzebNm+nSpQufZIzgi6hduzbz58/nzDPPDCmZiCSyQ4diuw3dcAOYhZdHpDRRUSAiMT788EOSk5PZtm1bTPtpp51GSkoK9evXDyeYiCS85cthx45g+YQTIJe5E0XkKGigsYjEGDt2bLaCoE2bNrz11lsqCEQkVCkpmcvdukEZfYoRKTT67yQiMR5//HHOOOOMw+s9e/Zk0aJF1KhRI8RUIiJwxx3w2mswbJjuNiRS2NR9SERiVKtWjfnz59OyZUsuvvhiJkyYQBmdjhORYqB8eTj//OAhIoVLRYGIZHPCCSfw/vvv6+qAiIhIgtDpP5EEtWfPHt55551ct6sgEJHiYN8+mD497BQipZ+KApEEtH37di644ALat2/P+++/H3YcEZEcpaUFYwf69YPhwyE9PexEIqWXigKRBLNu3TpatWrFO++8w//+9z9+97vfsXHjxrBjiYjEcIdrr4V//ztYnzABFi4MN5NIaaaiQCSBrFixglatWrF+/frDbV999RUjR44MMZWISHYvvABTpmSuDx8Oycnh5REp7VQUiCSIOXPmcP755/Ptt9/GtHfs2JEp0X95RURC5g733JO53q8fjBsXXh6RRKCiQCQBTJw4kYsvvpi9e/fGtA8YMID58+dTvXr1kJKJiGT3+uvwwQfBcqVK8Pe/a6IykXjTfzGRUszdueWWW7j22mtJzzJCb9SoUUydOpXy5cuHlE5EJGf33Ze5/Kc/Qc2a4WURSRSap0CklDpw4ACDBg3imWeeiWkvU+z+kR4AACAASURBVKYMEydO5MorrwwpmYhI7j7+OHNAsRnceGO4eUQShYoCkVJo165d9OzZk1dffTWmvUqVKjz//PN07do1pGQiInmbOTNzuUcPaNQovCwiiURFgUgps3XrVpKTk1m1alVMe61atZg/fz6//e1vQ0omIpK/Zcsyly++OLwcIolGRYFIKXPvvfdmKwgaNWrEwoULadiwYUipREQKJvoGaa1bh5dDJNFooLFIKXPvvfdy3nnnHV5v0aIFy5YtU0EgIiXCypXw1VcwZw40aBB2GpHEoSsFIqVMxYoV+fe//02bNm1o1KgR06dPp0qVKmHHEhEpsOOPh27dwk4hklh0pUCkFEpKSuL1119n5syZKgikWDKzzmb2qZl9bmZ/yWWfXmb2iZl9bGbTizqjiEgiUVEgUkIdOnSIlStX5rq9Vq1alC1btggTiRSMmZUF/gl0AZoAfcysSZZ9TgFuAVq7e1PghiIPKiKSQFQUiJRAe/fupVevXrRq1YrXX3897DgiR+ps4HN33+DuB4DngIuy7DMY+Ke77wRw9+1FnFGK0JYtMGsWXHIJLF8O+/eHnUgk8cStKDCzSmb2jpn9N3Lp96+R9gZmtsLM1pnZ82ZWIV4ZREqjHTt20KFDB2bNmsXBgwe5+OKLWb16ddixRI5EHWBL1HpqpC3aL4FfmtlbZva2mXUusnRSpGbPhoYNoWfPYI6CVq3gwgvDTiWSeOJ5pWA/cIG7NwN+A3Q2sxbAOOABdz8F2AkMimMGkVLliy++oFWrViyLupH3rl27GDx4MO4eYjKRI2I5tGV9A5cDTgHaAX2Ax8wsKdsTmV1pZivNbOU333xT6EElvt58Ey67DA4ejG1v2zacPCKJLG5FgQf+F1ktH3k4cAHwYqT9SaBHvDKIlCYrV66kRYsWfPbZZzHtTZs25YUXXsAsp89ZIsVSKlA3av0k4Msc9pnt7gfd/QvgU4IiIYa7T3b35u7evFatWnELLIXv88+he/fMrkK1awfFwODBMHx4uNlEElFcb0kaGUz2HtCIYFDZeuB7d0+L7JLTJeOMY68ErgSoV69ePGOKFHspKSlceuml7NmzJ6a9Xbt2vPTSSyQlZTuBKlKcvQucYmYNgK1Ab6Bvln3+TXCF4Akzq0nQnWhDkaaUuLrzTti1K1iuXRveegt+8YtwM4kksrgONHb3Q+7+G4KzQGcDjXPaLZdjdfZHBHj88cfp1q1btoKgd+/eLFy4UAWBlDiRE0NDgJeBNcAMd//YzMaYWffIbi8DO8zsE+B14GZ33xFOYilsqakwPeomsy+8oIJAJGxFMnmZu39vZouBFkCSmZWL/FHI6ZKxiADuzh133MGYMWOybRsxYgR33303ZcroBmJSMrn7AmBBlrbRUcsO3Bh5SCnzj39AWqTPwLnnQps24eYRkfjefahWxqAwM6sMdCA4I/Q6cElktwHA7HhlECmpDh48yKBBg7IVBGbGQw89xLhx41QQiEiJ9MMPMGlS5rrGD4gUD/G8UnAC8GRkXEEZgsvD8yKXgp8zs78BHwCPxzGDSImze/duLr30Ul5++eWY9kqVKjF9+nQuvvjikJKJiPx0P/sZPP443HdfMKaga9ewE4kIxLEocPdVwBk5tG8gGF8gIjmYOHFitoKgRo0azJ07l1atWoWUSkSkcJQtG0xS1rMn7NgBuugpUjzov6JIMXPTTTfRo0fmnXobNGjAsmXLVBCISKliBjVrhp1CRDIUyUBjESm4smXLMn36dNq3b8+BAweYN28etWvXDjuWiIiIlGIqCkSKocqVKzN37lwqVqxI1apVw44jIvKT7NsHycnQsWMwYVmTJsGVAhEpPtR9SCREq1atynXbscceq4JAREqFxYvh9dfh1lshqnekiBQjKgpEQpCens5NN93Eb37zG1566aWw44iIxNWcOZnL3bvrKoFIcaSiQKSI7du3jz59+nD//ffj7vTt25dly5aFHUtEJC7cYe7czPVu3cLLIiK5U1EgUoR27txJp06dmDFjxuG2ffv20bdvXw4cOBBiMhGR+PjwQ0hNDZZ//nNo3TrcPCKSMw00FikimzdvpnPnzqxZsyamvXbt2syaNYsKFSqElEzk6JlZGaAZcCKwF/jY3b8ON5UUFz/+CEOGZK4nJ0P58uHlEZHcqSgQKQIffvghycnJbNu2Lab91FNPZeHChdSvXz+cYCJHycwaAiOBDsA64BugEvBLM/sRmAQ86e7p4aWUMB08CL16QXTvyCuuCC+PiORNRYFInC1atIiePXuye/fumPbWrVszZ84catSoEVIykZ/kb8BE4Cp39+gNZnYc0Bf4I/BkCNkkZO4weDDMn5/Z9uCD0K5daJFEJB8qCkTi6KmnnmLQoEGkpaXFtPfs2ZNnnnmGSpUqhZRM5Kdx9z55bNsO/L0I40gxc8st8OSTsevXXx9eHhHJnwYai8SBuzN27FgGDBiQrSAYOnQozz//vAoCKVXMrJGZPWNmM82sZdh5JFxnnpk5dmDQILjrrnDziEj+dKVApJClpaUxZMgQJk2alG3b/fffz7Bhw0JIJVK4zKySu++LaroTuB1w4AXgN6EEk1Ds2QNpaXDMMcF6r15w7LEwdSo8+qjmJRApCXSlQKSQvfDCC9kKggoVKvD888+rIJDSZK6Z/TFq/SBQP/I4FEYgCUdKChx3HGQ9D9K+PTzzDJTT6UeREkFFgUgh6927N4MGDTq8npSUxKJFi+jVq1eIqUQKXWfgGDNbaGbnAsOBtkAXoF+oyaTILFsGPXsGtx4VkZJN9btIITMzJk6cyJdffsnHH39MSkoKTZo0CTuWSKFy90PAw2b2NDAaOAH4P3dfH24yKSrr1kHXrrB3b7C+enVw1yF1FRIpmVQUiMRB+fLlmTFjBj/88AMnnnhi2HFECp2ZnQPcDBwAxhJMXHaXmaUCd7r7rjDzSfzddhvs3Bks16oFo0erIBApydR9SOQn+Oyzz3LdVrVqVRUEUpo9SjB52Thgkruvd/fewFxgRqjJJO42bICZMzPXZ8+GU04JL4+I/HR5Xikwsxvz2u7u9xduHJGSY+LEiQwZMoTJkyfHjCEQSRCHCAYVVyG4WgCAu78BvBFSJikiDzwA6ZG5qi+8EFrqJrQiJV5+VwqqRR7NgWuAOpHH1YA6SUtCSk9P55ZbbuHaa68lPT2dq666ipSUlLBjiRS1vkAy0AroH3IWKUI7dsC//pW5Pnx4eFlEpPDkeaXA3f8KYGavAGe6++7I+h0E96EWSSgHDhzg8ssvZ9q0aYfbDh06RO/evdm4cSM///nPQ0wnUqTWuftNee1gZubuXlSBpGg8/njm3YaaNYMOHcLNIyKFo6BjCuoRdXk4sly/0NOIFGO7du0iOTk5piAAqFKlCtOmTVNBIInmdTP7s5nVi240swpmdoGZPQkMCCmbxFH0WIKhQzW4WKS0KOjdh54G3jGzlwhmq7wYeCpuqUSKmdTUVJKTk1m9enVMe61atZg3bx5nn312SMlEQtMZuBx41swaAN8DlYCywCvAA+7+YYj5JA62bYN33gmWy5aFHj3CzSMihadARYG732VmKcC5kaY/ufsH8YslUnx89NFHdOnShdTU1Jj2Ro0asXDhQho2bBhSMpHwuPs+4BHgETMrD9QE9rr79+Emk3iqWBHGj4c5c6BCBdAFUpHS40jmKagC/ODuU82slpk1cPcv4hVMpDhYvHgxPXr0YNeu2Fuun3POOcydO5datWqFlEyk+HD3g8C2sHNI/NWoATfdFDwOHQo7jYgUpgKNKTCz2wnuR31LpKk88Ey8QokUB88++yydOnXKVhB0796d1157TQWBiCS0smXDTiAihamgA40vBroDewDc/UuCW5WKlDruzn333Uffvn05cOBAzLZrrrmGWbNmUaVKlZDSiYgUja1bITkZTjst8/Gf/4SdSkTipaDdhw64u5uZA5jZz+KYSSRUS5cuZcSIEdna7777bkaOHInpVhsih5nZEGCau+8MO4sUrv/7P8g6Bcv//hdOFhGJv4JeKZhhZpOAJDMbDPwHeCx+sUTCc+655zJy5MjD6+XLl+fpp5/mL3/5iwoCkexqA++a2Qwz62z6T1IqHDoEc+eGnUJEilJB7z403sw6Aj8ApwKj3X1RXJOJhGjs2LFs2bKFefPmMWvWLNq3bx92JJFiyd1Hmdn/ARcCfwIeNrMZwOPuvj7cdHK03n4bvv02WD7hBHj11WA+gjp1ws0lIvFToKLAzMa5+0hgUQ5tIqVOmTJl+Ne//sXGjRs59dRTw44jUqxFupd+BXwFpAE/B140s0Xunr0vnhR7c+ZkLnftCo0bh5dFRIpGQbsPdcyhrUthBhEJw6ZNm3D3HLdVrFhRBYFIPszsejN7D7gXeAv4tbtfA5wF9Aw1nBy16KKge/fwcohI0cmzKDCza8xsNXCama2KenwBrM7rWJHibsGCBTRt2pT77rsv7CgiJVlN4Pfu3sndX4jMWYC7pwNdw40mR2PdOli7NliuXBnUe1IkMeTXfWg6kALcDfwlqn23u38Xt1Qicfb4449z1VVXcejQIUaOHMlJJ51E3759w44lUhItAA7/PTCzakATd1/h7mvCiyVHK/oqQceOQWEgIqVfnlcK3H2Xu28EHgS+c/dN7r4JOGhm5xRFQJHC5O7cfvvtXHHFFRyKmo5z4MCBrF+vMZEiR2EiEH2jyj2RNimhqleHjLkZ1XVIJHEUdEyBfulLiXfw4EEGDRrEmDFjYtrNjPvvv5+GDRuGlEykRDOPGpgT6TZU0DlwpBgaPBg2bYLJk6FXr7DTiEhRKegv7my/9M1Mv/SlxNi9ezeXXnopL7/8ckx7pUqVmD59OhdffHFIyURKvA1mdj2ZJ4quBTaEmEcKQeXKQXEgIomjoFcKNkTuMFE+8hiKfulLCbFt2zbOO++8bAVBjRo1ePXVV1UQiPw0VwOtgK1AKnAOcGWoiURE5IgV9Gz/1cA/gFGAA6+iX/pSAqxZs4YuXbqwadOmmPYGDRqQkpKiW46K/ETuvh3oHXYO+em2btXkZCKJrKAzGuuXvpQ4S5cupXv37uzcuTOm/ayzzmL+/Pkcf/zxISUTKT3MrBIwCGgKVMpod/fLQwslR2zdumCCsi5d4OaboW3bsBOJSFHLsygwsxHufq+ZPURwhSCGu18ft2QiP8HMmTPp168f+/fvj2lPTk7m+eefp2rVqiElEyl1ngbWAp2AMUA/QLciLWEeeAAOHYJ58yA9XUWBSCLK70pBxi/2lfEOIlJYPv/8c3r16kV6enpM++DBg3nkkUcoV05j5EUKUSN3v9TMLnL3J81sOvByvkdJ6NLTYfly+PxzmDo1s3348PAyiUh48vx05O5zI1+fLJo4Ij9do0aNGDduHDfffPPhtjFjxjBq1CjMLMRkIqXSwcjX783sV8BXQP3w4khBPfooXHddbNuZZ0K7dqHEEZGQ5dd9aC45dBvK4O6a1kSKpZtuuoktW7bwyCOPMGXKFAYOHBh2JJHSarKZ/ZzgRhRzgKrA/4UbSQrimWeyt40cCTp3IpKY8utHMT7y9fdAbSDjV0gfYGOcMon8ZBkTkvXv35+zzjor7DgipZKZlQF+cPedwJvAL0KOJAW0bx+8917mep8+cO65cOml4WUSkXDl133oDQAzu9Pdo4cdzTWzN/M61szqAk8RFBPpwGR3f9DMagDPE1xe3gj0ivxBETliX3/9NbVq1aJMmexTbpQtW1YFgUgcRSayHALMCDuLHJn334cDB4LlU06B6dPDzSMi4Svo5GW1zOzwGSAzawDUyueYNOAmd28MtACuM7MmwF+AV939FIL5Dv5y5LFF4MMPP+SMM85g5MiRYUcRSWSLzGy4mdU1sxoZj7BDSd6WLctcbtUqvBwiUnwU9DYsw4DFZpYxi3F94Kq8DnD3bcC2yPJuM1sD1AEuAtpFdnsSWAzoU50ckUWLFtGzZ092797N+PHjqVu3LtdfrzvkioQgYz6C6CGrjroSFWtly0K9erB5M7RuHXYaESkOzD3XccSxO5pVBE6LrK519/157Z/l2PoE/U1/BWx296SobTvd/ec5HHMlkVmT69Wrd1bWGWklcT311FMMGjSItLS0w21mxpIlS2itv26SoMzsPXdvHnaOMDVv3txXrtQdtI9EaipUrQpJSfnvKyIlX15/Kwp0pcDMqgA3Aie7+2AzO8XMTnX3eQU4tiowE7jB3X8o6C0h3X0yMBmCX/QFOkhKNXdn7NixjBo1Ktu2oUOH0rJlyxBSiSQ2M+ufU7u7P1XUWeTInXRS2AlEpLgoaPehqcB7QManrlTgBSDPosDMyhMUBNPcfVak+WszO8Hdt5nZCcD2I48tiSYtLY0hQ4YwadKkbNvuv/9+hg0bFkIqEQF+G7VcCWgPvE9wowkRESkhCloUNHT3y8ysD4C777V8TvlHtj8OrHH3+6M2zQEGAPdEvs4+8tiSSPbs2UPv3r2ZNy+2Bq1QoQJPP/00vXr1CimZiLj7n6PXzewY4OmQ4oiIyFEqaFFwwMwqE5nIzMwaAvmNKWgN/BFYbWYfRtpuJSgGZpjZIGAzoLsiS662b99O165deffdd2Pak5KSmD17Nm3bts3lSBEJyY/AKWGHkNzdfDPUrRvcdeiMM4JBxyIiBS0KbgcWAnXNbBrBB/6BeR3g7kuB3K4mtC9oQElc69ato0uXLqxfvz6mvV69eqSkpNCkSZOQkolIBoud+b4M0IQCzFtgZp2BB4GywGPufk8u+11C0F31t+6uUcQ/0Q8/wIQJ4A5lysD330O1amGnEpHiIN+iININaC3BrMYtCD7oD3X3b+OcTRLY22+/Tbdu3fj229i3WbNmzViwYAEnnnhiSMlEJIvxUctpwCZ3T83rADMrC/wT6EgwRu1dM5vj7p9k2a8acD2wonAjJ64VK4KCAKBZMxUEIpIp36LA3d3M/u3uZwHziyCTJLjvvvuOTp068cMPP8S0d+zYkRdffJHq1auHlExEcrAZ2Obu+wDMrLKZ1Xf3jXkcczbwubtviBzzHMEcNp9k2e9O4F5geKGnTlBvvZW5rEnLRCRaQWc0ftvMfpv/biI/XY0aNXjwwQdj2vr378/8+fNVEIgUPy8A6VHrhyJteakDbIlaT420HWZmZwB1C3Lrayk4zWQsIrkpaFFwPkFhsN7MVpnZajNbFc9gktgGDhzImDFjABg1ahRPPPEE5cuXDzmViOSgnLsfyFiJLFfI55icxpsdno/GzMoADwA35ffiZnalma00s5XffPNNASMnpkOH4O23M9dVFIhItIIONO4S1xQiORg1ahRt27blvPPOCzuKiOTuGzPr7u5zAMzsIiC/MWepQN2o9ZOAL6PWqwG/AhZH7n5dG5gTeZ2Ywcaa6LLgPv4Ydu8Olk84AU4+Odw8IlK85FkUmFkl4GqgEbAaeNzd04oimCSGH374gcqVK+d4FcDMVBCIFH9XA9PM7OHIeiqQ4yzHUd4FTjGzBsBWoDfQN2Oju+8Camasm9liYLjuPvTTRHcdat0a8p5tSEQSTX7dh54EmhMUBF2ACXFPJAlj69atnHvuuQwePBh3neATKYncfb27tyC4FWlTd2/l7p/nc0waMAR4GVgDzHD3j81sjJl1j3/qxLRkSeayug6JSFb5dR9q4u6/BjCzx4F34h9JEsFHH31Ely5dSE1NZdWqVdSrV+/wGAIRKTnMbCxwr7t/H1n/OXCTu4/K6zh3XwAsyNI2Opd92xVO2sSVlgYpKZnruggrIlnld6XgYMaCug1JYVm8eDFt2rQhNTXzVuZ33nknM2fODDGViBylLhkFAYC77wSSQ8wjOVi+HHbuDJbr1AlmMhYRiZbflYJmZpZxs3gDKkfWjWAKA90fUo7Ic889x4ABAzhw4EBM+0UXXUSXLhrPLlIClTWziu6+H4J5CoCKIWeSLFq1CsYUzJkTTFim8QQiklWeRYG7ly2qIFK6uTsTJkzg5ptvzrbtmmuu4aGHHqJsWb3dREqgZ4BXzWwqwW1FLweeCjeSZFW2LLRsGTxERHJS0FuSihy1Q4cOMWzYMB566KFs2+655x5GjBiB6bSVSInk7vdG5q3pQHAV+U53fznkWCIicoRUFJQGa6bBkttg92aoVg/OvQsa9ws7FQB79+7lD3/4A7NmzYppL1++PFOnTqVfv+KRU0SOnrsvBBYCmFlrM/unu18XciyJOHgQNPejiOSnoDMaS3G1Zhq8ciXs3gR48PWVK4P2kO3YsYMOHTpkKwiqV69OSkqKCgKRUsLMfmNm48xsI/A3YG3IkSTitdegceNg4jIRkbyoKCjpltwGaT/GtqX9GLSH6IsvvqBVq1Ysi54tB6hTpw5Lliyhffv2ISUTkcJgZr80s9FmtgZ4mGDSMnP38909e19BKXIffQQ9esD69XDuubBiRdiJRKQ4U/ehkm735iNrLwIHDhygQ4cObNiwIaa9adOmpKSkULdu3ZCSiUghWgssAbplTFZmZsPCjSTR7rwTdu8OlqtUgdq1w80jIsVb6b1SsGYaTK4PE8oEX4tBd5q4qFbvyNqLQIUKFXjggQcoUybz7dWuXTuWLl2qgkCk9OgJfAW8bmZTzKw9wUBjKSbeiZpu9Nln4eSTw8siIsVf6SwKinE/+0J37l1QrkpsW7kqQXuIunfvzj//+U8A+vTpw8KFC0lKSgo1k4gUHnd/yd0vA04DFgPDgOPNbKKZXRhqOGHXLti4MVguVw7OOSfUOCJSApTOoqCY9rOPi8b94MLJUO1kwIKvF04uFncfuvrqq0lJSeGZZ56hYkXNZSRSGrn7Hnef5u5dgZOAD4G/hBwr4a1enbncuDFUqBBeFhEpGUrnmIKj6WdfjG/rma/G/ULLevDgQdLS0qhcuXKO2zt37lzEiUQkLO7+HTAp8pAQ/fe/mcvNmoWXQ0RKjtJ5peBI+9knUnejQrR79266devGZZddRlpaWthxREQkYtWqzGUVBSJSEKWzKDjSfvaJ1N2okHz11Ve0a9eOl19+mblz53L99dfj7mHHEhERYq8UnH56eDlEpOQonUXBkfazL4a39SzO1q5dS8uWLXn//fcPt02cOJFHHnkkxFQiIgJw6FDsmAJdKRCRgiidYwrgyPrZV6sX6TqUQ7vEWLp0Kd27d2fnzp0x7WeddRaXXHJJSKlERCTDpk3wY+Ti93HHwfHHh5tHREqG0nml4EgV09t6FjczZ86kQ4cO2QqC5ORkFi9ezPH6yyMiErpf/AK+/x6WLIHJk8NOIyIlRem9UnAkMq4olNS7DxWBf/zjH9xwww3Zxg1cccUVTJw4kXLl9FYSESkujjkG2rQJO4WIlCT6JJchxNt6Fmfp6emMGDGCCRMmZNs2ZswYRo0ahZkmMRUREREpyVQUSK7279/PgAEDeP7552Pay5Urx5QpUxg4cGA4wURERESkUKkokBzt3LmTiy++mDfeeCOmvWrVqrz44ot06tQppGQiIpKblSth3Tro0gWSksJOIyIliQYaSzbuTvfu3bMVBLVr1+aNN95QQSAiUkxNmgR9+0KtWhpkLCJHpvQWBWumweT6MKFM8FWzExeYmXHXXXdRoUKFw22nnXYay5cv58wzzwwxmYiI5CY9HebODZbT0uBXvwo3j4iULKWzKFgzDV65MjL3gAdfX7lShcERaNu2LU8//TQAbdq04a233qJ+/frhhhIRkVy9+y58/XWwXKsWnHNOuHlEpGQpnWMKltwGaT/GtqX9GLTrDkMF1qtXLypXrkzHjh2pVKlS2HFERCQPGVcJALp2hbJlw8siIiVP6bxSsHvzkbUnMHdn7969uW7v1q2bCoJSYMuWLbRv357GjRvTtGlTRowYkW3OCREp2ebMyVzu1i28HCJSMpXOoqBavSNrT1BpaWlcc801dO7cmX379oUdR+KoXLlyjBs3jjVr1vDBBx+wYsUKZs2aFXYsESkkGzfC6tXBcsWK0LFjqHFEpAQqnUXBuXdBuSqxbeWqBO0CwJ49e/j973/PpEmTePPNNxkwYADp6elhx5KjMGzYMMzs8OP444+nf//+7Nix4/A+J5xwAs2bNwegQoUKnH766WzZsiWsyAX2ySef0L59e6pUqcKJJ57I6NGjOXToUKEd99xzz3HmmWdStWpV6tSpQ//+/fnyyy8Pb3/xxRdp1aoVxx57LJUqVeLUU0/lb3/7GwcOHDji5xKJp+iuQ+3bQ9Wq4WURkZKpdBYFjfvBhZOh2smABV8vnKzxBBHbt2/nggsuYG7UX5EZM2YwevToEFPJ0Vq9ejUtWrRg+fLlvPXWW9x0001MmzaN6667Lsf9d+zYwb///e9if2vZnTt30qFDB8yM2bNnM3r0aCZMmMDtt99eKMfNmTOHPn360KpVK2bPns24ceN488036dq16+ECeceOHZx//vk89thjpKSkcPnll3PXXXdx4403HvFzicSTug6JyE/m7sX+cdZZZ7kUjnXr1nnDhg0diHnUq1fPP/7447DjyVE47rjj/Lrrrotp69atm9eoUSPbvvv27fN27dr5+PHjiyreURs7dqwnJSX5rl27DreNGzfOK1euHNN2tMdddtllfuaZZ8YcO3v2bAf8k08+yfX5b731Vj/mmGM8PT39Jz9XYQJWejH4fR3mI1H/Vnz/vXu5cu4QPLZsCTuRiBRXef2tKJ1XCiRHK1asoGXLlqxfvz6mvVmzZixfvpwmTZqElEyO1vbt29m+fTuNG/9/e3ceH1V59n/8c5EQIOwCVhQEVATUVkCqoOIumyzFFUotdcPqQ4vKo62FonVpUatPBevCU0VaeAIu5QdSEBGwoqBCtYIsCiqbWkFBREEgyfX7YybJTDKTzCSZzGTm+3695pUz95xz5jpJZs59nXMvXcLKDz/8cLJKDT1SUFDAiBEj6NatSL0AowAAIABJREFUG2PHjq3JMCtlwYIF9O3blyZNmhSXDRs2jP3795eZWK8y2x06dIimTZuGbdssOAWsl9MJu0WLFmWaD1V2XyLVYeHCwLwEAKecAm3aJDceEamdlBQUSfPJzubOncu5557LF198EVZ+4YUX8uqrr3LkkUcmKTKpitWrVwOByeWKFBYWsmLFCgaVakNw/fXX07hxYx588MFy9+nu5OfnV/hItA0bNoQdF8DRRx9Nbm4uGzZsqPJ2V199NcuWLeOvf/0rX3/9NR988AHjx4/n3HPPLZMgFxQUsG/fPl577TUmTZrEDTfcgJlVal8i1a1vX8jLg+HDA7MZi4hUhpICSPvJzh577DGGDh1aZujRkSNH8o9//CPsiqrULkVJwXHHHUd+fj5bt27l+uuvp0mTJtx///3F673++us8+eSTrFq1im7dutG1a1cmTZoUcZ/Tpk2jbt26FT4Sbffu3cVX20M1b96c3bt3V3m7iy66iKeffppRo0bRtGlTOnXqREFBQcRRmRo2bEjDhg3p3bs3Z599Ng888EDY6/HsS6S6NW0Kw4bB//0flOruIiISs/ScvCxeaTrZmbvzm9/8hokTJ5Z5bfz48dx1111hVzul9lkTHIMwdLbpNm3asHLlSlq0aFFcdsYZZ8TcjGXQoEGsXLmyWuLbs2cPn332WYXrlb6yXyTS/6e7V/h/G8t2S5cu5ec//zljxoyhf//+fP7559x5550MHTqUl19+Oaz51fLly9m3bx9vvfUWd911F6NHj+bRRx+t1L5ERERSkZICSMvJzg4ePMg111zD9OnTw8rr1KnDY489xqhRo5IUmVSnopGHJk+ezKFDh1i0aBF33HEH99xzD4888kil9nnYYYeVaR9fWc8++yzXXXddhetFSliaN2/OV199VaZ8z549Ee8ExLvd2LFjGTx4MPfdd19xWdeuXencuTNz5szh4osvLi7v3r07AGeeeSYtW7Zk5MiRjB07lmOPPTbufYmIiKQiNR+CtJzsLFJCkJuby5w5c5QQpInCwkLWrVvH6aefTo8ePejVqxcTJkygX79+PPvss5UeCrM6mw9de+21MY0aE0nnzp3L9B3Ytm0b3377bdQ7C/Fst2HDBrp27Rq2XqdOnWjQoEGZzvihihKEjz/+uMr7EqkKd1i+vKSTsYhIVSQsKTCzp8xsh5m9F1J2mJktMrONwZ/NE/X+cUnDyc7Gjh1Lo5DZa1q1asUrr7zCwIEDkxiVVKeNGzeyf/9+unXrFlb+4x//mB07dvDWW29Var9FzYcqeiRa//79WbhwIXv37i0umzVrFg0aNODss8+u8nbt2rXj7bffDtt2/fr17N+/P6w5Vmmvv/46AB06dKjyvkSq4v334Ywz4PDD4ec/T3Y0IlLrxXIVrzIP4CygO/BeSNn9wK+Dy78G7otlXzUy9vS66e5PtHP/owV+rpue+PdMsJdeesmzs7P9uOOO802bNiU7HKlmzz77rAO+Zs2asPKdO3d6nTp1/I477khOYNVk165dfsQRR/gFF1zgixYt8ieeeMIbNmzo48aNC1tv2rRpnpWV5Zs3b45ruz/96U9uZn7LLbf4okWLfPr06X788cd7+/bt/ZtvvnF39759+/oDDzzg8+fP94ULF/qECRO8YcOGfsUVV8S9r0RD8xRk3DwF999fMjfB0KHJjkZEaoPyzhUJ/YIG2pdKCt4HWgeXWwPvx7KfSn3Rp2ElvzJeeOEF37FjR7LDkASYMGGC169f3w8dOlTmtV69evlpp52WhKiq19q1a/3cc8/1+vXr+xFHHOHjx4/3/Pz8sHWmTp3qgH/88cdxbVdYWOiPPvqof//73/fc3Fw/8sgj/fLLL/cPP/yweJ3x48f7iSee6A0bNvSmTZt6t27dfNKkSX7w4MG495VoSgoyLyk499ySpGDq1GRHIyK1QXnnCgu8nhhm1h6Y5+4nBZ9/5e7NQl7f7e4VNiHq0aOHr1q1KvY3LhpiNHREoexc6DOl/NGE1s8IjDi0d2ugP0HvexM7+lA1vd/BgwfJyclJQIAiUluY2b/cvUey40imuM8VtdihQ9CsGewLnua2bdOkZSJSsfLOFSnb0djMRpnZKjNbtXPnzvg2Lm+I0Whqeq6Cang/d+ePf/wjPXv25Ouvv05MnCIiknJWry5JCI4+WgmBiFRdTScFn5tZa4Dgzx3RVnT3Ke7ew917tGrVKr53qcwQo5VJJKqiiu9XUFDAmDFjuPXWW3nnnXe45JJLOHjwYAICFRGRVLN8ecny6acnLw4RSR81nRTMBUYGl0cCcxLyLpUZYrSm5yqowvvt37+fyy67jMmTJxeXvfzyy9x4443VFZ2IiKQwJQUiUt0SOSRpHrAC6GRm283sGmAicKGZbQQuDD6vfpUZYrTeYfGVV1Ul50b48ssvOf/885k9e3ZYeZMmTRg+fHh1RSciIilMSYGIVLeEzWjs7tFqqOcn6j2LdRkBn7wOq6eAF4BlwYkjy+/Ea3GWV1XveyN3hi4ncfn444/p168fH3zwQVj5UUcdxfz58/nBD36QoGBFRCRVbN8OW4M3lXNz4eSTkxuPiKSHlO1oXCXrZ8DaaYGEAAI/104rvxPvd7viK6+qLiMCoyE1bgdY4Gc5oyOtWrWKnj17lkkITjzxRFasWKGEQEQkQ4TeJTjtNMhO2OU9Eckk6flVUl4n3mh3CxofHRwJKEJ5onQZEdMQpAsWLOCyyy7j22+/DSs/55xzmD17Ns2aNYuypYiIpJvcXDjrLHjrLTUdEpHqk553CiJV7ssrh8j9EOrkwMFv4ME6MKV94oYnLceTTz7JoEGDyiQEw4cP58UXX1RCICKSYQYOhH/+E/bsgdtuS3Y0IpIu0jMpsKwoLwQr95Eq+aWb89RvEZgo8sCX1Mi8BaW4O3fccQfXXnstBQUFYa/ddtttTJ8+nXr16lW8o/Uzoh+ziIjUWjk50KRJsqMQkXSRnkmBF0R5obD8ycK6jIBRm2FsIdRtBH4ofPNEzltQyoQJE7jrrrvCysyMyZMnc99991GnTgx/upqekE1EREREaqX0TAoat4ttvfIq+VHnEdgS+5X3KlylHzlyJC1btix+Xr9+fZ5//nlGjx4d8z5qfEI2ERFJmM2bYe/eZEchIukqPZOCug1jX3fvlsiV9agdjC22K+/xXqUvlUAcd+hN5s2bR4MGDWjRogVLlixh6NChsR8XlJ/YiIhIrXLTTdC2LfzqV7BjR7KjEZF0k55Jwa518a0fqbIeqeMxBnh4UbQr7/FcpY+SQJzWZBOzZ89m+fLl9OrVK75jgvITGzUhEhGpNfbtg4ULA52L778fdiVotGwRyVzpmRTEK1JlPdI8AqUTgiKRrshHvUpftjz/ld9ETSD69u3L8ccfX+EhRNT7XiLPvuZqQiQiUou8/DJ8911g+YQToHPn5MYjIulHSUGRSJX40I7HozZH76sQ6Yp8tKv0pcqff/55uv5uKzsitRONlljEqssI4kpkREQkJc2dW7I8eHDy4hCR9JV5SUE8FfvSIjUpys4NXpGPsK7VDS+zumHrPvzww1x22WWs/RwGPgXfHqhETBWpyvGKVJaGwhWpNoWFMG9eyfNBg5IXi4ikr8xLCuKp2JcWqUlRnynRZyU2i/i8sLCQsWPHctNNN+EeuJK/chsMmx6YGiGwbt3YYqpIVY5XpDI0FK7EwMz6mdn7ZrbJzH4d4fVbzGydma02s8VmFuOwculn5Ur4/PPAcqtWcNppyY1HRNJTdrIDqHFFHYAtKzCfQeN2gQpytIp9aV1GxLbusnFQeDC8rPAg3y35DSPvnMszzzwT9lJWHbj4+yF5ROmEorKKYl02LtBkqPHR8R2vlFg/Q7/HWJTXyV6/LwHMLAv4M3AhsB1YaWZz3T10lIh3gB7uvs/MbgDuB66o+WiTL7Tp0MCBkBVtfk4RkSrIvKSgaDhOLyi5Yp6IikqENvu798GPHt3Kqx+Fv9aonvHcT52+nUIKCw9WXyUq1kRGoiu6+l1U2S26+g363ZYWRyd7yVinApvc/SMAM5sJDAGKkwJ3Xxqy/hvAT2o0whSRnw8zQm6yqemQiCRK5jUfCpXIibxKtdnfuhvOeARe/Sh8tSOOOIJ/3lAqIShSGytRoW3JH2kJf26ZHu3KNRFc7GLsZC8Z7ShgW8jz7cGyaK4BFiQ0ohT13HOwJXgtq2VL6NcvufGISPrK7KQAElfxDmnL/+9PoOdkWF9qsplOnTqxYsUKundJk87ApduSH/gSvvuStGhXHuvVb3WwVT8WiUWUsZIjrGj2E6AH8ECU10eZ2SozW7Vz585qDDH53OGBkKP+r/+CBg2SF4+IpDclBfFWvGOt9AU7JS/adjhnPQqffR3+8hlnnMHy5ctp3759+lSiIl1ND1Wbr6zHcvVbHWwD4u2QL5loO9A25Hkb4NPSK5nZBcA4YLC7lx6fDQB3n+LuPdy9R6tWrRISbDKNGwc9e0L9+oGkQEQkUdI0KYi1k67FV/GOs9L315UFDPjzLvaWOpVdcsklvPzyyxx22GGBgspWolLtqnQsd11qY5MoiC1xUxOjEqXn+FBCIOFWAh3NrIOZ5QDDgLmhK5hZN+AJAgnBjgj7SHtmcPHFsHw5rFsXGHlIRCRR0rSjcZQJuyKtF09lJY5RVaZOncrVV19dZhdjxozhwQcfJKv08BHxdgZOxY6vjY8u6chd3jq1USyjOKmDbfk0epMEuXu+mY0GFgJZwFPuvtbM7gJWuftcAs2FGgHPWmA0tq3unpHTdplBhw7JjkJE0l2aJgUxKj2xV0WVljgqff369aNdu3Zs2VJSSX7ooYe4+eabqyPy1Bz2sfe9MP9KoiZltbFJVKiKErdoSVFtTYSqUyomsZJU7j4fmF+qbELI8gU1HlSKKCyEOml6H19EUlfmfu2UrqDG0jQojlFVWrduzYIFC2jevDk5OTnMmjUrkBBUV5OfVLwq3WUE5d6lSfd25enSNyQR1LRKJCaffhq4K/D738Pu3cmORkQySQYmBVHa7MdSaYmz0telSxfmzp3LSy+9xOWXX169HVFTddjH0ndfQsvTOSEAdbAtTyomsSIpaNIk2Lo10MH48suTHY2IZJLMaz40tjByeSyVlijtygs7DY+aXZ155pklT6qzyU/ve8ObY0BqXJVO1bhqiiaKi0xNq0QqtHcvPP54yXONNiQiNSkD7xREEeuV91Kjqryx51i6du3K5s2bK36PqInHlvjvFqTqVelUjUviU90jW6lplUiFnnwS9uwJLHfsqNmLRaRmZd6dgmgqcYV77ty5DBs2jP3799O/f39ef/31kmFGIylvdJ75P4GF10PBvthHZkm1q9KlO2oP+FtqxZep4h31JxGdgmMZvSmTaWQmAaZMKVm+5RYoPUidiEgiZd6dggomG4v1Cvdjjz3G0KFD2b9/PwAbNmxgyJAhHDhwoOR9Sl9pjXS1NFTBt9TaSa/i7S+RanMspKvK9GOJp1NwPH9HzV0QmSa9E2DjRli/PrDcoAFceWVy4xGRzJN5SUF5J9sYKi2FhYXcfvvt3HjjjRQWhvdPOOecc8jJyYl+kodAohGL2jYyS7wVSVWCakZlRv2JtVOw/o7VQyMzCfDCCyXLF1wADRsmLxYRyUyZlxTk74MFIyt1hfrgwYP89Kc/ZeLEiWHlderU4fHHH+fuu+/GzCruUFyvRWxvWJ0jsyT6ynw8o8uoElRzKjPqT6z9a/R3DKjqZ0sjMwnhScHgjJyiTUSSLfOSAgAvoPjK5vwr4UGr8GS+Z88eBgwYwIwZ4evk5uYyZ84crr/++pLCik7yFmOc1TUyS01c0Y1niFRVgmpOZYaujbVTsP6O1fPZStXhhaXG7N4Ny5aVPL/oouTFIiKZKzOTgjDBybZCT+alrvxtXzqJ3r17s3jx4rAtW7VqxdKlSxk4cGD4Lis6yX+3q+KwqnNklpq4ohvP6DKqBIVL5F2cyoz6E2v/Gv0dq+ezpZGZMt6CBVBQEFg+9VRo3Tq58YhIZlJSECp/HywZE3bl770PttDr4ptYs2ZN2KrHHXccK1as4NRTTy27n4pO8lErTXVIyDCeNXFFN56O2qoElUj0XZzKDhEbS6dg/R2r57OlYXwz3sGDcHTwtKCmQyKSLBqStLTvvixefGUT/Ohp2POdh61y2mmn8cILL9CqVavI+6ho+MVow58mqiJQUxNHxTpEqoanLFGdE9pFk6iha/V3rL7PVqoNLyw16mc/g5EjYfVqiHZaERFJNCUFUSz6AAY+CQcLwssHDx5MXl4eubnlDC0K5Z/kY61MVdfY5ak4y7AqQQG1vV1+pv8dU/GzJbWSGZx8crKjEJFMpqQgih+2hY4tYe3nJWU33HADkydPJivSjDLRKvDRyiuqTFXnBFK6opu6auoujiSGPlsiIpImzN0rXivJevTo4atWrYp9gwdjHd4nkjpAYP6Brbuh12T49Gv4w+BcfvX/vgkMOVpa6Qo8BK4WnjgS1k6Lv5nQ+hmBYVO9oOxrjdsF2nhLeoj2v1PVpmSaITfjmNm/3L1HsuNIprjPFUn0t7/B11/DVVdBRTeeRUSqS3nnigzvaFyqgp+dS1FCAHB0c5h/Lcz4Mfz6rP2REwKI3i589ZT4RyYpqiRGSgig9jQrkdgkopOpJhUTSWmHDsH48TB6dKCDcS3JY0QkzWVu86GiK/kfzce/3oI1aRe4mrpsXFhzjpOPDDzKbc4RraJemYp9pAQjlJqVpJ/qbpcfLUldPEZ3C0SSaMsWePxx+PBD2Bo8DdSpAyeemNy4REQgU5MCyyq+GvvRS//D8Ot/zVOXbuHEb4qa7BjF8xdAxR0H6x0GB74sW25ZUZoAVSLBiCUOEYj+P3Tgy8DdAiUGIjVu2zY480zYvj28fPRoaNAgOTGJiITKvOZD2bnQfxp0GcGqZ++m1yW38Nbmg/T/C3zyVVEF3iluWlRRc471M+DQ3rLlVhd+MCr+cdyjJQwhiYxIucpLOqtzwjoRicmXX0KfPmUTghYt4MYbkxOTiEhpmZcUBCvW8+fP55wr72DHN4HibV/BgL/A198VreglnXrLq4gvGweFB8uW12sCFzwaf3vxaBNCBRMZkQqVl3SqT4pIjbvuOtiwIbBcty7cfTc88gi88Qa0bJnc2EQqy93p0KEDZsamTZvKvH7nnXfSMso/+H//93/Tvn37MuWvvPIKAwcOpGXLluTk5NC+fXt++ctfsnVrzZy75syZw/e//33q16/PCSecwKxZs2LaLj8/n4kTJ9KxY0fq1atHmzZtuPnmm8PWad++PWYW9jjiiCPC1nn66afLrGNmPP7448XrvPLKKxHXMTP69u1bpePPvOZDXUbw5JNPcv3111NQED7y0gnfg3qhv5HSFahII7pEq2R9t6v4/eKqzGuIQ6mqLiMCM3N/F6FJm/qkiNSoXbtg7tyS53/7G1xxRfLiEakuK1asYPPmzQDMnDmT8ePHV2l/kyZN4qabbuKSSy7hiSeeoFWrVnz44Yc89dRTDBkyhHfeeacaoo7utdde45JLLuHGG29k0qRJzJ8/n+HDh9O8eXP69OlT7rZXXXUVixcv5o477qBz585s27aNdevWlVnvxz/+Mb/4xS+Kn+fk5ETc35IlS2gQ0q7wmGOOKV7u3r07K1asCFt/69atXHHFFfTv3z+mY40mo5ICd7jzjju46667yrx22znwhwGBTl/FiipQ62cEOmmG9hsoGtElWn+CqlS+Mn1CKKm68x7WpFoiKWDBAigItkw97TQlBJI+8vLyaNiwISeddBJ5eXlVSgreeecdbrnlFsaPHx9WRzvrrLO46qqrmDdvXnWEXK67776bs846i0mTJgFw7rnnsnbtWu66665yk4IXX3yRmTNn8u6773LCCSeU+x6tW7emZ8+eFcbywx/+kEaNGkV8rUmTJmX2sWzZMurUqcPll19e4b7Lk57Nh+qW/UUeKoBrnssukxCYweQfwX0DSyUERRWoouEdI1X88/cFuh7E22+gIutnwJT28GCdwM9UGUoyVeOSshIx1KmIxO2442DkyEAzocGDkx2NSPUoKCjg2WefZfDgwVx99dWsW7eO1atXV3p/kydPpmXLlvz2t7+N+PrAgQMrve9YHDhwgKVLl5apVA8bNowVK1awZ8+eqNs+9dRTnHfeeRUmBImUl5fH2WefzZFHHlml/aRnUnDh44TOQbD3Oxj0FEx9Mz9stfrZ8Pw1jRh9QYtAgQVnKg6tQFU0ROh3u6q38pWqY8ynalwSXZcRgT4xYwsr7hsjIglx2mnw9NPwn//ATTclOxqR6rFkyRI+//xzhg0bxqWXXkrdunXJy8ur9P7++c9/cv7551O3bt1KbZ+fn1/ho7zJej/88EMOHTpE586dw8q7dOlCYWEhH3zwQdRt33zzTY4//nhGjx5NkyZNyM3N5eKLL+bTTz8ts+5TTz1FTk4OTZs25dJLL2XLli0R9gjHHnss2dnZdOrUiSeeeKLcY9+4cSPvvPMOw4cPL3e9WKRx86E6QAGffQ0XPQnvfBL+aotceOFq6NX+G8gvhAHTI1eaKuqY2fjo6m3uE22M+WXjklupS9W4RERqgawszVws6SMvL49mzZrRr18/cnJyuPDCC5k5cya///3vo0/0Wo5PPvmEo4+ufLPrWJKJqVOn8rOf/Szia7t37wagWbNmYeXNmzcPez2S//znPzz99NOcfPLJzJw5k71793LbbbcxdOhQ3njjjeLfx5AhQ+jZsydt2rRh/fr1/O53v6N3796sWbOGpk2bAoHmRXfffTennnoqBQUF5OXl8fOf/5x9+/aV6bhcJC8vj7p163LJJZdU+DuoSFKSAjPrBzwMZAF/cfeJ1foGi8cABaz/HPr/BbaU+lt2OAxevA6ObxUsKK9yG63PACSmjXa0JCTZo8akalwiIinmtdfgoYdg/354/nklA5JeDhw4wOzZsxk6dGhxR9nhw4dz5ZVX8sYbb9CrV69K7bcyyUSRlStXVrhOhw4d4o6h6O5CebG5O+7OnDlzaNEi0PKkdevWnH322SxZsoTzzz8fgIcffrh4m969e3P66afTtWtXpk6dyk3B24h9+/YNG0Gof//+HDhwgHvuuYcxY8ZQp07ZBj4zZ86kT58+HHbYYRUeX0VqPCkwsyzgz8CFwHZgpZnNdfey3bQr68CXrPkMzn4Udu8Pf+mUNvCPa+B7jUttE61yG/X/oE5i2mg3PjpsRuWw8mRK1bhERFLIm29C376wL3hjtSDKxPYitdWCBQv46quvGDBgAF999RUA55xzDvXq1SMvL684KcjOzqYgygegoKCA7OySKuhRRx1VpWFHu3btWuE6WVlZUV8ruiNQdDxFip6XvoNQettjjjmmOCEAOPPMM8nJyWHdunXFSUFpJ510Ep06deLtt98uN+5LL72UZ555hs2bN4eNQgTw7rvvsn79esaNq545iJLRp+BUYJO7f+TuB4GZwJDqfpPjWkLnw8PLBnSGV25tWzYhgOiV26KhRcvwxDSbiTZPQbJHjUnVuEREkmz7dvjDH2DCBBgwoCQhEElHRX0HLrvsMpo3b07z5s1p27YtBw4c4JlnnilOBFq1asXXX3/NvggfiM8++4zDDy+ppJ1zzjksXryY/Pz8MuvGom7duhU+pk2bFnX7Y489lrp167KhaEKRoA0bNlCnTh2OP/74qNt26dIlYrm7R7yyX1qsd0girTdz5kwaNGjAkCHVU41ORlJwFLAt5Pn2YFkYMxtlZqvMbNXOnTvjfpMGdWHu1dAxOG/GdafBnKug0YV/iK9yGy1ZSNQV8lQdNSZV4xIRSaJt26BXL/jNbwKTku0KXkdq2BAmT4aQocZFar1vvvmGefPmMXz4cJYuXRr2eOihh/j8889ZunQpEGgiU1hYWGY40W+//ZbFixfTu3fv4rJf/OIX7Ny5k3vvjVwXmz9/frlxrVy5ssLHoEGDom5fr149zj33XJ599tmw8lmzZtGrV6/iNv+RDBw4kNWrV/PFF18Ul7366qscOnSIk08+Oep27733Hu+//z6nnHJKucf2/PPP07JlS9q1a1fmtVmzZjFo0KCow5fGragtVE09gMsI9CMoen4lMLm8bU455RSPyx8pfnx4O/7AQLzwgWCZu/u66e5PtHP/owV+rpsefV/rprv/KTdsn/6n3PK3ERGpYcAqr+Hv81R7xH2uqKIvvnDv3Nk9MAtOyaNhQ/e33qrRUERqxPTp0x3wN954o8xrBw8e9BYtWvjVV19dXHbFFVd4o0aNfOLEib5o0SKfMWOGd+/e3Vu0aOHbt28P2/7hhx92M/PLL7/c//73v/urr77q06ZN87PPPtu7du2a8GNbtmyZZ2Vl+ZgxY3zp0qV+6623upn5woULi9fZvHmzZ2Vl+bRp04rL9uzZ423btvWePXv63LlzfcaMGd6mTRu/4IILiteZN2+eDxs2zKdPn+5LlizxRx991I888kjv0KGD79mzp3i9iy++2CdOnOjz58/3F154wX/yk5844JMmTSoT74oVKxzw2bNnx3Wc5Z0rkpEU9AIWhjy/Hbi9vG2qkhSUeVRGPEmEiEgSKCmo2aTgm2/cTzutJBGoW9f9V79yv+ce9zVraiwMkRp10UUXeceOHaO+fsMNN3izZs38u+++c3f3AwcO+G9/+1s/5phjPDs725s3b+5Dhw719evXR9x+yZIlPmDAAD/ssMM8Ozvb27Vr56NGjfKNGzcm5HhKmz17tp944omek5PjnTp18ry8vLDXP/74Ywd86tSpYeUbN270/v37e25urjdr1sxHjhzpu3btKn793Xff9fPOO89btmzp2dnZ/r3vfc9Hjhzpn3zySdh+br/9dj/++OO9QYMGXr9+fe/evbv/9a9/jRjrmDFjvGnTpsW/61iVd64fAPpZAAAPVklEQVSwwOs1x8yygQ+A84FPgJXAj919bbRtevTo4atWrYr9TR4sp33W2Jo9XhGRmmBm/3L3HsmOI5niPldUwZQpcP31gWUzmDkTqjiZqIhIwpV3rqjxPgXung+MBhYC64FnyksIREREUs1118E99wSWJ01SQiAitV9S5ilw9/lA+b1GEuXBOoFOwr3vVSdZERGpFDMYNw769YMK+gmKiNQKaTyjcTQeGG//pVGBp0oMRESkAhMmQOPGMGoUhA5EooRARNJFMoYkTTyLPkFFsaJZjEVERMrx+edw//1w223Qti18/HGyIxIRqX7pmRT8YFRs60WbxVhERCToz3+GAwcCy506Qfv2SQ1HRCQh0rP50AWPBn6ungJezhzziZqATERE0sK+fYGkoMittwb6E4iIpJv0vFMAgcTglvzAEKQDpsc3i7GIiAgwdWrJLMXt28PFFyc1HBGRhEnPOwUAL99YcqfAsqDNOfDVpkCTIY0+JCIiFSgogIceKnl+882Qnb5nTRHJcOn59fbyjfDuYyXPvQC2LYaTbyhpWiQiIlKO2bPho48Cy82bw9VXJzceEZFESs/mQ6unxFcuIiISwh0eeKDk+Q03QKNGyYtHRCTR0jMpiNa5uLxOxyIiIkGvvw5vvRVYzsmB0aOTG4+ISKKlZ1IQbZ6CWOYvEBGRjBd6l+DKK6F16+TFIiJSE9IzKYg2T0Gs8xeIiEjGcoczzihJBG65JbnxiIjUhPRMCi54NNCpuOjOgGWpk7GIiMTELDB78ccfw8KFcMIJyY5IRCTx0nP0IQgkAEoCRESkkurVgz59kh2FiEjNSM87BSIiIiIiEjMlBSIiIkH5+cmOQEQkOZQUiIiIAIWF0LEj9O8Pjz0GBw4kOyIRkZqTvn0KRERE4rBqFWzeHHj8618wSgPWiUgG0Z0CERERYO7ckuWLLoIsTW0jIhlESYGIiAjhScHgwcmLQ0QkGZQUiIhIxtu8GdasCSzn5MCFFyY1HBGRGqekQEREMt4LL5Qsn38+NGqUvFhERJJBSYGIiGS80KRg0KDkxSEikixKCkREJKPt2QOvvFLyXEmBiGQiJQWxWD8DprSHB+sEfq6fkeyIRERqNTPrZ2bvm9kmM/t1hNfrmdms4Otvmln7RMWycCEcOhRY7t4d2rRJ1DuJiKQuJQUVWT8DXhoFe7cAHvj50iglBiIilWRmWcCfgf7ACcBwMzuh1GrXALvd/Tjgf4D7EhWPmg6JiCgpqNiycZC/L7wsf1+gXEREKuNUYJO7f+TuB4GZwJBS6wwBpgWXnwPONzNLRDD160PDhoFlDUUqIplKSUFF9m6Nr1xERCpyFLAt5Pn2YFnEddw9H9gDtCi9IzMbZWarzGzVzp07KxXM//4vfPFFoBlRt26V2oWISK2npKAijY+Or1xERCoS6Yq/V2Id3H2Ku/dw9x6tWrWqdED160OfPpCYexEiIqlPSUFFet8L2bnhZdm5gXIREamM7UDbkOdtgE+jrWNm2UBTYFeNRCcikoGUFFSkywjoMwUatwMs8LPPlEC5iIhUxkqgo5l1MLMcYBgwt9Q6c4GRweVLgSXuXuZOgYiIVI/sZAdQK3QZoSRARKSauHu+mY0GFgJZwFPuvtbM7gJWuftc4Engb2a2icAdgmHJi1hEJP0pKRARkRrn7vOB+aXKJoQsfwdcVtNxiYhkKjUfEhERERHJcEoKREREREQynJICEREREZEMp6RARERERCTDKSkQEREREclwSgpERERERDKckgIRERERkQynpEBEREREJMMpKRARERERyXDm7smOoUJmthPYUsnNWwJfVGM4tUGmHXOmHS/omDNBvMfbzt1bJSqY2iADzhWpHqPiq5pUjw9SP0bFV7Go54pakRRUhZmtcvceyY6jJmXaMWfa8YKOORNk2vEmW234fad6jIqvalI9Pkj9GBVf1aj5kIiIiIhIhlNSICIiIiKS4TIhKZiS7ACSINOOOdOOF3TMmSDTjjfZasPvO9VjVHxVk+rxQerHqPiqIO37FIiIiIiISPky4U6BiIiIiIiUI22TAjPrZ2bvm9kmM/t1suNJBDN7ysx2mNl7IWWHmdkiM9sY/Nk8mTFWNzNra2ZLzWy9ma01szHB8rQ8bjOrb2Zvmdm7weP9XbC8g5m9GTzeWWaWk+xYq5uZZZnZO2Y2L/g8rY/ZzDab2Roz+7eZrQqWpeX/dapJ5fNFtO+8VFP685pqzKyZmT1nZhuCv8teyY4plJndHPz7vmdmeWZWP8nxpHz9IkqMDwT/xqvNbLaZNUul+EJe+28zczNrmYzYoknLpMDMsoA/A/2BE4DhZnZCcqNKiKeBfqXKfg0sdveOwOLg83SSD4x19y5AT+C/gn/bdD3uA8B57n4y0BXoZ2Y9gfuA/wke727gmiTGmChjgPUhzzPhmM91964hQ9al6/91yqgF54to33mppvTnNdU8DLzo7p2Bk0mhWM3sKOCXQA93PwnIAoYlN6paUb94mrIxLgJOcvcfAB8At9d0UCGepmx8mFlb4EJga00HVJG0TAqAU4FN7v6Rux8EZgJDkhxTtXP3V4FdpYqHANOCy9OAH9VoUAnm7p+5+9vB5b0EvtiPIk2P2wO+CT6tG3w4cB7wXLA8bY63iJm1AS4C/hJ8bqT5MUeRlv/XKSalzxflfOeljNKf11RjZk2As4AnAdz9oLt/ldyoysgGGphZNpALfJrMYGpD/SJSjO7+krvnB5++AbSp8cBKYon0OwT4H+A2AufylJKuScFRwLaQ59tJsS/RBPqeu38GgZMJcHiS40kYM2sPdAPeJI2PO3hb/t/ADgJXQT4Evgr54kvH/+8/EfjSLAw+b0H6H7MDL5nZv8xsVLAsbf+vU0itOV+U+s5LJaU/r6nmGGAnMDXYxOkvZtYw2UEVcfdPgD8SuHL8GbDH3V9KblQR1bbvo6uBBckOIpSZDQY+cfd3kx1LJOmaFFiEspTLyKTyzKwR8Dxwk7t/nex4EsndC9y9K4ErHqcCXSKtVrNRJY6ZDQR2uPu/QosjrJo2xxx0hrt3J9CM5b/M7KxkB5QhasX/Vqp+50X5vKaabKA78Ji7dwO+JflNX4oF2+YPAToARwINzewnyY2qdjOzcQSa3s1IdixFzCwXGAdMSHYs0aRrUrAdaBvyvA1JvhVXgz43s9YAwZ87khxPtTOzugROjjPc/e/B4rQ/7uDt7lcItCtuFrzNDOn3/30GMNjMNhNoynEegSuR6XzMuPunwZ87gNkEEsC0/79OASl/vojynZcqynxezWx6ckMqYzuw3d2L7rA8RyBJSBUXAB+7+053PwT8HTg9yTFFUiu+j8xsJDAQGOGpNe7+sQQSv3eDn5c2wNtmdkRSowqRrknBSqBjcLSSHAIdduYmOaaaMhcYGVweCcxJYizVLti2/Elgvbs/FPJSWh63mbUqGj3BzBoQOHmsB5YClwZXS5vjBXD32929jbu3J/DZXeLuI0jjYzazhmbWuGgZ6AO8R5r+X6eYlD5flPOdlxKifF5T6iq3u/8H2GZmnYJF5wPrkhhSaVuBnmaWG/x7n08KdYQOkfLfR2bWD/gVMNjd9yU7nlDuvsbdD3f39sHPy3age/D/MyWkZVIQbHc8GlhI4IP1jLuvTW5U1c/M8oAVQCcz225m1wATgQvNbCOB3u0TkxljApwBXEngatS/g48BpO9xtwaWmtlqApWXRe4+j8CX3i1mtolAe/snkxhjTUnnY/4e8JqZvQu8BfzD3V8kff+vU0YtOF9E+86T+PwCmBH8Lu0K/D7J8RQL3sF4DngbWEOgbpbUmW9rQ/0iSoyPAI2BRcHPyuMpFl9K04zGIiIiIiIZLi3vFIiIiIiISOyUFIiIiIiIZDglBSIiIiIiGU5JgYiIiIhIhlNSICIiIiKS4ZQUSMYxsxYhQ/v9x8w+CXmek6D37B4cP1lERFKImRUEv//Xmtm7ZnaLmVWqfmRmPcxsUgXrdA0dVtbMBptZysywLJlLQ5JKRjOzO4Fv3P2PcWyT5e4Fcb7PtcBJ7n5TnCGKiEgCmdk37t4ouHw48H/A6+5+R4Le72dAD3cfnYj9i1SW7hSIhDCzF8zsX8ErRtcGy7LN7Cszu8fM3gJODV7Zed/MlpnZZDP7f8F1G5nZ02b2lpm9Y2aDgjMRTwBGBK9GXVpOCCIikiTuvgMYBYy2gCwze8DMVprZajO7HsDMZpW62v+0mV1iZueY2bxg2almtjx4LlhuZp2Cd6PvAq4Ing+uMLOfmdkjwW3amdni4HstNrOjQ/Y/Kbifj3QekURQUiASbqS7nwL8kMDsuc2D5U2Bt939VOBd4FGgD3AWcETI9hOAF4PrnQc8CDiBk8AMd+/q7s/VzKGIiEi83P0jAvWjw4FrgD3u/kMC54XrzKwDMBO4AiBY0T8fmF9qVxuAs9y9G4Fzw+/d/WBweVbwfDCr1DaPAH919x8AM4DQpkitgTOBgWh2c0mA7GQHIJJibjazwcHlNsCxwL+Bg8DsYPkJwPvuvgWKpzL/afC1PkD/kPah9YGjayJwERGpNhb82Qf4QciV+aZAR2ABMMnM6gH9gFfdfb+Zhe6jKTDNzDoSuDhUN4b37QVcHFz+G3B/yGv/z90LgXVm9r1KHJNIuZQUiASZ2QUErvz3DH65v0agUg+w30s64FjEHZS89iN3/7DUvs+q9oBFRKTamdkxQAGwg8B3+i/cfWGE9V4B+hK4Y5AXYVd3A0vdfaiZtQdeqUQ4oR0/D4S+fSX2JVIuNR8SKdEU2BVMCE4kcKs4krVAJzNra4HLQleEvLYQ+GXREzPrFlzcCzROQMwiIlJNzKwV8DjwSPBC0ELgBjOrG3z9eDNrGFx9JnAV0Du4XmlNgU+Cyz8LKS/vfLAcGBZcHgG8VrkjEYmfkgKREv8Acs3sXQJtPt+MtJK77wNGAy8Dy4BPgT3Bl38X3McaM1sL3BksXwKcHOxwpg5iIiKpo0HRkKQEvtdfIvBdDvAXYB3wtpm9BzxBSSuLlwjcXX452FegtPuBP5jZ60BWSPlS4ISijsaltvklcJWZrQauBMZU/fBEYqMhSUUqwcwaufs3wTsFTwBr3H1ysuMSERERqQzdKRCpnBvM7N8EriA1AP43yfGIiIiIVJruFIiIiIiIZDjdKRARERERyXBKCkREREREMpySAhERERGRDKekQEREREQkwykpEBERERHJcEoKREREREQy3P8HpAEqOBK3YY4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 936x432 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "import time\n",
    "import numpy as np\n",
    "from sklearn.model_selection import StratifiedShuffleSplit\n",
    "from sklearn.model_selection import cross_val_score\n",
    "from sklearn import linear_model\n",
    "from sklearn.linear_model import LinearRegression\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from sklearn.linear_model import Lasso\n",
    "from sklearn.metrics import mean_absolute_error\n",
    "from sklearn.metrics import explained_variance_score\n",
    "from sklearn.metrics import r2_score\n",
    "from scipy.integrate import simps\n",
    "from sklearn.pipeline import Pipeline\n",
    "from sklearn.preprocessing import PolynomialFeatures\n",
    "%matplotlib inline\n",
    "# make stable output\n",
    "np.random.seed(123)\n",
    "\n",
    "# start the execution time\n",
    "start_time=time.time()\n",
    "\n",
    "# store data into a pandas.DataFrame\n",
    "try:\n",
    "    data = pd.read_csv('forestfires.csv')\n",
    "    print(\"File loaded locally! \")\n",
    "except FileNotFoundError:\n",
    "    try:\n",
    "        data = pd.read_csv('./data/forestfires.csv')\n",
    "        print(\"File loaded locally! \")\n",
    "    except FileNotFoundError:\n",
    "        url = 'https://raw.githubusercontent.com/haichenzhu/Forest-Fire-Area-Prediction/master/forestfires.csv'\n",
    "        data = pd.read_csv(url)\n",
    "        print(\"File loaded online! From: \", url)\n",
    "except:\n",
    "    print(\"No file named forestfires.csv in the ./data directory or no Internet access! Mission Failed! \")\n",
    "\n",
    "# remove these month due to there are few record\n",
    "remove_month = ['nov', 'jan', 'may']\n",
    "fire_data = data.drop(data[data.month.isin(remove_month)].index, axis=0)\n",
    "\n",
    "# classifier area values into 6 groups\n",
    "fire_data['area_bins'] = pd.cut(fire_data['area'], bins=[0, 5, 10, 50, 100, 150, 1100], \n",
    "                                include_lowest=True,\n",
    "                                labels=['0-5', '5-10', '10-50', '50-100', '100-150', '>150'])\n",
    "\n",
    "# delete the record when the area is larger than 50, due there is few records\n",
    "fire_data.drop(fire_data[fire_data.area>50].index, axis=0, inplace=True)\n",
    "\n",
    "# make these two object type feature's value to numerical values\n",
    "fire_data_num = fire_data\n",
    "fire_data_num.month.replace(('jan', 'feb', 'mar', 'apr', 'may', 'jun',\n",
    "                             'jul', 'aug', 'sep', 'oct', 'nov', 'dec'),\n",
    "                            (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12), inplace=True)\n",
    "fire_data_num.day.replace(('mon', 'tue', 'wed', 'thu', 'fri', 'sat', 'sun'),\n",
    "                          (1, 2, 3, 4, 5, 6, 7), inplace=True)\n",
    "fire=fire_data_num\n",
    "\n",
    "#select predictor variables and target variable as separate data sets  \n",
    "predvar = fire[[\"X\", \"Y\", \"month\", \"day\", \"FFMC\", \"DMC\", \"DC\", \n",
    "                \"ISI\", \"temp\", \"RH\", \"wind\", \"rain\"]]\n",
    "target = fire.area\n",
    "\n",
    "# split data into train and test sets\n",
    "spliter = StratifiedShuffleSplit(random_state=34, n_splits=1, test_size=0.3)\n",
    "\n",
    "train_index = fire.index\n",
    "test_index = fire.index\n",
    "for tr, te in spliter.split(fire.values, fire.area_bins.values):\n",
    "    train_index = tr\n",
    "    test_index = te\n",
    "pred_train = predvar.iloc[train_index]\n",
    "pred_test = predvar.iloc[test_index]\n",
    "tar_train = target.iloc[train_index]\n",
    "tar_test = target.iloc[test_index]\n",
    "\n",
    "\n",
    "def rec(y_true , y_predict):\n",
    "    '''\n",
    "    Regression Error Characteristic curve with X-axis represents deviation, and y-axis shows \n",
    "    the current accuracy under deviation, return the range of deviation, the accuracy of every\n",
    "    deviation, the area under the curve\n",
    "    '''\n",
    "    accuracy = []\n",
    "    begin = 0\n",
    "    end = 15\n",
    "    interval = 0.1\n",
    "    \n",
    "    epsilon = np.arange(begin, end, interval)\n",
    "    \n",
    "    # loop to find the accuracy for every deviation\n",
    "    for i in range(len(epsilon)):\n",
    "        count = 0\n",
    "        for j in range(len(y_true)):\n",
    "            if np.abs(y_true[j] - y_predict[j]) < epsilon[i]:\n",
    "                count = count + 1\n",
    "        \n",
    "        accuracy.append(count/len(y_true))\n",
    "    \n",
    "    #  using Simpson's rule to calculating area-under-curve\n",
    "    area_under_curve = simps(accuracy, epsilon) / end\n",
    "   \n",
    "    return epsilon, accuracy, area_under_curve\n",
    "\n",
    "\n",
    "def NLL(y_true, y_predict):\n",
    "    '''\n",
    "    defined as : −log𝑝(𝑦|𝐷,𝑥), \n",
    "    where y is the true value, D is the dataset, x represents a query point\n",
    "    '''\n",
    "    y_predict_mean = np.mean(y_predict)\n",
    "    y_predict_variance = np.var(y_predict)\n",
    "    \n",
    "    return np.log(2 * np.pi * y_predict_variance) / 2 + \\\n",
    "            np.sum((y_true - y_predict_mean) ** 2) / (2 * y_predict_variance)\n",
    "\n",
    "\n",
    "def mean_absolute_deviation(y_predict):\n",
    "    '''\n",
    "    the mean absolue deviation of every prediction to the prediction mean, defined as:\n",
    "    mean(abs(p - y_predict.mean() for p in y_predict))\n",
    "    '''\n",
    "    prediction_mean = np.mean(y_predict)\n",
    "    distance = np.abs(y_predict - prediction_mean)\n",
    "    return np.sum(distance) / len(y_predict)\n",
    "\n",
    "\n",
    "# make a pipline to do polynomial features and lasso regession sequantly\n",
    "model = Pipeline([('poly', PolynomialFeatures(degree=2)),\n",
    "                  ('lasso', Lasso(fit_intercept=False, max_iter=800000))])\n",
    "\n",
    "print(\"Start training Lasso regressor...\")\n",
    "model.fit(pred_train, tar_train)\n",
    "\n",
    "# predicting\n",
    "print('Training finished! Start predicting...')\n",
    "y_predict = model.predict(pred_test)\n",
    "y_true = tar_test\n",
    "\n",
    "# print out performance\n",
    "print('The Performance of Our Lasso regressor:')\n",
    "mse = mean_squared_error(y_true, y_predict)\n",
    "rmse = np.sqrt(mse)\n",
    "print('MSE:  ', mse.round(3))\n",
    "print('RMSE: ', rmse.round(3))\n",
    "print('MAD:  ', mean_absolute_deviation(y_predict).round(3))\n",
    "print('MAE:  ', mean_absolute_error(y_true, y_predict).round(3))\n",
    "print('NLL:  ', NLL(y_true, y_predict).round(3))\n",
    "\n",
    "rr = r2_score(y_true, y_predict)\n",
    "print('RR:   ', rr.round(3))\n",
    "\n",
    "deviation, accuracy, area_under_curve = rec(y_true.to_numpy(), y_predict)\n",
    "print('AUC:  ', area_under_curve.round(3))\n",
    "\n",
    "print(\"Explained variance regression score:\", explained_variance_score(y_true, y_predict).round(3))\n",
    "\n",
    "# record the running time \n",
    "elapsed_time = time.time() - start_time\n",
    "print(\"Time taken to run Logistic Regression =\", elapsed_time, \"seconds\")\n",
    "\n",
    "# Plotting\n",
    "plt.figure(figsize=(13 , 6))\n",
    "\n",
    "plt.subplot(1, 2, 1)\n",
    "plt.scatter(y_true, y_predict, color=\"darkorange\")\n",
    "plt.xlabel(\"Target\")\n",
    "plt.ylabel(\"Predicted\")\n",
    "plt.plot([y_true.min(), y_true.max()], [y_true.min(), y_true.max()], 'k--', lw=4)\n",
    "plt.text(30, 10, r\"$R^2 = %0.4f$\" %rr, fontsize=15)\n",
    "\n",
    "plt.subplot(1, 2, 2)\n",
    "plt.title(\"Regression Error Characteristic (REC)\")\n",
    "plt.plot(deviation, accuracy, \"--b\", lw =3)\n",
    "plt.xlabel(\"Deviation\")\n",
    "plt.ylabel(\"Accuracy (%)\")\n",
    "plt.text(10, 0.09, \"AUC = %0.4f\" %area_under_curve, fontsize=15)\n",
    "\n",
    "plt.show()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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": 2
}