2Bases_model.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Expansion modification system - 2 bases"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### DNA Sequence generation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 93,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy import stats\n",
    "from scipy.sparse import csr_matrix\n",
    "import subprocess\n",
    "import numpy as np\n",
    "from Bio import SeqIO\n",
    "\n",
    "import tensorflow as tf\n",
    "\n",
    "from tensorflow.keras.preprocessing.text import Tokenizer\n",
    "from tensorflow.keras.preprocessing.sequence import pad_sequences\n",
    "from sklearn.feature_extraction.text import CountVectorizer\n",
    "from sklearn.feature_extraction.text import TfidfTransformer\n",
    "\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "from tensorflow.keras.models import Sequential\n",
    "from tensorflow.keras.callbacks import ModelCheckpoint\n",
    "from tensorflow.keras.optimizers import SGD, Adam, Adadelta, RMSprop\n",
    "from tensorflow.keras.layers import Conv1D, Dense, MaxPooling1D, Flatten, Dropout\n",
    "from tensorflow.keras.layers import Embedding, GlobalAveragePooling1D, LSTM, SimpleRNN, GRU, Bidirectional\n",
    "\n",
    "from sklearn.metrics import confusion_matrix\n",
    "import itertools"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def DNA_gen(N0, N, Pmut, g, M = 0):\n",
    "    print('Pmut: ', Pmut)\n",
    "    DNA = list(np.random.randint(0,2, size = N0)) \n",
    "    limit = [0,g,1]\n",
    "    base = [1,0]\n",
    "    print(limit)\n",
    "    while (len(DNA) != N):\n",
    "        j = np.random.randint(0,len(DNA))\n",
    "        nuc = DNA[j]\n",
    "        p = np.random.random()\n",
    "        if(p>Pmut):\n",
    "            DNA.insert(j+1, nuc)\n",
    "        else:\n",
    "            if(False):\n",
    "                pg = np.random.random()\n",
    "                v = np.array(limit)-pg\n",
    "                v[v<0] = 1\n",
    "                idx = np.argmin(v)\n",
    "                DNA[j] = base[idx-1]\n",
    "            if(nuc==1):\n",
    "                DNA[j] = 0\n",
    "            else: DNA[j] = 1\n",
    "    count = 0\n",
    "    while(count != M):\n",
    "        j = np.random.randint(0,len(DNA))\n",
    "        pg = np.random.random()\n",
    "        v = np.array(limit)-pg\n",
    "        v[v<0] = 1\n",
    "        idx = np.argmin(v)\n",
    "        DNA[j] = base[idx-1]\n",
    "        count += 1\n",
    "\n",
    "    return DNA"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Corr(s,r,alphabet):\n",
    "    bases = ['A','C']\n",
    "    C = 0\n",
    "    freq = 0\n",
    "    s1 = s[:-r]\n",
    "    s2 = s[r:]\n",
    "    long_dist = np.equal(s1,s2)\n",
    "    long_dist = np.sum([int(long_dist[i]) for i in range(len(long_dist))])/len(long_dist)\n",
    "    for i in bases:\n",
    "        freq += (alphabet[i][0]/len(s))**2\n",
    "    C = long_dist -freq\n",
    "    return C"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "def smooth(y, box_pts):\n",
    "        box = np.ones(box_pts)/box_pts\n",
    "        y_smooth = np.convolve(y, box, mode='same')\n",
    "        return y_smooth\n",
    "def scaling(x,a,c):\n",
    "    return c*x**(-a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Sequence Generation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "**Human**\n",
      "Pmut:  0.1\n",
      "[0, 0.4, 1]\n",
      "**Neanderthal**\n",
      "Pmut:  0.12\n",
      "[0, 0.4, 1]\n"
     ]
    }
   ],
   "source": [
    "N = 5000000\n",
    "N0 = 1\n",
    "\n",
    "print('**Human**')\n",
    "human2 = DNA_gen(N0,N,0.10,0.4, M = int(N/3))\n",
    "\n",
    "print('**Neanderthal**')\n",
    "neanderthal2 = DNA_gen(N0,N,0.12,0.4, M = int(N/6))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "change = {0:['A'], 1:['C']}\n",
    "N_DNA2 = []\n",
    "H_DNA2 = []\n",
    "for i in range(N):\n",
    "    N_DNA2.append(change[human2[i]][0])\n",
    "    H_DNA2.append(change[neanderthal2[i]][0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "**Human**\n",
      "{'A': [2482016], 'C': [2517984]}\n",
      "GC-content:  0.5035968\n",
      "**Neanderthal**\n",
      "{'A': [2618255], 'C': [2381745]}\n",
      "GC-content:  0.476349\n"
     ]
    }
   ],
   "source": [
    "human_alphabet2 = {'A':[0], 'C':[0]}\n",
    "neanderthal_alphabet2 = {'A':[0], 'C':[0]}\n",
    "for (i,j) in zip(H_DNA2,N_DNA2):\n",
    "    human_alphabet2[i][0] += 1\n",
    "    neanderthal_alphabet2[j][0] += 1\n",
    "    \n",
    "print('**Human**')\n",
    "print(human_alphabet2)\n",
    "print('GC-content: ',(human_alphabet2['C'][0]#+human_alphabet['G'][0]\n",
    "                     )/len(H_DNA2))\n",
    "\n",
    "print('**Neanderthal**')\n",
    "print(neanderthal_alphabet2)\n",
    "print('GC-content: ',(neanderthal_alphabet2['C'][0]\n",
    "                      #+neanderthal_alphabet['G'][0]\n",
    "                     )/len(N_DNA2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 108,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "C[H]: 0.2272081990841874\n",
      "G[H]: 0.43267001061430865\n",
      "C[N]: 0.34988264488625825\n",
      "G[N]: 0.5832441197389489\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAD8CAYAAACPWyg8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deZyNdf/H8ddnjJ3Idkej7MNgzDAGkT1kJ4VkzzKWSL+E7tByxy2Kyb4Xhey7yi5hDMY2Y802SkJNKrJ9f3985dY0mDNm5jpz5vN8PObRnHOu5TOXHvOe73V9FzHGoJRSSrnCy+kClFJKpTwaHkoppVym4aGUUsplGh5KKaVcpuGhlFLKZRoeSimlXObtdAFJLVeuXKZAgQJOl6GUUinGrl27Lhhjct9vG48PjwIFChAeHu50GUoplWKIyKkHbaO3rZRSSrlMw0MppZTLNDyUUkq5zOOfeSil3NP169eJjo7m6tWrTpeSamXIkAEfHx/Spk3r8r4aHkopR0RHR5M1a1YKFCiAiDhdTqpjjOHixYtER0dTsGBBl/fX21ZKKUdcvXqVnDlzanA4RETImTNnglt+Gh73sOLICs79ds7pMpTyaBocznqY66/hEYffr/1Ou8Xt8B3ry9iwsdy8ddPpkpRSSSBLlix/ez1z5kx69erlUDUpi4ZHHDKny8yOl3dQ4fEK9F7dm+CpwYSdDXO6LKWUchsaHvdQNGdRvnzpS+a1mMcPl3+g4tSKhKwI4ecrPztdmlIqGXTo0IEFCxbcef1XK2Xjxo1Uq1aNF154gWLFijFgwAA+++wzgoODKV26NMePHwdg+fLlVKhQgcDAQGrXrs2PP/4IwNChQ+nUqRPVq1enUKFChIaGJv8Plwi0t9W9vP02EhDAC41aUK9IPYZsGEJoWCgLoxYyss5I2vq31fu1SiWWvn0hIiJxjxkQAKNH33eTK1euEBAQcOf1pUuXaNy48QMPvXfvXqKiosiRIweFChXi5ZdfJiwsjDFjxvDxxx8zevRoqlSpwvbt2xERpk6dyogRIxg1ahQAhw4dYsOGDVy+fBlfX19CQkIS1F3WSdryiMuVKzB7NjRtCv7+PPLFUj6qNYJdXXdROEdh2i9pT/VPqnPw/EGnK1VKPYSMGTMSERFx5+udd96J137ly5cnb968pE+fnsKFC1OnTh0ASpcuzcmTJwHbFblu3bqULl2aDz74gIMH//f7okGDBqRPn55cuXKRJ0+eO62SlERbHnHJmBGiouCLL2DYMGjXDgYPJuD119naYS3Toj7njbVvEDApgH4V+zGk+hAypc3kdNVKpVwPaCE4wdvbm1u3bgF2TMS1a9fufJY+ffo733t5ed157eXlxY0bNwDo3bs3/fr1o3HjxmzcuJGhQ4fGuX+aNGnu7JOSaMvjXry94cUXYe9eWLYM8uaFnj3xKliILl9f5HC7nbT1b8uIb0dQbnI5Is4lcpNbKeWoAgUKsGvXLgCWLl3K9evXXdo/JiaGxx9/HIBPPvkk0etzmobHg3h5QaNGsHUrbNwIgYEwcCC5S5Rjelhevm70BTFXY6gwtQIfbfuIW+aW0xUrpRJBly5d2LRpE8HBwezYsYPMmTO7tP/QoUN5/vnnefrpp8mVK1cSVekcMcY4XUOSCgoKMom+nseuXTB8OCxcCOnTc6FvVzr7HWPZd6uoW7guM5vO5LEsjyXuOZXyMFFRUZQoUcLpMlK9uP4dRGSXMSbofvt5bMtDRBqJyOSYmJjEP3i5cjB/vn0u0qIFuYaHsuT/djL+kRfZdGoT/hP8WXV0VeKfVyml3ITHhocxZrkxpmu2bNmS7iS+vjBrFoSFIcV8Cen3OeGr8vMYWWjweQP6rO7D1Rs6Y6hSyvN4bHgkq/LlYfNmmD+fkj/cIGzQCV75sQChYaEETwnWLr1KKY+j4ZFYRKBFC4iKIsP7Ixjz2SVWzvHi3LljBE0OYmzYWH2YrpTyGBoeiS19enj9dTh2jPo1u7Hvo6tUP36D3qt7U2tGDU78fMLpCpVS6qFpeCSV3Llh/Hge27afVT/WZsoy2HV8C6VDizNhW6i2QpRSKZqGR1IrWRJZtZqXx2zhwO6KVDp+jR5f9aHOcD9O/XTU6eqUStVEhNdee+3O65EjR/5tJHhSSMi07xEREaxa9b8enEOHDmXkyJEuHaNAgQJcuHDBpX3uR8MjuVSpwhOrv+Wr9muZdLAQOy4fptQYXyaHtse4OHJVKZU40qdPz6JFixL1l2piu3Hjxj/Cwx1oeCQzqVWLrvOOsb/CTIJjstDt50+p2ysbp2eOgZu66JRSycnb25uuXbvy0Ucf/eOzn376ieeee47y5ctTvnx5tm7dCkBYWBhPPfUUgYGBPPXUUxw+fBiwLYrmzZtTr149ihYtSv/+/e8ca8aMGRQrVoxq1ardOc79zjF06FC6du1KnTp1aNeuHYMHD2bevHkEBAQwb948ACIjI+Oc1r1p06aUK1eOkiVLMnny5MS/aLfpxIhOEKFAk/Z83eglJs3oxes3JlPqSF9GNR7Gy60+QFq3tnNrKZVK9F3TN9Hnhwt4LIDR9R484WLPnj3x9/f/2y97gD59+vDqq69SpUoVTp8+Td26dYmKiqJ48eJs3rwZb29v1q5dy6BBg1i4cCFgby/t2bOH9OnT4+vrS+/evfH29mbIkCHs2rWLbNmyUaNGDQIDA+97DoBdu3bxzTffkDFjRmbOnEl4eDhjx44FbLjca1r36dOnkyNHDq5cuUL58uV57rnnyJkzZ2JeWkDDw1FeXmkI6TyBehf/j84zmtA1+CDz1rZjyqi3KNh3KLRpAylsjn+lUppHHnmEdu3aERoaSsaMGe+8v3btWiIjI++8/vXXX7l8+TIxMTG0b9+eo0ePIiJ/mzCxVq1a/DUw2c/Pj1OnTnHhwgWqV69O7ty5AWjZsiVHjhy57zkAGjdu/Ld6YvtrWvf06dPfmdbdx8eH0NBQFi9eDMCZM2c4evSohoenKpizMGtf28eU8Em8vuY1Sj9xhuETOtLjnbfxGvSmnRI+XTqny1QqycSnhZCU+vbtS9myZenYseOd927dusW2bdv+8Qu8d+/e1KhRg8WLF3Py5EmqV69+57N7TbV+r4Xj7nUO4IETMcZ1ro0bN7J27Vq2bdtGpkyZqF69OlevJs0sF/rMw014iRfdyodw4JVDVPF9ht71oVrDCxwZ0AWKFYNJk+DPP50uUymPlCNHDl544QWmTZt25706dercuU0E9pYU/H2q9ZkzZz7w2BUqVGDjxo1cvHiR69evM3/+/AeeI7asWbPeaZHcT0xMDI8++iiZMmXi0KFDbN++/YH7JJSGh5t5ItsTrG6zmplNZnLgcW/KvJKOkZUMN0O6Q5EiMGGChohSSeC11177W6+r0NBQwsPD8ff3x8/Pj4kTJwLQv39/Bg4cSOXKlbkZj04uefPmZejQoVSqVInatWtTtmzZB54jtho1ahAZGfm3B+ZxqVevHjdu3MDf35+33nqLihUrxvfHd5lOye7Gfrj8AyErQ1h6eCnBWXyZ9lUGSn21F3x8YOBA6NzZjmhXKgXSKdndg07J7oHyZs3L4paLmfvcXL67dZGyVSIZPLUNfxbIDz172pbI+PHaElFKJTsNDzcnIrQs1ZLIHpG0LNWSd6M/o0zLS2xZ+BE8+aQNkcKFYdw4SKIHY0opFZuGRwqRO3NuZjWbxZo2a/jz5p9U3f8q3QeUJGbNEihYEHr1gkKFYORIiMeDNaWUehgaHilM3SJ1ORBygH4V+zFl91RKHAhh0aS+sG4d+PnZGX2ffBKGDIGLF50uV6n78vRnru7uYa6/hkcKlDldZkbVHcWOl3fwryz/4rn5LWj208ecXfQJ7NgB1arBO+/YEOnXD86edbpkpf4hQ4YMXLx4UQPEIcYYLl68SIYMGRK0v/a2SuGu37zOR9s/YsjGIaRLk45htYbRPag7XpFR8N//wuefg5cXtG8P//d/dulcpdzA9evXiY6OTrJBbOrBMmTIgI+PD2ljzWQRn95WGh4e4vil43Rf2Z21362lkk8lJjeaTKk8peDECfscZNo0uHYNmjaF/v0hCft/K6VSNu2qm4oUzlGYr176ik+bfsrRS0cJnBTIm+ve5IrPY7Yn1qlT8OabsHEjVKoEVavC8uVwSxelUkq5TsPDg4gIbcu0JapnFG1Kt+H9b97Hf6I/60+sh3/9C959F06fhtGjbZg0bgylSsGMGTpWRCnlkhQZHiLSVESmiMhSEanjdD3uJlemXMxsOpO1bdcCUOvTWnRY0oELf1yALFmgTx84dgw++8xOuNipk+3uO3w4/PKLw9UrpVKCeIWHiGQXkQUickhEokSkUkJOJiLTReS8iByI47N6InJYRI6JyID7HccYs8QY0wXoALRMSC2pQa1CtdjXfR+Dqgzis/2f4TvWl+l7ptv109OmhRdfhD174MsvbQtk4EDInx9efdW2TJRS6h7i2/IYA6wxxhQHygBRd38oInlEJGus94rEcZyZQL3Yb4pIGmAc8CzgB7QWET8RKS0iK2J95blr13/f3k/dQ8a0GflPrf8Q0S2CkrlL0nlZZ6rNrMbB8wftBiJQpw589ZUNkqZNYexYO2r9xRdh925nfwCllFt6YHiIyCNAVWAagDHmmjEm9r2NasBSEclwe58uQGisbTDGbAYuxXGaYOCYMeY7Y8w1YC7QxBiz3xjTMNbXebH+C6w2xsT5201EGonI5JiYmAf9iKlCyTwl2dRhE9MbTyfqpygCJgUwcO1A/rj+x/82CgiAWbPgu++gb19YsQLKlYNnnrGDED28Z55SKv7i0/IoBPwEzBCRPSIyVUT+tkqJMWY+sAaYKyJtgE7ACy7U8Thw5q7X0bffu5feQG2ghYh0j2sDY8xyY0zXv1b1UvaBesfAjhzqdYh2/u0YvnU4JceXZOWRlX/fMH9+2733zBk7VuTAAahdG4KDYeFCXWtdKRWv8PAGygITjDGBwO/AP55JGGNGAFeBCUBjY8xvLtQR1zJb9/wz1xgTaowpZ4zpboyJewJ8dU+5MuViWpNpbOqwiUxpM9FwTkOaz2vO6ZjTf98wWzY7JuTECZg82T5Mb9ECSpSAqVO1h5ZSqVh8wiMaiDbG7Lj9egE2TP5GRJ4GSgGLgSEu1hEN5L/rtQ/wvYvHUC6q+mRV9nTbw7Baw1hzbA0lxpXgg60fcP3m9b9vmCEDdOkChw7BF19A1qz2dcGCMGIE6K1BpVKdB4aHMeYccEZE/prXohYQefc2IhIITAGaAB2BHCLyngt17ASKikhBEUkHtAKWubC/SqB0adIxoMoAIntGUrtQbfqv7U/gpEC2nNryz43TpIHnn4fwcPuAvWRJeOMNe5urf3+dQ0upVCS+va16A5+JyD4gAHg/1ueZgOeNMceNMbeA9sA/+nqKyBxgG+ArItEi0hnAGHMD6AV8ie3J9YUx5mBCfiCVMAWyF2Bpq6UsbbWU3679RtWZVemwpAM//f7TPzcWsQ/Rv/4adu2CBg1g1CjbEunUCSIj/7mPUsqj6NxW6h9+v/Y7721+j5HbRpI1XVaG1RpGl3Jd8JL7/K1x4oQNkOnT4coV+2zkrbfA3z/5CldKJQqd20olSOZ0mRlWexh7u+/F/1/+dF/ZnUrTKrHnhz333qlgQTs+5K85tL78EsqUgebNISIi+YpXSiULDQ91T365/djQfgOzms3i5C8nCZoSxCurXyHm6n0ekOfODe+9BydPwuDBsH49BAbawYc64FApj6Hhoe5LRHjJ/yUO9zpMSFAIY8PGUnxccebsn3P/RXxy5IC337Yh8vbbsGmTHXD47LP2ew+/XaqUp9PwUPGSPUN2xtYfS1iXMHwe8eHFRS9Se1ZtDl049IAds9sWyMmT8J//2NZH9epQuTIsW6ZTwiuVQml4KJcE5Qtie+ftjK8/nt0/7MZ/gj8D1w7k92u/33/HbNlg0CAbIuPGwQ8/QJMmULo0fPop3LiRLPUrpRKHhodyWRqvNISUD+Fwr8O08W/D8K3DKTGuBIuiFj14PeqMGaFHDzh6FGbP/t8SuSVLwty52hJRKoXQ8FAJlidzHmY0mcGWjlt4NOOjPPfFc9T/vD7HLh178M7e3tCmDezbB4sW2SniW7e2kzMuXarPRJRycxoe6qFVeaIKu7ruYnTd0Ww9vZWS40syeMNgrly/8uCdRaBZM9i7Fz7/HK5etT2zKlSANWs0RJRyUxoeKlF4e3nTp2IfDvc6TAu/Fry7+V38xvux/PDy+B0gTRrb8oiMtAMNz5+3PbPKl4fFi/V2llJuRsNDJaq8WfPyWfPP2NB+A5nSZqLx3MY0ntOYEz+fiN8BvL2hY0c4csTO3BsTYwca+vvbZXP1wbpSbkHDQyWJ6gWqE9Etgg+e+YD1J9bjN96Pdza9w9UbV+N3gHTpoHNniIqyt7NE4KWXwNfXTg+v08Er5SgND5Vk0qZJy/899X8c6nWIJr5NGLJxCH7j/FhyaMmDe2X9xdvb3s7auxeWLIGcOaFbNyhUyM6l9Zsry8YopRKLhodKcj6P+DC3xVzWtVtHprSZaDavGXVn1yXqp6j4H8TLy44L2bED1q6F4sXh//4PnnzSjmC/FNfqxkqppKLhoZJNzYI1iegewcfPfszO73fiP9GfV9e8yi9Xf4n/QUSgVi27pvq2bVClCgwdaidmfOcd+PXXJKtfKfU/Gh4qWXl7edMruBdHeh2hc2BnxuwYg+9YXxZGLnT9YBUr2jEhe/dCzZowZIi9nfXBB/DHH4lfvFLqDg0P5YjcmXMzseFEdnXdRf5H8tNifgtaLWjFhT8uuH4wf3/bnTcsDIKC7KqGhQtDaKhdW0Qpleg0PJSjAvMGsq3zNt6r8R6LohZRcnxJFkctTtjBype3Aws3b4ZixaBPHxsiY8ZoiCiVyDQ8lOPSpknLm1XfJLxrOPmy5qP5F81ps6gNF/+4mLADPv20nfZ9wwYbIn372ttZo0driCiVSDQ8lNvw/5c/YS+HMbTaUL44+AV+4/34fP/n8e/WG1v16rBxo/0qUQJefdWGyJgxdhoUpVSCaXgot5I2TVqGVB/Czi47eTLbk7RZ1IY6s+vEb7LFe6lWza5ouGmTDZG+fe3trHHjdLChUgmk4aHcUsBjAWzrvI2xz44l7GwYpcaX4t1N7/LnjYf4ZV+1qg2R9ettePTqBUWKwIQJGiJKuUjDQ7mtNF5p6Bnck6ieUTQp3oTBGwcTMCmATSc3PdyBa9SwrZCvv4YnnrDrixQpAuPHa4goFU8aHsrt5cuaj3kt5rG6zWr+vPEn1T+pTselHRPWrfcvIlC7NnzzjQ2RJ5+Enj1ti2TsWH2wrtQDaHioFKNekXoc6HGAgVUGMnvfbIqPLc4nEZ8k/IE6/C9Etmyx054ULAi9e9v/jhihI9aVugcND5WiZEqbifdrvc+ebnvwzeVLh6UdqPlpTQ5fOPxwB/5r2pPNm20XX39/eOMN2yIZPBguPEQrRykPpOGhUqRSeUqxpeMWJjWcRMS5CPwn+vPW+rfit3rh/YjYLr5ffWVHrNeoAe++C/nzQ5cudtlcpZTnhoeINBKRyTExMU6XopKIl3jRtVxXDvU8RMuSLXlvy3uUHF+SVUdXJc4Jype366sfPAjt29vFqMqUseGycKEuTKVSNY8ND2PMcmNM12zZsjldikpi/8ryLz5t9ikb2m8gg3cGGnzegObzmnMm5kzinMDPDyZOhOhoO+niyZPQooUdcDhsmN7SUqmSx4aHSn2qF6hORPcI3q/5PmuOraHEuBJ8sPUDrt+8njgnyJHDriFy/LhdmKpYMRg0CHx87NK5ektLpSIaHsqjpEuTjoFPDySyZyS1CtWi/9r+BE4KZPOpzYl3kjRp7MJUa9faW1qdOsH8+faWVpMm9lmJUh5Ow0N5pALZC7C01VKWtlrKb9d+o9rMarRb3I4ff/sxcU/k52cHF54+bRel2rIFKlSAZ56x3yvloTQ8lEdr7NuYyJ6RvPn0m8w9MJcS40owfc/0hxsbEpccOexiVKdO2fEh+/fb6VCeecaueKiUh9HwUB4vU9pMvFfzPfaF7KNUnlJ0XtaZ2rNqP9xki/eSNSu8/jqcOAEffmifgzz1FDz7LGzfnvjnU8ohGh4q1SieqzgbO2xkYoOJhH8fTukJpRn+zfDEe6B+t4wZ7RTw331nWyI7d0KlSlC5su3me/Nm4p9TqWSk4aFSFS/xoltQN6J6RlG/aH0GrhtI0JQgws4m0UPuzJltS+TkSbss7rlztptv0aJ2XZHLl5PmvEolMQ0PlSrly5qPhS8sZHHLxVz44wIVp1ak75q+XP4ziX6ZZ8li58w6csS2PPLmteuK+PjYcDl9OmnOq1QS0fBQqVrT4k2J7BFJSFAIoTtCKTm+JCuPrEy6E6ZJA82bw9at9hnIs8/CRx/ZAYetWkF4eNKdW6lEpOGhUr1sGbIxrsE4vun0DVnTZ6XhnIa0WtAq8bv1xlahAsyda5+L9O0Lq1bZKVGqVYNly+DWraQ9v1IPQcNDqdueyv8Ue7rt4e3qb7P40GJKjCvBzIiZid+tN7YnnoCRI+30J6NG2ecjTZrYJXMnTdK1RZRb0vBQ6i7p0qRjcLXBRHSLwC+3Hx2XduSZWc9w/NLxpD/5I49Av352+pM5c+zr7t1tuAwZAt9/n/Q1KBVPGh5KxaFE7hJs7riZ8fXH2zXUJ5RKum69sXl72+cfYWF2udxKleCdd2yIPPecXflQb2kph2l4KHUPXuJFSPkQIntG3unWW3ZyWb49823yFCBiR6kvWwZHj9pWyaZNUKcOFC8O48bB778nTy1KxaLhodQD+Dziw8IXFrK01VJirsZQeXplui3vxi9Xf0m+IooUsYMNz56F2bPtdCi9etmuvm+8obe0VLLT8FAqnv6aJ6tfxX5M3TMVv3F+LIpalLxFpE8PbdrYbr7ffmvXXx850q65HhJiH7YrlQw0PJRyQZZ0WRhVdxRhL4fxWJbHeO6L52g2rxlnfz2b/MVUqmSngj961K4nMn26baF06GAfuiuVhDQ8lEqAcvnKEdYljBG1R7Dm2Br8xvsxYecEbhkHHmQXKmRXOvzuO3sra948+0ykSxc7y69SSUDDQ6kE8vby5vXKr3Mg5ADl85Wnx6oeVJlehQPnDzhT0OOPw+jRttUREgKffmrn0GrfHnbvdqYm5bE0PJR6SIVzFObrtl/zSdNPOHLxCIGTAhm0bhBXrjs0uC9fPjsJ47FjdpzIwoVQrtz/Rq4n9aBHlSpoeCiVCESEdmXacajXIdqUbsOwb4YROCmQXd/vcq6o/PltiPw1cv3UKTtyvUwZe2tLp4VXD0HDQ6lElCtTLmY2ncnXbb/mt2u/UXFaRYZtGcbNWw7+os6e3Y4ROXbM3sq6ft0OQixVChYv1paIShAND6WSQO1CtdkXso/mJZozaP0ganxSg5O/nHS2KG9vaNsWDhyAL76w7zVvDlWq2Fl+lXKBhodSSSRHxhzMfW4unzT9hIhzEZSeUJqPd3zsbCsE7LTwzz9v11mfPNkumVulCtSsCV99pS0RFS8aHkolob+ehewL2Ufl/JV5Zc0rVJnhYI+su3l72+68x47ZZyJHjkDduhAUBF9+6XR1ys1peCiVDApkL8DqNquZ3Ww2xy4do+yksgzeMJg/b/zpdGmQKdP/ZvOdNg1++QXq1bNBsm+f09UpN6XhoVQyERHa+LchqmcULUu15N3N7xI0JYjw791k9cD06aFTJ4iKsqsbhodDQIAdsa6DDVUsGh5KJbNcmXIxq9ksVrRewaUrl6g4tSKD1g3i6o2rTpdmpUtnVzY8dgxee82udlismG2dXLrkdHXKTWh4KOWQBsUacLDHQdqXac+wb4YRNDmIPT/scbqs/3n0UfjgAzt3Vps2MGaMnTvr449td1+Vqml4KOWg7BmyM63JNFa3Wc2lK5eoMLUCw78Z7nyPrLvlz28nXYyIgLJl4ZVX7EBDHa2eqml4KOUG6hWpx/6Q/TQp3oSB6wZSbWa15Fn61hWlS9tVDJcsgRs37Gj1ypXtAlUq1dHwUMpN5MyUky9afMGsZrPYf34//hP9Cd0R6sxMvfciYkPj4EGYNMk+SK9eHRo3hkOHnK5OJSMND6XciIjwkv9LHAg5QLUnq9FnTR+qzqjKkYtHnC7t79Kmha5d7UP1YcNg40Y73UnPnna1Q+XxNDyUckP5s+Vn5YsrmdlkJgd/OkiZiWUYtmUY12+62YPqjBlhwAAbIt262dZIoUJ2Nt8TJ5yuTiWhFBkeItJURKaIyFIRqeN0PUolBRGhfUB7IntEUr9ofQatH0S5yeXYHr3d6dL+KU8eGDfOhkinTjBjBvj62ofrP/3kdHUqCcQ7PEQkjYjsEZEVCT2ZiEwXkfMi8o+5GUSknogcFpFjIjLgfscxxiwxxnQBOgAtE1qPUilB3qx5WfjCQpa2WsrPV3/mqWlP0XtVb/64/ofTpf1TgQIwYYJd1bBTJxg/HgoXhrffhl9/dbo6lYhcaXn0AaLi+kBE8ohI1ljvFYlj05lAvTj2TwOMA54F/IDWIuInIqVFZEWsrzx37frv2/sp5fEa+zYmskckvYJ7MXbnWAInBRJ2NszpsuL2+ON2adwDB6B2bRg6FAoWhOHD4bffnK5OJYJ4hYeI+AANgKn32KQasFREMtzevgsQGnsjY8xmIK4hqsHAMWPMd8aYa8BcoIkxZr8xpmGsr/Ni/RdYbYyJc31NEWkkIpNjYmLi8yMqlSJkTZ+V0GdDWdduHVeuX+GpaU/x9sa3uXHrhtOlxa14cVi0yE51UrEiDBwITzxhw+TiRaerUw8hvi2P0UB/IM4+g8aY+cAaYK6ItAE6AS+4UMfjwJm7Xkfffu9eegO1gRYi0v0eNS03xnTNli2bC2UolTLULFiTfSH7aFWqFUM3DaXC1ArsPbfX6bLurVw5WLkSduywy+G+/ba9xTVsGFx1k2lZlEseGB4i0hA4b4y573qaxpgRwFVgAtDYGONK21TiOuR9zhVqjClnjOlujJnownmU8hjZM2RndvPZLHh+AdG/RhM0JVKCu/UAABafSURBVMh9Zuq9l+Bgu3rh/v1QqxYMGgR+fvDZZzrlSQoTn5ZHZaCxiJzE3k6qKSKzY28kIk8DpYDFwBAX64gG8t/12gf43sVjKJUqPef3HJE9ImldqjXvbn6XspPLsiN6h9Nl3V+pUnak+tq1kCULvPSSfbD+wQfw449OV6fi4YHhYYwZaIzxMcYUAFoB640xL929jYgEAlOAJkBHIIeIvOdCHTuBoiJSUETS3T7PMhf2VypVy5kpJ582+5SVL67k1z9/5anpT/Hal6+5Z4+su9WqZefMWrHCTrrYvz/ky2fXEpk1C65ccbpCdQ+JNc4jE/C8Mea4MeYW0B74xwIAIjIH2Ab4iki0iHQGMMbcAHoBX2J7dH1hjDmYSLUplWrUL1qfgz0O0qVsFz7c/iH+E/zZdNLN557y8oIGDWD9ets7a+BAO5Nvu3a219arr9r3lVsR4+GzYgYFBZnwcDdZbEepZLTx5EY6L+vMdz9/R58KfXi/1vtkSpvJ6bLixxg75cmkSba31vXrdnncDh2gVSvImdPpCj2aiOwyxgTdb5sUOcJcKfVg1QtUZ1/3ffQs35MxO8YQMDGArae3Ol1W/IhAjRp2IaqzZ2H0aLh2DXr1grx5oXlzGy4e/sevO9PwUMqDZU6XmbH1x7Ku3Tqu3bzG0zOepufKnvz6Zwoa7Z07N/TpA3v32ucjvXvDN9/YcAkOhjlz4E837mHmoTQ8lEoFahasyYEeB3ilwitMCJ+A3zg/VhxJ8ExDzilTBkaNslPBT5wIMTHw4ot2waoBA+xU8SpZaHgolUpkSZeF0fVGs63zNnJkzEGjOY147cvX3G+m3vjImNHO4nvoEKxZYxelGjnSdgEODLTB8vvvTlfp0TQ8lEplKvhUYGeXnfQs35MPt39IjU9qcPbXFLoGh5eX7da7eDF8/z2Ehtr3QkLAx8d2/T192ukqPZKGh1KpUHrv9IytP5Y5z80h4lwEAZMCUuZtrLvlyWOfh4SH22cizzwDH35o1xdp1QoiI52u0KNoeCiVirUq1YrwruHky5qPRnMa0WNlD/cfWPggIvY21hdf2Knh+/Wz82qVKgVt2ugiVYlEw0OpVK54ruKEvRzGa5VeY0L4BMpNLse3Z751uqzE8cQTMGKEDYz+/e3treLF4fXX4VJcE3yr+NLwUEqR3js9I+uMZG3btfx+7XcqT69Mt+Xd+PnKz06Xljhy5bJriRw5YntnjRoFTz5pJ2bUqeETRMNDKXVHrUK1iOwZSb+K/Zi2ZxrFxxVn1t5ZeMxMFD4+doncffugfn0bKAULwpAh8MsvTleXomh4KKX+Jku6LIyqO4rwruEUzF6QdkvaUXVmVfb9uM/p0hJPqVIwb56dGr5OHXjnHXuLq39/O6JdPZCGh1IqTgGPBfBt52+Z2mgqUT9FUXZSWd74+g2uXPegmW5LloQFC2DPHjs546hRtiXSvr1tnah70vBQSt2Tl3jRuWxnjvQ+QseAjoz4dgQBk1LQHFnxFRBgpzk5etSOEVm40I5mr1kTVq92ujq3pOGhlHqgHBlzMKXxFL5u+/WdObLaL2nP6RgPG4BXqBCMGWMHFv73v3DsmH020ry5HYSo7tDwUErFW+1Ctdkfsp/Xn3qdeQfmUezjYrz+1eue0yvrLzly2Ocfx4/bh+qrV0PRoraL77lzTlfnFjQ8lFIuyZIuC/995r8c6X2E1qVbM2rbKIp8XIQx28dw7eY1p8tLXGnTwhtv2MWomjWzI9YLFIDOne3D9lRMw0MplSBPZHuCGU1msKfbHsrmLUvfL/tScnxJlh5a6jlde/9SuDDMnm0nYuzY0T4f8feHKlXgs8/sYlWpjIaHUuqhlHmsDF+99BUrX1xJWq+0NJ3XlGdmPeNZXXv/UrQoTJgAZ87YWXx//BFeesn20Bo2DM6fd7rCZKPhoZR6aCJC/aL12dt9Lx8/+zF7zu0hcFIgLy97me8ve+CD5pw54bXX4PBhWLXKTnkyaJAdhNi6tV2P/dYtp6tMUhoeSqlEkzZNWnoF9+Jo76P0qdCHT/d+StGPi/Le5vc873kI2Onfn30W1q61s/aGhNj1RWrVgtKl4dNPPfaWloaHUirR5ciYgw/rfkhUzyjqF63PWxveotzkcoSdDXO6tKRTooTt5vv99zY00qSxgw3/WuXQw7r6angopZJM4RyFmf/8fJa1WsbPV36m0rRK9F7Vm1+uevA8UhkzQtu2ds31VaugYkX7fKRwYTs9vIcslavhoZRKco18G3Gwx0FCgkIYHz4e37G+fBLxief1yrqbiL2ltWSJHbneqpVtmZQqZW9pjR8Ply87XWWCaXgopZJFtgzZGFt/LOFdwin8aGE6LO1As3nNuPDHBadLS3oFC9rZfM+ehbFjIX166NkTHn8cevWyXYBTGA0PpVSyCswbyDedvuHDOh+y6ugq/Cf48/Xxr50uK3k89pgNjZ07Yft2aNoUpkyxz0vq14d16yCFtMY0PJRSyc5LvHi10quEdQkjW4Zs1Jldh16revH7td+dLi15iECFCvbB+pkzdkr43buhdm37/oYNTlf4QBoeSinHBDwWwO6uu3m14quM3zmeMhPLsPPsTqfLSl558sBbb8HJkzBpkh14WLMmNGzo1g/XNTyUUo7KmDYjH9b9kI0dNnLj1g2qzKjCuLBxnv0wPS4ZMkDXrvb5x/DhsGWLnQKlWTNYscLtBh1qeCil3ELVJ6uyu9tunin0DL1W96LlgpaeN1tvfGTMaCdjPH7czuy7bRs0agRBQbbn1s2bTlcIaHgopdxIjow5WNZ6GcNrDWfxocX4T/Rnwwn3v/+fJHLlsvNlnTkDn3wCMTG2FeLrC1OnOj5yXcNDKeVWvMSLN6q8wbbO28iUNhO1Pq1F/6/7e+b0JvGRNi20a2fn0Zo3Dx59FLp0sfNpzZgBN244UpaGh1LKLQXlC2J31910K9eND779gErTKnH4wmGny3KOtze88AKEhcHy5ZA9O3TqZFsi06cne4hoeCil3FbmdJmZ0HACS1ou4dQvpyg3uRyf7//c6bKcJWJ7YoWHw7JltiXSubMdub50abKNE9HwUEq5vSbFm7AvZB9l85alzaI2dFvejSvXrzhdlrNE7IP0nTvtg3QRO+iwbl3YujXJQ0TDQymVIuTLmo/17dczoPIAJu+eTPkp5dn/Y+peChawodGkiV0WNzTUtkiqVIGyZeG335LstBoeSqkUw9vLm2G1h7GmzRou/HGB8lPK8/GOj1PfmJC4eHtD7962d9akSRAcDFmyJNnpxNMvelBQkAkPD3e6DKVUIjv/+3k6Le3EyqMrqVO4DjOazCBf1nxOl+URRGSXMSbofttoy0MplSLlyZyH5a2XM6HBBLac2kLpCaWZvW+2tkKSiYaHUirFEhG6B3VnT7c9+Ob0pe3itjSc05Czv551ujSPp+GhlErxfHP5sqXjFj6q+xEbTmyg9ITSzDswz+myPJqGh1LKI6TxSkPfin2J6B5B0ZxFabWwFS0XtOTH3350ujSPlCLDQ0SaisgUEVkqInWcrkcp5T6K5SzG1k5bea/Geyw5tITi44ozfc90fRaSyB4YHiKSQUTCRGSviBwUkbcTejIRmS4i50XkQByf1RORwyJyTEQG3O84xpglxpguQAegZULrUUp5Jm8vb96s+ib7uu/D/1/+dF7Wmfqf1yf612inS/MY8Wl5/AnUNMaUAQKAeiJS8e4NRCSPiGSN9V6ROI41E6gX+00RSQOMA54F/IDWIuInIqVFZEWsrzx37frv2/sppdQ/+ObyZUP7DXz87MdsPrWZkuNLMmLrCB2dnggeGB7G+muYYtrbX7Hbf9WApSKSAUBEugChcRxrM3ApjtMEA8eMMd8ZY64Bc4Emxpj9xpiGsb7Oi/VfYLUxZnd8f1ilVOrjJV70Cu7F3u57qfJEFd5Y+wbFxhZj/sH5eivrIcTrmYeIpBGRCOA88LUxZsfdnxtj5gNrgLki0gboBLzgQh2PA2fueh19+7176Q3UBlqISPd71NxIRCbHxMS4UIZSylMVyVGElS+uZEP7DeTOlJsXFrxAg88bcOLnE06XliLFKzyMMTeNMQGADxAsIqXi2GYEcBWYADS+q7USHxLXae9TT6gxppwxprsxZuI9tllujOmaLVs2F8pQSnm66gWqE9YljNF1R7Pl9BZKji/J+1veT73rhSSQS72tjDG/ABuJ+7nF00ApYDEwxMU6ooH8d732Ab538RhKKRUv3l7e9KnYh6ieUdQvWp83179JwMQAvjn9jdOlpRjx6W2VW0Sy3/4+I/Z20aFY2wQCU4AmQEcgh4i850IdO4GiIlJQRNIBrYBlLuyvlFIu83nEhwUvLGBF6xX8cf0Pnp7xNN1XdOfyn5edLs3txaflkRfYICL7sL/kvzbGrIi1TSbgeWPMcWPMLaA9cCr2gURkDrAN8BWRaBHpDGCMuQH0Ar4EooAvjDEHE/pDKaWUKxoUa8DBHgfpV7Efk3dNpszEMmw+tdnpstyazqqrlFJ3+eb0N7Rf0p7vfv6OnuV7MqzWMLKmz/rgHT2IzqqrlFIuqvJEFfZ238srwa8wfud4/Cf6sz16u9NluR0ND6WUiiVLuiyMeXYMWzpuAeDpGU8zbMswbty64XBl7kPDQyml7qHyE5XZ020PzYo3Y9D6QQRPCWb3DzouGTQ8lFLqvrJnyM68FvOY//x8fvjtB4KnBPPv9f/mzxt/Ol2aozQ8lFLqAUSEFn4tiOoZRdsybfnPlv9Qfkp5Is5FOF2aYzQ8lFIqnrJnyM6MJjNY1moZP/3xE8FTgnl/y/vcvHXT6dKSnYaHUkq5qJFvIw6EHKBZiWa8uf5Nnpn1DN9fTl2TYmh4KKVUAuTMlJO5z81leuPp7Di7g5LjSzJh54RU0wrR8FBKqQQSEToGdmR3190EPhZIj1U9CJ4azI7oHQ/eOYXT8FBKqYfkm8uXde3WMee5OZz77RyVplWiy7IuXPjjgtOlJRkND6WUSgQiQqtSrTjU8xCvVnyVGREzKPZxMcbvHO+Rt7I0PJRSKhFlTZ+VUXVHsbf7XgIeC6Dnqp6Um1yOXd/vcrq0RKXhoZRSSaBknpKsa7eO+c/P58IfF6g0rRKjt4/mlrnldGmJQsNDKaWSyF+DC/d238uzRZ/l1S9fpe7supyOOe10aQ9Nw0MppZJYzkw5WdJyCRMbTGTbmW34jfPj3+v/zS9Xf3G6tATT8FBKqWQgInQL6sb+kP00KNaA/2z5D4VDCzMubFyKnK1Xw0MppZJRwUcLMq/FPPZ020OZf5Wh1+peBEwM4MtjXzpdmks0PJRSygEBjwWwrt06Fr2wiKs3rlLvs3o0n9ec6F+jnS4tXjQ8lFLKISJCsxLNONjjIMNqDWPNsTX4jfNj9r7ZTpf2QBoeSinlsPTe6RlQZQAHehygzGNlaLu4La0XtubnKz87Xdo9aXgopZSbKPRoITa038C7Nd5lQeQCSk8ozbLDy5wuK04aHkop5Ua8vbz5d9V/s63zNrJnyE6TuU1o+HlDTv5y0unS/kbDQyml3FBQviD2dNvDqDqj2HRqEyXHl+TDbR9y/eZ1p0sDNDyUUsptpU2Tln6V+nGwx0FqFKjBa1+9RuCkQLZHb3e6NA0PpZRyd09ke4LlrZezpOUSLl+7TOXplRm8YbCjs/VqeCilVAogIjQp3oT9Iftp69+Wdze/y/Pzn+eP6384Uo+Gh1JKpSCPpH+EmU1nMrruaJYcWkL1mdX54fIPyV6HhodSSqVAfSr2YXHLxRz86SDBU4OJ/CkyWc+v4aGUUilUk+JN2NppKzdu3aDazGpEnItItnNreCilVAoW8FgAmztsJqN3RqrNrMbqo6uT5bwaHkoplcIVzVmUbzt/S+FHC9NwTkM+2vYRxpgkPaeGh1JKeQCfR3zY0nELTYs3pd9X/Xh52cv8eePPJDufhodSSnmIzOkyM//5+bxV9S12n9vN9VtJNxpdkrpp47SgoCATHh7udBlKKZWsrt64SgbvDAnaV0R2GWOC7reNtjyUUsoDJTQ44kvDQymllMs0PJRSSrlMw0MppZTLNDyUUkq5TMNDKaWUyzQ8lFJKuUzDQymllMs8fpCgiPwEnHK6DjeQDYhxugg3kFKvgzvW7URNyXHOpDpHYh43MY6VC7hwj8+eNMbkvt/OHh8eyhKRycaYrk7X4bSUeh3csW4nakqOcybVORLzuIlxLBEJf9Ao8vvR21apx3KnC3ATKfU6uGPdTtSUHOdMqnMk5nEd//9BWx5KKZUKactDKaVUQkx+mJ215aGUUspl2vJQSinlMg0PpZRSLtPwUPclIoVEZJqILHC6Fqek5GuQkmtPLHoNkoaGRwogIvlFZIOIRInIQRHp8xDHmi4i50XkQByf1RORwyJyTEQGABhjvjPGdH6Y+hODiGQQkTAR2Xv7Grz9EMdy5BqISBoR2SMiKx7iGCny3w9ARLKLyAIROXT7/+VKCTxOir0G7kpEmorIFBFZKiJ14rWTMUa/3PwLyAuUvf19VuAI4BdrmzxA1ljvFYnjWFWBssCBWO+nAY4DhYB0wN67zwEscPgaCJDl9vdpgR1AxZR0DYB+wOfAijg+c+vaE+nf8BPg5dvfpwOyp7ZrkMzXezpwPo7rVA84DBwDBsT67FFgWnyOry2PFMAY84MxZvft7y8DUcDjsTarBiwVkQwAItIFCI3jWJuBS3GcJhg4ZuxfadeAuUCTxPspHo6xfrv9Mu3tr9hdBd32GoiID9AAmHqPTdy29sQgIo9gf+lPAzDGXDPG/BJrM4++Bg6YiQ2KO0QkDTAOeBbwA1qLiN9dm/z79ucPpOGRwohIASAQ+5f3HcaY+cAaYK6ItAE6AS+4cOjHgTN3vY4GHheRnCIyEQgUkYEPUfpDu33bJwL719TXxpiUdA1GA/2BW3F96Oa1J4ZCwE/AjNu37qaKSOa7N0gF1yBZ3SNk4wxYsf4LrP7rD9UH8U7cclVSEpEswEKgrzHm19ifG2NGiMhcYAJQ+K6/1ON1+DjeM8aYi0D3BBWcyIwxN4EAEckOLBaRUsaYA7G2cbtrICINgfPGmF0iUv1e27lj7YnIG3urqbcxZoeIjAEGAG/dvZGHXwN3EFfAVgB6A7WBbCJSxBgz8UEH0pZHCiEiabHB8ZkxZtE9tnkaKAUsBoa4eIpoIP9dr32A7xNQapK7fbtjI7Ga5OC216Ay0FhETmL/0qspIrNjb+SmtSeWaCD6rtbiAmyY/I2HXwN3cK+ADTXGlDPGdI9PcICGR4ogIoK9VxxljPnwHtsEAlOw93g7AjlE5D0XTrMTKCoiBUUkHdAKWPZwlSceEcl9u8WBiGTE/pV0KNY2bnkNjDEDjTE+xpgCt4+53hjzUkqoPbEYY84BZ0TE9/ZbtYDIu7fx9GvgJhIvYJ3uEaBf8eo1UQX7cHgfEHH7q36sbSoDpe96nRboEsex5gA/ANdv/4/U+a7P6mN7ch0H3nT6545Vtz+w5/Y1OAAMjmMbt78GQHXi7m3l9rUnws8eAITf/jdcAjya2q6BA9e8AHf1tsLePvwOKMj/eqSVTMixdW4rpZTyQCIyB/vHSi7gR2CIMWaaiNTHduBIA0w3xvwnQcfX8FBKKeUqfeahlFLKZRoeSimlXKbhoZRSymUaHkoppVym4aGUUsplGh5KKaVcpuGhlFLKZRoeSimlXKbhoZRSymX/D/Rp71skngtVAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "r = np.linspace(10,100,100)\n",
    "n = 5\n",
    "m = 5\n",
    "\n",
    "cr_H2 = [Corr(human2,int(ri),human_alphabet2) for ri in r]\n",
    "cr_H2 = smooth(cr_H2,10)\n",
    "\n",
    "cr_N2 = [Corr(neanderthal2,int(ri),neanderthal_alphabet2) for ri in r]\n",
    "cr_N2 = smooth(cr_N2,10)\n",
    " \n",
    "\n",
    "gamma_H2 = scipy.optimize.curve_fit(scaling,r[n:-m],cr_H2[n:-m])[0]\n",
    "gamma_N2 = scipy.optimize.curve_fit(scaling,r[n:-m],cr_N2[n:-m])[0]\n",
    "        \n",
    "plt.plot(r[n:-m],cr_H2[n:-m], c = 'red', label = 'Human')\n",
    "#plt.plot(r[n:-m],scaling(r[n:-m],gamma_H2[0],gamma_H2[1]), c = 'black')\n",
    "\n",
    "plt.plot(r[n:-m],cr_N2[n:-m], c = 'green', label = 'Neanderthal')\n",
    "#plt.plot(r[n:-m],scaling(r[n:-m],gamma_N2[0],gamma_N2[1]), c = 'black')\n",
    "\n",
    "print('C[H]:',gamma_H2[1])\n",
    "print('G[H]:',gamma_H2[0])\n",
    "\n",
    "print('C[N]:',gamma_N2[1])\n",
    "print('G[N]:',gamma_N2[0])\n",
    "\n",
    "plt.yscale('log')\n",
    "plt.xscale('log')\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 121,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAqIAAAHFCAYAAAApGJuMAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd1yX5f7H8dcFiICCey8ER+6FmBvLMld5NBtq5TZbdrSysmGd0tM5TTPNsjLLzI6pZWpppoi591ZQ0VREHIkDBeH+/XEjPwciKnDz/fJ+Ph48hO/3Hh/I6s11X9fnMpZlISIiIiKS0zycLkBERERE8iYFURERERFxhIKoiIiIiDhCQVREREREHKEgKiIiIiKOUBAVEREREUd4OV2A3LjixYtbgYGBTpchIiIicl1r1649allWifTeUxB1QYGBgaxZs8bpMkRERESuyxiz71rv6dG8iIiIiDhCQVREREREHKEgKiIiIiKOUBAVEREREUcoiIqIiIiII7RqXkRERNLEx8dz5MgRkpKSnC5FXEC+fPkoWbIkAQEBN3W+gqiIiIgAdgiNjY2lXLly+Pr6YoxxuiTJxSzLIiEhgYMHDwLcVBjVo3kREREB4MiRI5QrVw4/Pz+FULkuYwx+fn6UK1eOI0eO3NQ1FERFREQEgKSkJHx9fZ0uQ1yMr6/vTU/lUBAVERGRNBoJlRt1K39nFERFRERExBEKoiIiIiJyTcYYvv3222y5toKoiIiIuIXevXvTtm3bTB//1ltvERgYmH0FyXUpiEr6zp1zugIREZFcITEx0aWvn5spiMrVNm6EwED4+muwLKerERERuWEXR0c/++wzKlWqREBAAPfddx9xcXEATJo0iVdffZV9+/ZhjMEYw8iRIwEIDAzklVde4YknnqBYsWI0b94cgJiYGB566CEKFy6Mr68vYWFhrFmz5rL7/v7779SpUwcfHx/q1q1LeHj4ZY+2o6OjMcYwZcoUOnToQIECBXj55ZcBiIqKolu3bhQuXJgiRYpw9913s3nz5rRrx8fH06dPH0qXLk3+/PmpUKECQ4cOTXt/6dKlNG/eHH9/f/z9/alXrx6//fZb2vuxsbH07t2bEiVK4O/vT/PmzVmyZMll9S9atIi6deum1b9o0aIs+ieSPgVRuZqvL1SpAr17Q7t2sHev0xWJiIjcsNWrV7No0SLmzJnDr7/+yoYNG3juuecAePDBBxk+fDjly5cnJiaGmJiYtPcAxowZQ8mSJVm+fDlff/01lmXRpUsXduzYwS+//MKqVasoVaoUd911F0ePHgXg4MGD3HvvvTRp0oR169bxwQcfXBYULzV8+HB69OjB5s2befLJJ4mNjaVFixaULFmSiIgIVqxYQfXq1QkLC0sLz6+88grr1q3jp59+IjIykmnTplGjRg0AkpOTL7v3unXrGDlyJH5+fgAkJCTQpk0bTp06xbx581i/fj0dOnTgrrvuYvv27QAcOnSITp060ahRI9atW8d7773HkCFDsucfTirtrCRXq1YNliyBCRNg+HCoXRv+9S945hnw0l8ZEZE85dlnYcOGnL9v/frw4Ye3dAlvb28mTZpE/vz5ARg8eDAfffQRYPe+LFiwIJ6enpQuXfqqcxs3bpw2QgqwcOFCVq1axdatW6lZsyYAkydPJjAwkHHjxvHaa68xbtw4SpYsyYQJE/D09KRmzZq8/fbbtG/f/qrrDxo0iF69eqV9PXLkSAIDAxk/fnzaa2PGjGHu3LlMmTKFZ599ln379tGgQQOaNGkCQMWKFWnWrBlgj5aeOHGCe++9l6pVqwKk/Qkwbdo04uPjmTZtGl6p/y8fMWIECxcuZMKECXz44YeMGzeO4sWL8/nnn+Pl5UXNmjUZNWoUnTt3vvEffiZpRFTS5+EBgwfD1q1wxx0wbBg0bWo/thcREXEBNWrUSAuhAOXKlSM2NjZT54aGhl729datWylWrFhaCAXInz8/TZo0YevWrQBs27aNxo0b4+npmXZM06ZNM3X91atXs3btWgoWLJj24e/vT3R0NJGRkQA88cQTTJ8+ndq1azNkyBDmzZtHSkoKAEWKFKF///60a9eO9u3b8+9//5udO3dedv3Dhw9TuHDhy+4RERGRdv1t27YRGhqaFlQBWrRokamf183S8JZkrEIF+Pln+N//4OmnISQEXngBXn0VfHycrk5ERLLbLY5KOsnb2/uyr40xWJlc+1CgQIGrXkuvcbtlWZe9fuUx12r2fuX1U1JSuPPOOxk7duxVxxYqVAiAdu3asX//fn777TcWL15Mr169qFOnDgsXLsTT05PPP/+cIUOGMH/+fBYsWMCrr77K2LFjGTRoECkpKdSoUYOZM2dedf2Lj++v/F4yqj+raERUrs8YeOAB2L4devWCUaOgXj0ID3e6MhERkZvm7e1NcnJypo6tVasWR48eZdu2bWmvnT9/nlWrVlGrVi0AatasyerVqy+75vLlyzN1/ZCQELZu3Uq5cuWoUqXKZR8lSpRIO65o0aI8/PDDTJgwgTlz5hAeHn5ZTbVr12bo0KHMmzePfv368dlnn6Vdf8+ePQQEBFx1/bJly6Z9jytXrrys/qVLl2aq/pulICpXSU5JYfHWQySnDvenKVoUvvoKFiyApCQIC4NBg+DkSUfqFBERuRWVK1fm8OHDLF++nKNHj3L27NlrHnvHHXcQGhpKjx49+PPPP9myZQuPPvoo586dY/DgwYD96Dw2NpbBgwezfft2Fi1axIgRI4Drjyw+9dRTJCcn06VLFyIiIoiOjmbp0qWMGDGCZcuWAfaczhkzZrBz504iIyOZMmUKBQsWpGLFikRFRTF8+HCWLl3Kvn37WL58OREREWlTCXr27EnlypXp2LEj8+fPJzo6mpUrVzJ69GhmzZoF2HNo4+LiGDhwINu3b2fhwoVp9WcXBVG5yopdRxg9Yz39xoUzf+NfVwfStm1h82Z47jmYOBFq1IDUv8QiIiKuokuXLnTv3p2OHTtSokQJ/vOf/1zzWGMMs2bN4rbbbqNjx440btyYw4cPs2DBAooXLw7Yc1B//vlnli1bRv369RkyZAhvvfUWAD7Xmc5WqlQpli9fTvHixenatSvVq1enZ8+e7Nu3jzJlyqRd47XXXqNRo0aEhISwadMm5s2bR6FChShQoACRkZE89NBDVKtWjW7dutGsWbO0R/0+Pj6Eh4cTEhJCnz59qFatGl27dmXVqlVUqlQprf7Zs2ezatWqtPrff//9W/45Z8Rkdq6E5B4hISHWlX3LspJlWSzfFcuUJZFEHY6nTBE/Hm5RhTvrlMPL84rfXdasgX79YNMm6NYNPv4YUv+FERER17J9+/a0dkCSNZYsWULr1q3ZtGkTderUcbqcbJPR3x1jzFrLskLSe08jonIVYwzNqpdmbP8WvPFgCAXye/H+7E30G7eYX9fv50LyJSOkISF2GB01Cn75xR4dnThRjfBFRCRPGj9+PMuWLSM6Opq5c+cyYMAAmjRp4tYh9FYoiMo1GWO4vVqptEAa4OvNB79spu+4xcxdt5+ki4E0Xz546SV7VLRePRgwwG75tHu3s9+AiIhIDtu3bx8PPfQQ1atXZ/DgwbRs2ZI5c+Y4XVaupUfzLii7H81fi2VZrI6K45slu9h16CSlCvnyUIsq3FWvPPkuPrJPSYEvvoDnn7c///JLuP/+HK9VRERunB7Ny83So3nJdsYYQquWZEzf5rz1cGOKFMzPR3M20/eTxfyydh+JF5LtRvgDBtiN72vUgO7d7f6j5887Xb6IiIjkMgqicsOMMTSuUpIP+zTj7R6hFPPPz8dzt9D3k8XMXpMaSCtVgogIe2u4sWOhRQvtWS8iIiKXURCVm2aMISS4BB/0bsaonqGUCPBl7Lwt9PlkMT+vjibRwxM++ABmzIDISGjQQG2eREREJI2CqAsxxnQ2xnx2Mpc1kDfG0CioBO/3bsronk0oVciXT37dysBPl7BiVyz84x+wbh1UqWJ/PnQoJCY6XbaIiIg4TEHUhViWNduyrIEX95zNbYwxNAwqznuPNWVUj1C8PAyvT1vDq9+v5mDhUvDnn/DUU/YoaatWsG+f0yWLiIiIgxREJcsZY2gUXILxg1oxoG0NNu87xqBPlzBpWTTn3vsAfvgBtm2zH9X/8ovT5YqIiIhDFEQl2+Tz9OD+pkF88UQYrWqWYerSKPqPDyeiVgustWvtBU2dO8MLL9h714uIiLiQ3r1707ZtW0drMMbw7bffOlrDrVAQlWxXzN+HF7rU573HmuLv681bP67jxZXH2Dd7ATz+OPz3vxAWBgcOOF2qiIi4sN69e2OMwRiDl5cXRYsWpWnTprzxxhscP37c6fIAaNu2Lb1793a6jFxDQVRyTO2KRRnbvzlPta9FVMxJBn+9igldn+bMt9/ZuzLVrw/z5jldpoiIuLCWLVsSExPD/v37iYiIYMCAAUydOpVatWqxa9cup8uTKyiISo7y9PCgc0ggXzwRxt31yjNzxV76HynB79P/wCpbFjp0gJdfhgsXnC5VRERckLe3N6VLl6Zs2bLUqlWLvn37smrVKvz8/Hj88cfTjvv++++pX78+Pj4+BAYGMnToUM6cOZP2flhYGH379uXFF1+kePHiBAQE0L9/fxISEjK8f0bX7d27NwsXLuTrr79OG7ldvHgxALGxsfTu3ZsSJUrg7+9P8+bNWbJkyWXXXrRoEXXr1sXHx4e6deuyaNGiLPqpOUdBVBxRuEB+nu1Ul4/6NadEgC//XXWEYc+OI2rwMBg92t6rfv9+p8sUERE3EBAQwODBg1m8eDFxcXFMmjSJwYMHM2zYMLZt28bkyZP5/fffLwuqANOnT+fYsWNEREQwZcoUfv75Z4YPH37N+1zvuh999BEtW7bkgQceICYmhpiYGJo1a0ZCQgJt2rTh1KlTzJs3j/Xr19OhQwfuuusutm/fDsChQ4fo1KkTjRo1Yt26dbz33nsMGTIk+35oOcTL6QIkb6tetjAf9m3Ggo0H+GLhDp4u24YOHzTlsbefJKBuXRg/Hh5+2OkyRUTyrPG/bWVPbHyO3zeoVACD29XKsuvVrl0by7LYu3cvI0eOZPTo0TzyyCP2vYKCGDt2LK1bt2bMmDEUKVIEgKJFi/Lpp5/i6elJjRo1eOutt3j66acZPXo0BQoUuOoembmut7c3vr6+lC5dOu28SZMmER8fz7Rp0/DysqPZiBEjWLhwIRMmTODDDz9k3LhxFC9enM8//xwvLy9q1qzJqFGj6Ny5c5b9jJygEVFxnIcxtKtfgS+eCKNzSCBzT/vQ97mv+LFDbxIffQx69oS//3a6TBERcWGWZQH2KvN9+/YxdOhQChYsmPbRvn17AKKiotLOCQ0NxdPTM+3r5s2bk5iYyO7du6+6flxcXKave6XVq1dz+PBhChcufNm5ERERREZGArBt2zZCQ0PTgipAixYtbuEnkjtoRFRyDX/ffDxxTy3uaVCBib9v57PEu5j1WisenTmeO+rVw3PyZGjd2ukyRUTylKwclXTSli1bMMZQsWJFwH5M3qZNm6uOK1++/DWvcTHMpiclJeWmr5uSkkKNGjWYOXPmVe/5+fml3dsYc9l7V37tihREJdcJKhXAqJ5NWLfnKF/+sYN3//Es008cou/A5wn9xx2YN98Eb2+nyxQRERcRHx/P+PHjufPOOylVqhQVKlRg586dDBgwIMPzVq9eTXJyctqo6PLly/H29iY4OPiqYzN7XW9vb5KTky97LSQkhMmTJxMQEEDJkiXTPa9WrVp88803l9WzdOnSDOt3BXo0L7lWw6DijOnXnJe7NuB8YBCv9Xqd5+MKs/2ebvbOTCIiIldITEzk8OHDxMTEsG3bNr788ktCQ0M5f/4848ePB+Dtt99mzJgxvPXWW2zZsoWdO3cya9YsBg0adNm1jh07xpNPPsn27duZM2cOr776KgMGDEh3fmhmr1u5cmXWrl3L7t27OXr0KElJSfTs2ZPKlSvTsWNH5s+fT3R0NCtXrmT06NHMmjULgMGDBxMXF8fAgQPZvn07CxcuZMSIEdn0U8w5GhGVXM3DGFrXKkvz20ozb/1+piyweDawDs3fnkqf0PJUeGYguMGjCRERyRoRERGUKVMGT09P/P39qV69Oj169ODpp59OW4T0yCOP4O/vzzvvvMOoUaPw8vIiKCiIrl27Xnat+++/H39/f1q0aEFiYiLdu3fnP//5zzXvnZnrDhs2jM2bN1OvXj3OnDnDokWLCAsLIzw8nFdeeYU+ffoQFxdHiRIlCA0N5Z577gGgXLlyzJ49m2effZb69etTtWpVxowZw5133pkNP8WcYzKa7yC5U0hIiLVmzRqny3BEQuIFZvy+mf+tiuY8ntwTt51eQx+mWHBFp0sTEXF527dvp0aNGk6XkSuEhYVRpUoVJk6c6HQpLiGjvzvGmLWWZYWk954ezYtL8fX2omeHBkx6rj33FjjD/BLV6fP1Gr789GdOn9N+9SIiIq5EQVRcUuGCPgx+/mEmdgyi+eEdTIvzpPe/f+HH8O0kXki+/gVERETEcZojKi6tTOO6DK93G93e+C9f/pXMZ57ezFgZTc+2tbi7Xnm8PPW7loiI3LiLW29K9tL/pcX1eXtT5e0RjOofxjtzP6RE1FY+mrOZ/uMW8/umAySnaB60iIhIbqQgKu6jVSvq/zqdD5K28q9vRlIgejf//Wkjj09YQsS2GFK0ME9ERCRX0aN5cS+FCmG++orQWbMIGTiIP8vUYHL3Z3jrx3UElwrg0bBqNKla0i12oxARyQ7p7eAjkpFb6cCkEVFxT1264LFlMy0rF+bT17rzwtY5nD17jtenreGfXy1j/d6jTlcoIpLr5MuXj4SEBKfLEBeTkJBAvnz5bupcBVFxXyVLwsyZeH75BXfO/ZaJIx9gSOF44uLP8eK3K3nhmxVs/eu401WKiOQaJUuW5ODBg5w9e/aWRrkkb7Asi7Nnz3Lw4MFrbk16PWpo74LyckP7mxYdDb17Q3g4iV26Mu/pkUzdEMuJM+dpXKUEj4VVp2qZQk5XKSLiuPj4eI4cOUJSknozy/Xly5ePkiVLEhAQcM1jMmporyDqghREb1JKCnz4Ibz8MgQEcG7C5/xcqjY/LNvNqYQkmt9WmsfCqlGphL/TlYqIiLgNBVE3oyB6i7ZuhUcegfXroXdvzrzzLjO2H2PGir2cS7rAnXXK06t1VUoX9nO6UhEREZenLT5FLlWrFqxYAa+8ApMnU6BxQx6xDjHp6TZ0vT2I8G2H6PfJYsbO28KxU+ecrlZERMRtaUTUBWlENAutWAGPPgqRkTBkCIwaxdELHkyJiOS3DX/h5WHoElqZ+5sFEeDr7XS1IiIiLkeP5t2MgmgWO3MGXnwRxo6FKlXgq6+gRQsOHj/Dt+G7WLTlEH75vbi/aRD/aFIZX2+13xUREcksPZoXyUiBAvDxx/DHH3DhArRqBcOGUc7Xg+H/aMD4gS2pW6kYXy/exWMfL2Lmyr0kXkh2umoRERGXpxFRF6QR0Wx0+jQMHw7jxkG1ajBpEjRtCsD2Ayf4atFONkYfo0SAD4+0rkbbuuXw9NDvcyIiIteiEVGRzCpYED75BBYuhPPnoUULeP55SEigRvki/OeR2/l3ryYULejD+7M30X98OAs2HiA5JcXpykVERFyORkRdkEZEc8ipU3YInTABqle3R0dvvx2wd5NYviuWb8Mj2R0bT5kifvRoWYU762iEVERE5FJarORmFERz2O+/Q79+cOAAPPccvPEG+PgAdiBdsesI3y7ZRdRhBVIREZErKYi6GQVRB8TH2yH088+hRg17dDQ0NO1ty7JYGXmEb8L/P5A+3MIOpF6eCqQiIpJ3KYi6GQVRB/32G/TvD4cOwQsvwMiRkD9/2tsXA+m3SyKJjDlJ6cK+9GhZVYFURETyLAVRN6Mg6rCTJ2HYMPjiC6hZ0x4dbdz4skMsy2JV1BG+DY9kV2ogfbhFFdrWLa9AKiIieYqCqJtREM0lfv0VBgywR0eff94eHU2dO3qRZVmsjorjmyW72HXoJKUK+/Jgs2Duqlceby9PZ+oWERHJQQqibkZBNBc5edKeOzpxItx2m70rU+rK+ktdDKRTIiLZcfBvivnnp3vTYNo3rIhPPgVSERFxXwqibkZBNBeaP9+eO3rwIAwdCm++Cb6+Vx1mWRYboo/xXUQkm/Ydp5CfN91ur0ynkEoUyJ/PgcJFRESyl4Kom1EQzaXi4+0FTBMm2LsyffUVNGt2zcO3/nWcqUujWB0VR0EfL+5rXJkuTQIJ8PXOwaJFRESyl4Kom1EQzeV+/90eHd2/H559Ft56C/z8rnl4ZMxJpkZE8ufOWHy9PenUqBLdbg+iSMH81zxHRETEVSiI5mLGmArAJKAskALMAYZbGfyDURB1AadOwYsv2nvWV6kCX34JLVtmeEr0kVN8/2cU4VsP4eXpQfsGFbm/aRAlC139iF9ERMRVKIjmYsaYMkA5y7LWGGO8gQXAGMuyfrzWOQqiLmTRIntXpuhoePppGDUKChTI8JSDx84wbVkUv286iAHa1i3PA82CKVcs4/NERERyo4yCqEs1NDTGeBljXjTGRBpjzhtjDhhjPsiG+1Qxxkwwxmw0xiQbYxZncGxNY8xCY8xZY8whY8ybxphML4O2LCvGsqw1qZ8nApuACrf8TUju0KYNbNoETz4JY8ZA3bp2OM1AuWIFGNq5Hl89GUb7hhVZuPkg/ccvZvSM9eyJjc+hwkVERLKfSwVR4CvgGeBd4G7gRSAhG+5TC+gA7Er9SJcxpgjwO2AB9wFvAsOAN27mpsaYYkAX4LebOV9yqYIF4eOPYfFiMAbuuAMGDbJbP2WgVGE/nmpfm8nPtOH+psGsijzC4M8ieO371Ww7cCJnahcREclGLvNo3hhzDzAbqGdZ1rZMntMAiLMs60A673UAfrMsKzmd9zwsy0pJ/Xw6UNyyrLB0jnsJeAGoZFlWfOprLwAjgdKXvLYUKJ9OiQsty+qXekx+4FfgF8uy3svo+9KjeRd29qzd+P6996B0afj0U+jcOVOnnkpIYvaaaGau3Et8QhL1AovxUPMqNKhcDGNM9tYtIiJyk9xijqgx5gegkGVZ7W7gnJ+AakBLy7KOXvL6o9ijq/dYlrXgOtfIKIguAQ5ZlvXQJa9VBPYB91qWNTuTdXoC04D9lmUNvd7xCqJuYM0a6NsXNm+Ghx6Cjz6CkiUzdWpC4gXmrdvP9BV7OHbqPNXLFuahFsHcXq0UHgqkIiKSy7jLHNEmwC5jzFhjTHzqnMwZxpiyGZzTG0gCfjXGBAAYY+4DvgBevl4IzYTbgB2XvmBZ1n7gbOp7mTUBOIX9WP+ajDGdjTGfnbzOI11xASEhdhh980348Ud7z/opUyATvxj6envR9fYgJj3VhiEd6xCfkMgbP6zl8QlL+GPzQZJTXOOXSxEREVcKoqWxg2V94CGgD9AImGmu8VzSsqwTQDugCPCzMaY99sjj+5ZlvZMFNRUB/k7n9ROp712XMaY50A8IAdYbYzYYY55J71jLsmZbljWwUKFCN1uv5Cbe3vDqq7B+PVStCr162Y/p//orc6d7edKhYUW+eKI1w7vUx2B4Z9YGBo4PZ8HGAySnpGTzNyAiInJrvJwu4AaY1I/7LMs6BmCMiQHCgTuAhemdZFlWjDGmLfAnMBeYaFnW8CysK73hJ3ON19Or78/U4yWvqlULli6FsWPh5Zftr//zHxg4EDyu/7uip4cHd9QpR1jtsizbcZgpEVG8+/NGpkRE8nCLKtxZpxxenq70O6eIiOQVrvR/pxPA5oshNNVSIBGoeZ1zfQFv7IbxBY0xWfV9nwAKp/N6IdIfKRVJn6cnDBlizxkNDYXBg+3WT5GRmb6EhzG0qFGGcQNaMPKBEAr65OP92Zvo+8li5qzdR+KFq9bliYiIOMqVguj2a7xusANm+m8aEwjMB1YCrYBOwCdZVNMOrpgLmrpTUgGumDsqkilBQbBgAXzxBWzcaPcdfecdSErK9CWMMTStXoqP+zXnrYcbU6RgfsbM3UKfTxbz8+poBVIREck1XCmI/gLUNcYUv+S1VkA+YGN6JxhjSmHvVLQbuD/1Mfh9QB9jzKgsqGke0M4Y43/Jaw9i9zYNz4LrS15kjL2ifts2aN/e3io0NBTWrbvByxgaVynJh32aMapnKKUK+fLJr1t57ONFzFi5l3NJCqQiIuIsV2rfFABsAQ4CowB/4B1gh2VZd13jnF+wFzndcbGnZ+rr9wHTgU6WZV3VPN4Y44fd0B7slewBwOupX8+1LOts6nFFgG2pdb0DBAHvAx9alvXKLX3DGVD7pjxmxgx7Z6a4OBg61O5D6ud3w5exLIuN+44xZUkkm/Ydp5CfN91uD6JzSCX88rvSdHEREXElbtFHFOytN4ExQGvsuaE/Af9MXR2f3vHVgOOX9hC95L1mwHIrnR9A6uP8vdcoo7JlWdGXHFsTGAs0xZ4XOhEYmV6j/KyiIJoH/f03vPACfP45BAfDZ5/ZOzTdpM37jjF1aRRr9xyloE8+/hEayH2hlfH3zZeFRYuIiLhREBWbgmgetmiRvZo+Ksp+fP/uu1AkU53C0rXj4N9MXRrFil2x+Hl70blxJbo2qUzhAvmzsGgREcnLFETdjIJoHpeQYDfC/+9/oXhxex/7+++355bepN2H4/n+zygitsXg7eVBx0aVuL9pEMX8fbKwcBERyYsURN2MgqgAsGED9OtnL2K67z745BMoV+6WLrn/6Gmm/RnFH5sP4elhaFe/PA80C6ZU4RufkyoiIgIKom5HQVTSXLgAH34Ir70G+fLZrZ4y2Qg/IzEnzvLDst3M3/AXFtCufgV6tqxK8QCNkIqIyI1REHUzCqJyld27YdAgWLgQWra0FzVVr37Ll42LT2Dan7uZt24/xhg6h1TiwebBmkMqIiKZpiDqZhREJV2WBV9/bbd4OnPGHiV9/nl7T/tbdPjEWb6NiGThpgN4e3nStUllujUNoqCPVtmLiEjGFETdjIKoZCg21t4udNo0qFMHJk60G+Jngf1HT/NN+C6WbIuhoI8X3ZsGc19oIL7e6kMqInsyR6MAACAASURBVCLpyyiIutLOSiKSGaVKwfffw08/wfHjcPvt8M9/wunTt3zpisULMqJbQ8YNaEGtCkX5atFOeo+1d2rS1qEiInKjNCLqgjQiKpkWHw8vvQTjxkGlSvDpp3DPPVl2+W0HTvD1op1siD5G8QAfHmoeTLv6FfD28syye4iIiGvTo3k3oyAqN2zpUhgwAHbsgF694IMP7B6kWWTD3qN8vXgX2w6coLi/Dw80C6J9w4oKpCIioiDqbhRE5aacPw+jRsHo0VCokB1Ge/a8pUb4l7Isiw3Rx/h2SSRb9h+naMH8dG8WTIeGFfHJp0AqIpJXKYi6GQVRuSVbtkD//rByJbRrB+PHQ+XKWXqLjdHHmBIRycboYxQu4M39TYPo3KgSPlrUJCKS5yiIuhkFUbllycl2AH3pJfvzN9+EZ58Fr6wNipv3H2fKkkjW7z1KIT9vut0eROeQSvjlVyAVEckrFETdjIKoZJm//oInn4TZs6FBA7vVU8OGWX6brX8dZ0pEFGt3xxHgm4+HW1ShU0glzSEVEckD1L5JRNJXoYLd5ul//4OYGGjcGJ57zm6In4VqVSjKqB6hfNS3GVXKFGLCgu30Hx/OH5sPkqJfhkVE8iyNiLogjYhKtjhxAoYPt7cHDQy0Wz21a5ctt1q7J44vF+4g6nA8waUC6Nf2NhoFlciWe4mIiLM0Iioi11ekCHz2GYSHQ/78dr/RXr0gLi7Lb9UoqAQf92/B8C71OX0+iZenrOKlKSuJijmZ5fcSEZHcSyOiLkgjopLtzp2z2zyNHg3+/vD++/Doo1nW6ulSiReS+WXtfqZGRBKfkESb2mXpHVad0kX8svxeIiKS87RYyc0oiEqO2boVBg6EZcvgjjtgwgSoUiVbbnXmXBLTlu1m1sq9pFjQoWFFHmweTDF/n2y5n4iI5AwFUTejICo5KiXFfmQ/fDgkJsJrr9kLmvLly5bbHY0/xzdLdjF/wwG8PE1aIC1aUIFURMQVKYi6GQVRccShQ/DMM/Djj1Cnjh1Ob789224Xc+Is30VE8vumg3h5GjqFVOKBpsEUKZg/2+4pIiJZT0HUzSiIiqN++snuPXrokP3n229DQEC23e7g8TNMjYhi4eYD5PP0oHPjQLo3DaJwAQVSERFXoCDqZhRExXHx8fDKKzB2LJQtC598Avfdl623PHjsDFMiIlm05SD5vDy5N6QS3ZsFU8jPO1vvKyIit0ZB1M0oiEqusXIlDBgAmzdD164wZgyUK5ett/zr6GmmRESyeMsh8ufz5L7GgXRrGqRAKiKSSymIuhkFUclVkpLgvffgjTfsver/9S946qks37f+SvvjTjElIorwrQqkIiK5mYKom1EQlVxpzx47gM6bB/Xrw/jx2bqY6aL9caf4bmlU2gjpvY0D6XZ7Zc0hFRHJJRRE3YyCqORalgUzZtir62Ni7B6ko0fbuzZlMwVSEZHcSUHUzSiISq536hS8/jp89BEUK2Y/uu/VK1t2ZrrSlYG0S2ggDzQLpoBP9vQ9FRGRjCmIuhkFUXEZGzbA44/bi5rCwmDcOKhRI0duvf/oaaYsiWTx1kME+Obj4ZZV6dSoIt5enjlyfxERsWUURD1yuhgRyUPq17e3B/30UzuU1qsHI0bA+fPZfuuKxQvyUtcGjO3fguDShZgwfxsDxofzx+aDpOgXcBGRXEEjoi5II6Liko4cgeefh8mToW5dmDIFatfOsduv3RPHF7/vYHdsPFVKB9Dvzho0DCqeY/cXEcmrNCIqIs4rWRK+/hpmz7YXMoWE2HNIU1Jy5PaNgkowdkALhnepz6lzSbw0ZSUvT1lJVMzJHLm/iIhcTSOiLkgjouLyYmOhXz+YMwfuugsmTbJ3aMohiReS+WXNPr5bGsWphCTCapXlkdZVKV+sYI7VICKSV2ixkptREBW3YFkwYQIMHQq+vvD55/buTDno9Lkk/rdsN7NWRZN4IYW765WnZ6uqlCzkm6N1iIi4MwVRN6MgKm5l507o2RPWroU+fezH9f7+OVrCidPn+f7PKOas3Q9Ax0YVebhFFfUgFRHJAgqibkZBVNxOUpK9Rejo0RAYCN9+C02b5ngZR04mMGVJJPM3HsDby4N/NKnM/U2DKKgepCIiN01B1M0oiIrbWroUHnkE9u+HF1+E116D/Dk/KvnX0dN8E76L8G0xFPTx4oFmwdwXWhmffOpBKiJyoxRE3YyCqLi1kydhyBB7hX2NGvDllzmyZ316dh8+yaTFu1gVeYSiBfPTs1VV7qlfAS9PNRwREckstW8SEddRqJC9in7uXDh9Gpo1g3/+E86cyfFSgksX4l8PNebdx5pSpogfH8/dwoBPw1m85ZCa4ouIZAGNiLogjYhKnhEfbz+iHz8egoLslfV33OFIKZZlsSrqCF/9sZO9R05RpXQAfe64jUZBxTHGOFKTiIgr0IioiLimgAB7f/rFi8HDA+68EwYOtB/f5zBjDE2qlmLcwJYM71Kf0+eSGPHdKl74ZgXbDpzI8XpERNyBRkRdkEZEJU86exZefx3efx/KlLH3r+/UybFykpJTmLduP99FRHHizHmaVC3JY2HVCS4d4FhNIiK5kRYruRkFUcnTVq2Cvn1h61a7/+hHH0GxYo6Vcy7xAj+tjuaHZbs5fe6CdmkSEbmCgqibURCVPO/8eRg1yv4oWhQ++QTuv9/Rkk6fS2L68j3MXLlXuzSJiFxCQdTNKIiKpNq40R4dXbfO3h70k0+gdGlHS0pvl6YHmwdTtKCPo3WJiDhFQdTNKIiKXOLCBXj3XRg5Evz84MMP7ab4Dq9kP3IygSkRkczfcAAvT0OnRpXo3ixIgVRE8hwFUTejICqSjh07oF8/WLYM2reHCROgQgWnq+Lg8TNMjYhi4eaD5PM0dAypRPemCqQikncoiLoZBVGRa0hOth/Pv/QSeHraI6UDBjg+OgqXBtID5PP0UCAVkTxDQdTNKIiKXMeePXYA/eMPCAuzG+FXqeJ0VQAcPHaGqUv/P5B2aFSJf4QGUqqwn9OliYhkCwVRN2GM6Qx0rlKlyoDIyEinyxHJ3SwLJk6E55+3V9mPHAnDhoGXl9OVAXYg/W5pJH9sPgRAq5pl6HZ7ZaqVLexwZSIiWUtB1M1oRFTkBhw6BE89BTNnQoMGdjht2NDpqtIcOZnArFV7mbfuL84mXqBOxaLc3zSI0Kol8cgFUwpERG6VgqibURAVuQkzZsCTT0JcHAwd+v+r7HOJM+eT+HX9X8xcuZe4+HOUL1aArk0q07ZuefLn83S6PBGRm6Yg6mYUREVu0okT8MIL9qhocDB89hnccYfTVV3mQnIKEdtjmL58D1GH4ynk503HRhXpHFJJC5tExCUpiLoZBVGRW/THHzBwIOzebTfEf/ddKFLE6aouY1kWm/YdZ8bKvazcFYuXpwdtapel2+1BBJb0d7o8EZFMUxB1MwqiIlng7Fn78fz770Px4vDxx/Y2oblwXuaBY6eZtSqa+Rv+4vyFFBoFFafr7UE0CiqOyYX1iohcSkHUzSiIimShdeugf39Yvx7uvdfuQ1q+vNNVpSv+bCJz1+3np9XRHD99niqlA3ikdTWaVC2pQCoiuZaCqJtREBXJYhcu2FuDvvaa3d7p3/+Gxx8HDw+nK0tXUnIKf2w+yNSlUcScOEu1MoV4pHU1GlcpoUAqIrmOgqibURAVySZ79sCgQfD779C8ub2YqWZNp6u6pgvJKSzcfJApEZHE/p3AbeUK80jranpkLyK5ioKom1EQFclGlgXffAP//CecOgUjRthbhnp7O13ZNSUlp7Bg4wGmLo3iyMkEapQvTO821akfWNzp0kREFETdjYKoSA44cgSefRamToU6dWDSpFzVCD89Sckp/LbhL6YujeJo/DkaBZegb5vqVClTyOnSRCQPyyiI5s4JUCIiTitZEr77DmbPhqNHoUkTeP11SEx0urJryufpQadGlfjyiTAGtK3BrkN/8+TEpYyesZ5Dx884XZ6IyFU0IuqCNCIqksNOnIAhQ+xH9vXq2aOj9es7XdV1nT6XxP+W7WbmqmguJKfQoWFFerSsosb4IpKj9GjezSiIijjk55/txUxHj8Irr8DLL0O+fE5XdV3HTp3ju4hI5q3/Cy9PDzo2rEjX2ytTIsDX6dJEJA9QEHUzCqIiDjp+HJ55BqZMgQYN7NHRunWdripTDh4/w5QlkSzacggPA23qlOOBpkFULKGdmkQk+yiIuhkFUZFcYNYse3T0xAl49VV48UWXGB0FiP37LD+u2Muv6/dz/kIKzaqX4oFmwdQon7u2ORUR96Ag6mYUREVyiWPH4Omn7ZX1DRvCV1+5zOgowN9nzvPz6n38tDqa0+eSqFupKA80CyYkWI3xRSTrKIi6GQVRkVxmxgwYPNglR0cBEhIvMG/dfn5cuZej8ecIKhXAA82CaFWzDJ65dHcpEXEdCqJuRkFUJBc6etSeOzp1qsvNHb0oKTmFRVsO8r9le9h/9DSlC/tyf9Mg7q5Xgfz5PJ0uT0RclIKom1EQFcnFZs6096l30dFRgBTLYsWuWH74czfbD/5NIT9vuoQGcm/jQAr6uNb3IiLOUxB1MwqiIrncpXNHXXR0FMCyLLbsP860ZbtZHRWHn7cXHRvZrZ/Ui1REMktB1M0oiIq4iEtHR19/HYYPBy8vp6u6KbsPx/PDst0s2XYITw8P7q5fnu5NgylTxM/p0kQkl1MQdTMKoiIu5OhRe3T0++8hJMQeHa1Vy+mqbtrB42eYvnwPCzYeIDklhda1yvJgs2AqlwpwujQRyaUURN2MgqiIC5o+3V5ZHx8Pb74Jw4a57Ogo2Ls1zVi5lzlr95GQmEyTqiV5sHkwtSoUdbo0EcllFETdjIKoiIs6cgSeeAJ+/BGaNIGvv4bq1Z2u6pbEJyQye/U+Zq3aS3xCEnUqFuXB5upFKiL/T0HUzSiIirgwy4Jp0+DJJ+HsWXjrLXj2WfB07fZI5xIvMHf9X/y4Yg9H488RXCqA+0IDaVO7LN5erv29icitURB1MwqiIm7g8GF7IdNPP9mjo19+CTVrOl3VLUtKTuGPzQf5ccUe9sWdppCfN/c0qECnRpUoWcjX6fJExAEKom5GQVTETViWvYjp6afh1Cl47TV44QWX6zuaHsuy2Bh9jJ9WR7NiVyxgaH5bKe4LrUztCkX02F4kD1EQdTMKoiJu5sgRO4z+8APUr2/vWV+/vtNVZZnDf5/llzX7mLf+L06fS+K2coV5oFkwTauXwkOBVMTtKYi6GQVRETc1c6a9sv7YMXtHpldegfz5na4qy5xLSmbBxr+YvnwPh/9OoHzRAtzfLIg765TTPFIRN6Yg6mYUREXc2PHjMHSovaK+Zk2772jjxk5XlaWSU1KI2HaY/y3fTdTheIoWzE/XJpXpFFIJX2/XbWklIulTEHUzCqIiecC8eTBwIMTE2DsyvfaaW42Ogj2PdN3eo/ywbDcb9h6jkJ83DzQLplNIJXzyaYRUxF0oiLoZBVGRPOLkSXt09MsvoXZte5S0YUOnq8oW2w6c4JvwXazbc5QiBfLzQPNgOjasSH4FUhGXpyDqZhRERfKYuXNhwACIjYWXX7bnjnp7O11Vttiy/zjfhO9iQ/QxihbMz8MtqtC+YUXyeXo4XZqI3CQFUTejICqSB504YTe+nzwZ6tWzR0fr1XO6qmyzad8xJi/exeb9xylTxI9HW1cjrHZZrbIXcUEZBVH9iiki4gqKFLHD508/2c3wQ0LgjTcgMdHpyrJF3UrF+O+jt/PWw43x9fbinVkbePLzpayOOoIGUETch0ZEXZBGREXyuGPH4Jln4LvvoG5de2V9gwZOV5VtUiyLxVsO8fXinRz+O4G6lYrySOtq1KlYVI3xRVyAHs27GQVREQHs0dHHH4e4OHjpJbfrO3qlpOQU5q7bz3cRkfx9JpHqZQvTvVkQzaqXxtNDgVQkt8qWIGqMKQZYlmUdv5Xi5MYpiIpImuPH4Z//tOeO1qpl78rkZn1Hr3Q+KZkFmw4wffkeYk6cpWxRP7rdHsRddctrlb1ILpQlQdQYEwo8BLQBagEX/21PAbYC4cBUy7JW3HLFkiEFURG5ypw5MGiQ3Xf0+edh5Ejw8XG6qmyVnGKxbMdhfli+m12HTlK4gN2HtGMj9SEVyU1uKYgaY+4DXgcuXZ555TOQSy+yERhpWdbPN1GrZIKCqIik6++/4bnn4Isv7BX1P/wA1ao5XVW2syyLTfuO893SSDbstds+Pdg8mA4NK2rrUJFc4KaDqDFmMdDy4peXvHUBOJ76WhHg0j3ZLl4wwrKssJsrWTKiICoiGZozBx57DM6dg08/hV69nK4ox2zed4zJ4bvYtO84xfzz81DzKtzToIICqYiDbiWIpqR+egr4CZgDrLQsK/qK44KAUKAjcB9QEHv+qP7NzwYKoiJyXQcOQI8eEBEBffrAxx9DgQJOV5VjNkQf5ZvwSLak9iEd2LYGTauX0ip7EQfcShDdDbwDTLYs61wmb+YD9AWesywr6CbqletQEBWRTLlwwe41+vbbcNttMG0a1KnjdFU5xrIs1uyO47MF29l/9DR1KxXl8btrEly6kNOlieQptxJEPS3LSr7Jm970uZIxBVERuSELF0LPnvbe9R99ZG8XmodGBpNT7LZPkxfv4lRCEu0aVOCxsGoULejei7lEcousWjXvA9yd+uVuy7K2ZlF9coMUREXkhsXGwiOPwIIFcO+9MGEClC7tdFU56lRCEt8tjeSnVdF4eRg6hlSie9MgBVKRbJZVQdQASdgLlB6xLOu7rCtRboSCqIjclJQU+PBDGDECfH3teaM9euSp0VGAg8fOMCUikkVbDuLl6UGHhhV5oFkwxfwVSEWyQ5Y1tDfGHARKA+0ty5qfRfXJDVIQFZFbsnOnvYBp+XLo0gXGj89zo6NgB9Kpf0axcNNBPD0MHRpWpGerqhTy83a6NBG3klEQ9bjBa03HHhG985arEhERZ1Svbq+m/+9/Yd48e0emqVMhj235XK5YAZ67tx5fPhnGnXXKMXvNPvp+spifVkeTnJJy/QuIyC270RHROthhNAj4AJgBxHB5Q3ssy9qfhTXKFTQiKiJZZscO6N0bVq6Erl3tvqMlSjhdlSOij5xi/PytbNh7jMAS/gxuV5P6lYs7XZaIy8vKR/Mp2KHTcEX4vIRlWZbXNd6TLKAgKiJZKjkZ3nsPXn0VihSxd2bq2NHpqhxhWRbLdsYyYcE2Yv9OoPltpRnQtgZlivg5XZqIy8rKR/MXXQyh5pIPrvhcRERcgacnvPACrFkDpUpBp072vvWnTztdWY4zxtD8ttJMHNyax8KqsWZ3HAPGhzPx9+2cOZfkdHkibudmguiV4fPS10VExFXVqQOrVtmh9PPPoX59e0FTHuTt5UmPllX56skw2tQuy/Tle+jzyWJmr9H8UZGsdENB1LIsj0x8aFtPERFXlT8/vPMOLF5s78zUogW88gokJjpdmSOK+fsw7N56jO3fgkolCjJ23lYenxDB+r1HnS5NxC3c0BxRyR00R1REckR8PAwZApMmQePG8P33EJR3d262LIsVu44wYcE2Yk6cpVXNMgy8qwYlAnydLk0kV8uOOaIiIuLuAgLgq69g+nSIjLQf1U+d6nRVjjHG0LR6KT57vBWPtq7Gil2x9B8Xzv+W7SYpWY/rRW5GhkHUGLPYGNP2Ri9qjGlnjAm/+bJERCTX6NYNNmyw55D26AF9+8KZM05X5RhvL096tqrK54+3pkHl4kxcuIPBE5awKvIIesoocmMyfDR/SbumKGAqMAfYYFlW0hXH5QMaAB2Bh4FgAM0XzR56NC8ijrhwAUaOhFGjoFo1+1F9/fpOV+W4lZGxTJi/nYPHz1C/cjEGtq1BcOlCTpclkmvcdB9RY8xuoDKX9wxNBPYDx7BXyhcFKgH5Lj0V2G1ZVtVbK13SoyAqIo764w/o1QuOHbN3Z3rqKfDI2zO9kpJTmLt2H98uieRUQhJt65XnsbBqmj8qwq0F0XzAM8ALwKVbbVx50qWtm+KAfwNjrxw5layhICoijouLs/ernzMH7rrLnktarpzTVTnu9Lkkvl8axaxV0XgYeLB5Fbo3C8LbSw8IJe+65Z2VjDHewAPAg0AL4MpnDieBCOB7YLplWXmzz0cOURAVkVzBsmDCBBg2DLy9Ydw4ePhhp6vKFQ6fOMvnv29n6Y7DlC7sy6C7a9K0WimMUcttyXuybIvP1IsZoCL2CKkBjgD7Lc3QzjEKoiKSq0RGwqOPwooV8OCDdiAtWtTpqnKF9XuPMu7Xrew/eppGQcV5/O6aVCzh73RZIjkqK0ZEA1I/PWtZ1oV03vcC/AAsy4q/hVolExRERSTXuXDBboQ/ciSULAlffgnt2jldVa5wITmF2Wv28U34LhISL9Cmdjl6tqxKuWIFnC5NJEfcUh9RY0wn4AQQCwRe47BA7JHR46nHSyYZYyoYYxYaY7YbY7YaY/5j9OxGRFyNlxeMGAErV0LhwnDPPTBwIJw86XRljvPy9OAfTSrz5ZNhdLs9iKXbY+g/Ppz3ft5IzImzTpcn4qjMLHN8EPsR/EzLsqLSOyD19emp13sw68rLEy4Awy3LqoHdAqsJ0NXZkkREblLDhrBmjb1f/RdfQO3a8OuvTleVKxQukJ/+bWsw6ek2dAkNZPHWQ/Qbt5ix87Zw8qyWVkjelJkg2gh7lfy86xx38f1Gt1RRJhhjyhljThtjLGNMwWy4fhVjzARjzEZjTLIxZnEGx9ZMHdE8a4w5ZIx50xiT6eWRlmXFWJa1JvXzRGATUOGWvwkREaf4+tqP6Zcts3dnat/eboL/999OV5YrFC3ow6C7azLpqTa0b1CBOWv30/eTRcxYuVc7NEmek5kgerEfx1/XOe7gFcdnp/8Cp7Px+rWADsCu1I90GWOKAL9jB/X7gDeBYcAbN3NTY0wxoAvw282cLyKSqzRpAuvWwcsvw+TJUKsW/PKL01XlGsX8fXi6Qx3GD2xJ9bKFmTB/G4M+XcKKXbHaoUnyjMwEUa/UP6+3TcTFBU1eGR51i4wxLYF7gHczcWwDY0z5a7zXIYORy9mWZVWwLKs7sDWDWzwO+AJdLctaYFnWp9ghdOglC7wwxiw1xkSn8/HFJcfkx57e8KFlWduv972JiLiE/Pnh7bftuaPFikHnzvDkk5CQ4HRluUZgSX/e7hHKvx5qjIeB16et4bnJK9h24ITTpYlku8wE0SOpf3a8znEXFykdyfCoW5AaHD/GHnk8molTRgILjDHFr7jOo8Bs4I70TrIsK7PPRtoDv13RKeB77HDa+pLrtbAsKzCdj36XfF9TgPWWZb2XyXuLiLiORo3suaPPPWe3d2rSBLbrd+6LjDGEVi3Jp4Na8VT72hw8doZ/frWMkdPWsC/ulNPliWSbzATRVdiLlfoYY7qnd4AxphvQB/sR9aqsK+8qjwM+wCeZPL43kAT8enGE0hhzH/AF8LJlWQtusZ7bgB2XvmBZ1n7gbOp7mTUBOIX9WF9ExD15e9tbgs6dC4cP2+F04kS7Mb4A9gr7ziGV+OqpMB4Lq8bG6GM8PmEJ7/60kYPHzzhdnkiWy0wQ/e6SY783xiwwxjxvjOmb+udvwA/AxcfcU7Kj0NT5k/8ChmZ261DLsk4A7YAiwM/GmPbANOB9y7LeyYKyigDpzb4/kfredRljmgP9gBBgvTFmgzHmmWsc29kY89lJtUMREVfWvj1s3AjNmsGAAfZuTPrv2mV8vb3o0bJq6gr7yoRvO0T/ceEKpOJ2rtvQPrWn5WKgZepL6Z1gUl9fYllWm6ws8JI6PgUqWZbVPvXr3sBXgL9lWRkuXDLGVAb+BMoAEy3LGnAD950OFLcsKyyd95KA5yzL+uiK1w8CkyzLGpHZ+9wINbQXEbeQkmKvrn/1VahQAb7+Glq1crqqXOnYqXNMX76HX9bu40KyxZ11y9GnTXWK+fs4XZrIdd1SQ/vUrTu7Y7cVSvf6qX9uIpt6iBpjagF9gTeMMYWNMYVJ3ckJKGSM8b3OJXwBbyAFKGiMycxIcGacAAqn83oh0h8pFRGRizw84KWXICICPD0hLMyeQ3runNOV5TrF/O2WT19f7EG6xR4h/XHFHi6o5ZO4sEwFMsuyjgBNgRHY7YzMJR87gZeApqnHZYeqQD5gOXb4O8H/zxM9gL2AKV3GmEBgPrASaIW9qCqzc0yvZwdXzAU1xlQACnDF3FEREbmGpk1hwwYYNAjee8+eO7pundNV5UoXe5BOeLwVtSoW4bMF23ni8wg27TvmdGkiNyVTe81fdZIxftgjgX9blpXt+5OlrnqvfcXL9wDDsft97rEsa2c655UClgKHgHssy0owxtwBzMWeJ/pyJu6d0aP5l4DnsacMnEp97TnsVf2lr1hNn2X0aF5E3Navv0K/fnDkCLz2mj1i6pWtXQFdlmVZLN8Vy6e/bSP2ZAJNqpbkwebB1KpQ1OnSRC6T0aP5mwqiuUFm5ogaY34BSgN3XBoKU1fOTwc6WZZ1VfP41KDdIfXLYdg9Ul9P/XruxfCd2tB+G7AFeAcIAt7H7gX6yq1+j9eiICoibu34cXjqKZg6FRo3tpvh33YjjUjylnNJycxYsYeZK/cSn5BE7YpFeah5MCHBJbCXeYg4Ky8H0WrAccuyruo5aoxpBiy30vkBpD7O33uNW1e2LCv6kmNrAmOxpy78DUwERlqWlXwD384NURAVkTzhhx9g8GA4exb+/W94+ml7Xqmk61ziBeat/4vpK/ZwNP4cFYoV4K565bmjTjlKBFxvKYVI9nHLIPp/7d13nFTV/f/x12Ep0gQEQQWkCCqgKIqFIqigQELERFCJSayJNeYbNBqjiaDGn9EYYzSxxBobREFFLBG7CCgIqIiA0kGCihTp7HJ+f9wlIbi7LLDL3Zl5PR+P+2CZhT8QLgAAIABJREFUmTvzmbmU95xz7+fkMoOopJyxeDGcd17Se/S44+DBB6FZs7SrqtA2Fmzijamf8+Lk+Xy8YBkBOKRFfU5o34RubfemauXiFhWUyodBNMsYRCXllBjh/vvhl7+EEOD22+Gss5KfVaLPv17Nqx8t4tWPFrF42Rrq1KhK7w5N6Xt4MxrWcZRUu4ZBNMsYRCXlpDlzkgD61ltw8snw979Dgwbb3E3JhU1T5i5l5IS5jJ+5BIBO+zfiu4c3o0PLBlQy1KscGUSzjEFUUs7atAluuw1+8xvYYw946CHo1SvtqjLKkuVreP79+bw4eT4r125kr7rV6XVoU3od2tQG+SoXBtEsYxCVlPM++AB++EOYNg1+8YvkYqbdDFHbY0N+AWOnL+HFyfOZMncplULgyNYN6XVoE45s1ZDKeV4YprJhEM0yBlFJAtauhSuvhDvugIMOgsceg/bt064qIy1aupoXJ8/nlQ8XsWz1evaoVY0eBzfmpCOaey6pdppBNMsYRCVpCy++CGefDcuWwZAhyTKhNsHfIfkFm5jw2Ze8NGUB7336BSFAz/aNOa1zKxrXr5l2ecpQBtEsYxCVpK18+WXSc3T4cDjyyOTc0TZt0q4qo32xYi1PjpvFS5MXkF+wiW5t96F/p5a03rtO2qUpwxhEs4xBVJKKEGPSBP/ii2HVKrj+ehg0CPLsm7kzlq1az/Dxsxn1/jzWbiigXdN6nHxkC7oc2Ig8FxhQKRhEs4xBVJJKsGRJMjr69NNw9NHwyCPQqlXaVWW81es28q8PFjJywlwWL1tDg91348RDmtCzfRMa7+G0vYpnEM0yBlFJ2oYYk7XqL7kk+XnoUNs8lZGCTZEJn33ByInzmDTrSyLQrmk9TjykCccd1JhqVRyB1v8yiGYZg6gkldKcOUnz+6lT4cYb4YorXJGpDH25ci2vfbSI0R8sZMHS1dSpUZWTOjajb8dm1K1ZLe3yVEEYRLOMQVSStsPq1XDuuTBsGJx6KjzwANR0KrksxRj5aP7XDB83m/GffkHVypU44ZAmnHJ0S6ftVWIQtb+FJCm71ayZTNMfdhj8+tcwfTo88wy0aJF2ZVkjhED7ZvVp36w+87/8huHj5/DylIW88P58urbZiwGd9+OAfeqmXaYqIEdEM5AjopK0g156CQYOTH5+8MFk2l7lYuk363h2wlxGTZzH6vX5HNqiPuccf6CBNAc5NZ9lDKKStBNmzYLTToP334dLL4Wbb4Zqns9YXtasz+eFSfP559hZrFizgWPa7MWZxx5A0wa10i5Nu4hBNMsYRCVpJ61fnywPevvtyZT9sGG2eCpna9bnM3z8bIaPn836jZv47uH78pNj92f36lXTLk3lzCCaZQyiklRGnnkmWR60oADuvRdOPz3tirLe8tXrefStT3n+/XnU3K0KZx57AN85bF/yKtnNIFuVFERdEkGSlLtOPhmmTIF27ZJzRwcMSBriq9zUrVmNS/ocxN9+egwtGtbmzhencsl9Y5jw2Rc4OJZ7DKKSpNzWrBm89Rb8/vcwciS0bQuPPpo0wle5adFod27+8dFcc8phrF6/kWuemMBlD4/jw3lL0y5Nu5BT8xnIqXlJKifTpiU9R8ePh+98B+6+G5o2TbuqrLexYBMvTV7A429/yter1nN4ywacc/yBtNq7TtqlqQx4jmiWMYhKUjkqKIA77oCrr4a8PLjpJrjgAqjkJGJ5W7+xgJET5zLsnVl8s3Yjxx+0D2ceewB71auRdmnaCQbRLGMQlaRdYPZsOP98eOUV6NQpuZjpoIPSrionrFq3kX+OncUz785hU4QTD23CaZ33o1FdA2kmMohmGYOoJO0iMSbniw4aBMuXJ2vVX3MNVK+edmU54auV63h8zKe8PGUhm2KkZ/vGnNallcuGZhiDaJYxiErSLvbVV3D55fDww0m/0QcfhK5d064qZ3y5ci1Pjp3Ni5Pnk1+wiW5t9+HUzvux3167p12aSsEgmmUMopKUkldfhZ/+FObOTUZHhwxxVaZd6OtV6xgxfg7Pvz+fNRvyObLVngzovB8H77sHIdiHtKIyiGYZg6gkpeibb+Cyy+Dvf4eDD4ZHHoFDDkm7qpzyzdqNPDdxLs+8N5cVazaw/z51GNBpP7ocuJeN8Ssgg2iWMYhKUgXw/PNw3nmwdGkyMnr55VClStpV5ZR1GwsY/cFCho+fzeJla2i2Zy0uO+kQDtinbtqlaQsG0SxjEJWkCmLpUrjwQnjySTjgALj11qT/qNPEu1TBpsg70//NPS9P4+tV6zm1c0vO6NaaqpXz0i5NuMSnJEnlo359GDYMnnsuucK+b1/o3Rs+/jjtynJKXqVAt7Z7c88F3ejZvjFD35nFxX8fw7ufLnHZ0ArOICpJ0s4IIQmgU6fCn/8M772XnDN68cXJiKl2mVq7VeGykw7h+tOPYEN+Ab8bOpFfPjiWSbO/MpBWUE7NZyCn5iWpAlu6FAYPhrvugt13hxtugJ/9DCpXTruynJJfsImXP1jIY29/ylcr19G2ST0GdG7J0fs3opKnTuxSniOaZQyikpQBpk6FX/wCXnsN2reHv/wFundPu6qcsyG/gJcmL+Cp8bNZsnwtTevX5JROLTnuoMbsVsVzSHcFg2iWMYhKUoaIEUaMSFZmmj8fTj01uaCpSZO0K8s5BZs28fa0f/PPsbOYtWQltatXofehTflex2YuHVrODKJZxiAqSRlm7Vq4+Wa46SbIy4Nrr01GS6tWTbuynBNj5KP5X/Pse3MZO2MJEOneLlmpqWUjV2oqDwbRLGMQlaQMNWcO/N//wciR0KYN3HknHH982lXlrC9WrOWZ9+bwwqT5rN1QwJGt9uQHR7fk0Ob1XampDBlEs4xBVJIy3KhRcOmlSTAdMABuuQWaNUu7qpy19UpNzfasRb8jmtPj4MbsVtWLzHaWQTTLGEQlKQusXZsE0JtuSs4lveIKuPJKqOH5imnZkF/Amx8v5pn35vDZv1dSs1plTjikCd89vBn7NqiVdnkZyyCaZQyikpRF5s9PQuiwYdC0Kfzxj8koqVPDqYkx8vGCZYx6fx5vT1tM/qZIx/325IxurWnbpF7a5WUcg2iWMYhKUhZ6++1kun7KFOjXD+65Bxo1SruqnLds1XpenDz/P9P2h7VswMCurTh43z08j7SUDKJZxiAqSVmqoCBZnenqq6FWraQp/oABaVclYO2GfEZNnMeT42azYs0GWjSsTd+OzehxcGOqex5piQyiWcYgKklZ7pNP4MwzYcIEOO205Or6Bg3SrkrAuo0FvD51Ec9NmPeffqQ/7NqKvh2bUbWyDfKLYhDNMgZRScoB+flJ79HBg6FOHbjtNjjjDM8drSBijExbuIxH3/qUSbO/olHd6pxxTGu6td3bEdKtGESzjEFUknLIRx8la9WPHw8nnJBM1++3X9pVaQvvz/6S+1+ZzqwlK6lWJY9O+zeiZ/vGHL7fnq5rj0E06xhEJSnHbNoEd98NV10FGzYko6SXX56s0qQKYVOMTFuwjNemLuKtaYv5Zu1GmuxRk5OPas4J7ZvkdD9Sg2iWMYhKUo5atCi5sn7EiGRFpsceg732SrsqbWVjwSbGfLKYEe/OYebnK6i1W2X6dNiXk45oTsM61dMub5cziGYZg6gk5biHHoKLLoLdd4fHH3eZ0Apq83mkT787h3em/xsIdGu7N/07taT13nXSLm+XMYhmGYOoJImpU5PWTjNmwLXXwjXXOFVfgS1ZvoZnJ8zlxUkLWLMhn0Nb1OfUzvtxWIsGWd+P1CCaZQyikiQAVq2Ciy+Gf/wDuneHRx5JVmdShbV63UZemDSfp9+bw9Jv1nNg47oM7NqKo1o3zNpAahDNMgZRSdL/ePhhuOQSqFwZ7r3XJvgZYEN+Aa98uIih73zGkuVrOWCfulzUuy0HNs6+JUQNolnGICpJ+pbPPkv6jL73Hpx9Ntx+O9SunXZV2ob8gk28+tEiHnp9Bl+vWk+Pgxtzbo8DqV97t7RLKzMG0SxjEJUkFWnjRhgyBG68EVq3hmeegTZt0q5KpbBmfT5Dx3zGiHfnkFcpcHrXVpxydIusWK3JIJplDKKSpBK9+SaceiqsXZucN9qvX9oVqZQ+/3o1f3/lE8bOWEKjutU55agW9GjfhFq7VUm7tB1mEM0yBlFJ0jYtWAA/+AFMnAi/+11yZX2lSmlXpVKaMucrHnhtBjM+X061ypU49qB9GNi1NXvXq5F2advNIJplDKKSpFJZtw4uvDDpO9q3b3JR0x57pF2VtsOni1fwwqT5vPrhQgo2Rfp2bMbArq2oW7Na2qWVmkE0yxhEJUmlFiP89a8waFCyCtMTT0CXLmlXpe209Jt1PPLmTP41ZQFVK+fRu0NTTjm6ZUas1GQQzTIGUUnSdps4EU47DebNg+uug1//2qn6DDT/y28YNnYWr0/9HICe7Rvzo277V+hAahDNMgZRSdIOWbkSzj8fhg6Fnj3h0UehUaO0q9IO+GLFWoaPn83z788H4KQjmnF6l1bsXqNqypV9m0E0yxhEJUk7LEZ44IGkAX7dusla9ccdl3ZV2kFLlq/hkTc/5ZUPF1K9amVO6dSSHxzVghrVKqdd2n8YRLOMQVSStNM++ihp8TRzZnJF/dVXu1Z9Bpv7xTc8/MYMxs5YQp0aVTmty370PbwZ1aqkf0wNolnGICpJKhOrViVX1T/6KPTokYyONmyYdlXaCdMXLeOh12cyec5XNKi9G2d0a02vQ5uSVym9dewNolnGICpJKjMxwoMPwsUXQ4MGMGIEHHFE2lVpJ30wdykPvj6dTxYup0XD2lzQqy2HNm+QSi0lBVEvl5MkKZeFAOecA++8k0zNH3NM0ndUGe2Q5vW57azOXH3KYaxZn8+Vj7zLDU9N4quV69Iu7X8YRCVJEhx2WNLiqWtXOPvsZIR0w4a0q9JOCCHQre3e/P3C7vy4+/6Mn7mEn971Jk+Om8Wa9flplwc4NZ+RnJqXJJWb/Hy46ir44x/hqKNg2DBo1iztqlQGPv96NX996WMmzvqSGtUq07tDU/od0Zy96pbvsqFOzUuSpNKpXBluuQWefBKmTYMOHWDUqLSrUhnYZ4+a/P6HR3L7OZ05slVDnnl3LqMmzku1JkdEM5AjopKkXeKzz2DAAJgyBS6/HG68EapUSbsqlZEvVqylauVK5b5uvSOikiRp+7VqBePGwQUXJFP13brB3LlpV6Uy0rBO9XIPodtiEJUkScXbbTe4665kWdDNU/XDh6ddlbKEQVSSJG3baafB5MnQujX07w8XXQRr16ZdlTKcQVSSJJVOy5YwZkxyvuhdd/235ZO0gwyikiSp9KpWTa6qHz06WSL06KNhyBDYuDHtypSBDKKSJGn79ewJH30EAwfC4MHQuTPMnJl2VcowBlFJkrRj6taFRx6Bp56C2bOTqfp//CNZv14qBYOoJEnaOaecAh98AB07wplnwo9/DCtXpl2VMoBBVJIk7bwmTeDVV+G66+CJJ5I2T2+9lXZVquAMoikLITQNIbwaQvgkhPBxCOHmEEJIuy5JkrZbXh789rf/DaDdu8MvfgGrV6dblyosg2j68oErY4xtgA7AUcAP0i1JkqSd0KULfPgh/Pzn8Je/wCGHwNtvp12VKqCMCaIhhP4hhLEhhKUhhHUhhBkhhGtCCFXL4bVahRDuCSF8EEIoCCG8UcJj2xaOaK4JIXweQrguhJBX2teKMS6OMU4s/HkD8CHQdKffhCRJaapZMwmhb7yRXLzUvTtcdRVs2JB2ZapAMiaIAvWB14HzgD7AA8DVwJ/K4bXaAd8BZhZuRQoh1ANeASLQD7gOuAwYsiMvGkKoD5wM/GtH9pckqcLp3j25kOncc+Gmm5I2TzNmpF2VKogQM7jFQgjh98DFQL1YxBsJIXQAvowxLizivu8A/4oxFhRxX6UY46bCn58CGsQYjy3icVcBVwDNYowrC2+7AhgM7LXFbWOAJkW8hVdjjOcWPqYa8BIwKsZ4a0nvu2PHjnGiK1lIkjLN00/DeefBunVJE/xLL00a5CurhRDejzF2LOq+TBoRLcpSoKQ/wYOB0SGEBlveGEL4CfAccHxRO20OoaXQhyTMbtmjYihQHei+xfN1jTE2L2LbHELzgMeAydsKoZIkZazvfz9pgn/ccfCrX0G7dvDss/YdzWEZF0RDCHkhhBohhK7ApcBdRY2GFjoL2Ai8FELYvXD/fsD9wG9ijKN3spwDgelb3hBjnA+sKbyvtO4BviGZ1pckKXvtsw+MGgUvvABVqsDJJyerNLkqU07KuCAKrC7c3gbeBH5V3ANjjMuAXkA9YGQIoQ8wDPhTjPEPZVBLPWB5EbcvK7xvm0IIXYBzgY7A5BDClBDCpcU89nshhHtXrFixo/VKklQx9OmTnDt6550waVJyZf0tt0B+ftqVaRfKxCDaGTiGZPSwH3BnSQ+OMS4GegL7Ay8Aj8QYryzDeooajQ3F3P7tnWN8J8YYYowHxxgPLdz+Usxjn4sx/qxOnTo7U68kSRVDlSpw8cUwbRr07g1XXAGdOiXT98oJGRdEY4yTYoxjYox/IpmavzCEsN82dqtOci7pJqBWCKGs3vcyoG4Rt9eh6JFSSZK0tb33hhEjYNgwmDcvWbP+2mth/fq0K1M5y7ggupVJhb+2KO4BIYTmwMvAu0A3oC/w1zJ6/elsdS5oCKEpUJOtzh2VJEklCAFOPTUZHT399GSp0A4dYNy4tCtTOcr0INql8Nc5Rd0ZQmgEjAZmAf1jjO+QTOefHUK4sQxe/0WgVwih9ha3nQasJTl/VZIkbY8GDeCRR5KLmVatSlZp+uUvXSY0S2VMEA0hvBRCuDyE0CeEcGIIYQhwKzAsxjirmN3uB1YA34sxrgWIMb5GEhZ/FULoVcxr1Shcyak/0BjYc/PvQwg1tnjo3cB6YEQIoWcI4WckLaP+tFVLJ0mStD369IGPP4YLL4Q//zm5mOmNN9KuSmUsYxrahxCuB74PNCdZn3028CBwd4xxYzH77A98HWP8qoj7OgPjimmE35xiRlmBFjHGuVs8ti3JBVOdSM4LvQ8YXFSj/LJiQ3tJUk55881kZaZZs+Cii5Kr62vU2PZ+qhBKamifMUFU/2UQlSTlnNWr4eqr4fbboU0beOKJZJRUFV42r6wkSZJyQc2ayRT9yy/DsmVw5JFJKHVALaMZRCVJUuY44QT48EM48UT4v/+DXr1g4cK0q9IOMohKkqTMsueeMHIk3HUXvPMOHHQQPPaYo6MZyCAqSZIyTwhwwQXJMqHt2sGPfpT0H11p05pMYhCVJEmZq1UreOstuPFGGD4cjjjCJUIziEFUkiRltrw8uOoqeO21ZET0qKOSpviq8AyikiQpO3TrBpMnJ1fU/+QnyVT9V99qJa4KxCAqSZKyx157wSuvwA03wIgR0LZtMmWvCskgKkmSskvlyknz+/ffh6ZNoX9/OOMMWL487cq0FYOoJEnKTgcfDOPHw3XXwbBhyUpMb72VdlXagkFUkiRlrypV4Le/TfqNVq0Kxx6bXNi0cWPalQmDqCRJygVHHZVcyHTOOXDTTdC5M3z6adpV5TyDqCRJyg21asF998FTT8GsWdChAzzwgCsypcggKkmScssppyTr1R95JJx7Lpx1FqxZk3ZVOckgKkmSck+TJjB6NAwenDS/79TJqfoUGEQlSVJuysuDa6+FF16AhQuhY0d46CGn6nchg6gkScptvXvDpElJu6ezz4YePWDmzLSrygkGUUmSpGbNkh6jd9+dhNL27eH662HDhrQry2oGUUmSJIBKleD88+GTT6BfP/jd7+Cww2DcuLQry1oGUUmSpC3tvXeyEtNzz8GKFdClC1x8cfKzypRBVJIkqSh9+8K0afDzn8Ndd0HbtjB8uBczlSGDqCRJUnFq14bbb4d334WGDaF/fzjpJJg3L+3KsoJBVJIkaVuOOAImTIBbb4XXX09GR//wB9es30kGUUmSpNKoXBkGDUouZurVC37966T3qK2edphBVJIkaXs0bQojRsCzz8KiRUkYHTEi7aoykkFUkiRpR5x0UtJztE2bZP36yy6z7+h2MohKkiTtqH33TRrhX3QR/OlPyZr106enXVXGMIhKkiTtjGrV4K9/Tabn585NmuD/7W9QUJB2ZRWeQVSSJKksfP/78NFH/22Af+ih8Pzz9h0tgUFUkiSprOyzD7z8Mvzzn7BuXdIUv0ePpDG+vsUgKkmSVJZCgAEDkvB5xx0wZQoccghccQWsWpV2dRWKQVSSJKk8VKkCl1yS9Bk980y45RY48EAYOtTp+kIGUUmSpPLUoAHcdx+MHQuNGsHAgXDccTB1atqVpc4gKkmStCt06gTvvQd3351c1HTYYTBkSE73HjWISpIk7Sp5eXD++TBjBpx6KgwenKzMNGlS2pWlwiAqSZK0qzVoAI8+Cs89B0uXJqOld9yRc+eOGkQlSZLS0rcvfPghnHACXHop9O8Py5enXdUuYxCVJElKU/36MHJkclX9yJFJI/yxY9OuapcwiEqSJKWtUiW4/HJ4++3k527d4IYbsn6ZUIOoJElSRXH00TB5cnIh029/m/w+iy9kMohKkiRVJHXqwGOPwRNPwIIFcMQRMGgQrF6ddmVlziAqSZJU0YQAp58On3wC550Ht92WnDs6ZkzalZUpg6gkSVJFVa8e3HMPvPYa5Ocn544OGpQ1a9YbRCVJkiq6445L2jxdcEEyOtquHTz/fNpV7TSDqCRJUiaoXRv+9rfkyvpatZIepD/6EaxcmXZlO8wgKkmSlEm6dk2urB88GIYOhQ4dkjXsM5BBVJIkKdNUrQrXXgtvvpmcO9qlC9x/f9pVbTeDqCRJUqbq0gWmTIEePZKr62+6KaPWqzeISpIkZbJ69ZKlQQcOhKuugssug02b0q6qVCqnXYAkSZJ2UtWq8Oij0KBBclX955/Dww9DtWppV1Yig6gkSVI2qFQJbr8dmjSBK6+Ef/8bnn46GTGtoJyalyRJyhYhwBVXwOOPw7hxcOyxsHRp2lUVyyAqSZKUbQYOhFGjYOZM6NkTvv467YqKZBCVJEnKRiecAM8+m6xX37NnhRwZNYhKkiRlqxNPhGeegWnTklZPs2enXdH/MIhKkiRls9694ZVX4MsvoVMnmDAh7Yr+wyAqSZKU7bp2hbFjoUaN5AKmkSPTrggwiEqSJOWGAw6A8eOhXTv4/vfhzjvTrsggKkmSlDMaNYLXX4e+feHnP4frr0+1HIOoJElSLqlZE0aMgMsvhz59Ui3FlZUkSZJyTV4e3HJL2lU4IipJkqR0GEQlSZKUCoOoJEmSUmEQlSRJUioMopIkSUqFQVSSJEmpMIhKkiQpFQZRSZIkpcIgKkmSpFQYRCVJkpQKg6gkSZJSYRCVJElSKgyikiRJSoVBVJIkSakwiEqSJCkVIcaYdg3aTiGEL4F5adehHVIHWJF2EdpuHrftl2ufWba830x4HxWtxrTr2dWv3wD4ajv3aRZj3LOoOwyi0i4UQrg3xviztOvQ9vG4bb9c+8yy5f1mwvuoaDWmXc+ufv0QwsQYY8eyej6n5qVd67m0C9AO8bhtv1z7zLLl/WbC+6hoNaZdT9qvv1McEZUkSVKpOCIqSZKktNxblk/miKgkSZJS4YioJEmSUmEQlbJACOGuEMKiEIJTHBnCY7Zj/Nwyj8dMJTGIStnhCeCwtIvQdvGY7Rg/t8zjMcsRIYSmIYRXQwifhBA+DiHcHEIIJe1jEJUKhRAGhBBGFn5zXxVCeD+EMLCcXqtVCOGeEMIHIYSCEMIbJTy2beFf7DUhhM9DCNeFEPK2fEyM8a0Y45LyqLUiCyH0DyGMDSEsDSGsCyHMCCFcE0KoWg6vlZXHLITQuPDPewwh1CqH58/Kzy0NIYTKIYRfhxA+DSGsDyEsDCHcVg6v4zHTjsoHrowxtgE6AEcBPyhph8q7oiopQwwC5gC/JFk14jvA4yGEBjHGO8r4tdoVPv94oNjQFEKoB7wCTAP6AfsBt5J8ibymjGvKRPWB14FbgOXAkcBgYC/gkjJ+rWw9ZrcAq4Ca5fT82fq5peFBoAcwBJgONAXalsPreMxyTAihFfAr4GjgIODtGOOxRTyuLXAH0Ink39z7gCExxgKAGONiYHHhzxtCCB+S/DktXozRzc0tRoAGRdz2ODCnhH06AE2Kue87QF4x91Xa4uengDeKedxVwDJg9y1uuwJYs+VtW9wX0/4c096A3xf+Axk8Ztv8rI4BvgYuByJQq4TH+rmluAG9gY1A2+3Yx2PmVto/K/2ABcCTwCdFHXOgHvA5yZePE4ALgNXADcU8Z/3C52xT0ms7NS8VijEWtXbuZKBhCbsNBkaHEBpseWMI4Sckq10cX8xrbSplWX2Af8UYV25x21CgOtC9lM+Ra5ZSwigOHjMACqdP7wCuo3TrRg/Gzy1N5wCvxRinbcc+g/GYqXSeizE2jTEOAD4u5jEXkBzbH8QYR8cY7yYZnR8UQth9yweGEKqRfIn5c4zxk5Je2CAqlawzybRTcc4iGaV4afNfxBBCP+B+4DcxxtE7+foHkkzB/UeMcT7JiMOBO/ncWSOEkBdCqBFC6ApcCtwVC7+SF+EsPGaQ/KeyG/DXUj7+LPzc0nQUMDOEcGcIYWXhOZkjQgj7lLDPWXjMVAql/PJRqi8ehV9yHwMmxxhv3daTGkSlYoQQepBMVxT7H3WMcRnQi2TKYmQIoQ8wDPhTjPEPZVBGPZJp5q0tK7xvc633hRAWFv68MIRwXxm8diZZXbi9DbxJcq5TkTxmEEKoD1wPDIoxbizNPn5uqduLJFgeCpwOnA0cDjxd3FXJHjOVsdJ+8bgH+Aa4rDRP6sVKUhFCCM1Jzg99Nsb4UEmPjTEuDiH0BN4BXgDuizFeWYblFDWyF7a8PcZ4Xhm+XibqDNR9NvMNAAAEN0lEQVQguVjpd8CdwEXFPdhjxu+Bd2OML2zPTn5uqQqFW78Y41KAEMJiki9exwOvFrWTx0xlaJtfPEIIXYBzganA5MLvSA/EGP9S3JMaRKWthBD2AF4E5gM/KuVu1UnOS9wE1AohVNqO86xKsgyoW8TtdSj6H4ScFGOcVPjjmBDCV8DDIYRbY4yzStgtJ49ZCKEdyfmG3UIIm+usUfhrnRBCQYxxbQlPkZOfWwWwDJi9OYQWGgNsILlyvsggWshjprJS4hePGOM7hb8vNafmpS2EEGoAo0j+0f5ujHF1KfZpDrwMvAt0A/pS+vPutmU6W51rFUJoStJqZ3qRe2hzKG1R3ANy/Ji1BqoA40iCxDL++94XklzAVKQc/9zSVtwFH4EkYBZ9p8dMZadcvng4IioVCiFUJmld0RroEmP8ohT7NAJGA7OA/jHGtYUXA7wQQlgWY/zNTpb1IvCrEELtGOM3hbedBqwlmZLTt3Up/HVOUXd6zBgDHLfVbb2BK0la+swuaic/t9SNAoYU9jXe3OWgG8mXig+K2sFjpjJWPl880u5d5eZWUTbgXpLphUtJmvpuuVUrZp9RwES26plHcpHTRqBXMfvVAPoXbuNI2mVs/n2NLR5Xj6Q58GigJ/AzkubjRfZty7UNeImkB2Yf4ESSViKrgKEl7OMx+/Z7O4tt9xH1c0v3GO1OcrrQOOB7wA9JejSO9pi5lfGftSJ7x5L0jf0aqL3FbZdTTN/YUr9e2m/Yza2ibMDcwv+Mi9qaF7PP/hTRCL/wvs4U31S9eWlfi+T8r9dIRhkWk1ztXGQj6lzbCj+LqYX/+S0nmZb/OVClhH08Zt9+b2ex7SDq55b+cWpFcsHRapJp0oeAeh4ztzL4s7XNLx+U0xePUPjkkiRJykGF5xIXeToT0CLGOLfwcW1JupJsucTn4Fi4xOcOvbZBVJIkSWnwqnlJkiSlwiAqSZKkVBhEJUmSlAqDqCRJklJhEJUkSVIqDKKSJElKhUFUkiRJqTCISpIkKRUGUUmSJKXCICpJkqRUGEQlSZKUCoOoJEmSUmEQlSRJUioMopIkSUqFQVSSJEmpMIhKkiQpFQZRSZIkpcIgKkmSpFQYRCVJkpQKg6gkSZJSYRCVJElSKgyikqQihRBqpV2DpOxmEJWkHBZCGBxCiIXbQyGEXiGEsSGE1cCYtOuTlN0qp12AJKnCOAb4MQ5SSNpFDKKSpM1aAtOA/wd8ATRNtxxJ2c4gKknabA1wfIxxSdqFSMoNTr9IkjZ7xxAqaVcyiEqSNlucdgGScotBVJK0WUy7AEm5xSAqSZKkVBhEJUmSlAqDqCRJklIRYvSUIEmSJO16johKkiQpFQZRSZIkpcIgKkmSpFQYRCVJkpQKg6gkSZJSYRCVJElSKgyikiRJSoVBVJIkSakwiEqSJCkVBlFJkiSl4v8Dq6z11y1JXd8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize = (10,7))\n",
    "plt.plot(r[n:-m],cr_N2[n:-m], '-', color = 'red', markersize = 10, label = 'Introgressed')\n",
    "plt.plot(r[n:-m],cr_H2[n:-m], '-', color = 'steelblue', markersize = 10, label = 'Depleted')\n",
    "plt.xlabel('r', fontsize = 18, labelpad = 15, weight = 'semibold')\n",
    "plt.ylabel('C(r)', fontsize = 18, labelpad = 15, weight = 'semibold')\n",
    "plt.yscale('log')\n",
    "plt.xscale('log')\n",
    "plt.tick_params(labelsize=15,  which='both')\n",
    "plt.legend(fontsize = 'x-large')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Write DNA Sequence"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "string = ''\n",
    "with open('human.txt', mode='w') as f:\n",
    "    f.write(string.join(H_DNA2))\n",
    "    f.close()\n",
    "string = ''\n",
    "with open('neanderthal.txt', mode='w') as f:\n",
    "    f.write(string.join(N_DNA2))\n",
    "    f.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Divide in Sentences and then in KMers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('human.txt', mode='r') as f:\n",
    "    H_DNA = f.read().rstrip('\\n')\n",
    "    f.close()\n",
    "with open('neanderthal.txt', mode='r') as f:\n",
    "    N_DNA = f.read().rstrip('\\n')\n",
    "    f.close()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [],
   "source": [
    "H_DNA = [i for i in H_DNA]\n",
    "N_DNA = [i for i in N_DNA]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 86,
   "metadata": {},
   "outputs": [],
   "source": [
    "sentence_length= 200 \n",
    "\n",
    "H_DNA = np.array(H_DNA[0:(len(H_DNA)//sentence_length)*sentence_length])\n",
    "N_DNA = np.array(N_DNA[0:(len(N_DNA)//sentence_length)*sentence_length])\n",
    "\n",
    "#divide in sentences\n",
    "H_DNA_s= np.array(np.split(H_DNA, len(H_DNA)//sentence_length))\n",
    "N_DNA_s= np.array(np.split(N_DNA, len(N_DNA)//sentence_length))\n",
    "\n",
    "H_DNA_s = [''.join(i) for i in H_DNA_s]\n",
    "N_DNA_s = [''.join(i) for i in N_DNA_s]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 89,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "KMER:  10\n"
     ]
    }
   ],
   "source": [
    "def getKmers(sequence, size):\n",
    "    return [sequence[x:x+size] for x in range(len(sequence) - size + 1)]\n",
    "\n",
    "kmer = [10]\n",
    "for k in kmer:\n",
    "    print('KMER: ', k)\n",
    "    intr_texts = [' '.join(map(str, (getKmers(i, k)))) for i in N_DNA_s]\n",
    "    depl_texts = [' '.join(map(str, (getKmers(i, k)))) for i in H_DNA_s]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 90,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "50000\n",
      "50000\n"
     ]
    }
   ],
   "source": [
    "merge_texts = intr_texts + depl_texts\n",
    "\n",
    "tokenizer = Tokenizer()\n",
    "tokenizer.fit_on_texts(merge_texts)\n",
    "#X = tokenizer.texts_to_matrix(merge_texts, mode = 'freq')\n",
    "\n",
    "encoded_docs = tokenizer.texts_to_sequences(merge_texts)\n",
    "max_length = max([len(s.split()) for s in merge_texts])\n",
    "X = pad_sequences(encoded_docs, maxlen = max_length, padding = 'post')\n",
    "labels = list(np.ones(len(intr_texts))) + list(np.zeros(len(depl_texts)))\n",
    "print(len(X))\n",
    "print(len(labels))\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, labels, test_size = 0.20, random_state = 42)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 96,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1025"
      ]
     },
     "execution_count": 96,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "y_train = np.asarray(y_train)\n",
    "y_test = np.asarray(y_test)\n",
    "vocab_size = len(tokenizer.word_index) + 1\n",
    "vocab_size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([491, 457, 260, 264, 455, 408, 182, 370, 180, 399, 613, 549, 229,\n",
       "       214, 175,  70,  10,   7,  42,   8,   9,  12,  14,  13,  11,  10,\n",
       "         7,   3,   1,   1,   1,   4,   8,  55, 102, 116, 126, 125, 115,\n",
       "       103,  53,   7,   3,   1,   4,  24,  37,  43,  45,  49,  48,  44,\n",
       "        38,  23,   3,   1,   1,   4,   8,   9,  12,  14,  13,  78, 184,\n",
       "       224, 158,  36,  34,  30,  26,  39,  16,  18,  20,  21,  22,  19,\n",
       "        63,  87,  64,  59, 164, 275, 619, 768, 789, 666, 263, 247, 120,\n",
       "       246,  87,  64,  59, 164, 275, 619, 695, 958, 888, 909, 969, 703,\n",
       "       625, 267, 162,  60,  17,  15,   5,   2,   2,   2,   2,   2,   2,\n",
       "         2,   2,   2,   2,   2,   6,  16,  56,  84,  89,  90,  91, 334,\n",
       "       423, 308,  99,  68,  22,  19,  17,  15,  39,  16,  56,  84,  89,\n",
       "       333, 576, 812, 831, 688, 287, 226,  46,  41,  26,   5,   2,   2,\n",
       "         2,   2,   2,   2,   2,   2,   2,   2,   2,   2,   2,   2,   2,\n",
       "         6,  16,  18,  59, 164, 197, 427, 519, 513, 724, 523, 454, 533,\n",
       "       593, 748, 422, 465, 912, 905, 471, 441, 305], dtype=int32)"
      ]
     },
     "execution_count": 97,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X_train[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Model: \"sequential_5\"\n",
      "_________________________________________________________________\n",
      "Layer (type)                 Output Shape              Param #   \n",
      "=================================================================\n",
      "embedding_5 (Embedding)      (None, None, 10)          10250     \n",
      "_________________________________________________________________\n",
      "bidirectional_4 (Bidirection (None, 10)                640       \n",
      "_________________________________________________________________\n",
      "dense_8 (Dense)              (None, 10)                110       \n",
      "_________________________________________________________________\n",
      "dense_9 (Dense)              (None, 1)                 11        \n",
      "=================================================================\n",
      "Total params: 11,011\n",
      "Trainable params: 11,011\n",
      "Non-trainable params: 0\n",
      "_________________________________________________________________\n",
      "None\n",
      "Train on 32000 samples, validate on 8000 samples\n",
      "Epoch 1/5\n",
      "32000/32000 - 36s - loss: 0.5759 - accuracy: 0.7017 - val_loss: 0.5158 - val_accuracy: 0.7477\n",
      "Epoch 2/5\n",
      "32000/32000 - 31s - loss: 0.5094 - accuracy: 0.7579 - val_loss: 0.5250 - val_accuracy: 0.7491\n",
      "Epoch 3/5\n",
      "32000/32000 - 31s - loss: 0.5010 - accuracy: 0.7600 - val_loss: 0.5220 - val_accuracy: 0.7499\n",
      "Epoch 4/5\n",
      "32000/32000 - 29s - loss: 0.4908 - accuracy: 0.7680 - val_loss: 0.5138 - val_accuracy: 0.7533\n",
      "Epoch 5/5\n",
      "32000/32000 - 29s - loss: 0.4849 - accuracy: 0.7692 - val_loss: 0.5058 - val_accuracy: 0.7515\n",
      "Accuracy: 75.06%\n"
     ]
    }
   ],
   "source": [
    "d = 10\n",
    "model = Sequential()\n",
    "model.add(Embedding(vocab_size, d)) #dropout = 0.2 #input_length = max_length\n",
    "model.add(Bidirectional(LSTM(int(d/2)))) #dropout = 0.2, recurrent_dropout = 0.2\n",
    "model.add(Dense(d, activation = 'relu'))\n",
    "model.add(Dense(1, activation = 'sigmoid'))\n",
    "\n",
    "epochs = 5\n",
    "model.compile(loss = 'binary_crossentropy', optimizer = 'rmsprop', metrics = ['accuracy'])\n",
    "#checkpoint = ModelCheckpoint(\"weights.best.hdf5\", monitor = 'val_acc', verbose = 1, \n",
    "#                             save_best_only = True, mode = 'max')\n",
    "\n",
    "print(model.summary())\n",
    "\n",
    "history = model.fit(X_train, y_train, \n",
    "                    epochs = epochs, verbose = 2, validation_split = 0.2, batch_size = 32, shuffle = True, \n",
    "#                    callbacks = [checkpoint]\n",
    "                   )\n",
    "\n",
    "predicted_labels = model.predict(X_test)\n",
    "cm = confusion_matrix(y_test, [np.round(i[0]) for i in predicted_labels])\n",
    "cm = cm.astype('float') / cm.sum(axis = 1)[:, np.newaxis]\n",
    "\n",
    "scores = model.evaluate(X_test, y_test, verbose = 0)\n",
    "print(\"Accuracy: %.2f%%\" % (scores[1]*100))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 102,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVIAAAEdCAYAAABe2+9WAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dd5xU1fnH8c8XpAsCIkWkKiB2I2JXbIi9R7Bh5acGjRpj78aIxtiiRohijw0bxoKAghUF7KJ0BKQX6VLk+f1x7uDsMLN7d3d29y48b1/z2plbzj13Fp8991SZGc4550quSkVnwDnnKjsPpM45V0oeSJ1zrpQ8kDrnXCl5IHXOuVLyQOqcc6XkgRSQdLMkkzQoy74BkoZVQLZKJLqP3mmfh0kaUE7Xnifp5vK4VhyS6kh6XtL86Hs5K0/pto7SOyof6VU0SZ2L83uT1CW6/x3KMFuVyiYVnYGE6SppdzMbWdEZyaOLgNUVnYkKciFwNHAm8DMwMU/pzgT2An7MU3oVrTNwE3BzzOO/INx/vr7PSs8D6e8WANOB64Dj8p24pFpmtiLf6RbFzMaU9zUTZFtgrJm9nM9EzWwlMCKfaVYGkgTUMLPFbIT3Xxh/tP+dAX8HjpG0Y2EHStpF0lBJyyUtlPSspCZp+1OPfqdJekrSL8Ab0b4pku6WdLWkmZIWSfqngiMkfS9piaTXJDVIS7OOpAcljY2uO1nSQ5LqFZHXAo/2Ub6yvc5KO2ZfScOj68yX9B9JdTPS3V/S15J+lTRa0t5xvmRJVSVdI2mcpJWSpkt6IuOY3pLGR/snSLosY//NUTXCrpJGRPn8UtJ+acdMAc4Fdk3dY7T9CUmjMtJb71Fd0rnR72JFdK3hkrYv5PiqUb6mRvn+XtKpGdd5QtIoSYdK+kbSMkkfpdIt5DtLPUofLOn16LzxkrpG1/1HlMefJV2ece5ekgZKmhGd95Wk09L2nwX8K3qf+rcwLON73lfSSOBX4GRlPNpLOlnSWkkHZ3yniyX9rbB721B4IC3oJWAcoVSalaQtgGFAbeBU4GLgAGCwpOoZh98NLAFOJgTplO6Ex6mzgbuAy4F7gNuAG4ALojTvSDunNlA1ytvh0XEHRXkujr0yXv8i/BGZEN3fPsBQYBZwEnApcATweNp3sCXwNqEUfxLQF3g2ymNR+gK3AC8CRwF/AeqkpX1+lKeBhMfyl4B/Sro6I53awJNReicCK4FXJaXycDzwFuHxO3WvsUjaH3gEeIbwXZ8DfAJsVshptxJ+N/2AY4CPgWcl9cg4riXwD+B2oAfQGHhRkmJkrS/wEeHefgIGAA8CdQn/FgcQvqs9085pFeXlPML3+TLweFq+3gT+Gb1PfU8XpZ2f+p4fBboBn2dmysxeAl4A+kuqF91Lf2By9L1s+Mxso38R6obmRe/PAn4D2kefBwDD0o7tA/wC1Evb1pkQjHpEn1tHn1/Ncq0phKBVNW3b58AaoE3atruA2YXkeRNgn+g6LdO2G9A77fMwYECONP4ArABuStv2IfB+xnEHRenukJa3+UDttGNOi465uZA8bxsdc0mO/VUIdZmPZ2x/GFgE1Ez7fRlwUNoxu0TbuqVtewIYlZFWtm2p39dR0ecrgNGF3Efm8Q2BZenfY7T9LULVQvq11wDt0rYdF6W1bSHX6xIdk/572i7a9l7G9zcLuDNHOor+3fTNOK83YDn+vzDg2Bz52SFtW0NgBvAYcAmwCtg5X/+PJv3lJdL1PQNMBa7Jsb8z8K6FeiIAzOxzQoDcN+PYN3OkMczMfkv7PAGYYmaTM7ZtkV7KlXRG9Ai7lNCA9FG0q33ht7S+qGT9KjCEqNQQleb2IpSQNkm9ouusBnaLTu8MDDaz5WlJvhLjsgdGP5/IsX8rYEvWL2W/ANQD0qtcVhP+SKSk6oK3ipGPonxFqBK4N6rCyHzSyLQDoeSWLd/tJTVO2zbFzManfS5OvoemvZ8Q/XwvtcHM1gKTgOapbZIaSHpA0k+E72w10Iv4/2aM8PRR+EFmC4DzCaX3fwC3mNnXMa9R6XkgzWBmawglrtMltcpySDNgdpbtswl/lTO3ZfNLxudVObYJqA4g6XjgKeBTQlXBnoRHPICaOa6TVRQcX4yucbpFRQqgAaH64GF+/59uNeGxuRrQIjquKTAnPU0LDWlLi7j05sCy9D9CGZpFPzO/t9Tn9O93cRQ4UtdfFb0t1neRjZkNIVS77E8I1vMkPSypTo5Tisp3g7Rt2X7PEC/f685Nu99s6aWn9QRwCiG4dQV2Jzx2x/2eFqZdqyjvEe65CvCfmOdsELzVPrv+wPXAVVn2zSTUa2VqAozO2JbPOQpPBj4zs3X1V5IOKGFadxP+h9rTzBalbf+F6PGc8FiaaUb0cxYZ34GkWsCmRVx3PlBHUr0cwXRm9DPz+0015C0oIv04fiX645Qm8w8gZvYk8GRUcj8BuBdYDGTW1ULBfM9P257PfBebpJrAkYSqnkfSthenAFWcf8N9CH+IZwH3EeptNwpeIs3CQveWuwmPKc0ydn8GHKa0VmxJuxPqzT6i7NQilAzTnZbtwMJIOgP4M3CumX2Xvs/MlhG6tXQws1FZXqlAOhI4NK1hB0KwKUrqMfTMHPunE4L1yRnb/0gIYt/GuEZRpgOtoyCTcmiug81srpn1JdQdb5fjsO+A5WTP9zgzm1uK/JZGDUJgW/fvJvp3e0zGcauifSUuzUd/1C8m9N09F+gh6cSSplfZeIk0t77AtcDewPC07fcQ/rEMknQnoRTWh/A/eV77K2YYDDwk6TpCMD8COLjwUwqStDWhVflt4KeM1t2J0f/wVwJDJa0lNLQtIbQ0HwlcZ2bjCKWNPwH/k3QPoV7zGkLDVU5mNlZSP0LLcmPgA6A+cJKZdTeztQojbPpKmh/d8wGE7/taM/u1OPebw2uEOuFHFbpd7Up4jF9H0i2EUuowYF50zAFkL41iZgsk3QdcL2kNMIrwh+UIQst8hTCzRVG3pRslLQbWEu5hEaHOOSU1sODPkt4jVJuMjXsdSZsSenW8YGYDom19gX9L+qAC/5CUGy+R5hA1pNybZftcQqPJr8BzwEOE0sqhxahLKom+hG4qfyY07LSi+I9OLQh1Y4cT6lrTX0cCmNlHhLrBLYCnCf1frwSmEdX5mdnPhCDRiPDH4yLgdEKprCgXEbo/nU6oPriPtABsZv8htPoeD/yPEIj+YmZ9inmvWUWl8HMIjWoDCQHynIzDRhJKn48AgwiB/Gbg/kKSvpHQXe3CKN/7E+qfn89HvkvhVEI3pKcI+X85ep/uQ0Id6p8Jf6T7FvMa/yQ8MfVO23YFoc78kaxnbGD0ezuDc865kvASqXPOlVJiAqmkhpIGR0PfBitteGTGcb9Fw9y+kjQwbXsbSZ9F578Qo++fc87lRWICKaESfKiZtSN0PM5asQ+sMLNdold66+OdwL3R+QsJLYfOOVfmElNHKmks0MXMZkpqRhj90yHLcUvNbNOMbQLmAk3NbI2kvQhDFQ8rl8w75zZqSer+1MTMZgJEwTRbp3eAmgqz96wB+pjZa4QRM79Eo5Ig9BVsnuN8JPUiDJODqjV2q1Ivs6uoS7KOLbLW+riE+nnaVBYumBdnUpacqtZrZbYm3iyUtmLuIDPrVprrFVe5BlJJQwjDCzPlnG0pi5ZmNkNSW+A9Sd8SOmtnylnUNrN+hP6UVG3Yxmp3vbkYl3cV7ZX7Mvu9uyQ74bDMKSiKz9asoEaHP8Y69tevHmpU6gsWU7kGUjM7JNc+SbMlNUt7tJ+T7bjU6BozmxTNm7groW9cfUmbRKXSrfh9OKNzrtITFGtka/lKUs4GAj2j9z2B1zMPiGayqRG9b0SYRm5MNOnG+4S5MXOe75yrpARUqRrvVQGSFEj7EMZvjyeMfe4DIKmTpEejYzoCoyR9TQicfez3pTSuAi6XNIFQZ/pYuebeOVe2pHivCpCYxiYzm0+WseNmNoowuzdm9gkF56RMP24SYZ5M59wGJ9mP9okJpM45V6gKKm3G4YHUOZd8wkukzjlXOhVX/xmHB1LnXOVQQS3ycXggdc5VAt7Y5JxzpSP80d4550rNS6TOOVca/mjvnHOlI6CqNzY551zpeB2pc86Vhj/aO+dc6XmJ1DnnSslLpM45VwoVOEVeHB5InXOVgw8Rdc650vDGJuecKz1/tHfOuVJI+HykicmZpIaSBksaH/1cb/FySbtI+lTS95K+kXRK2r4nJE2W9FX02qV878A5V3aiR/s4rwqQmEAKXA0MNbN2wNDoc6blwJlmtj3QDbhPUv20/X81s12i11dln2XnXLnxVURjORZ4Mnr/JHBc5gFmNs7MxkfvZwBzgC3KLYfOuYqT4FVEkxRIm5jZTIDoZ+PCDpbUGagOTEzbfHv0yH+vpBpll1XnXLlSsh/ty7WxSdIQoGmWXdcVM51mwNNATzNbG22+BphFCK79COvc35rj/F5ALwDV3rw4l3bOVRRvtQ/M7JBc+yTNltTMzGZGgXJOjuPqAW8C15vZiLS0Z0ZvV0p6HLiikHz0IwRbqjZsY8W/E+dceVOCA2mSHu0HAj2j9z2B1zMPkFQdeBV4ysxeytjXLPopQv3qd2WaW+dcuQkrjSjWqyIkKZD2AQ6VNB44NPqMpE6SHo2O+SOwP3BWlm5Oz0r6FvgWaAT8rXyz75wrMxKqEu9VERLTId/M5gMHZ9k+Cjgvev8M8EyO8w8q0ww65ypUPkubkroB9wNVgUfNrE/G/nuBA6OPtYHGZlafHBITSJ1zrjD5CqSSqgIPEZ58pwMjJQ00szGpY8zssrTjLwZ2LSzNJD3aO+dcTnmsI+0MTDCzSWa2Cnie0I89lx7Ac4Ul6IHUOZd8KsYLGkkalfbqlZFac2Ba2ufp0bb1Lyu1AtoA7xWWPX+0d84lnihWi/w8M+tUaHLry9UNsjswwMx+K+yCHkidc5VClSp5e4CeDrRI+7wVMCPHsd2BPxWVoD/aO+cqhTzWkY4E2klqE/VN707ox555vQ5AA+DTohL0QOqcS77i1ZEWyszWAL2BQcAPwItm9r2kWyUdk3ZoD+B5Myty9KM/2jvnKoV89iM1s7eAtzK23Zjx+ea46Xkgdc4lXjEbm8qdB1LnXKVQUcM/4/BA6pxLPiV79icPpM65SsEDqXPOlZIHUuecKwVvbHKxbdmwNtefvAuH7NychnVrMOuXFfxv5FT6vPw1vyxbVeT5+27XhLdv7Fbkcdv+6SV+nr983eczD9yG3bZuxI6tG7J9iwbUrrEJd73yDbe9+GWp7mdjMGvGz9z/j9v48P3B/LJwAY0bN+XgbkfR+y/Xsln99VYUX8/y5csY8vYbDB/6Dt9/8zWzZkxHVarQZut2HHX8yZx+zoVUr159vfNWr17N0489zOsDnmPKxAlUqVqVDh23p0fP8zn2pB5lcasVL7lx1ANpUrRpUpchtxxO4/q1+N/IqYybsYjdtm7En47YjkN3bs6hN73NgqUrC01j6tyl/H1A9lWot2/RgGP3aMWYaQsLBFGA20/fnfp1qrNg6UpmLlzO1k3r5e2+NmRTp0yi+9EHMX/eXA4+7CjatmvPN1+O4qlHH+bD94fw3MAhNGhY+Jpgo0Z8zF97n0v9Bg3ZY+/9OeTwo1i0cCHvD36LO2+5lnffGsiTL75JjZo1152zatUqzj/tOEZ8NJzmLVpx/CmnA/DB0EFcefF5jPn2K6655c4yvfdyp7wOEc07D6QJcc85e9C4fi2uePwz+g76cd32O87oRO8jt+fGU3bl0sdGFJICTJ27jDsGfJ11X/+L9wfg8aHj1tt39gPDGfvzIqbNW8ZpB2zNIxfuW4o72XjccvWlzJ83l+v/djdnnHvhuu133HQVT/R7kHv73MKtdz1QaBpbNG7CPx58jG5Hn1Cg5Ll06RLOPKEbX44cwbNP9OWcC/68bt9/n+jHiI+Gs2unPej/whvUrl0HCKXbnicdwRP9HuSgw45kj733z/MdV6wkP9onN8RvRFo33pRDdm7OlDlL6PfujwX23f7SVyz9dTXd92tL7Rol+7vXcNMaHL17S5avXMPzH05ab/+Qr2cwbd6yEqW9sZr202Q+Gj6U5i1acdrZ/1dg38V/vZ7ateswcMBzLF9e+PfacYedOebE7us9vm+6aV3OvuASAD7/5MMC+wa/FYaFX/DnK9cFUYDatetw0WVXA/DMY4+U7MaSLE9DRMuCB9IEOGD7ZgC8980MMkf1Lv11DSPGzqFOzWrs3q5RidI/7YCtqVm9Kq+NmBKrrtUVbcRHwwHY94CD13vk3HTTuvxh9z1ZsWI5X4/+vMTX2KRaNQCqVi34B3Te3NkAtGjZer1zUts+/WhYia+bVL74XTFI6iZprKQJkq7Osr+GpBei/Z9Jap2275po+1hJh5Vnvkuj3ZahTnLCzMVZ90+ctQSAbZptVqL0ex7UDoD+WR7rXclMmhi+y9Zbb5N1f6u2YfvkSRNKfI2Xn3sKgP0OPLTA9lS96/RpU9Y7Z9rUsG3J4kXMnTOrxNdOmrhB1AMpBdZSORzYDughabuMw84FFprZNsC9wJ3RudsRpsPaHugGPByll3j1aofHusXLV2fdv3h5KEXWr12t2Gnv07EJHZrXZ8y0hXw2bm7JM+kKWLo4/NGrWzf7H7e6dcMfxyWLfilR+s/0f4QP3x9Mxx124sQeZxbY1+WQwwF45P67+XXFinXbVyxfziP337Xu8+JfSnbtpEpyIE1aY9O6tVQAJKXWUhmTdsyxwM3R+wHAg9Fa9scSprxaCUyWNCFKr8i5BJMu9W+j6Mm81nf2we0B6D/ES6PlKTXzWkn+x373zdf5+41XskXjJvzr0f9SrVrBP6Bnnn8R7775Gl+M/JQju3TigIMPw8wYPnQQy5YupXHTZsyZNZMqVStFOSK2SjnWXtIRxUkompaqtLKtpbJHrmPMbI2kRcDm0fYRGefmWoelF9ALQLUL755SHlIlzno5Spx1a4US66IV2UusuTSoU51jO7eKGpkmli6TroBN60UlziWLsu5funRJdFzxqmOGvP0Gl1/Yk4aNtuCpAW/TolWb9Y6pXbsOz7z2Lv3+9U/eeeNVXnz2cWrXrsOe+3bhL9fewmnHhaqAhpuXrE49qZLcal9YifR/hHVM4uTeCOtDl1actVRyHRN7HRYz6wf0A6jasE0Jynn5NX5GeEzcpln2/ptbN60LwISZ2f+nzeXUA7ahZvWqPDt8AotyVBu4kmm7dSjpT5mYvQ70p6hutE3b7HWo2bz9xitccdHZNGrchCdfeovWhZxbu3YdLr3qRi69qsAUmkybOoW5c2bTqs3WsQYEVBqVeNKS9f8Ulr04a6mkjpkuaRNgM2BBzHMT6YMxoVHgoJ22RCr4CL9pzU3Ys0Njlq9cw8jx84qV7llRI1O2vqOudPbYJ/TR/Gj4UNauXVug5X7p0iV8MXIENWvWYufdOsdK741XXuCqS86nSdMteerl7CXROF565nEAjj7hlBKdn1Ti9yquJMrZ2GRmPxXnlaf8xFlLZSDQM3p/EvBetBTAQKB71KrfBmgHlLzvSTmaPHsJQ77+mdaN69Kr67YF9l138i5sWrMaz30wkeUr16zb3n7LerTfMvcIpL23bcy2W9Xn+6neyFQWWrZuy74HHMzP037i2cf7Ftj3r3/8jeXLl3HsyacW6Oc5cfxYJo4fu15ar774DFdefB7NmrfgmVffjRVEly5Zv4fH8KHv0L/vAzRptiVnnndRCe4qyZLdah+7sUlSDeAcoBOh5PcnMxsv6RTgGzP7obSZieo8U2upVAX6p9ZSAUaZ2UDgMeDpqDFpASHYEh33IqFhak2Uv0KXUE2Sy/t/xpBbDufus/egyw7NGPvzIjpt04gDdmjG+BmLuPWFguPeR99zPAB1uz+ZNb2zokamOKXRnge2Y69tGwPQtkmoRjh8t61ovnltAMb9vIh7Bn5XshvbgN3U5z66H30Qf7v+Cj79cBhbt+vA11+O5LOPP6D11u247OqbChx/xP5/AGDszN876Y/4eDjXXnYha9euZY999ueVF55e7zp1623GWb16F9jWbb9d6dBxB9pu057qNWrw3ddfMOKj4TTcvBGPPPnShvVYH6lSGRub0klqDwwmPEaPBroAdaPd+wFHAmdmPbmYilpLxcx+BU7Oce7twO35yEd5mzx7Cftf9z+uP3lXDtl5S7ru2pxZC1fw8Ntj6DPgaxYWoyN9/TrVOW6P+I1Me23bmNMOKFgft2OrhuzYqiEAH46Z5YE0i5at2/LyOx/xQDRpyQfvDWKLxk0549wL6f2Xa6nfoGGRacyYPpW1a9cCv/cbzdR8q5brBdKjTziFj94fzJejPmPNmtU0a96Ccy64hPP/dDkNG21R+ptLGiX70V4xFshD0jtAHeBoYCmwCuhkZl9IOhm408zalmlOy0jVhm2sdtebKzobrhhG35f176hLqBMO25fvvv6iVGGwVrP21ubsB2Md+8Mdh402s06luV5xxX203w842cx+ydLJfTbQLL/Zcs65gpJcIo0bSH8FauXY1xzYsIZQOOcSJ8ndn+IOER0MXCspvXexRQ1QF5NRp+mcc3kV1ZHGeVWEuCXSvwIfAxMIQdWAGwnj2qsDJ5RJ7pxzjrDUSJIndo6VMzObBuwMPAK0BiYS6kVfAnYzsw1nmhnnXCJtCCVSzGwhcEP0cs65cpXkOtJizf4kqT6wA6E0OgP43sy8ock5V7YS3o801qO9pE0k3UkYz/4B8ALwIWG8+12Sij9RpnPOxRTG2udviGhRE8hHx/xR0hhJ30v6b2HpxS2R3kOYdu5W4BVgDtAYOJHwqF8TuCRmWs45V2z5KpGmTSB/KKFwOFLSQDMbk3ZMO+AaYB8zWyipcWFpxg2kZwDXmtk9adsWALdL+hW4Hg+kzrkylMex9nEmkD8feChqG8LM5hSat5gXXgt8n2Pfd+SY99M55/JCxXq0byRpVNqrV0Zq2SaQz5wEvj3QXtLHkkZI6lZY9uKWSJ8GziPMypTpfOCZmOk451yxFXM+0nlFjLWPMwn8JoSpOLsQ5jb+UNIOuRrXC1tqJH1CwynASZK+J8z7maojPZYwC9TdhWTaOedKKa9zjcadQH6Ema0mrAE3lhBYR2ZLsLASabapVrYEOmbZfg9wfyFpOedcqeSx+9O6CeSBnwlzGp+accxrQA/gCUmNCI/6k3IlmDOQmllyx2M55zYuyl9jU8wJ5AcBXSWNAX4D/mpm83OlmbTlmJ1zbj2pfqT5EmMCeQMuj15FKu7Ipq0IRdyaOTLmnHNlotIPEZVUF3gR6JraFP1Mb+nKx3LMRN0M7o/Se9TM+mTsv5zQg2ANMBc4J7X4nqTfgG+jQ6ea2TH5yJNzruIlOI7G7kd6B9CSMFO+gOMJ3QIeAyYDe+YjM2kjDg4HtgN6SNou47AvCcuc7AQMAO5K27fCzHaJXh5EnduAJHkV0biB9AjConKfRZ9nmNkHZtYLeJ0wX2k+rBtxYGargNSIg3XM7H0zWx59HEHouuCc25AlfGLnuIG0CTAtWt54GZC+POJb/P7IX1pxRhykOxd4O+1zzWgkwwhJx+U6SVKv1KgHW7mkdDl2zpW5MLFzvFdFiNvYNA1oFL0fDxzF76Oc9iCs6ZQPcUYchAOl04FOwAFpm1ua2QxJbYH3JH1rZuutR2xm/YB+EFYRLX22nXNlrUqCK0njBtLBwCHAq8C9wJOSdgNWAvsD/8xTfuKMOEDSIcB1wAFmtjK13cxmRD8nSRoG7EqYzd85V8klOI7GDqRXAbUBzOxpSUuBkwgri/YG+uYpP0WOOJC0a3S9bukzskhqACw3s5XRSIR9KNgQ5ZyrpKQNoPtT1LizPO3zq4TSaV7FHHHwD2BT4KXoi011c+oI9JW0llD32yd9fkHnXOVWQdWfsSRuZFOMEQeH5DjvE2DHss2dc66iVFRDUhyFzf40l2LMM2pmhc4g7ZxzJSVCy31SFVYifQifsNk5lxAJLpAWOvvTzeWYD+ecy60CRy3Fkbg6UuecyybBcdQDqXMu+cSG0SHfOecqVKVstXfOuaSoyAlJ4vBA6pyrFPzR3jnnSim5YbTwDvlrKV6H/LzMkO+cc9lU1u5Pl/B7IK0G/AVYSpjIeQ5hjtJjgTrkb/Yn55xbT2i1r+hc5FZYh/x169pLuocwO/7J0ep6qe1XAy8Bbcoyk865jZwqbtLmOOLOkH8m8J/0IArrliz9D3B6vjPmnHPpNoQ1m6oSpqnLZvtipOOcc8WWerSP86oIcVvtnwX+LmkTYCChjrQxoY70VsJqos45V2Yqa2NTusuB1YSgeWfa9pWE2eqvzHO+nHOugOSG0fgz5K8CLpN0G7ATocV+FvCtmS0ow/w55xwSVE1wY1OxOuRHQXNY2WTFOedyS/KjfexGIkk7SXpB0kRJKyX9Idp+u6TD85UhSd0kjZU0Iepelbn/LElzJX0Vvc5L29dT0vjo1TNfeXLOVbzUePuiXvHSKnmcySZWII0C5WigKfAUoYN+ykrg4njZL/I6VQkz8x8ObAf0kLRdlkNfMLNdotej0bkNgZuAPYDOwE3RyqLOuUpOiCqK9yoyrVLEmVzilkjvAJ4wswOA2zP2fQXsEjOdonQGJpjZpKhe9nlCz4A4DgMGm9kCM1sIDAa65SlfzrmKFLM0GrNEWpo4k1XcOtJtgSui95nj7xcDDUuTiTTNgWlpn6cTSpiZTpS0PzAOuMzMpuU4t3m2i0jqBfQCaNGyJeOe8VqAyqTB7r0rOguuGFaOn56XdIpRR9pI0qi0z/3MrF/a59LEmazilkjnAG1z7NsemBoznaJk+6YyA/cbQGsz2wkYAjxZjHPDRrN+ZtbJzDpt0WiLEmfWOVc+BFSVYr2Aean/v6NXvyzJZYobZ7KKG0ifB26VtG/6hSW1B64idNjPh+lAi7TPWwEz0g8ws/lmtjL6+B9gt7jnOucqrzyObCpNnMmet5j3cAMwChjO76XP14HvgG+Av8dMpygjgXaS2kiqDnQnjKRaR1KztI/HAD9E7wcBXSU1iBqZukbbnHMbgDwG0tLEmazidshfCaRP9nIAABwdSURBVBwl6WDgYKARsAAYamaDY2U93nXWSOpNCIBVgf5m9r2kW4FRZjYQuETSMcCaKA9nRecuiAYMjIySu9UHCzi3YQgNSfnpR1qaOJNLrEAqqSUw08yGAkMz9m0CbGlmeaknNbO3gLcytt2Y9v4a4Joc5/YH+ucjH865ZMnnwKbSxJmseYt53GRg1xz7do72O+dcmclnh/x8i9v9qbDs1SR0ynfOuTIhYJMEDxEtbM2mnSjY0f4ISdtmHFYT+COhn5VzzpWZBMfRQkukxxOGXELoY3VjjuMmA/+Xz0w551w6xRz+WVEKC6R/B+4mlKoXAwfxe4t4yiozW11GeXPOuXUSHEcLXfxuNWEyZ/ClRJxzFSzB05HG7v50CaGLU7bppu4Afk5fddQ55/JJJHti57glzYuACTn2jYv2O+dc2Yg5qinpi9+1IncgnQy0zktunHMuByV41aa4JdKFQIcc+zoQGqOcc65MJH055riB9A3gZkk7pm+UtAOhi9Tr+c6Yc86lS3Igjftofw2wN/ClpC+BmUAzwrDR74D1GqGccy6fKv3id9EsSrsDfwImArWinxcCe0RLezjnXJkIyzHHe1WE2Msxm9mvQN/o5Zxz5aqyjmxyzrlESDU2JVVhk5bMAQ4zsy8lzSXH+kcpZtY435lzzrmUBBdICy2RPgTMTntfaCB1zrmyI6okuB9pYWPtb0l7f3O55MY557IQlbdEWiEkdQPuJ6yl8qiZ9cnYfy9wYPSxNtDYzOpH+34Dvo32TTWzY8on1865MiXYJMGVpIXVkb5XnITM7KDSZkZSVUI1wqGEJVNHShpoZmPSrnNZ2vEXU3AJlBVmlj4ZtXNuA5D0Emlhva7mZ7zaA/sRSoFLo5/7Au2AeXnKT2dggplNMrNVwPPAsYUc3wN4Lk/Xds4lWJVocueiXhWhsDrSk1PvJZ1LGFO/d/pqodHqov8D8rUkc3NgWtrn6cAe2Q6U1ApoA6SXnGtKGkVYQrWPmb2W49xeQC+AFi1b5iHbzrmyVllLpOmuA27MXHI5+nwTcG2e8pPtq8rVW6A7MMDMfkvb1tLMOgGnAvdJ2jrbiWbWz8w6mVmnLRptUbocO+fKnAjBKs6rIsS9blOgRo59NYB89SGdDrRI+7wVMCPHsd3JeKw3sxnRz0nAMHIvIe2cq0yU7Ef7uIF0GHCnpE7pGyXtDtwJDM9TfkYC7SS1kVSdECwHZh4kqQPQAPg0bVsDSTWi942AfYAxmec65yqfMLIpuYE0bvenXoSA9pmk2cAcQim0CfBNtL/UzGyNpN7AIEL3p/5m9r2kW4FRZpYKqj2A580s/bG/I9BX0lrCH4g+6a39zrnKLcFVpPECqZlNB/4g6QjCLFBNgVnASDN7K58ZitJ7K2PbjRmfb85y3ifAjpnbnXMbhiQ3NhWrQ362IOecc2VPlX8+UgBJNSRdKOkxSYMktYu2nyKpY9ll0Tm3sct3q72kbpLGSpogKefE9JJOkmSZ7UOZ4i7H3J7QV3QzYDTQBagb7d4POBI4M05azjlXEvlqSIozgjI6ri5wCfBZkXmLee0HgKmE1UIPo2C973DCCCfnnCsbCkuNxHnFEHcE5W3AXcCvRSUYN5DuB9xhZr+wfgf52YT1m5xzrkwU89G+kaRRaa/MXkXZRlA2L3A9aVeghZn9L07+4jY2/UpYpymb5sAvMdNxzrkSKUZj07xohGPOpLJsW1dAlFQFuBc4K+4F45ZIBwPXStos/cJRB/iL8ZZ851wZU8xXDEWNoKwL7AAMkzQF2BMYWFiDU9wS6V+Bj4EJhKBqwI3A9kB14ISY6TjnXLEJqJq/7k/rRlACPxNGUJ6a2mlmi4BG664tDQOuMLNRuRKMuxzzNGBn4BFCg9NEQr3oS8BuZjarmDfinHPFIsV7FcXM1gCpEZQ/AC+mRlBKKtFk8EWWSCVVI7RyTTazG4AbSnIh55wrOaE8DhKNM4IybXuXotKLUyL9jTDnp3e6d85VmHyVSMtCkSVSM1sraTxhghLnnCt3oftT5R8ieh1woySfFMQ5V/5ilkYTWyKNXA9sDnwl6WdCJ/wCHfPNrHOe8+acc+tU1FyjccQNpN8D35VlRhxMnz6d226+kXfffYcF8+fTtFkzjj7mOK674SYaNGhQ5PnLli1j4Ouv8c7bb/LVl18wfdo0qlSpQvv2HTi5ew8u+tPFVK9ePeu5P4wZw99uu5kPhw9j8eLFtGzVipP/2J0rrryaWrVyjcVwzRvX54YLj6TrPtvRcLPazJq3mDfe/4bb+77FL0tWxE5n713acmnPQ9ipfXOabF6PuQuW8P3EmTz83DAGf/JD1nO2bduU6//vCPbr1I56dWoydeYCXho0mrsfH8yvK1fn6xYTIUzsXNG5yE0F50be+Oy2Wyf7+LOc3cPKzaSJEzlw/72ZM2cORx1zLB06bMuokZ8zfNj7tO/QgfeGf8zmm29eaBrvDnqHY486nIYNG7J/lwPZeuttWLhgAW+9+QazZs1iz7325u13h1KzZs0C533+2Wcc3vUgVq9ezfEnnsRWW7Vg2Pvv8cXoUey19z68/e5QatTItdJM+Wuwe++KzgIAbbZqxPtPXE6TzevxxvtfM3bKbDpt34ounTswdvIsDjr7XhYsWlZkOuefvC8PXNudpctXMvC9r/l5zi80b1yfYw/emTq1anDTg29w12ODCpyz+w6teLvfJVTbpCqvDvmK6bMW0qVze3bbvhWffDmRw//vX6xavaasbr1YVo59kbXL55QqDHbYYRf794ChsY49uGOj0UWMbMq7QkukkmoBRxD6js4EhprZ7HLI10bnzxdfxJw5c/jnvQ9wUe+L122/8orL+df993LzDdfxr4cfKTSNJk2a0v/JZzjxpJMLlDyXLFnCYQd3YcSnn/DIvx/i0sv+sm7fb7/9xv+dfzbLly/npVde56ijQze6tWvXclqPP/LaKy/zwP338tcrc840ttG6/5pTaLJ5PS6/8yX+/fzvq+3c+ZcTuOT0g7i599FccvvzhaaxySZVuPXiY1jx6yr2PvVOxv80Z92+u/o3YcRzV3PVuYdx31ND1wXGKlVE31tOp06tGpx0aV/eHP4tEIZQPnvXORx/yK5ccvqB3P14vhb3TYYEP9nnbmyS1JbwSP8S8A/gGWCspK7llLeNxuRJkxgy+F1atW7NBRf9qcC+G266hTp16vDfZ59m2bLCSzc777ILPU49bb3H97p163JJFDw/HD6swL4PPxjOjz/8wL777b8uiAJUqVKFv99xFwCP9nuEjf3JJVPr5ptz6N4dmfLzPB554YMC+27795ssXb6SU4/cndo1s1elpDSsV4f6dWszfuqcAkEUYOzk2Yz/aQ61a1Vn09q/PxHst1s7OrZtxoejx68LogBmxrX3hRXIzztpw5uQTTH/qwiFtdrfBawlzPxUmzAc9Eugb1lmSFJ/SXMkZa2TVfBANCHrN5L+kLavp6Tx0atnWeYzn4a9/x4AhxzSlSpVCv5K6taty15778Py5cv5/LMRJb5GtWrVAKi6ScGHkNS1ux7Wbb1z2rRtS7v27Zn6009MnjSpxNfeEHXp3B6AIZ/+uN4fmaXLV/LpV5OoU6sGnXdqXWg6cxYsYc6CJbRr2ZitWxZcGnyblo3ZpuUWfPXjtAJVBF12D9d+N0vd6ZSf5zNuymxabbk5bbZqtN7+yipVRxrnVREKC6R7Adeb2cdm9quZ/QD8H9BSUllOm/cEsP7/1b87HGgXvXoB/waQ1BC4CdiDMBLrJklFt9AkwLhxYwHYpn37rPu33qYdAOPHjSvxNZ56oj8AXbsW/GrXXbtd9mtvk7r2+JJfe0PUvlXoVj1h6pys+ydG29u1Knql8svueJEqVcQnz17Jf249g1svPoZHbzuDT/57JWMmzeK0Kx8reO3WIc0JP2W/9oSpc2Nfu9KIuYJoElcRbQZkFkMmEv44NCXUmeadmX0gqXUhhxwLPBWtIDpCUv0osHcBBpvZAgBJgwkB+bmyyGc+LV60CIDN6m2Wdf9mm4XtixaVbLbCfz/0IO8Oeoedd96Fnmefk/3am2W/dr0oT7/84jMlpqu3aejJsGhp9pb5RUvDXMCb1a1dZFqvDPmSmXMX8eQdZ3H60Xus2z5r3mKefn0Ek6fPL9a1F0fb69fdsHpbJLiKtMgO+UmsGMs1KWuRk7WmSOqVmvR17ry5ZZbRfFn36FiCv7avvfoKf/3LpTRt2pTnXnx53SN+7GtH/wSSvPBYEq37umLULXc/YnfefKQ3H385kZ2Pv40Ge17GzsffxrDPx3LfNX/k6T5nF/PainvpSiPp69oXFUgHRfWVcyTN4fdS6ND07dG+8pJrUtZCJ2stsNGsn5l1MrNOWzTaItsh5apeqsS5eFHW/YsXLwZyl1hzGfj6a5x5Wne2aNyYQUOG0aZt29zXXpT92ktS185RYt1YpUp9m22avdRXr07oYpar1JiyTcvG9L35NMZMmsU51z/FuCmz+XXlasZNmc051z/F6DFTObHrH9hvt3axr1035rUrmzzOR5p3hT3a31JuuSieXJOyTic83qdvH1ZuuSqF9u07ADAhRx3oxAnjAWiXow41m5cHvMRZZ5xKk6ZNeefd99imXbusx627do460Ampa+eoQ91Yjfsp9ALcpmX2esito+2ZLfGZDtlrW6pX24SPRo9fr9HKzPj4iwnstl1L/rBdCz4cHX4X46aENLfJUQe6TdRoVdS1K50EPxTlLJGa2S3FeZVjngcCZ0at93sCi8xsJmFuwa6SGkSNTF2jbYl3QJcDARgy5F3Wrl1bYN+SJUv49JOPqVWrFp332DNWes8/9196nt6DZltuyeChw3MGUYAuBx4EhM78mSZPmsT4ceNo2apV1tLsxmz4yBDUDtlr2/WqPTatXYO9dmnL8hWr+PybKYWmU6N6KMs0arBp1v2N6oftq1b/tm7bsJHhj17XvdefkK11881p37oJP82Yz+Tp8+LdTCVRmR/ty52k54BPgQ6Spks6V9IFki6IDnmL0Ag2AfgPcBFA1Mh0G2H265HAramGp6Rru/XWHHJoV36aMoVHHn6owL7bbrmJZcuWcdrpZ1KnTp1128f++CNjf/xxvbSeeepJzj3rDFq0bMng9z4oMgDut/8BbNuxIx99+AH/e2Pguu1r167lumuvAuC8Xhd4HWmGydPnMfiTH2jdvBEXnLJ/gX03XHgkm9auwbNvfs7yX1et296+dRPaty44idrHX0wE4PiDd2WHdlsW2LdT++Ycf8gurF27lmGfj123/cPR4/lh0kz2260dRx7w+zxCkrj9z2ExzEcHfJSfG02QJD/a+xDRhA4R3Xbbjoz8/DOGD3ufdu3b8/4HnxQYIlqrWvgns2L177+/4cPe54jDDmHt2rX0POsctmrRYr3rbLZZfS7+86UFtmUOEW3RoiXvvzfUh4gWIXOI6I+TZ7P7DmGI6LgpsznwrHsK9P9c8eWDANTatWD+H7npNHoetxcrV61m4PvfMHXGAlpt2ZCjD9yJGtWr8a9n3uPKf75S4JzMIaLTZi3gwM4dNtghoh133NWeGjgs1rGd29Yv9yGiHkgTEkgBpk2bxm0338jgd99hfsakJQ0bNixwbLZA+vSTT9DrvMJbeFu2asXYCVPW2/7DmDHcdutNfDDsfZYsWULLVq344yk9EjlpSVICKcBWTepzw4VHcejeHdm8fp1o0pKvub3v2yxcvLzAsbkCKcDpR+/BGcfsyY7tm1O3dg0WL/uVr3+czuOvfsJLg0Znvfa2bZtywwVHsH+n9tStU4OpMxfy4jujEjdpST4C6XY77mpPDRxe9IHA7m0380Ba3pIUSF08SQqkrmh5CaQ77WpPxwykndqUfyCNO42ec85VqCTX0nsgdc5VAkp0g6cHUudcpZDgOOqB1DmXfBXZtSkOD6TOucohwZHUA6lzrlKoqEmb4/BA6pyrFLyO1DnnSqMC16yPwwOpc65SSPKjfeImLXHOuUwilEjjvGKlJ3WTNDZa+229JXKjiZK+lfSVpI8kbVdYeh5InXOVQr5mf5JUFXiIsP7bdkCPLIHyv2a2o5ntQlgI9J7C0vRA6pyrHPI3j15nYIKZTTKzVcDzhLXg1jGzxWkf61DEskteR+qcqxSKMWlzI0npMxH1M7N+aZ+zre+2Bxkk/Qm4HKgOHFTYBT2QOucqhWI0Nc0rYvanWOu7mdlDwEOSTgWuB3rmStAf7Z1zlUP+Hu1zrfuWy/PAcYUl6IHUOZd4IUbG+y+GkUA7SW0kVQe6E9aC+/16UvpCZ0cC4wtLMHGP9pL6A0cBc8xshyz7TwOuij4uBS40s6+jfVOAJcBvwJryntzVOVdG8tgh38zWSOpNWByzKtDfzL6XdCswyswGAr0lHQKsBhZSyGM9JDCQAk8ADwJP5dg/GTjAzBZKOhzoR8GK4gPNbMNaPtE5l9fu+Gb2FmEhzfRtN6a9/3Nx0ktcIDWzDyS1LmT/J2kfRxDqN5xzG7RkT+xc2etIzwXeTvtswLuSRkvqVUF5cs6VgXyObMq3xJVI45J0ICGQ7pu2eR8zmyGpMTBY0o9m9kGWc3sBvQBatGxZLvl1zpVc0id2rpQlUkk7AY8Cx5rZ/NR2M5sR/ZwDvEoYwbAeM+tnZp3MrNMWjbYojyw750orf92f8q7SBVJJLYFXgDPMbFza9jqS6qbeA12B7yoml865fMtj96e8S9yjvaTngC6EYV7TgZuAagBm9ghwI7A58HBU+Zzq5tQEeDXatglh0oF3yv0GnHNlIsFtTckLpGbWo4j95wHnZdk+Cdi5rPLlnKtAgioeSJ1zrrSSG0k9kDrnEi81sXNSeSB1zlUKCY6jHkidc5WDl0idc66UkjxE1AOpc65SSG4Y9UDqnKsEKnIcfRweSJ1zlUKS17X3QOqcqxySG0c9kDrnKocEx1EPpM65ykDFWY653Hkgdc4lXtJHNlW6afSccy5pvETqnKsUklwi9UDqnKsUvPuTc86VhnfId8650kl6Y5MHUudcpeCP9s45V0pJLpEmrvuTpP6S5kjKugKopC6SFkn6KnrdmLavm6SxkiZIurr8cu2cK2sJXo05kSXSJ4AHgacKOeZDMzsqfYOkqsBDwKHAdGCkpIFmNqasMuqcK0deIo3PzD4AFpTg1M7ABDObZGargOeBY/OaOedchRBQRYr1qghJLJHGsZekr4EZwBVm9j3QHJiWdsx0YI9sJ0vqBfSKPq6sVS17NUIl1wiYV9GZKCMb6r1tqPfVobQJfPHF6EG1qqlRzMPL/TusjIH0C6CVmS2VdATwGtCO7AV/y5aAmfUD+gFIGmVmncoqsxVlQ70v2HDvbUO+r9KmYWbd8pGXspK4R/uimNliM1savX8LqCapEaEE2iLt0K0IJVbnnCtTlS6QSmqqaBUsSZ0J9zAfGAm0k9RGUnWgOzCw4nLqnNtYJO7RXtJzQBegkaTpwE1ANQAzewQ4CbhQ0hpgBdDdzAxYI6k3MAioCvSP6k6L0i//d5EIG+p9wYZ7b35flZRCDHLOOVdSle7R3jnnksYDqXPOldJGF0glNZQ0WNL46GeDHMf9ljYMNbGNVkUNi5VUQ9IL0f7PJLUu/1wWX4z7OkvS3LTf0XkVkc/iijEEWpIeiO77G0l/KO88lkRphnZvCDa6QApcDQw1s3bA0OhzNivMbJfodUz5ZS++tGGxhwPbAT0kbZdx2LnAQjPbBrgXuLN8c1l8Me8L4IW039Gj5ZrJknsCKKxP5OGEftHtCING/l0OecqHJyj8viAM7U79vm4thzyVm40xkB4LPBm9fxI4rgLzUlpxhsWm3+8A4OBU97EE22CH+8YYAn0s8JQFI4D6kpqVT+5KrhRDuzcIG2MgbWJmMwGin41zHFdT0ihJIyQlNdhmGxbbPNcxZrYGWARsXi65K7k49wVwYvT4O0BSiyz7K6O4914Z7SXpa0lvS9q+ojOTT4nrR5oPkoYATbPsuq4YybQ0sxmS2gLvSfrWzCbmJ4d5E2dYbOyhswkSJ89vAM+Z2UpJFxBK3QeVec7KXmX8fcWRa2j3BmGDDKRmdkiufZJmS2pmZjOjR6Y5OdKYEf2cJGkYsCuQtEAaZ1hs6pjpkjYBNiP5j2BF3peZzU/7+B8qQd1vTBvkUGczW5z2/i1JD0tqZGYbxCQtG+Oj/UCgZ/S+J/B65gGSGkiqEb1vBOwDJHFe0zjDYtPv9yTgPUv+KIwi7yuj3vAY4IdyzF9ZGgicGbXe7wksSlVFVWaFDO3eIGyQJdIi9AFelHQuMBU4GUBSJ+ACMzsP6Aj0lbSW8Avvk8QJos0s67BYSbcCo8xsIPAY8LSkCYSSaPeKy3E8Me/rEknHAGsI93VWhWW4GGIMgX4LOAKYACwHzq6YnBZPKYZ2bxB8iKhzzpXSxvho75xzeeWB1DnnSskDqXPOlZIHUuecKyUPpM45V0oeSF1OkizGq0sF5u95SR8V85xto3znHLRRjLS6RWltU9q0XOW2MfYjdfHtlfa+FvAe8DfgzbTtietf61x580DqcopmHwJA0qbR24np23ORVNPMfi2zzDmXIP5o70pN0gXRI+4fJH0oaQVwca5H32hGrWcyth0o6SNJKyTNk/RvSbWLmY8Wkp6UNDlKZ6ykmyRVy3J4g6hqYKmkWZKuyZLeLpLekbQkmpT4OUlbFCdPbuPggdTl0wvAy4Qhju/GPUnSQdHxU4ATgCuA4yn+6pONgVnApYRJhu8FLgTuznLsfcA84ETCzFF/j4YNp/LUEfgw+ngacB6wG/BqMfPkNgL+aO/y6W4z65v6UIwJie8EhpjZ6WnnzgHekHSLmY2Pk4iZjQZGR+cL+BhYBdwv6XIz+y3t8NFm1jt6P0jSlsD1hLkJAG4lBPajonlckTQG+FbSwWY2NOa9uY2Al0hdPr1Z9CEFSapPKOm9KGmT1AsYHh0Se80iSVUk/VXSj4SJMVYTAuOmQGZQzyxZvgK0lpSa6PuQaBtpeRoLzAQ6xb5Bt1HwQOryaXYJztmcMJlxf0LgS72WEv59Fmfm+6uAvxOqGI4mLFlyWbSvZsaxmfPQpj43U1gzqj5wY0aeVgNbFjNPbiPgj/YunzKnEku12lfP2N4w7f3C6Oc1wJAsaU4vxvVPBp41s5tSG5R7Fc7MJWZSn2ea2W+SFhMWdHs6y7lZJwN3Gy8PpK4spYJgR6L+ppK2BtoCnwOY2QJJXwLtzKxPKa9XC1iZse20HMceDzye9vkE4CczSwXJocD2ZjaqlHlyGwEPpK7MmNkESd8Cd0QT+lYHrmX9mdH/CrwtqQqhXnIZ0Bo4CrjMzH6KecnBwLmSvgB+IqwMsFWOY3eT9C/C2k8HA6cTlj9OuQEYIWkgoWS6IErrMODfZvZJzDy5jYDXkbqydgqh7vS/wC2EBQgnpx8QtYAfSAhUzxKW2/gLMIniLUdxPSEQ94nSWUToSpXNpYQFEl8hzK5/g5n9Jy1P3wN7Ar8RGqzeJsz6viQz/875DPnOOVdKXiJ1zrlS8kDqnHOl5IHUOedKyQOpc86VkgdS55wrJQ+kzjlXSh5InXOulDyQOudcKf0/Nax3KZxnksAAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn.metrics import confusion_matrix\n",
    "import matplotlib.pyplot as plt\n",
    "import itertools\n",
    "\n",
    "predicted_labels = model.predict(X_test)\n",
    "cm = confusion_matrix(y_test, [np.round(i[0]) for i in predicted_labels])\n",
    "cm = cm.astype('float') / cm.sum(axis = 1)[:, np.newaxis]\n",
    "\n",
    "plt.imshow(cm, cmap = plt.cm.Blues)\n",
    "plt.title('Normalized confusion matrix', fontsize = 15)\n",
    "plt.colorbar()\n",
    "plt.xlabel('True label', fontsize = 15)\n",
    "plt.ylabel('Predicted label', fontsize = 15)\n",
    "for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):\n",
    "    plt.text(j, i, format(cm[i, j], '.2f'),\n",
    "             horizontalalignment = 'center', verticalalignment = 'center', fontsize = 20,\n",
    "             color='white' if cm[i, j] > 0.5 else 'black')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(20, 25000, 10)\n",
      "(20, 25000, 10)\n",
      "words,sentences,letters\n",
      "50000\n",
      "50000\n",
      "(40000, 20)\n",
      "(10000, 20)\n",
      "1025\n"
     ]
    }
   ],
   "source": [
    "sentence_length= 200\n",
    "word_length = 10\n",
    "\n",
    "H_DNA = np.array(H_DNA[0:(len(H_DNA)//sentence_length)*sentence_length])\n",
    "N_DNA = np.array(N_DNA[0:(len(N_DNA)//sentence_length)*sentence_length])\n",
    "\n",
    "#divide in sentences\n",
    "H_DNA_s = np.array(np.split(H_DNA, len(H_DNA)//sentence_length))\n",
    "N_DNA_s = np.array(np.split(N_DNA, len(N_DNA)//sentence_length))\n",
    "\n",
    "#divide in words\n",
    "H_DNA_w = np.array(np.split(H_DNA_s, sentence_length//word_length, axis=1  ))\n",
    "N_DNA_w = np.array(np.split(N_DNA_s, sentence_length//word_length , axis=1 ))\n",
    "\n",
    "print(H_DNA_w.shape)\n",
    "print(N_DNA_w.shape)\n",
    "print(\"words,sentences,letters\")\n",
    "\n",
    "#make words strings.\n",
    "Human_DNA = np.array([[''.join(H_DNA_w[i,j]) for i in range(H_DNA_w.shape[0]) ]for \n",
    "                  j in range(H_DNA_w.shape[1])])\n",
    "#make words strings.\n",
    "Neanderthal_DNA = np.array([[''.join(N_DNA_w[i,j]) for i in range(N_DNA_w.shape[0]) ]for \n",
    "                  j in range(N_DNA_w.shape[1])])\n",
    "\n",
    "depl_texts = [' '.join(Human_DNA[i,:]) for i in range(Human_DNA.shape[0])]\n",
    "intr_texts = [' '.join(Neanderthal_DNA[i,:]) for i in range(Neanderthal_DNA.shape[0])]\n",
    "\n",
    "merge_texts = intr_texts + depl_texts\n",
    "print(len(merge_texts))\n",
    "labels = list(np.ones(len(intr_texts))) + list(np.zeros(len(depl_texts)))\n",
    "print(len(labels))\n",
    "\n",
    "tokenizer = Tokenizer()\n",
    "tokenizer.fit_on_texts(merge_texts)\n",
    "\n",
    "encoded_docs = tokenizer.texts_to_sequences(merge_texts)\n",
    "max_length = max([len(s.split()) for s in merge_texts])\n",
    "X = pad_sequences(encoded_docs, maxlen = max_length, padding = 'post')\n",
    "\n",
    "from sklearn.model_selection import train_test_split\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, labels, test_size = 0.20, shuffle = True)\n",
    "\n",
    "print(X_train.shape)\n",
    "print(X_test.shape)\n",
    "y_train = np.asarray(y_train)\n",
    "y_test = np.asarray(y_test)\n",
    "\n",
    "vocab_size = len(tokenizer.word_index) + 1\n",
    "print(vocab_size)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 106,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Train on 32000 samples, validate on 8000 samples\n",
      "Epoch 1/5\n",
      "32000/32000 - 13s - loss: 0.5304 - accuracy: 0.7367 - val_loss: 0.5189 - val_accuracy: 0.7419\n",
      "Epoch 2/5\n",
      "32000/32000 - 9s - loss: 0.4952 - accuracy: 0.7592 - val_loss: 0.5150 - val_accuracy: 0.7441\n",
      "Epoch 3/5\n",
      "32000/32000 - 9s - loss: 0.4875 - accuracy: 0.7629 - val_loss: 0.5240 - val_accuracy: 0.7473\n",
      "Epoch 4/5\n",
      "32000/32000 - 9s - loss: 0.4822 - accuracy: 0.7689 - val_loss: 0.5150 - val_accuracy: 0.7477\n",
      "Epoch 5/5\n",
      "32000/32000 - 8s - loss: 0.4791 - accuracy: 0.7688 - val_loss: 0.5129 - val_accuracy: 0.7429\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<tensorflow.python.keras.callbacks.History at 0x7fd5b013ae50>"
      ]
     },
     "execution_count": 106,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.fit(X_train, y_train, \n",
    "                    epochs = epochs, verbose = 2, validation_split = 0.2, shuffle = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 107,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVIAAAEdCAYAAABe2+9WAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dd5gUVdbH8e8PJAcBiSJRAUFWURBFUTFhBrOYWVFMiDmLoq6K7prWVwVU1ixmwRUlCkZUXDGgklGRJEGyxPP+UTXQ03TPVE/3zPQM5+PTz3RXuHWrB8/culFmhnPOuYIrU9wZcM65ks4DqXPOpckDqXPOpckDqXPOpckDqXPOpckDqXPOpckDKSCpvySTNDLBvjckjS+GbBVIeB99Yj6Pl/RGEV17saT+RXGtKCRVkTRU0pLwe+mZoXSbhukdn4n0ipukjqn83iR1Ce+/bSFmq0TZobgzkGW6StrXzL4q7oxk0GXAhuLORDG5FDgBOA/4HZiZoXTnA52AnzOUXnHrCNwB9I94/P8I7j9T32eJ54F0q6XAXOBW4MRMJy6pkpmtzXS6+TGzH4v6mllkd2Cqmb2ZyUTNbB0wMZNplgSSBFQwsxVsh/efF3+038qAe4Fukv6W14GS2kkaK2mNpGWSXpJUL2Z/zqPf2ZKel/Qn8G64b46kf0m6SdJ8ScslPajAsZKmSFop6R1JNWPSrCLp/yRNDa87W9Ljkqrnk9dcj/ZhvhK9esYc01nShPA6SyQ9JalaXLoHS/pW0l+SvpZ0QJQvWVJZSTdLmiZpnaS5kp6NO6aPpOnh/hmSro7b3z+sRthb0sQwn99IOijmmDlAL2DvnHsMtz8raVJcets8qkvqFf4u1obXmiBpjzyOLxvm69cw31MknRV3nWclTZJ0pKTvJK2W9ElOunl8ZzmP0odLGhaeN11S1/C6/wzz+Luka+LO7SRpuKR54XmTJZ0ds78n8Fj4Puffwvi477mzpK+Av4DTFPdoL+k0SZslHR73na6Q9I+87q208ECa2+vANIJSaUKS6gDjgcrAWcAVwCHAaEnl4w7/F7ASOI0gSOfoQfA49XfgAeAa4CHgbqAfcEmY5n0x51QGyoZ5OyY87rAwz6noFPd6jOCPyIzw/g4ExgILgFOBq4Bjgf/EfAc7A+8TlOJPBQYBL4V5zM8g4E7gNeB44FqgSkzaF4V5Gk7wWP468KCkm+LSqQw8F6Z3CrAOeFtSTh5OAkYQPH7n3Gskkg4GBgIvEnzXFwCfATvmcdpdBL+bwUA34FPgJUlnxh3XGPgncA9wJlAXeE2SImRtEPAJwb39ArwB/B9QjeDf4hsE39X+Mec0CfNyIcH3+Sbwn5h8vQc8GL7P+Z4uizk/53t+Gjga+DI+U2b2OvAqMERS9fBehgCzw++l9DOz7f5FUDe0OHzfE9gEtAw/vwGMjzl2APAnUD1mW0eCYHRm+Llp+PntBNeaQxC0ysZs+xLYCDSL2fYAsDCPPO8AHBhep3HMdgP6xHweD7yRJI19gLXAHTHbPgY+jDvusDDdtjF5WwJUjjnm7PCY/nnkeffwmL5J9pchqMv8T9z2J4DlQMWY35cBh8Uc0y7cdnTMtmeBSXFpJdqW8/s6Pvx8HfB1HvcRf3wtYHXs9xhuH0FQtRB77Y1Ai5htJ4Zp7Z7H9bqEx8T+ntqE28bFfX8LgPuTpKPw382guPP6AJbk/wsDuifJT9uYbbWAecAzQF9gPbBXpv4fzfaXl0i39SLwK3Bzkv0dgVEW1BMBYGZfEgTIznHHvpckjfFmtinm8wxgjpnNjttWJ7aUK+nc8BF2FUED0ifhrpZ539K2wpL128AYwlJDWJrrRFBC2iHnFV5nA9A+PL0jMNrM1sQk+VaEyx4a/nw2yf5dgJ3ZtpT9KlAdiK1y2UDwRyJHTl3wLhHykZ/JBFUCD4dVGPFPGvHaEpTcEuW7paS6MdvmmNn0mM+p5HtszPsZ4c9xORvMbDMwC2iYs01STUn/lvQLwXe2AehN9H8zRvD0kfdBZkuBiwhK7/8E7jSzbyNeo8TzQBrHzDYSlLjOkdQkwSENgIUJti8k+Kscvy2RP+M+r0+yTUB5AEknAc8DnxNUFexP8IgHUDHJdRIKg+Nr4TXOsbBIAdQkqD54gq3/020geGwuBzQKj6sPLIpN04KGtFX5XHonYHXsH6E4DcKf8d9bzufY73dFGDhyrr8+fJvSd5GImY0hqHY5mCBYL5b0hKQqSU7JL981Y7Yl+j1DtHxvOTfmfhOlF5vWs8AZBMGtK7AvwWN31O9pWcy18jOO4J7LAE9FPKdU8Fb7xIYAtwE3Jtg3n6BeK1494Ou4bZmco/A04Asz21J/JemQAqb1L4L/ofY3s+Ux2/8kfDwneCyNNy/8uYC470BSJaBqPtddAlSRVD1JMJ0f/oz/fnMa8pbmk34UfxH+cYoR/wcQM3sOeC4suZ8MPAysAOLraiF3vpfEbM9kvlMmqSJwHEFVz8CY7akUoFL5NzyA4A/xAuARgnrb7YKXSBOwoHvLvwgeUxrE7f4COEoxrdiS9iWoN/uEwlOJoGQY6+xEB+ZF0rnAlUAvM/shdp+ZrSbo1tLKzCYleOUE0q+AI2MadiAINvnJeQw9L8n+uQTB+rS47acTBLHvI1wjP3OBpmGQyXFksoPN7A8zG0RQd9wmyWE/AGtInO9pZvZHGvlNRwWCwLbl303477Zb3HHrw30FLs2Hf9SvIOi72ws4U9IpBU2vpPESaXKDgFuAA4AJMdsfIvjHMlLS/QSlsAEE/5NntL9inNHA45JuJQjmxwKH531KbpJ2JWhVfh/4Ja51d2b4P/wNwFhJmwka2lYStDQfB9xqZtMIShuXA/+V9BBBvebNBA1XSZnZVEmDCVqW6wIfATWAU82sh5ltVjDCZpCkJeE9H0Lwfd9iZn+lcr9JvENQJ/y0gm5XexM8xm8h6U6CUup4YHF4zCEkLo1iZkslPQLcJmkjMIngD8uxBC3zxcLMlofdlm6XtALYTHAPywnqnHPkDCy4UtI4gmqTqVGvI6kqQa+OV83sjXDbIOBJSR8V4x+SIuMl0iTChpSHE2z/g6DR5C/gFeBxgtLKkSnUJRXEIIJuKlcSNOw0IfVHp0YEdWPHENS1xr6OAzCzTwjqBusALxD0f70B+I2wzs/MficIErUJ/nhcBpxDUCrLz2UE3Z/OIag+eISYAGxmTxG0+p4E/JcgEF1rZgNSvNeEwlL4BQSNasMJAuQFcYd9RVD6HAiMJAjk/YFH80j6doLuapeG+T6YoP55aCbynYazCLohPU+Q/zfD97E+JqhDvZLgj/SgFK/xIMETU5+YbdcR1JkPTHhGKaOt7QzOOecKwkukzjmXpqwJpJJqSRodDn0brZjhkXHHbQqHuU2WNDxmezNJX4Tnvxqh759zzmVE1gRSgkrwsWbWgqDjccKKfWCtmbULX7Gtj/cDD4fnLyNoOXTOuUKXNXWkkqYCXcxsvqQGBKN/WiU4bpWZVY3bJuAPoL6ZbZTUiWCo4lFFknnn3HYtm7o/1TOz+QBhME3U6R2gooLZezYCA8zsHYIRM3+Go5Ig6CvYMMn5SOpNMEwOylZor6r1kh3qslDbprWLOwsuBXN/+4WlSxZHmZQlqbLVm5htjDYLpa39Y6SZHZ3O9VJVpIFU0hiC4YXxks62lEBjM5snqTkwTtL3BJ214yUtapvZYIL+lJSp0cQqHJJsWL3LRu8O6VncWXApOOHwA9NOwzaupUKr0yMd+9fkx4v8L22RBlIzOyLZPkkLJTWIebRflOi4nNE1ZjYrnDdxb4K+cTUk7RCWSndh63BG51yJJ0hpZGvRyqacDQfOD9+fDwyLPyCcyaZC+L42wTRyP4aTbnxIMDdm0vOdcyWUgDJlo72KQTYF0gEE47enE4x9HgAgqYOkp8NjWgOTJH1LEDgH2NalNG4ErpE0g6DO9Jkizb1zrnBJ0V7FIGsam8xsCQnGjpvZJILZvTGzz8g9J2XscbMI5sl0zpU62f1onzWB1Dnn8lRMpc0oPJA657Kf8BKpc86lp/jqP6PwQOqcKxmKqUU+Cg+kzrkSwBubnHMuPcIf7Z1zLm1eInXOuXT4o71zzqVHQFlvbHLOufR4HalzzqXDH+2dcy59XiJ1zrk0eYnUOefSUIxT5EXhgdQ5VzL4EFHnnEuHNzY551z6/NHeOefSkOXzkWZNziTVkjRa0vTwZ80Ex7ST9LmkKZK+k3RGzL5nJc2WNDl8tSvaO3DOFZ7w0T7KqxhkTSAFbgLGmlkLYGz4Od4a4Dwz2wM4GnhEUo2Y/debWbvwNbnws+ycKzK+imgk3YHnwvfPASfGH2Bm08xsevh+HrAIqFNkOXTOFZ8sXkU0mwJpPTObDxD+rJvXwZI6AuWBmTGb7wkf+R+WVKHwsuqcK1LK7kf7Im1skjQGqJ9g160pptMAeAE438w2h5tvBhYQBNfBBOvc35Xk/N5AbwAq1Url0s654uKt9gEzOyLZPkkLJTUws/lhoFyU5LjqwHvAbWY2MSbt+eHbdZL+A1yXRz4GEwRbytRoYqnfiXOuqCmDgVTS0cCjQFngaTMbELf/YeDQ8GNloK6Z1SCJbHq0Hw6cH74/HxgWf4Ck8sDbwPNm9nrcvgbhTxHUr/5QqLl1zhWZYKURRXrlm5ZUFngcOAZoA5wpqU3sMWZ2dU7DNfAY8FZeaWZTIB0AHClpOnBk+BlJHSQ9HR5zOnAw0DNBN6eXJH0PfA/UBv5RtNl3zhUaCZWJ9oqgIzDDzGaZ2XpgKEFjdzJnAq/klWDWdMg3syXA4Qm2TwIuDN+/CLyY5PzDCjWDzrlilcKjfW1Jk2I+Dw6r83I0BH6L+TwX2C/JNZsAzYBxeV0wawKpc87lJYVAutjMOuSVVIJtydpKegBvmNmmvC7ogdQ5VyJksLFpLtAo5vMuwLwkx/YALs8vwWyqI3XOucSUwit/XwEtJDULG7B7EDR2576k1AqoCXyeX4IeSJ1zWU9Ea7GPUmo1s41AH2Ak8BPwmplNkXSXpG4xh54JDDWzfLtI+qO9c65EKFMmc+U+MxsBjIjbdnvc5/5R0/NA6pwrETLZIT/TPJA657Jf9PrPYuGB1DlXIniJ1Dnn0pDT2JStPJA650qEiMM/i4UHUudc9pM/2jvnXNo8kDrnXJo8kDrnXBq8sclF1nCnKvQ7a1+67tOIWtUqsmDpGt79Yjb3DJ3En6vX53v+QW13ZtQ93fI9rkWvF5i7eDUAO9eqQvdOzTiqfWN236Um9WtVZtXaDUyetZin3p/CsImz076v0mz+vLk8NOBuJowdxZ/LllKnXn26HnsCV11/KzvW2GZF8W2sWb2aUSOGM270B/zw3WTmz5tLGZWh+W4t6XbyaZx/0WWUL19+m/M2bdrEu2+9xkvPPc3sWTNZtXIFDRo0pP1+neh9+VW03L1NgquVcNkbRz2QZotm9avz4f0nUq9GZd6dOJupv/9JhxZ16dNtT47cpxGH3fQOS1euyzONXxat5B+vTEq4r22TWpx4QHOm/LJ0SxAFuPT4tlx3yt7MXrCCCT/8zsJla2hcpxrdOzXj8Ha78O9h33LjkHznbNgu/TJ7FqcceyiL/1jEkcccz64tWvHt/ybxn0GPM2HsaN4cMY6atXbKM40vJ37KVZdeQI2atejU+WC6HnsCy5ctY8zIEdxzx8188N4wXnrrfSpWrJjrvL69z+e9YW/SYOeGHH1cd6pUrcrUn6bw5tAXGf7mqzw7dBgHHNylEO++iCmzQ0QzzQNplnj0koOoV6My1wz+hCff27pKyv0XdKJv973of05H+j75cZ5p/LpoJfcMTRxIn7s2mDN7yKgfc22fNG0RR94yjE+mzM+1vdUuNZjwwEn07b4XQydM55uZiwtyW6XabTdcyeI/FtH/vgfpedFlW7bffdsNPDPwMf55T3/uffCxPNOoU7cejzw5hGO7n5Kr5HnrypX06H4UX385kReeGchFl1+1Zd+3/5vEe8PepOXubRg26mMqVa68Zd9rLz/PDX0v5rGHBpSuQEp215Fmb4jfjjStV40j927EnIUrGDgi91JTd78yiVVrN3BWl5ZUrlCwv3u1qlWg2/7NWLNuAy9/OC3XvmETZ28TRAGmzv2TNz4JVro+uO3OBbpuafbrnNl8/OEYdmnchPN6XZJr39U39qNylSq8/frLrFm9OkkKgT3+thcnnnbmNo/vVatV46LL+gIw8dPcf0B//SWobjng4C65gihA12OOB2DpklL4hy9z0+hlnAfSLNBlz4YAjPlmLvETdq1au4HPf15AlYrl6NiqXoHSP/ewVlQsvwNvfTorUl1rjg2bgpWuN27yhVbjffbxeAAO6nLENo+cVatVo33HTqxds4Zvvv6ywNfYoVw5AMruUDbX9pz6z88/nsBfa9fm2jd21PsAHHjwoZQ2mZpGrzBkXSCVdLSkqZJmSLopwf4Kkl4N938hqWnMvpvD7VMlHVWU+U5Hy4bBKq8z5v2ZcP/MecsBaLHzjgVKv+eRrQF4ZuSP+Ry5VbVK5TixU3M2bzbGTP4t/xO2M7NmBCX75rvulnB/s+a7hsdNL/A1Xnv5eQAOOaxrru2tWu9Br0uuYOpPUzi80170u+EqBtx1G73OOoUbr7yEE046jetu6V/g62ajqEG0uAJpVtWRxiyTeiTBcgBfSRpuZrERoBewzMx2k9QDuB84I1xOtQewB7AzMEZSy/zWWskG1SsHj3XL1yQuLeZs37FKhZTT7rxHA3ZvVJMpvyxl4s8LI5/3ZJ8u1K9ZmYEjfmDq3MQBfnu2csUKAKpVT/zHLWf7ihXLC5T+c08/yYSxo2jzt704/ezzt9nf7x8P0Hy3ltzd7wZeGDJoy/a/7bUPp/Q4h8pVqhToutnM60iji7JManfgufD9G8Dh4Vr23Qlms15nZrOBGWF6Jd7Wfz6pP2L3Oip4DEylNHr/BZ04pfOufDJlHjc+81nK13SQM6l6Qf7f/+C/73DXrddTp259Bv7nFcqFj/ixafe/+Rpuv/Eq+l53M59/N50pc/7g9f+OQYKeZ3Tn+WcGZuI2skoGl2POuKQlUknHppJQOON0uqIsk7rlGDPbKGk5sFO4fWLcuQ0TXURSb6A3AJVqZSDb6VmRU+KsvG1/QYgpsaZQvwlQs2oFTuwUNjKNn5b/CcC9Pfenb/e9+PiHeZx09wjWb9yc0jW3F9WqVwdgZZIS56qVYYm1WmrVMSNHDOeKi85jp9p1eOWdkTRu2mybY94Y+iLPPvUkvS65gsuuvH7L9n33P5CnX3qTgzu04f67+nHKGedQpWrVlK6fzbK5RJrXo/1/CYpAUXJvQNl8j8pflGVSkx0TeYnVcI3rwQBlajQp9paUab8Hj8677Vwj4f5dw7rR6fNSe0w8J2xkemHs1EhB+IFeB3BFtz0Z/93vnHz3+6xdvzGl621Pmu/WEoBZM2ck3D971szwuBaR03xv2JtceXFP6tStx8tvf0CzJPWv40YFZZZOnQ/ZZl/devXZdbdWTPl+MrNmTONv7faJfP2sVoInLdn2T2Hhi7JMas4xcyXtAOwILI14blaa8H2QzSP23gWJXC33VSuVo9Pu9VmzbgNfTo1exwnw9667A9v2HU3k4Ys7c8mxbRnzzW+cdu8H/LU+66uWi1VOEPt4/Bg2b96cq+V+1cqVfP3l51SsVIm920erXRr2xlCuufxC6jfYOWlJNMf6dcEfxWRdnJYu+QOAcglGRJVUomDVJEUlaR2pmf2SyitD+YmyTOpwIKf2/VRgXLjK33CgR9iq3wxoARS870kRmr1gBaO/+Y2m9apzybFtc+3rd2YHqlYqx0sfTmPNuq0lxJYNa2xp7U/kwDb1ad2oFj/MWZJvI9Pjlx/CJce25YNJv3LqPR5Eo2jSrDkHHXoEc3/9ZZv6yIfvv5s1q1dz8uln52r0mTF9KjOmT90mrTeGvsjVl/Vi510a8eq7o/MMogD77n8AAE8/+e9tGrNe/M9TzJ/3O3Xq1qdFq9YFvb0sVEpa7SVVAC4AOhCU/C43s+mSzgC+M7Of0s1MWOeZs0xqWWBIzjKpwCQzGw48A7wgaQZBSbRHeO4USa8BPwIbw/yVmIhw5cCP+fD+E3mod2cO3bMhP89dxr4t69Flz4ZM+/1P+r+Y+2/Ct0/0AKBS98SNChd0DRqZhozK+9dyyxntuaBra9as28B3sxdz3Sl7b3PMd7MX8+4XcwpwV6XbPx54lFOOPZT+N1/Lpx99yG4td2fy11/x+ScTaL5rC66/tX+u44/o1A6AOYu39v387OMJ3ND3YjZv3kynzofw+ssvbHOd6jvuSK9Lrtjy+dwLLuadN1/l5ynfc9h+e3LEUcdRfccd+eG7yXz28XjKli3L3Q88TNmymahtyx5lSvrEzpJaAqMJHqO/BroA1cLdBwHHAedlIkP5LZNqZn8BpyU59x7gnkzko6jNXrCCzte+Sb+z9uXIvRtxVPvGLFi2hsff/Z57hk5i2aq8x9nHqlGlPCcd0DxSI1PTekGjSeUK5bjhtMT1aS+MneqBNIEmzZozfMwn4aQloxk/ZiR169WnZ+/LuOr6W6lRM/+GzN/n/srmzUGD3msvPZfwmIaNGucKpFWqVuXN98bx9JP/5oP3hjHsrVfZsH49tXaqzXHdT+aiy6+i3T77ZuYms4Wy+9FeFj+UJtFB0gdAFeAEYBWwHuhgZv+TdBpwv5k1L9ScFpIyNZpYhUNuLu5suBT8PKRncWfBpeCEww/ku8lfpxUGKzVoac3+/n+Rjv3pvqO+NrMO6VwvVVEf7Q8CTjOzP8NO87EWAg0ymy3nnMstm0ukUQPpX0ClJPsaAj70xTlXqLK5+1PUkU2jgVskxfYutrAB6gri6jSdcy6jwjrSKK/iELVEej3wKcGwy9EEHd1vJxjXXh44uVBy55xzBEuNZPPEzpFyZma/AXsBA4GmwEyCetHXgfZmtqCwMuicc5DdJdLIId7MlplZPzM7wMxamtn+ZnarmS0pzAw65xxkdj7S/KbrDI85XdKPkqZIejmv9FKaRk9SDaAtQWl0HjDFzLyhyTlXuDJY2owyXaekFsDNwIFmtkxS3bzSjNohfweCju6XA7FrG6yR9ARwq5ltSOlunHMuomCsfcae27dM10mQbs50nbGTUlwEPG5mywDMbFFeCUZ9tH8IuBK4F2gD1A5/3gf0BR6Mfg/OOZe6DNaRJpquM37KzZZAS0mfSpoo6ei8Eoz6aH8ucIuZPRSzbSlwj6S/gNsIAqpzzhWKFMba15YUu5zu4HDqzBxRptzcgWDioy4EM8l9LKltsqrMqIF0MzAlyb4fEmTCOecyJ7X5SBfnM0Q06nSdE8Mqy9mSphIE1q8SJRj10f4F4MIk+y4CXoyYjnPOpSxnPtIMPdpHma7zHeBQAEm1CR71ZyVLMK+lRi6L+TgHOFXSlPCCi4C6BBW01YB/Rcq+c84VSObmGo04XedIoKukH4FNwPV5dfXM69E+0VQrOwOJZot9CHg04n0451zKMtnZPsJ0nQZcE77ylTSQmln2jsdyzm1fVAomdnbOueKU4X6kGZfqyKZdCCpdK8bvy9ByzM45l1CJD6SSqgGvAV1zNoU/Y7s9ZWSBmLDj66Nhek+b2YC4/dcQ9CDYCPwBXJCz+J6kTcD34aG/mlm3TOTJOVf8sjiORu7+dB/QmGCmfAEnEXRUfQaYDeyficzEjIE9hmDk1JmS2sQd9g3BMid7Am8AD8TsW2tm7cKXB1HnSpFsXkU0aiA9lmCs/Rfh53lm9pGZ9QaGEcxXmglbxsCa2XogZwzsFmb2oZmtCT9OJOhM65wrzbJ8YueogbQe8Fu4vPFqIHZ5xBFsfeRPV5QxsLF6Ae/HfK4oaVI4NvbEZCdJ6h0eN8nWr0ovx865QhdM7BztVRyiNjb9RjBRCcB04HiCDqsA+xGs6ZQJUcbABgdK5wAdgENiNjc2s3mSmgPjJH1vZjO3STAYdzsYglVE08+2c66wlcniStKogXQ0cATwNvAw8Jyk9sA64GAyN/tTlDGwSDoCuBU4xMy2LPhuZvPCn7MkjQf2JpjN3zlXwmVxHI0cSG8knIfUzF6QtAo4lWBl0T7AoAzlZ8sYWOB3gjGwZ8UeIGnv8HpHx84RKKkmsMbM1oVjYw8kd0OUc66EUmqTlhS5SIE0bNxZE/P5bYLSaUZFHAP7T6Aq8Hr4xeZ0c2oNDJK0maDud0DsjNfOuZItiwc2Zd/IpghjYI9Ict5nwN8KN3fOueJSIoeISvqDFOYZNbM81zRxzrmCEkHLfbbKq0T6OD5hs3MuS2RxgTTP2Z/6F2E+nHMuuWIctRRF1tWROudcIlkcRz2QOueynygdHfKdc65YlchWe+ecyxbFOSFJFB5InXMlgj/aO+dcmrI3jObdIX8zqXXIz8gM+c45l0hJ7f7Ul62BtBxwLbCKYCLnRQRzlHYHqpC52Z+cc24bQat9ceciubw65G9Z117SQwSz458Wrvecs/0m4HWgWWFm0jm3nVPxTdocRdQZ8s8DnooNogDh56eAczKdMeeci1Ua1mwqSzBNXSJ7pJCOc86lLOfRPsqrOERttX8JuFfSDsBwgjrSugR1pHcRrCbqnHOFpqQ2NsW6BthAEDTvj9m+jmC2+hsynC/nnMsle8No9Bny1wNXS7ob2JOgxX4B8L2ZLS3E/DnnHBKUzeLGppQ65IdBc3zhZMU555LL5kf7yI1EkvaU9KqkmZLWSdon3H6PpGMylSFJR0uaKmlG2L0qfn9PSX9Imhy+LozZd76k6eHr/EzlyTlX/HLG2+f3ipZWweNMIpECaRgovwbqA88TdNDPsQ64Ilr2871OWYKZ+Y8B2gBnSmqT4NBXzaxd+Ho6PLcWcAewH9ARuCNcWdQ5V8IJUUbRXvmmlUacSSZqifQ+4FkzOwS4J27fZKBdxHTy0xGYYWazwnrZoQQ9A6I4ChhtZkvNbBkwGjg6Q/lyzhWniKXRiCXSdOJMQlHrSHcHrgvfx4+/XwHUSicTMRoCv8V8nktQwox3iqSDgWnA1Wb2W5JzGya6iKTeQG+ARo0bM+3NS0ofwhsAAByDSURBVDKQdVdUau7bp7iz4FKwbtpv+R8UQQp1pLUlTYr5PNjMBsd8TifOJBS1RLoIaJ5k3x7ArxHTyU+ibyo+cL8LNDWzPYExwHMpnBtsNBtsZh3MrEOd2nUKnFnnXNEQUFaK9AIW5/z/Hb4GJ0guXtQ4k1DUQDoUuEtS59gLS2oJ3EjQYT8T5gKNYj7vAsyLPcDMlpjZuvDjU0D7qOc650quDI5sSifOJM5bxHvoB0wCJrC19DkM+AH4Drg3Yjr5+QpoIamZpPJAD4KRVFtIahDzsRvwU/h+JNBVUs2wkalruM05VwpkMJCmE2cSitohfx1wvKTDgcOB2sBSYKyZjY6U9WjX2SipD0EALAsMMbMpku4CJpnZcKCvpG7AxjAPPcNzl4YDBr4Kk7vLBws4VzoEDUmZ6UeaTpxJJlIgldQYmG9mY4Gxcft2AHY2s4zUk5rZCGBE3LbbY97fDNyc5NwhwJBM5MM5l10yObApnTiTMG8Rj5sN7J1k317hfuecKzSZ7JCfaVG7P+WVvYoEnfKdc65QCNghi4eI5rVm057k7mh/rKTd4w6rCJxO0M/KOecKTRbH0TxLpCcRDLmEoI/V7UmOmw1cnMlMOedcLEUc/llc8gqk9wL/IihVrwAOY2uLeI71ZrahkPLmnHNbZHEczXPxuw0EkzmDLyXinCtmWTwdaeTuT30Jujglmm7qPuD32FVHnXMuk0R2T+wctaR5GTAjyb5p4X7nnCscEUc1Zfvid01IHkhnA00zkhvnnEtCWbxqU9QS6TKgVZJ9rQgao5xzrlBk+3LMUQPpu0B/SX+L3SipLUEXqWGZzphzzsXK5kAa9dH+ZuAA4BtJ3wDzgQYEw0Z/ALZphHLOuUwq8YvfhbMo7QtcDswEKoU/LwX2C5f2cM65QhEsxxztVRwiL8dsZn8Bg8KXc84VqZI6ssk557JCTmNTtspr0pJFwFFm9o2kP0iy/lEOM6ub6cw551yOLC6Q5lkifRxYGPM+z0DqnHOFR5TJ4n6keY21vzPmff8iyY1zziUgSm6JtFhIOhp4lGAtlafNbEDc/oeBQ8OPlYG6ZlYj3LcJ+D7c96uZdSuaXDvnCpVghyyuJM2rjnRcKgmZ2WHpZkZSWYJqhCMJlkz9StJwM/sx5jpXxxx/BbmXQFlrZrGTUTvnSoFsL5Hm1etqSdyrJXAQQSlwVfizM9ACWJyh/HQEZpjZLDNbDwwFuudx/JnAKxm6tnMui5UJJ3fO71Uc8qojPS3nvaReBGPqD4hdLTRcXfS/QKaWZG4I/BbzeS6wX6IDJTUBmgGxJeeKkiYRLKE6wMzeSXJub6A3QKPGjTOQbedcYSupJdJYtwK3xy+5HH6+A7glQ/lJ9FUl6y3QA3jDzDbFbGtsZh2As4BHJO2a6EQzG2xmHcysQ53addLLsXOu0IkgWEV5FYeo160PVEiyrwKQqT6kc4FGMZ93AeYlObYHcY/1ZjYv/DkLGE/yJaSdcyWJsvvRPmogHQ/cL6lD7EZJ+wL3AxMylJ+vgBaSmkkqTxAsh8cfJKkVUBP4PGZbTUkVwve1gQOBH+PPdc6VPMHIpuwNpFG7P/UmCGhfSFoILCIohdYDvgv3p83MNkrqA4wk6P40xMymSLoLmGRmOUH1TGComcU+9rcGBknaTPAHYkBsa79zrmTL4irSaIHUzOYC+0g6lmAWqPrAAuArMxuRyQyF6Y2I23Z73Of+Cc77DPhb/HbnXOmQzY1NKXXITxTknHOu8Knkz0cKIKmCpEslPSNppKQW4fYzJLUuvCw657Z3paLVXlJLgtVC7yNY6O4IoFq4+yCCGfSdc67QZLKxSdLRkqZKmiEp6Qofkk6VZPEN7dvkLeI9/Bv4lSCIHkXuet8JBCOcnHOucChYaiTKK9+ktg5FPwZoA5wpqU2C46oBfYEv8kszaiA9CLjPzP5k2w7yCwnWb3LOuUKR4Uf7qEPR7wYeAP7KL8GogfQvgnWaEmkI/BkxHeecK5BMlUhJPBS9Ydy19gYamdl/oyQYNZCOBm6RtGPMNgs7wF+Bt+Q75wqZIr6A2pImxbzi+7nnORRdUhngYeDaqHmL2v3peuBTYAZBUDXgdmAPoDxwctQLOudcqgSUjd79aXE450Yy+Q1Frwa0BcaHJdz6wHBJ3cxsUqIEoy7H/BuwFzCQoMFpJkG96OtAezNbECUd55wrKCnaK4I8h6Kb2XIzq21mTc2sKTARSBpEIUKJVFI5gsrZ2WbWD+gXKavOOZcxQhkaJJrCUPTIojzabyKY8/NYks/E5JxzhSqTA5uiDEWP2d4lv/TyDaRmtlnSdIIJSpxzrsgF3Z9K/hDRW4HbJfmkIM65ohexfrS4huNHbbW/DdgJmCzpd4JO+Lk65ptZxwznzTnntiiuuUajiBpIpwA/FGZGHMydO5e7+9/OqFEfsHTJEuo3aMAJ3U7k1n53ULNmzXzPX716NcOHvcMH77/H5G/+x9zffqNMmTK0bNmK03qcyWWXX0H58uVznfP7778z7O23GPnBCH7++ScWzJ9P1apVabf3Plx08aWceJL3bMtLw7o16HfpcXQ9sA21dqzMgsUrePfD77hn0Aj+XLk23/MPat+CUU9fme9xLY6+jbkLt457KVNGnH5Uey48tTO7Nq5L9SoV+X3Rn3w+eRaPPD+Gn2aVro40wcTOxZ2L5JR7buTtT/v2HezTL5L2aigys2bO5NCDD2DRokUc3607rVrtzqSvvmTC+A9p2aoV4yZ8yk477ZRnGqNGfkD344+hVq1aHNzlUHbddTeWLV3KiPfeZcGCBezf6QDeHzWWihUrbjnntltu4sF/3k/TZs046KBDqFe/Pr/++gvD3n6LdevWccWVV/PAvx4q7NtPSc19+xR3FgBotkttPnz2GurtVJ13P/yWqXMW0mGPJnTp2Iqpsxdw2N8fZuny1Xmm0bhBLc7ttn/CfW1b7MyJh7djyox5dDjt3lz7Xhjwd049qj1zFyzjvY++Z9XqdezRYme6HtCaDRs30b3Pk0z4alrG7jUd66a+xuY1i9IKg63atrMn3xgb6djDW9f+Op9+pBmXZ4lUUiWC1vqmwHxgrJktLIJ8bXeuvOIyFi1axIMP/5vL+lyxZfsN113DY48+TP9+t/LYEwPzTKNevfoMee5FTjn1tFwlz5UrV3LU4V2Y+PlnDHzyca66euuAjQ77dmTU2PEcdPAhudL6+aefOKTz/jz26MP0OPNs9mnfPkN3Wno8evMZ1NupOtfc/zpPDt262s79155M33MOo3+fE+h7z9A80/h1/lLuGZR4YOBz9/UEYMhbn+ba3r5NY049qj1TZszjoHP/ydq/NmzZd263/Rl85zncdOFRWRNIMyWLn+yTNzZJak7wSP868E/gRWCqpK5FlLftxuxZsxgzehRNmjblkssuz7Wv3x13UqVKFV5+6QVWr867dLNXu3acedbZ2zy+V6tWjb5h8Px4wvhc+0486eRtgijA7q1bc+ppZwDw0Ufjt9m/vWvacCeOPKA1c35fzMBXP8q17+4n32PVmnWcddy+VK5YPkkKeau1YxW6HboXa9au5+X/fplrX7NdagMw/supuYIowH/HfwdA7ZrVKG0U8b/ikFer/QPAZoKZnyoTDAf9BhhUmBmSNETSIkkJ62QV+Hc4j+B3kvaJ2Xe+pOnh6/zCzGcmjf9wHABHHNGVMmVy/0qqVatGpwMOZM2aNXz5xcQCX6NcuXIAlN0h+qIIOefsUDalhRS2C106tgRgzOc/E189tmrNOj6fPIsqlSrQcc+mBUr/3G77UbFCOd4a8802da0/zpwPwCH7tqJihXK59h1zcFsAPvzi5wJdN1vl1JFGeRWHvAJpJ+A2M/vUzP4ys5+Ai4HGkgpz2rxngaPz2H8M0CJ89QaeBJBUC7gD2I9gJNYdkvJvockC06ZNBWC3li0T7t91txYATJ9W8Ee1558dAkDXrnl9tVutWLGCd95+E0kccaQ/hMRr2SToVj3j10UJ988Mt7doUrCVynueeAAAz7z5yTb7fpw5n3+/OI62LXbm27f78fBNp3N332688cjFDLz9bF77YBL9H480aVHJEXFS52xcRbQBMCtu20yCPw71CepMM87MPpLUNI9DugPPhyuITpRUIwzsXYDRZrYUQNJogoD8SmHkM5NWLF8OwI7Vd0y4f8cdg+3LlxdstsInH/8/Ro38gL32asf5f78g3+PNjEsvvpCFCxdy8SWXsXtrX0kmXvWqwaySy1clbplfviqYwnLHapVTTrtz+93YvXl9psyYx8RvZyc85sYH32LanIU8cO0pXHLGwVu2f/3jr7z47hes+Wt9ytfNdllcRZpvh/xsbNJPNpdgvnMM5pDUO2eKrT8W/1FoGc2ULY+OBfhr+87bb3H9tVdRv359XnntzS2P63m58fpreeuN1zmw80Hcn2Ut9iXFll9VAXrF9Dr5QACeefPTpMc8eMOpPHLT6dz31PvsdtRt1D7gGg7/+0NgxvDHL+fi0w9Oem5JlO3r2ucXSEeG9ZWLJC1iayl0bOz2cF9RSTaXYJ5zDObaaDbYzDqYWYc6tetkNHMFUT2nxLliecL9K1asAJKXWJMZPuwdzju7B3Xq1mXkmPE0a94833NuvvF6Hnv0YTofdDDvvDuCChUqpHTN7cWKsCS6Y9XE851XrxJ0MUtWYk2mZvXKnHh4u4SNTDnOOWE/LjuzC08MncC//jOa3xf9yeq16/ls8ixOuXIga9au5+6+3ahSqWANXdkqhflIi1xej/Z3FlkuUpNsLsG5BI/3sdvHF1mu0tCyZSsAZiSpA505YzoALZLUoSby5huv0/Pcs6hXvz4fjBrHbi1a5HvO9ddezf/9+xEO6XIobw37L5Urp/5Yur2Y9kvQC3C3xonrQHcNt0//JbUyxjknBI1MLwyfmDQIHxs2KCXq3rRwyUqmzlnI3q0b0bJpPb756bdtjimxsvjZPmkgNbNsDaTDgT6ShhI0LC03s/mSRgL3xjQwdaWErG56SJdDARgzZhSbN2/O1XK/cuVKPv/sUypVqkTH/RJ33I439JWXufDv57Fzw4aMHP1hviVRM+Pqvn0YNPAJDj/iSF5/axiVKiVbWcYBTPgq+ON2RKfdkZSr5b5q5Qp0atecNWvX8+V3c1JK9+8nB41M8X1HY5UvF/xvW7tm1YT764Tb12/YlNK1s102DxEtrmWgk5L0CvA50ErSXEm9JF0i6ZLwkBEEjWAzgKeAywDCRqa7CSZt/Qq4K6fhKds133VXjjiyK7/MmcPAJx7Pte/uO+9g9erVnH3OeVSpUmXL9qk//8zUn7ft4vLi88/Rq+e5NGrcmNHjPooURC+/pDeDBj7BUUcfwxtvD/cgGsHsuYsZ/dlPNG1YO1djD0C/S4+jauUKvPTel7kafVo2rUfLpsknUTtw711p3bwBP0xP3sgE8Nk3MwHoe85hVK9aMde+C0/tzC71azL/j+X8NKtQ2oOLTTY/2vsQ0SwdIrr77q356ssvmDD+Q1q0bMmHH32Wa4hopXLBP5m1G7b+/iaM/5BjjzqCzZs3c37PC9ilUaNtrrPjjjW44sqrtny+5+47+cdd/alUqRJ9+l61TWd+gD33ake37idm8nbTkq1DRH+evZB92wZDRKfNWcihPR/KNUR07Tf/B0ClvRPn/5m7z+Os4ztuM1IqXpVK5Rn37DXs2XIXFi5ZwXsTvmf5yrW0270Rh+7Xio0bN3HOjUMYNu7bzN5wAWViiGjrv+1tzw8fH+nYjs1rZNcQUVd0mu+6K59MnMTd/W9n9KgPGPn+COo3aMBlffpya787qFWrVr5p/PrLL2zevBmA58J+o/EaN2mSK5DOmR2UfNauXcs/778v4TnnnHt+VgXSbDF77mI6n/0A/S49niMPaM1RnfdgweIVPP7yh9wz6H2WrVgTOa0a1Spx0hF5NzLlWL12PYf1fIi+5x5G98PaccYxHSi/ww4s/nMVb476H488P5ZJU35J9/aySlDazN5Hey+RZkmJ1EWXLSVSF00mSqRt9tzbXhievJQeq0OzHb1E6pxziWRvedQDqXOuRBDK4lZ7D6TOuRIhi+OoB1LnXPYrzq5NUXggdc6VDFkcST2QOudKhGzu/uSB1DlXIngdqXPOpaMY16yPwgOpc65E8Ed755xLg/ASqXPOpS2L46gHUudcCZHFkTTr5iN1zrlEMrlmk6SjJU0Nl3W/KcH+SyR9L2mypE8ktckzbwW8J+ecK1KZmthZUlngcYKl3dsAZyYIlC+b2d/MrB3wAJDnKpAeSJ1zJUPmpsjvCMwws1lmth4YSrDM+xZmtiLmYxXyWVHZ60idc1kvxYmda0uKnWR4sJkNjvmcaOn2/ba5pnQ5cA1QHjgsrwtmXSCVNAQ4HlhkZm0T7D8buDH8uAq41My+DffNAVYCm4CNRT25q3OukKTWIX9xPv/vR1q63cweBx6XdBZwG3B+sgSz8dH+WeDoPPbPBg4xsz0JFrsbHLf/UDNr50HUudIlg4vfJVvSPZmhQJ5r7WRdIDWzj4Ckq3+a2Wdmtiz8OJHgS3DOlWrBxM5RXhF8BbSQ1ExSeaAHwTLvW68mtYj5eBwwPa8Es+7RPkW9gPdjPhswSpIBg+LqRZxzJVimRjaZ2UZJfYCRQFlgiJlNkXQXMMnMhgN9JB0BbACWkcdjPZTgQCrpUIJA2jlm84FmNk9SXWC0pJ/DEm78ub2B3gCNGjcukvw65wou0xM7m9kIYETctttj3l+ZSnpZ92gfhaQ9gaeB7ma2JGe7mc0Lfy4C3ibo5rANMxtsZh3MrEOd2nWKIsvOuXRlsJI000pcIJXUGHgLONfMpsVsryKpWs57oCvwQ/Hk0jmXaYr4X3HIukd7Sa8AXQj6gs0F7gDKAZjZQOB2YCfgibBiOaebUz3g7XDbDgQjEz4o8htwzhUKn/0pBWZ2Zj77LwQuTLB9FrBXYeXLOVeMBGU8kDrnXLqyN5J6IHXOZT2f2Nk55zIgi+OoB1LnXMngJVLnnEtTxOGfxcIDqXOuRMjeMOqB1DlXAsjXtXfOufT5uvbOOZeu7I2jHkidcyVDFsdRD6TOuZIg+lLLxcEDqXMu62X7yKYSN42ec85lGy+ROudKhGwukXogdc6VCN79yTnn0uEd8p1zLj3Z3tjkgdQ5VyL4o71zzqUpm0ukWdf9SdIQSYskJVwBVFIXScslTQ5ft8fsO1rSVEkzJN1UdLl2zhW2LF6NOStLpM8C/wc8n8cxH5vZ8bEbJJUFHgeOBOYCX0kabmY/FlZGnXNFyEuk0ZnZR8DSApzaEZhhZrPMbD0wFOie0cw554qFgDJSpFdxyMYSaRSdJH0LzAOuM7MpQEPgt5hj5gL7JTpZUm+gd/hxXaVyiasRSrjawOLizkQhKa33Vlrvq1W6Cfzvf1+PrFROtSMeXuTfYUkMpP8DmpjZKknHAu8ALUhc8LdECZjZYGAwgKRJZtahsDJbXErrfUHpvbfSfF/ppmFmR2ciL4Ul6x7t82NmK8xsVfh+BFBOUm2CEmijmEN3ISixOudcoSpxgVRSfYWrYEnqSHAPS4CvgBaSmkkqD/QAhhdfTp1z24use7SX9ArQBagtaS5wB1AOwMwGAqcCl0raCKwFepiZARsl9QFGAmWBIWHdaX4GZ/4uskJpvS8ovffm91VCKYhBzjnnCqrEPdo751y28UDqnHNp2u4CqaRakkZLmh7+rJnkuE0xw1CzttEqv2GxkipIejXc/4WkpkWfy9RFuK+ekv6I+R1dWBz5TFWEIdCS9O/wvr+TtE9R57Eg0hnaXRpsd4EUuAkYa2YtgLHh50TWmlm78NWt6LIXXcyw2GOANsCZktrEHdYLWGZmuwEPA/cXbS5TF/G+AF6N+R09XaSZLLhngbz6RB5D0C+6BcGgkSeLIE+Z8Cx53xcEQ7tzfl93FUGeisz2GEi7A8+F758DTizGvKQryrDY2Pt9Azg8p/tYFiu1w30jDIHuDjxvgYlADUkNiiZ3BZfG0O5SYXsMpPXMbD5A+LNukuMqSpokaaKkbA22iYbFNkx2jJltBJYDOxVJ7gouyn0BnBI+/r4hqVGC/SVR1HsviTpJ+lbS+5L2KO7MZFLW9SPNBEljgPoJdt2aQjKNzWyepObAOEnfm9nMzOQwY6IMi408dDaLRMnzu8ArZrZO0iUEpe7DCj1nha8k/r6iSDa0u1QolYHUzI5Itk/SQkkNzGx++Mi0KEka88KfsySNB/YGsi2QRhkWm3PMXEk7ADuS/Y9g+d6XmS2J+fgUJaDuN6JSOdTZzFbEvB8h6QlJtc2sVEzSsj0+2g8Hzg/fnw8Miz9AUk1JFcL3tYEDgWyc1zTKsNjY+z0VGGfZPwoj3/uKqzfsBvxUhPkrTMOB88LW+/2B5TlVUSVZHkO7S4VSWSLNxwDgNUm9gF+B0wAkdQAuMbMLgdbAIEmbCX7hA7JxgmgzSzgsVtJdwCQzGw48A7wgaQZBSbRH8eU4moj31VdSN2AjwX31LLYMpyDCEOgRwLHADGAN8PfiyWlq0hjaXSr4EFHnnEvT9vho75xzGeWB1Dnn0uSB1Dnn0uSB1Dnn0uSB1Dnn0uSB1CUlySK8uhRj/oZK+iTFc3YP85100EYKaR0dprVbumm5km177EfqousU874SMA74B/BezPas61/rXFHzQOqSCmcfAkBS1fDtzNjtyUiqaGZ/FVrmnMsi/mjv0ibpkvARdx9JH0taC1yR7NE3nFHrxbhth0r6RNJaSYslPSmpcor5aCTpOUmzw3SmSrpDUrkEh9cMqwZWSVog6eYE6bWT9IGkleGkxK9IqpNKntz2wQOpy6RXgTcJhjiOinqSpMPC4+cAJwPXASeR+uqTdYEFwFUEkww/DFwK/CvBsY8Ai4FTCGaOujccNpyTp9bAx+HHs4ELgfbA2ynmyW0H/NHeZdK/zGxQzocUJiS+HxhjZufEnLsIeFfSnWY2PUoiZvY18HV4voBPgfXAo5KuMbNNMYd/bWZ9wvcjJe0M3EYwNwHAXQSB/fhwHlck/Qh8L+lwMxsb8d7cdsBLpC6T3sv/kNwk1SAo6b0maYecFzAhPCTymkWSyki6XtLPBBNjbCAIjFWB+KAeX7J8C2gqKWei7yPCbcTkaSowH+gQ+QbddsEDqcukhQU4ZyeCyYyHEAS+nNcqgn+fqcx8fyNwL0EVwwkES5ZcHe6rGHds/Dy0OZ8bKFgzqgZwe1yeNgA7p5gntx3wR3uXSfFTieW02peP214r5v2y8OfNwJgEac5N4fqnAS+Z2R05G5R8Fc74JWZyPs83s02SVhAs6PZCgnMTTgbutl8eSF1hygmCrQn7m0raFWgOfAlgZkslfQO0MLMBaV6vErAubtvZSY49CfhPzOeTgV/MLCdIjgX2MLNJaebJbQc8kLpCY2YzJH0P3BdO6FseuIVtZ0a/HnhfUhmCesnVQFPgeOBqM/sl4iVHA70k/Q/4hWBlgF2SHNte0mMEaz8dDpxDsPxxjn7AREnDCUqmS8O0jgKeNLPPIubJbQe8jtQVtjMI6k5fBu4kWIBwduwBYQv4oQSB6iWC5TauBWaR2nIUtxEE4gFhOssJulIlchXBAolvEcyu38/MnorJ0xRgf2ATQYPV+wSzvq+Mz79zPkO+c86lyUukzjmXJg+kzjmXJg+kzjmXJg+kzjmXJg+kzjmXJg+kzjmXJg+kzjmXJg+kzjmXpv8HCA2Ly6ShoC8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sklearn.metrics import confusion_matrix\n",
    "import matplotlib.pyplot as plt\n",
    "import itertools\n",
    "\n",
    "predicted_labels = model.predict(X_test)\n",
    "cm = confusion_matrix(y_test, [np.round(i[0]) for i in predicted_labels])\n",
    "cm = cm.astype('float') / cm.sum(axis = 1)[:, np.newaxis]\n",
    "\n",
    "plt.imshow(cm, cmap = plt.cm.Blues)\n",
    "plt.title('Normalized confusion matrix', fontsize = 15)\n",
    "plt.colorbar()\n",
    "plt.xlabel('True label', fontsize = 15)\n",
    "plt.ylabel('Predicted label', fontsize = 15)\n",
    "for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):\n",
    "    plt.text(j, i, format(cm[i, j], '.2f'),\n",
    "             horizontalalignment = 'center', verticalalignment = 'center', fontsize = 20,\n",
    "             color='white' if cm[i, j] > 0.5 else 'black')\n",
    "plt.show()"
   ]
  },
  {
   "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}