{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true,
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "# Comparing a-Priori SBM Fitting with Different Groupings #"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "It is often desirable to create models of data so we can generate new samples, learn something about the data's structure,\n",
    "or for use in some larger analysis. SBM models define different connection probabilities for different communities in\n",
    "a graph.\n",
    "\n",
    "Here we enforce a community structure a-priori based a some connectome attributes or combinations thereof, and fit an SBM.\n",
    "In this case, the connection probability between two blocks is simply the total number of observed connections between them,\n",
    "divided by the total number of possible connections between them.\n",
    "\n",
    "We fit these SBMs using different attributes or combinations of attributes from the connectome, and then compute the\n",
    "BIC (Bayesian Information Criterion) for each model. This gives us an idea of which attributes best explain the topological\n",
    "structure of the connectome.\n",
    "\n",
    "For consistency, we use the attributes `hemisphere` (the side of the brain a neuron is on), `cell_type` (the high-level\n",
    "classification of a neuron as interneuron, sensory, or motor) and their combination as groups for fitting. The analysis\n",
    "is only preformed on the Ciona Intestinalis Connectome, two C. Elegans Worm Wiring Connectomes, and the 8 C. Elegans\n",
    "Developmental Connectomes, because only these have both those attributes defined.\n",
    "\n",
    "*Note that the code to do the actual SBM fitting is located in `analysis/apriori_sbm_code.py`*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import sys\n",
    "import matplotlib.pyplot as plt\n",
    "sys.path.append(\"../\")\n",
    "import pandas as pd\n",
    "from graph import GraphIO\n",
    "from analysis.apriori_sbm_code import fit_apriori_sbm\n",
    "import seaborn as sns\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## Ciona Intestinalis ##\n",
    "We'll start with the single Ciona Connectome. First we load the data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "path = '../json_connectomes/ciona.json'\n",
    "ciona_connectome, _, _, _ = GraphIO.load(path)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Now we'll fit an SBM for each grouping"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "                                                    bic    likelihood  \\\nCiona Intestinalis ['Side']               -23033.922165 -10863.237119   \n                   ['Annotation']         -22820.202349 -11323.658336   \n                   ['Side', 'Annotation'] -26206.331665 -10245.149497   \n\n                                           n_params       estimator  \nCiona Intestinalis ['Side']                     121  SBMEstimator()  \n                   ['Annotation']                16  SBMEstimator()  \n                   ['Side', 'Annotation']       529  SBMEstimator()  ",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th></th>\n      <th>bic</th>\n      <th>likelihood</th>\n      <th>n_params</th>\n      <th>estimator</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Ciona Intestinalis</th>\n      <th>['Side']</th>\n      <td>-23033.922165</td>\n      <td>-10863.237119</td>\n      <td>121</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['Annotation']</th>\n      <td>-22820.202349</td>\n      <td>-11323.658336</td>\n      <td>16</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['Side', 'Annotation']</th>\n      <td>-26206.331665</td>\n      <td>-10245.149497</td>\n      <td>529</td>\n      <td>SBMEstimator()</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# In Ciona, Hemisphere is called side and cell_type is called annotation\n",
    "groups = [['Side'], ['Annotation'], ['Side', 'Annotation']]\n",
    "rows = {'Ciona Intestinalis': {}}\n",
    "for group in groups:\n",
    "    rows['Ciona Intestinalis'][str(group)] = fit_apriori_sbm(ciona_connectome, group, plot_sbms=False)\n",
    "ciona_results = pd.DataFrame.from_dict({(i,j): rows[i][j]\n",
    "                           for i in rows.keys()\n",
    "                           for j in rows[i].keys()},\n",
    "                       orient='index')\n",
    "ciona_results.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Great, now lets plot -BIC for different groupings. The lowest -BIC is set to 0, since the scale doesn't really matter here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "<Figure size 432x576 with 2 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa8AAAIZCAYAAAARELD0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAA6ZElEQVR4nO3df7yUdZ3//8dTRMT8EcbRLwKKGSnamslkalaUbUpZ6Gd1w36Ard9oXe137Wrblu6u+7HfbWtWVCaUZVia5EqFGFmbiueYqIgmqygI4SExMQkRX58/rvfRi+OcM3MOzJl5n/O8325zm2ve16/XnLnOPOe6rvdco4jAzMwsJzs1uwAzM7O+cniZmVl2HF5mZpYdh5eZmWXH4WVmZtlxeJmZWXYcXrZDSTpf0vf6Md/XJf1LGp4iaXU/lrHNfJKWSZqyPXX1o4YJkkLSzjtoeZ+Q9K0dsawc7cjnn16Xl6ThZ7c3y5PDKwOSVkraJOkJSRsk/bek8aXxl0n699LjXdKb9X2S/pzmv1TShKY8gefq6jGUIuLvI+LfduT6IuKwiFi8I5fZCJLeIak9vb5rJS2QdBxARPxHRPz/TahppaQ31jntYknbXWO17aNRz78R25sNLIdXPt4aEbsDY4B1wH/1Mu2PgLcB7wD2Al4OdADHN7pI6xtJHwG+DPwHsC+wP3AJMK2JZZm1PIdXZiLiLxThdGi18enT8l8D0yLi1oh4OiL+FBFfjYhv9zDPuZL+V9JGSXdLOqU07gxJv5H0+bTX94CkqaXxB0r6VZp3ITC6P8+r+95jt3EfSHWNkzQi1fKQpHXp8M/IHubrvvewi6S5qdZlkiqlaSelPYjH0ri3lcbtlebrlPSgpE9K2imNG5bqWS/pfuAtfXjOewH/CpwdEVdFxJ8jYktE/DQiPp6m2eZwp6S3pfoeS/VO6vZ8PybpDkl/kvRDSbumcaMkXZuew4Y0PK7OOnvcBiRdCLwGuDjtOV6c2g+RtFDSo5LulfS3peW9Ob2eGyU9nGp+AbAA2C8t5wlJ+5Wfv547JDszvf7rJf1zablHSbop/W3WSrpY0i49PKdntzdJo9Pf47FU76+7Xl9rXX6BMiNpN+DtwM09TPJGYElErOrDYv+X4g1oL+AC4HuSxpTGvwq4lyKYPgt8W5LSuO9T7NWNBv4NmNmH9dak4rzEGcDrImI18BngpcARwEuAscCn6lzc24ArgBcC84GuN9rhwE+BXwD7AO8HLpd0cJrvvyj+Ni8GXgfMAN6Txr0XOAl4BVABTu1W/7mSru2hnmOAXYGr6yle0kuBHwAfAtqA64CfdnuD/lvgROBA4HCKvx0U/+vfAQ6g2Lvb1PX861R1G4iIfwZ+DZwTEbtHxDkpiBZSbBv7AKcDl0g6LC3r28D7ImIP4GXADRHxZ2AqsCYtZ/eIWNNDLccBB1McSfhUKcC3Ah9ONR6Txv9DHc/to8Bqir/pvsAnAF83r8U5vPLxE0mPAY9T7Fl9rofpXgSs7cuCI+LKiFgTEc9ExA+B+4CjSpM8GBHfjIitwByKQ5f7StofeCXwLxGxOSJupAiBHUGSvgicALw+IjpTYL4X+HBEPBoRGykOt02vc5m/iYjr0vP4LsXhVICjgd2BiyLiqYi4AbgWOF3SMIoPC+dFxMaIWAl8AXh3mvdvgS9HxKqIeBT4v+UVRsRFEXFSD/W8CFgfEU/XWf/bgf+OiIURsQX4PDASOLY0zVfSa/koxWtxRKrjjxHx44h4Mv3dLqQI4npV3QZ6mPYkYGVEfCft+d8G/Jjngn0LcKikPSNiQxrfFxdExKaIWAosJb2OEdERETenda4EvlHnc9ySns8Bac/31+GLvrY8h1c+To6IFwIjgHOAX0n6/6pM90eKf8S6SZoh6fZ02OQxik/D5cN/f+gaiIgn0+DuwH7AhvSpucuDfVl3L14IzAL+b0T8KbW1AbsBHaVaf5ba6/GH0vCTwK4qegXuB6yKiGdK4x+k2KsbDezCts+raxxd83YbV68/AqNVf8/E/crLT/WuKtUCz3+Ou0Oxxy7pG+mw5+PAjcALUzjXo6dtoJoDgFd1vUbpdXon0LW9/g3wZuBBFYecj6mzhufVwrbP8aXp8N8f0nP8D+o7jP05YAXwC0n3Szq3j/VYEzi8MhMRWyPiKopDJMdVmeR64Kg+nM84APgmRSC+KAXkXYB6my9ZC4xKh4m67F/PeuuwgeIT/HckvTq1rac43HVYRLww3fZKHVm2xxpgfLfzHPsDD6d1bqF4Q+4+Doq/wfhu4+p1E/AX4OQ+1PlsHWlPdHyplt58lOJQ26siYk/gtV2LqbfYXnTfS1kF/Kr0Gr0wHQY8CyCdi51GcUjxJ8C8HpbTV18D7gEmpuf4Cep4fmmP+qMR8WLgrcBHJLlzU4tzeGVGhWnAKGB59/ERcT3F+YarJU2WtLOkPST9vaS/q7LIF1C8aXSm5b+HYs+rpoh4EGgHLlDRPf84in/+Ws9h1263qm8wqZv7O9NzeVXa0/gm8CVJ+6RljZV0Qj319uIW4M/AP0oaruK7YW8FrkiHyeYBF6a/4wHAR4CuThTzgA+o6EwyCqj7U3vao/wU8FVJJ6e9o+GSpkr6bJVZ5gFvkXR8Ok/3UWAz8Ns6VrcHRfA/Jmlv4NP11lmHdRTnA7tcC7xU0rvT8xku6ZUqOsXsIumdkvZKhz4fp/gg1rWcF6noyNIfe6TlPSHpEOCsemaSdJKkl6TtsKuerTVmsyZzeOXjp5KeoPjnuhCYGRHLepj2VIqT+T8E/kSxJ1Wh2CvbRkTcTXEO5yaKN4+/Av6nD3W9g+Jk/qMUb4hza0w/luJNtHw7qKeJI2IhReeI+ZImA/9EcYjn5nRo6HqKPYp+i4inKDpzTKXY07oEmBER96RJ3k8RbvcDv6HoiHBpGvdN4OcU515uA64qL1vFl2wX9LLuL1KE4ScpPkCsotgL/kmVae8F3kXRgWQ9RcC+NdVfy5cpzo+tp+js87M65qnXfwKnquiJ+JV0Tu1NFOci11Ac5vsMxSFvKM4Xrkyv399TPCfS3/sHwP3pcON+fazjYxTb40aK1+WHdc43kWI7eoLi/+CSHL4fONTJ5yXNzCw33vMyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOzs3asGSdgVuBEak9fwoIj4t6XzgvUBnmvQTEXFdmuc84ExgK/CBiPh5ap8MXAaMBK4DPhgR0dv6R48eHRMmTNjBz8rMbHDr6OhYHxFtza6jloaFF7AZeENEPCFpOPAbSQvSuC9FxOfLE0s6FJgOHAbsB1wv6aURsRX4GjALuJkivE4EFtCLCRMm0N7evkOfkJnZYCfpwWbXUI+GHTaMwhPp4fB0621vaRpwRURsjogHgBXAUZLGAHtGxE1pb2sucHKj6jYzs9bX0HNekoZJuh14BFgYEbekUedIukPSpZJGpbaxwKrS7KtT29g03L292vpmSWqX1N7Z2VltEjMzGwQaGl4RsTUijgDGUexFvYziEOBBwBHAWuALaXJVW0Qv7dXWNzsiKhFRaWtr+UO2ZmbWTwPS2zAiHgMWAydGxLoUas8A3wSOSpOtBsaXZhsHrEnt46q0m5nZENXI3oZtwJaIeEzSSOCNwGckjYmItWmyU4C70vB84PuSvkjRYWMisCQitkraKOlo4BZgBvBfjarbrFEWd6xi7oLlrN+widGjRjJj6iSmTB5fe0Yze55G9jYcA8yRNIxiD29eRFwr6buSjqA49LcSeB9ARCyTNA+4G3gaODv1NAQ4i+e6yi+gRk9Ds1azuGMVF1+5lM1bik26c8MmLr5yKYADzKwfVOPrUtmqVCrhrvLWKv7u339B54ZNz2tvGzWSSz/5piZUZFadpI6IqDS7jlp8hQ2zAbC+SnD11m5mvXN4mQ2A0aNG9qndzHrn8DIbADOmTmLE8GHbtI0YPowZUyc1qSKzvDWyw4aZJV2dMtzb0GzHcHiZDZApk8c7rMx2EB82NDOz7Di8zMwsOw4vMzPLjs95lfjyPWZmeXB4Jb58j5lZPnzYMJm7YPmzwdVl85atzF2wvEkVmZlZTxxeiS/fY2aWD4dX4sv3mJnlw+GV+PI9Zmb5cIeNxJfvMTPLh8OrxJfvMTPLgw8bmplZdhoWXpJ2lbRE0lJJyyRdkNr3lrRQ0n3pflRpnvMkrZB0r6QTSu2TJd2Zxn1FkhpVt5mZtb5G7nltBt4QES8HjgBOlHQ0cC6wKCImAovSYyQdCkwHDgNOBC6R1NWD4mvALGBiup3YwLrNzKzFNSy8ovBEejg83QKYBsxJ7XOAk9PwNOCKiNgcEQ8AK4CjJI0B9oyImyIigLmleczMbAhq6DkvScMk3Q48AiyMiFuAfSNiLUC63ydNPhZYVZp9dWobm4a7t1db3yxJ7ZLaOzs7d+hzMTOz1tHQ8IqIrRFxBDCOYi/qZb1MXu08VvTSXm19syOiEhGVtra2PtdrZmZ5GJDehhHxGLCY4lzVunQokHT/SJpsNVDupz4OWJPax1VpNzOzIaqRvQ3bJL0wDY8E3gjcA8wHZqbJZgLXpOH5wHRJIyQdSNExY0k6tLhR0tGpl+GM0jxmZjYENfJLymOAOanH4E7AvIi4VtJNwDxJZwIPAacBRMQySfOAu4GngbMjousy72cBlwEjgQXpZmZmQ5SKDnyDT6VSifb29maXYWaWFUkdEVFpdh21+AobZmaWHYeXmZllx+FlZmbZcXiZmVl2HF5mZpYdh5eZmWXH4WVmZtlxeJmZWXYcXmZmlh2Hl5mZZcfhZWZm2XF4mZlZdhxeZmaWHYeXmZllx+FlZmbZcXiZmVl2HF5mZpYdh5eZmWWnYeElabykX0paLmmZpA+m9vMlPSzp9nR7c2me8yStkHSvpBNK7ZMl3ZnGfUWSGlW3mZm1vp0buOyngY9GxG2S9gA6JC1M474UEZ8vTyzpUGA6cBiwH3C9pJdGxFbga8As4GbgOuBEYEEDazczsxbWsD2viFgbEbel4Y3AcmBsL7NMA66IiM0R8QCwAjhK0hhgz4i4KSICmAuc3Ki6zcys9Q3IOS9JE4BXALekpnMk3SHpUkmjUttYYFVpttWpbWwa7t5ebT2zJLVLau/s7NyRT8HMzFpIw8NL0u7Aj4EPRcTjFIcADwKOANYCX+iatMrs0Uv78xsjZkdEJSIqbW1t21u6mZm1qIaGl6ThFMF1eURcBRAR6yJia0Q8A3wTOCpNvhoYX5p9HLAmtY+r0m5mZkNUI3sbCvg2sDwivlhqH1Oa7BTgrjQ8H5guaYSkA4GJwJKIWAtslHR0WuYM4JpG1W1mZq2vkb0NXw28G7hT0u2p7RPA6ZKOoDj0txJ4H0BELJM0D7iboqfi2amnIcBZwGXASIpehu5paGY2hKnowDf4VCqVaG9vb3YZZmZZkdQREZVm11GLr7BhZmbZcXiZmVl2HF5mZpYdh5eZmWXH4WVmZtlxeJmZWXYcXmZmlh2Hl5mZZcfhZWZm2XF4mZlZdhxeZmaWHYeXmZllx+FlZmbZcXiZmVl2HF5mZpYdh5eZmWXH4WVmZtlpWHhJGi/pl5KWS1om6YOpfW9JCyXdl+5HleY5T9IKSfdKOqHUPlnSnWncVySpUXWbmVnra+Se19PARyNiEnA0cLakQ4FzgUURMRFYlB6Txk0HDgNOBC6RNCwt62vALGBiup3YwLrNzKzFNSy8ImJtRNyWhjcCy4GxwDRgTppsDnByGp4GXBERmyPiAWAFcJSkMcCeEXFTRAQwtzSPmZkNQQNyzkvSBOAVwC3AvhGxFoqAA/ZJk40FVpVmW53axqbh7u3V1jNLUruk9s7Ozh36HMzMrHU0PLwk7Q78GPhQRDze26RV2qKX9uc3RsyOiEpEVNra2vperJmZZaGh4SVpOEVwXR4RV6XmdelQIOn+kdS+Ghhfmn0csCa1j6vSbmZmQ1QjexsK+DawPCK+WBo1H5iZhmcC15Tap0saIelAio4ZS9KhxY2Sjk7LnFGax8zMhqCdG7jsVwPvBu6UdHtq+wRwETBP0pnAQ8BpABGxTNI84G6KnopnR8TWNN9ZwGXASGBBupmZ2RClogPf4FOpVKK9vb3ZZZiZZUVSR0RUml1HLb7ChpmZZcfhZWZm2XF4mZlZdhxeZmaWnV7DS9JLJL26SvtrJB3UuLLMzMx6VmvP68vAxirtm9I4MzOzAVcrvCZExB3dGyOiHZjQkIrMzMxqqBVeu/YybuSOLMTMzKxetcLrVknv7d6Yro7R0ZiSzMzMelfr8lAfAq6W9E6eC6sKsAtwSgPrMjMz61Gv4RUR64BjJb0eeFlq/u+IuKHhlZmZmfWg1/CStHcaXJpu27RHxKONK83MzKy6WocNO3juByG7ruDb9eOQAby4QXWZmZn1qNZhwwMHqhAzM7N61bw8lKSd049AImm8pFMlHdHwyszMzHpQ6/JQ7wUeAR5Mw4uAU4EfSvqnAajPzMzseerpKn8QsAewHDggItZL2g24FfhMY8szMzN7vlrh9VREbAA2SFoREesBIuJJSU81vjwzM7Pnq3XOa6SkV0iaDOySho9Mj3u7dBSSLpX0iKS7Sm3nS3pY0u3p9ubSuPMkrZB0r6QTSu2TJd2Zxn2l6/ybmZkNXbX2vNYCX0zDf0jDY1L7H2rMexlwMTC3W/uXIuLz5QZJhwLTgcOA/YDrJb00IrYCXwNmATcD1wEnAgtqrNvMzAaxWl3lX9+9TdLvqrVXmfdGSRPqrGMacEVEbAYekLQCOErSSmDPiLgprXsucDIOLzOzIa0Zv6R8jqQ70mHFUaltLLCqNM3q1DY2DXdvr0rSLEntkto7Ozt3dN1mZtYi+hNe39yO9X2NovfiERSHHr+Q2qudx4pe2quKiNkRUYmISltb23aUaWZmrazP4RURl/R3ZRGxLiK2RsQzFCF4VBq1GhhfmnQcsCa1j6vSbmZmQ9iAHjaUNKb08BSgqyfifGC6pBGSDgQmAksiYi2wUdLRqZfhDOCagazZzMxaT63ehv0m6QfAFGC0pNXAp4Ep6dJSAawE3gcQEcskzQPuBp4Gzk49DQHOoui5OJKio4Y7a5iZDXGK6PEUUtYqlUq0t7c3uwwzs6xI6oiISrPrqKUZvQ3NzMy2i8PLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8tOw8JL0qWSHpF0V6ltb0kLJd2X7keVxp0naYWkeyWdUGqfLOnONO4rktSoms3MLA+N3PO6DDixW9u5wKKImAgsSo+RdCgwHTgszXOJpGFpnq8Bs4CJ6dZ9mWZmNsQ0LLwi4kbg0W7N04A5aXgOcHKp/YqI2BwRDwArgKMkjQH2jIibIiKAuaV5zMxsiBroc177RsRagHS/T2ofC6wqTbc6tY1Nw93bzcxsCGuVDhvVzmNFL+3VFyLNktQuqb2zs3OHFWdmZq1loMNrXToUSLp/JLWvBsaXphsHrEnt46q0VxURsyOiEhGVtra2HVq4mZm1joEOr/nAzDQ8E7im1D5d0ghJB1J0zFiSDi1ulHR06mU4ozSPmZkNUTs3asGSfgBMAUZLWg18GrgImCfpTOAh4DSAiFgmaR5wN/A0cHZEbE2LOoui5+JIYEG6mZnZEKaiE9/gU6lUor29vdllmJllRVJHRFSaXUctrdJhw8zMrG4OLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLTlPCS9JKSXdKul1Se2rbW9JCSfel+1Gl6c+TtELSvZJOaEbNZmbWOpq55/X6iDgiIirp8bnAooiYCCxKj5F0KDAdOAw4EbhE0rBmFGxmZq2hlQ4bTgPmpOE5wMml9isiYnNEPACsAI4a+PLMzKxVNCu8AviFpA5Js1LbvhGxFiDd75PaxwKrSvOuTm3PI2mWpHZJ7Z2dnQ0q3czMmm3nJq331RGxRtI+wEJJ9/Qyraq0RbUJI2I2MBugUqlUncbMzPLXlD2viFiT7h8BrqY4DLhO0hiAdP9Imnw1ML40+zhgzcBVa2ZmrWbAw0vSCyTt0TUMvAm4C5gPzEyTzQSuScPzgemSRkg6EJgILBnYqs3MrJU047DhvsDVkrrW//2I+JmkW4F5ks4EHgJOA4iIZZLmAXcDTwNnR8TWJtRtZmYtYsDDKyLuB15epf2PwPE9zHMhcGGDSzMzs0y0Uld5MzOzuji8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsZBNekk6UdK+kFZLObXY9ZmbWPDs3u4B6SBoGfBX4a2A1cKuk+RFxd3MrMzNrDYs7VjF3wXLWb9jE6FEjmTF1ElMmj292WQ2Ty57XUcCKiLg/Ip4CrgCmNbkmM7OWsLhjFRdfuZTODZsIoHPDJi6+cimLO1Y1u7SGySW8xgLlV2F1ajMzG/LmLljO5i1bt2nbvGUrcxcsb1JFjZdLeKlKWzxvImmWpHZJ7Z2dnQNQlplZ863fsKlP7YNBLuG1GigfvB0HrOk+UUTMjohKRFTa2toGrDgzs2YaPWpkn9oHg1zC61ZgoqQDJe0CTAfmN7kmM7OWMGPqJEYMH7ZN24jhw5gxdVKTKmq8LHobRsTTks4Bfg4MAy6NiGVNLsvMrCV09SocSr0NFfG8U0eDQqVSifb29maXYWaWFUkdEVFpdh215HLY0MzM7FkOLzMzy47Dy8zMsuPwMjOz7AzaDhuSOoEH+zn7aGD9DizHrMzblzXS9m5fB0REy39RdtCG1/aQ1J5DbxvLk7cva6Shsn35sKGZmWXH4WVmZtlxeFU3u9kF2KDm7csaaUhsXz7nZWZm2fGel5mZZcfhZWZm2Rl04SVpgqRNkm5Pj1em+50kfUXSXZLulHSrpAPTuOskvbDKss6X9LEa6ztf0hlp+HOS/lBrHmuenraP0vgPS/qLpL0atP6TJR3a1+kk/aukN/ZznWdIOj8Nf1jSQ5Iu7s+yBqvu20VqW5nuG/7e0cs0z752pbalkn5Q51PrM0mf6M90kn67HetcLGlCGv6lpCck9drdf9CFV/K/EXFEt7a3A/sBh0fEXwGnAI8BRMSbI+Kx7V1pRHwc+Pr2Lscartr20eV0it+PO6VB6z4ZqBle3aeLiE9FxPXbu/KI+BLwqe1dziDV03bR8PeOekmaRPG+/VpJL2jQauoKr+7TRcSxO2LlEfF6oOZPggzW8CrrTPdjgLUR8QxARKyOiA1QfMKSNDoN/7OkeyVdDxzctRBJB0n6maQOSb+WdEga9QQweH9re/Dr2j6QdBCwO/BJihDraj9D0lXp9b9P0mdL456QdGH6NHyzpH1T+wGSFkm6I93vL+lY4G3A5yTdnrap96ZP8ksl/VjSbj1Md5mkU9Oyj5f0u7QXcKmkEal9paQLJN2WxnVto5sotlPrm1Z47+j+2r0D+C7wC4ptpGsdiyV9RtISSb+X9JrU3tu2e3raTu6S9JnUdhEwMm13l6e2n6Tal0ma1ct0T6R7qTgK1bWn+vbUPiXV+SNJ90i6XJJSOY8CW2v8LbYVEYPqBkwA7qrSPg5YCdwOfAF4RWncSopLqkwG7gR2A/YEVgAfS9MsAiam4VcBN/Sw/vO75vGt9W49bR9p3CeBf6H4ULcS2Ce1nwHcD+wF7Epx2bHxaVwAb03DnwU+mYZ/CsxMw38H/CQNXwacWlrni0rD/w68v4fpLgNOTetfBbw0tc8FPpSGV5bm/wfgWz08zzOAi5v9WrTSrcZ2MSDvHXXW+XvgAOBNwPxS+2LgC2n4zcD1vW27FHuSDwFtFD9KfANwcprniW7r3DvdjwTu6tpmq0z3RLr/G2AhxQ8H75vWMwaYAvwp/T13Am4CjuvheS4GKr39LYbCnhdQfFqi+DR0HvAMsEjS8d0mew1wdUQ8GRGPA/MBJO0OHAtcqeKY+DcoXgwbXKYDV0TxCfsq4LTSuEUR8aeI+AtwN8UbCMBTwLVpuIPiTRDgGOD7afi7wHE9rPNl6dP4ncA7gcNq1Hgw8EBE/D49ngO8tjT+qiq12HZolfcOSa8EOiPiQYpAPFLSqNIkPb321bbdVwKLI6IzIp4GLmfb7ajsA5KWAjdTBN/EGqUeB/wgIrZGxDrgV2l9AEui2HN9huLDwITqi6ht5/7OmKOI2AwsABZIWkdxXmFR98mqzLoT8Fj0fJ7EMifpcIp/yoXpSMYuFJ9Yv5om2VyafCvP/e9sifRRsVt7dz19ofIyik+8S1WcvJ9Sq9Qa47vq7K0W66MWee84HThEz3Uy2pNiL+db6XFPr321bbfWdgQUh/qANwLHRMSTkhZT7MH1Olsv43r6P+qzIbPnJelISful4Z2Aw3n+VedvBE6RNFLSHsBbAdInqQcknZbml6SXD1z1NgBOB86PiAnpth8wVtIBtWbswW8p9uSg2KP6TRreCOxRmm4PYK2k4Wk6epiuyz3ABEkvSY/fTfHJ1hpkIN47JJ0j6ZxeatiJ4kjA4V3bKDCN0rnZProFeJ2k0ZKGpeV0bUdb0vYIxeHGDSm4DgGOLi2jPF3ZjcDbJQ2T1EaxR7ekn3X2aMiEF7AP8FNJdwF3AE8D23QXjojbgB9S7M7+GPh1afQ7gTPT7vMyig3HBo/pwNXd2q7muQDqqw8A75F0B0XAfDC1XwF8PHW4OIjiHNstFOcI7inN3306ANKhn/dQHIa6k+Iwlnu4NtZAvHccAvyxlxpeCzwcEQ+X2m4EDpXU58OQEbGW4jDoL4GlwG0RcU0aPRu4I3XE+Bmwc9qO/43i0CFVpiu7muLvtJTiXNo/RsQf+lpjLYPu8lAqvitwbUS8rEnrP5/ixOXnm7F+612zt49WkA5PViKix0/6Q02ztwtJ1wL/JyKeasb6W006PPmxiOixy/xg3PPaCuyl0pcNB4qkzwHvAv480Ou2ujVt+2gFkj5M8Yn78WbX0mKaul1ExEkOroKkXwIvBrb0Ot1g2/MyM7PBbzDueZmZ2SDn8DIzs+w4vMz6SdK+kr4v6f50+ZybJDXqmohmVuLwMuuHdE22nwA3RsSLI2IyRbf6cd2m8xeFzRrAHTbM+iFdHuhTEfG6KuPOAN5CcSWCF1Bck/BSih5UTwKzIuKO7l+rSN8jOikt5mcU3/96BcX17GZExJONfE5mOfGel1n/HAbc1sv4YyguzPsG4ALgdxFxOMXPSMytY/kHA7PTPI9TXGjXzBKHl9kOIOmrKn7W5NbUtDAiHk3Dx1FcnJeIuAF4kWr/2OWqiPifNPw9er6wr9mQ5PAy659lwJFdDyLibOB4ip+YgG2/qF7tQqVBcZmh8v/grt3Gd5/ezBKHl1n/3ADsKumsUttuPUx7I+miu+kq3evTBVtXkgJQ0pHAgaV59pd0TBo+necu7GtmuMOGWb+lC6J+ieIHBjsp9ra+TvGjfc9eO1DS3sB3KMKp3GFjJHANxYVfb6U4NDg1Lf46itA7FrgPeLc7bJg9x+Fl1mKafZFYsxz4sKGZmWXHe15mZpYd73mZmVl2HF5mZpYdh5eZmWXH4WVmZtlxeJmZWXYcXmZmlh2Hl5mZZWfQ/lDe6NGjY8KECc0uw8wsKx0dHesjoq32lM01aMNrwoQJtLe3N7sMM7OsSHqw2TXUw4cNzcwsOw4vMzPLjsPLzMyy4/AyM7PsDNoOG2ZmQ8nijlXMXbCc9Rs2MXrUSGZMncSUyeObXVbDOLzMzDK3uGMVF1+5lM1btgLQuWETF1+5FGDQBpgPG5qZZW7uguXPBleXzVu2MnfB8iZV1HgOLzOzzK3fsKlP7YOBw8vMLHOjR43sU/tg4PAyM8vcjKmTGDF82DZtI4YPY8bUSU2qqPHcYcPMLHNdnTLc29DMzLIyZfL4QR1W3fmwoZmZZcfhZWZm2XF4mZlZdhxeZmaWHYeXmZllp2HhJelSSY9IuqvUtrekhZLuS/ejSuPOk7RC0r2STii1T5Z0Zxr3FUlqVM1mZpaHRu55XQac2K3tXGBRREwEFqXHSDoUmA4clua5RFLXN+6+BswCJqZb92WamdkQ07DwiogbgUe7NU8D5qThOcDJpfYrImJzRDwArACOkjQG2DMiboqIAOaW5jEzsyFqoM957RsRawHS/T6pfSywqjTd6tQ2Ng13b69K0ixJ7ZLaOzs7d2jhZmbWOlqlw0a181jRS3tVETE7IioRUWlra9thxZmZWWsZ6PBalw4Fku4fSe2rgfJ1TcYBa1L7uCrtZmY2hA10eM0HZqbhmcA1pfbpkkZIOpCiY8aSdGhxo6SjUy/DGaV5zMxsiGrYhXkl/QCYAoyWtBr4NHARME/SmcBDwGkAEbFM0jzgbuBp4OyI6PpZ0LMoei6OBBakm5mZDWEqOvENPpVKJdrb25tdhplZViR1RESl2XXU0iodNszMzOrm8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8tOU8JL0oclLZN0l6QfSNpV0t6SFkq6L92PKk1/nqQVku6VdEIzajYzs9Yx4OElaSzwAaASES8DhgHTgXOBRRExEViUHiPp0DT+MOBE4BJJwwa6bjMzax3NOmy4MzBS0s7AbsAaYBowJ42fA5ychqcBV0TE5oh4AFgBHDWw5ZqZWSsZ8PCKiIeBzwMPAWuBP0XEL4B9I2JtmmYtsE+aZSywqrSI1anNzMyGqGYcNhxFsTd1ILAf8AJJ7+ptlipt0cOyZ0lql9Te2dm5/cWamVlLasZhwzcCD0REZ0RsAa4CjgXWSRoDkO4fSdOvBsaX5h9HcZjxeSJidkRUIqLS1tbWsCdgZmbN1Yzwegg4WtJukgQcDywH5gMz0zQzgWvS8HxguqQRkg4EJgJLBrhmMzNrITsP9Aoj4hZJPwJuA54GfgfMBnYH5kk6kyLgTkvTL5M0D7g7TX92RGwd6LrNzKx1KKLq6aPsVSqVaG9vb3YZZmZZkdQREZVm11GLr7BhZmbZcXiZmVl2HF5mZpYdh5eZmWXH4WVmZtlxeJmZWXYcXmZmlh2Hl5mZZcfhZWZm2XF4mZlZdhxeZmaWHYeXmZllx+FlZmbZcXiZmVl2HF5mZpYdh5eZmWWn119SlvSR3sZHxBd3bDlmZma19RpewB7p/mDglcD89PitwI2NKsrMzKw3vYZXRFwAIOkXwJERsTE9Ph+4suHVmZmZVVHvOa/9gadKj58CJuzwaszMzOpQ67Bhl+8CSyRdnR6fDMxpSEVmZmY11BVeEXGhpAXAa4AA3hMRv2toZWZmZj2od88LYCvwDEV4PdOYcszMzGqr65yXpA8ClwOjgX2A70l6fyMLMzMz60m9e15nAq+KiD8DSPoMcBPwX40qzMzMrCf19jYUxWHDLltTW79IeqGkH0m6R9JyScdI2lvSQkn3pftRpenPk7RC0r2STujves3MbHCoN7y+A9wi6XxJFwA3A9/ejvX+J/CziDgEeDmwHDgXWBQRE4FF6TGSDgWmA4cBJwKXSBq2Hes2M7PM1RVe6TJQ7wEeBf5I0dvwy/1ZoaQ9gdeSwi8inoqIx4BpPNf9fg5Fd3xS+xURsTkiHgBWAEf1Z91mZjY49OXCvFt5rqfh9vQ2fDHQCXxH0u8kfUvSC4B9I2ItQLrfJ00/FlhVmn91anseSbMktUtq7+zs3I4SzcyslTWjt+HOwJHA1yLiFcCfSYcIe1p9lbaoNmFEzI6ISkRU2tra+lmemZm1umb0NlwNrI6IW9LjH1GE1zpJYyJiraQxwCOl6ceX5h8HrOnHes3MbJAY8N6GEfEHYJWkg1PT8cDdFFesn5naZgLXpOH5wHRJIyQdCEwElvRn3WZmNjjUu+fV1duwfG3D7elt+H7gckm7APdTdAbZCZgn6UzgIeA0gIhYJmkeRcA9DZwdEVurL9bMzIYCRVQ9ffT8CaXJwKsp9rhubPVrG1YqlWhvb292GWZmWZHUERGVZtdRS1+ubXg7sLZrHkn7R8RDjSjKzMysN3WFV+pZ+GlgHc+d7wrg8MaVZmZmVl29e14fBA6OiD82shgzM7N61NvbcBXwp0YWYmZmVq9e97wkfSQN3g8slvTfwOau8emyUWZmZgOq1mHDPdL9Q+m2S7qZmZk1Ta/hFREXDFQhZmZm9ap12PDLEfEhST+lyvUEI+JtDavMzMysB7UOG3433X++0YWYmZnVq9Zhw450/6uBKcfMzKy2WocN76T6z48IiIjwl5TNzGzA1TpseNKAVGFmZtYHtQ4bPtg1LOkAYGJEXC9pZK15zczMGqXeX1J+L8WPRn4jNY0DftKgmszMzHpV797T2cBRwC0AEXGfpH0aVlWTLO5YxdwFy1m/YROjR41kxtRJTJk8vvaMZmY2oOoNr80R8ZRU/HiypJ2p3pEjW4s7VnHxlUvZvKX4ncvODZu4+MqlAA4wM7MWU++FeX8l6RPASEl/DVwJ/LRxZQ28uQuWPxtcXTZv2crcBcubVJGZmfWk3vA6F+gE7gTeB1wXEf/csKqaYP2GTX1qNzOz5qk3vM6PiG9GxGkRcSpwqaTLG1nYQBs9amSf2s3MrHnqDa/9JZ0HIGkX4CrgvoZV1QQzpk5ixPBh27SNGD6MGVMnNakiMzPrSb0dNt4DXJ4C7PXAgoj4UuPKGnhdnTLc29DMrPXVujzUkaWH/0nxPa//oejAcWRE3NbI4gbalMnjHVZmZhmotef1hW6PNwCHpvYA3tCIoszMzHpT6/JQrx+oQszMzOpV67DhuyLie5I+Um18RHyxvyuWNAxoBx6OiJMk7Q38EJgArAT+NiI2pGnPA84EtgIfiIif93e9ZmaWv1q9DV+Q7veoctt9O9f9QaD8DeBzgUURMRFYlB4j6VBgOnAYcCJwSQo+MzMbomodNvxGur+g+zhJH+rvSiWNA94CXAh07dVNA6ak4TnAYuCfUvsVEbEZeEDSCorrLN7U3/WbmVne6v2eVzVVDyXW6cvAPwLPlNr2jYi1AOm+68K/Y4FVpelWp7bnkTRLUruk9s7Ozu0oz8zMWtn2hJf6NZN0EvBIRHRsx3qqXhQ4ImZHRCUiKm1tbf0pz8zMMrA9PyjZ36vKvxp4m6Q3A7sCe0r6HrBO0piIWCtpDPBImn41UP7y1ThgTX+LNjOz/PW65yVpo6THq9w2Avv1Z4URcV5EjIuICRQdMW6IiHcB84GZabKZwDVpeD4wXdIISQcCE4El/Vm3mZkNDrU6bOwxUIUAFwHzJJ0JPASclmpYJmkecDfwNHB2RGzteTFmZjbYKWJQ/abksyqVSrS3tze7DDOzrEjqiIhKs+uoZXs6bJiZmTWFw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Di8zM8uOw8vMzLLj8DIzs+w4vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsOw4vMzPLjsPLzMyy4/AyM7PsOLzMzCw7Ax5eksZL+qWk5ZKWSfpgat9b0kJJ96X7UaV5zpO0QtK9kk4Y6JrNzKy1NGPP62ngoxExCTgaOFvSocC5wKKImAgsSo9J46YDhwEnApdIGtaEus3MrEUMeHhFxNqIuC0NbwSWA2OBacCcNNkc4OQ0PA24IiI2R8QDwArgqAEt2szMWkpTz3lJmgC8ArgF2Dci1kIRcMA+abKxwKrSbKtTW7XlzZLULqm9s7OzYXWbmVlzNS28JO0O/Bj4UEQ83tukVdqi2oQRMTsiKhFRaWtr2xFlmplZC2pKeEkaThFcl0fEVal5naQxafwY4JHUvhoYX5p9HLBmoGo1M7PW04zehgK+DSyPiC+WRs0HZqbhmcA1pfbpkkZIOhCYCCwZqHrNzKz17NyEdb4aeDdwp6TbU9sngIuAeZLOBB4CTgOIiGWS5gF3U/RUPDsitg541WZm1jIGPLwi4jdUP48FcHwP81wIXNiwoszMLCu+woaZmWXH4WVmZtlxeJmZWXYcXmZmlh2Hl5mZZcfhZWZm2XF4mZlZdhxeZmaWHYeXmZllx+FlZmbZcXiZmVl2HF5mZpYdh5eZmWXH4WVmZtlxeJmZWXYcXmZmlh2Hl5mZZcfhZWZm2XF4mZlZdhxeZmaWHYeXmZllx+FlZmbZcXiZmVl2HF5mZpYdh5eZmWVn52YXUC9JJwL/CQwDvhURFzW5JLM+WdyxirkLlrN+wyZGjxrJjKmTmDJ5fLPLMstSFuElaRjwVeCvgdXArZLmR8Tdza3MrD6LO1Zx8ZVL2bxlKwCdGzZx8ZVLARxgZv2Qy2HDo4AVEXF/RDwFXAFMa3JNZnWbu2D5s8HVZfOWrcxdsLxJFZnlLZfwGgusKj1endq2IWmWpHZJ7Z2dnQNWnFkt6zds6lO7mfUul/BSlbZ4XkPE7IioRESlra1tAMoyq8/oUSP71G5mvcslvFYD5RMD44A1TarFrM9mTJ3EiOHDtmkbMXwYM6ZOalJFZnnLosMGcCswUdKBwMPAdOAdzS3JrH5dnTLc29Bsx8givCLiaUnnAD+n6Cp/aUQsa3JZZn0yZfJ4h5XZDpJFeAFExHXAdc2uw8zMmi+Xc15mZmbPcniZmVl2HF5mZpYdRTzv61KDgqRO4MF+zj4aWL8DyzEr8/ZljbS929cBEdHyX5QdtOG1PSS1R0Sl2XXY4OTtyxppqGxfPmxoZmbZcXiZmVl2HF7VzW52ATaoefuyRhoS25fPeZmZWXa852VmZtlxeJmZWXYcXmZmlp1BF16SJkjaJOn29Hhlut9J0lck3SXpTkm3pp9YQdJ1kl5YZVnnS/pYjfWdL+mMNPw5SX+oNY81T0/bR2n8hyX9RdJeDVr/yZIO7et0kv5V0hv7uc4zJJ2fhj8s6SFJF/dnWYNV9+0ita1M9w1/7+hlmmdfu1LbUkk/qPOp9ZmkT/RnOkm/3Y51LpY0IQ3/UtITknr9rtqgC6/kfyPiiG5tbwf2Aw6PiL8CTgEeA4iIN0fEY9u70oj4OPD17V2ONVy17aPL6RS/H3dKg9Z9MlAzvLpPFxGfiojrt3flEfEl4FPbu5xBqqftouHvHfWSNIniffu1kl7QoNXUFV7dp4uIY3fEyiPi9UB7rekGa3iVdab7McDaiHgGICJWR8QGKD5hSRqdhv9Z0r2SrgcO7lqIpIMk/UxSh6RfSzokjXoC2DRgz8Z2tK7tA0kHAbsDn6QIsa72MyRdlV7/+yR9tjTuCUkXpk/DN0vaN7UfIGmRpDvS/f6SjgXeBnxO0u1pm3pv+iS/VNKPJe3Ww3SXSTo1Lft4Sb9LewGXShqR2ldKukDSbWlc1za6iWI7tb5phfeO7q/dO4DvAr+g2Ea61rFY0mckLZH0e0mvSe29bbunp+3kLkmfSW0XASPTdnd5avtJqn2ZpFm9TPdEupeKo1Bde6pvT+1TUp0/knSPpMslKZXzKLC1xt9iWxExqG7ABOCuKu3jgJXA7cAXgFeUxq2kuB7YZOBOYDdgT2AF8LE0zSJgYhp+FXBDD+s/v2se31rv1tP2kcZ9EvgXig91K4F9UvsZwP3AXsCuFNfMHJ/GBfDWNPxZ4JNp+KfAzDT8d8BP0vBlwKmldb6oNPzvwPt7mO4y4NS0/lXAS1P7XOBDaXhlaf5/AL7Vw/M8A7i42a9FK91qbBcD8t5RZ52/Bw4A3gTML7UvBr6Qht8MXN/btkuxJ/kQ0Ebxu443ACeneZ7ots690/1I4K6ubbbKdE+k+78BFlL8cPC+aT1jgCnAn9LfcyfgJuC4Hp7nYqDS299iKOx5AcWnJYpPQ+cBzwCLJB3fbbLXAFdHxJMR8TgwH0DS7sCxwJUqjol/g+LFsMFlOnBFFJ+wrwJOK41bFBF/ioi/AHdTvIEAPAVcm4Y7KN4EAY4Bvp+Gvwsc18M6X5Y+jd8JvBM4rEaNBwMPRMTv0+M5wGtL46+qUotth1Z575D0SqAzIh6kCMQjJY0qTdLTa19t230lsDgiOiPiaeBytt2Oyj4gaSlwM0XwTaxR6nHADyJia0SsA36V1gewJIo912coPgxMqL6I2rL5JeUdISI2AwuABZLWUZxXWNR9siqz7gQ8Fj2fJ7HMSTqc4p9yYTqSsQvFJ9avpkk2lybfynP/O1sifVTs1t5dT1cDuIziE+9SFSfvp9Qqtcb4rjp7q8X6qEXeO04HDtFznYz2pNjL+VZ63NNrX23brbUdAcWhPuCNwDER8aSkxRR7cL3O1su4nv6P+mzI7HlJOlLSfml4J+Bwnv+TKTcCp0gaKWkP4K0A6ZPUA5JOS/NL0ssHrnobAKcD50fEhHTbDxgr6YBaM/bgtxR7clDsUf0mDW8E9ihNtwewVtLwNB09TNflHmCCpJekx++m+GRrDTIQ7x2SzpF0Ti817ERxJODwrm0UmEbp3Gwf3QK8TtJoScPScrq2oy1pe4TicOOGFFyHAEeXllGeruxG4O2Shklqo9ijW9LPOns0ZMIL2Af4qaS7gDuAp4FtugtHxG3ADyl2Z38M/Lo0+p3AmWn3eRnFhmODx3Tg6m5tV/NcAPXVB4D3SLqDImA+mNqvAD6eOlwcRHGO7RaKcwT3lObvPh0A6dDPeygOQ91JcRjLPVwbayDeOw4B/thLDa8FHo6Ih0ttNwKHSurzYciIWEtxGPSXwFLgtoi4Jo2eDdyROmL8DNg5bcf/RnHokCrTlV1N8XdaSnEu7R8j4g99rbGWQXdtQxXfFbg2Il7WpPWfT3Hi8vPNWL/1rtnbRytIhycrEdHjJ/2hptnbhaRrgf8TEU81Y/2tJh2e/FhE9NhlfjDueW0F9lLpy4YDRdLngHcBfx7odVvdmrZ9tAJJH6b4xP14s2tpMU3dLiLiJAdXQdIvgRcDW3qdbrDteZmZ2eA3GPe8zMxskHN4mZlZdhxeZv0kaV9J35d0f7p8zk2SGnVNRDMrcXiZ9UO6JttPgBsj4sURMZmiW/24btP5i8JmDeAOG2b9kC4P9KmIeF2VcWcAb6G4EsELKK5JeClFD6ongVkRcUf3r1Wk7xGdlBbzM4rvf72C4np2MyLiyUY+J7OceM/LrH8OA27rZfwxFBfmfQNwAfC7iDic4mck5tax/IOB2WmexykutGtmicPLbAeQ9FUVP2tya2paGBGPpuHjKC7OS0TcALxItX/sclVE/E8a/h49X9jXbEhyeJn1zzLgyK4HEXE2cDzFT0zAtl9Ur3ah0qC4zFD5f3DXbuO7T29micPLrH9uAHaVdFapbbcepr2RdNHddJXu9emCrStJASjpSODA0jz7SzomDZ/Ocxf2NTPcYcOs39IFUb9E8QODnRR7W1+n+NG+Z68dKGlv4DsU4VTusDESuIbiwq+3UhwanJoWfx1F6B0L3Ae82x02zJ7j8DJrMc2+SKxZDnzY0MzMsuM9LzMzy473vMzMLDsOLzMzy47Dy8zMsuPwMjOz7Di8zMwsO/8PK5KjQKVqP1wAAAAASUVORK5CYII=\n"
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "colors = sns.color_palette(\"deep\")[0]\n",
    "fig, ax = plt.subplots(2, 1)\n",
    "out_index = ciona_results.index.levels[0]\n",
    "for i, key in enumerate(out_index):\n",
    "    data = ciona_results.loc[key]\n",
    "    bic_data = data['bic'].to_numpy()\n",
    "    norm_bic = bic_data - np.min(bic_data)\n",
    "    like_data = data['likelihood'].to_numpy()\n",
    "    norm_like = like_data - np.min(like_data)\n",
    "    _ = ax[0].scatter(data['bic'].index,\n",
    "                  norm_bic,\n",
    "                  label=key,\n",
    "                  color=np.array(colors))\n",
    "    _ = ax[1].scatter(data['likelihood'].index,\n",
    "                  norm_like,\n",
    "                  label=key,\n",
    "                  color=np.array(colors))\n",
    "    _ = ax[0].set(xlabel='Group', ylabel='-BIC')\n",
    "    _ = ax[1].set(xlabel='Group', ylabel='Likelihood')\n",
    "_ = fig.suptitle('BIC and Likelihood: Ciona Intestinalis')\n",
    "_ = fig.set_size_inches(6, 8)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It appears annotation has the best BIC but the worst likelihood. This could be because annotation has much fewer parameters than the other groupings. (For Ciona\n",
    "'Side' is broken into 11 groupings instead of just 'Left' and 'Right'). Therefore, based on\n",
    "the BIC, annotation explains the structure better, and the other groupings are likely\n",
    "fitting noise."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## C. Elegans Developmental Conectomes (Witviliet) ##\n",
    "\n",
    "We will preform this analysis for all 8 developmental stages at once."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "# load some elegans Witvliet data\n",
    "base = '../json_connectomes/witvilet/'\n",
    "graph_files = sorted(os.listdir(base))\n",
    "wit_connectomes = []\n",
    "for f in graph_files:\n",
    "    # These contain multiple synapse types, and thus can have multiedges and need to be flattened.\n",
    "    # Alternatively, we could preform the analysis separately for each synapse type\n",
    "    wit_connectomes.append(GraphIO.flatten_multigraph(GraphIO.load(os.path.join(base, f))[0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "Now we will fit sbms for each grouping for each connectome."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "                                                         bic   likelihood  \\\nWitvliet 0 hours  ['hemisphere']                -8016.577339 -3961.208692   \n                  ['cell_type0']                -7651.218552 -3694.831560   \n                  ['hemisphere', 'cell_type0']  -8917.613900 -3574.749594   \nWitvliet 5 hours  ['hemisphere']                -9913.659242 -4909.418897   \n                  ['cell_type0']                -9401.135504 -4510.924858   \n                  ['hemisphere', 'cell_type0'] -10818.117288 -4376.558445   \nWitvliet 8 hours  ['hemisphere']               -10035.203502 -4970.007348   \n                  ['cell_type0']                -9387.116384 -4561.351516   \n                  ['hemisphere', 'cell_type0'] -10628.877074 -4420.721409   \nWitvliet 16 hours ['hemisphere']               -11716.860426 -5810.567133   \n                  ['cell_type0']               -10907.811713 -5320.952857   \n                  ['hemisphere', 'cell_type0'] -12142.430570 -5172.453006   \nWitvliet 23 hours ['hemisphere']               -14771.539210 -7337.602882   \n                  ['cell_type0']               -13806.911934 -6710.789074   \n                  ['hemisphere', 'cell_type0'] -15164.255483 -6533.163547   \nWitvliet 27 hours ['hemisphere']               -14586.909618 -7245.077303   \n                  ['cell_type0']               -13437.319264 -6584.277672   \n                  ['hemisphere', 'cell_type0'] -14596.318736 -6389.737317   \nWitvliet 45 hours ['hemisphere']               -19049.862125 -9476.429417   \n                  ['cell_type0']               -17546.406610 -8579.196723   \n                  ['hemisphere', 'cell_type0'] -19108.022424 -8341.470073   \n\n                                                n_params       estimator  \nWitvliet 0 hours  ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      25  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       169  SBMEstimator()  \nWitvliet 5 hours  ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      36  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       196  SBMEstimator()  \nWitvliet 8 hours  ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      25  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       169  SBMEstimator()  \nWitvliet 16 hours ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      25  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       169  SBMEstimator()  \nWitvliet 23 hours ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      36  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       196  SBMEstimator()  \nWitvliet 27 hours ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      25  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       169  SBMEstimator()  \nWitvliet 45 hours ['hemisphere']                       9  SBMEstimator()  \n                  ['cell_type0']                      36  SBMEstimator()  \n                  ['hemisphere', 'cell_type0']       225  SBMEstimator()  ",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th></th>\n      <th>bic</th>\n      <th>likelihood</th>\n      <th>n_params</th>\n      <th>estimator</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 0 hours</th>\n      <th>['hemisphere']</th>\n      <td>-8016.577339</td>\n      <td>-3961.208692</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-7651.218552</td>\n      <td>-3694.831560</td>\n      <td>25</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-8917.613900</td>\n      <td>-3574.749594</td>\n      <td>169</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 5 hours</th>\n      <th>['hemisphere']</th>\n      <td>-9913.659242</td>\n      <td>-4909.418897</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-9401.135504</td>\n      <td>-4510.924858</td>\n      <td>36</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-10818.117288</td>\n      <td>-4376.558445</td>\n      <td>196</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 8 hours</th>\n      <th>['hemisphere']</th>\n      <td>-10035.203502</td>\n      <td>-4970.007348</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-9387.116384</td>\n      <td>-4561.351516</td>\n      <td>25</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-10628.877074</td>\n      <td>-4420.721409</td>\n      <td>169</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 16 hours</th>\n      <th>['hemisphere']</th>\n      <td>-11716.860426</td>\n      <td>-5810.567133</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-10907.811713</td>\n      <td>-5320.952857</td>\n      <td>25</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-12142.430570</td>\n      <td>-5172.453006</td>\n      <td>169</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 23 hours</th>\n      <th>['hemisphere']</th>\n      <td>-14771.539210</td>\n      <td>-7337.602882</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-13806.911934</td>\n      <td>-6710.789074</td>\n      <td>36</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-15164.255483</td>\n      <td>-6533.163547</td>\n      <td>196</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 27 hours</th>\n      <th>['hemisphere']</th>\n      <td>-14586.909618</td>\n      <td>-7245.077303</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-13437.319264</td>\n      <td>-6584.277672</td>\n      <td>25</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-14596.318736</td>\n      <td>-6389.737317</td>\n      <td>169</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">Witvliet 45 hours</th>\n      <th>['hemisphere']</th>\n      <td>-19049.862125</td>\n      <td>-9476.429417</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-17546.406610</td>\n      <td>-8579.196723</td>\n      <td>36</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-19108.022424</td>\n      <td>-8341.470073</td>\n      <td>225</td>\n      <td>SBMEstimator()</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "groups = [[\"hemisphere\"], [\"cell_type0\"], [\"hemisphere\", \"cell_type0\"]]\n",
    "rows = {}\n",
    "ages = []\n",
    "for i in range(len(wit_connectomes)):\n",
    "    age = wit_connectomes[i].graph['age']\n",
    "    ages.append(age)\n",
    "    g_row = {}\n",
    "    for group in groups:\n",
    "        g_row[str(group)] = fit_apriori_sbm(wit_connectomes[i], group, plot_sbms=False)\n",
    "    rows['Witvliet ' + str(age) + ' hours'] = g_row\n",
    "wit_results = pd.DataFrame.from_dict({(i,j): rows[i][j]\n",
    "                           for i in rows.keys()\n",
    "                           for j in rows[i].keys()},\n",
    "                       orient='index')\n",
    "ages = np.array(ages)\n",
    "wit_results.head(24)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's plot -BIC for different groupings. Again the lowest -BIC for each developmental stage is set to 0, allowing for easy comparison."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "<Figure size 432x576 with 2 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcAAAAIaCAYAAAC3aJ2nAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABpuklEQVR4nO3deZwU1bn/8c8zC5uCgCwigw5EVFYREGKixiVxi1EwLiheNJIQvRiXaIzGm8Tr/WGMWVyi0RjjgqhAXAIxaFTU4I4DIgqIEEAZJGyC7EsPz++POjP0DD0zPUt3z0x/36/XvLr61Kmq0z3V9dSpOnWOuTsiIiLZJifTBRAREckEBUAREclKCoAiIpKVFABFRCQrKQCKiEhWUgAUEZGslLUB0MxuNrMJtVjufjP7eZg+3syKa7GOcsuZ2TwzO74u5apFGQrNzM0sL9XbSqIsl5jZG5kuR0NiZsea2cK494eZ2ftmtsnMrjSzlmb2dzP70sz+msmyNjTp+g1lUm2PPVJevQZAM1tmZtvMbLOZrTezf5hZt7j5j5jZ/4t73yzsrIvMbEtY/iEzK6zPctVUVTuXu1/m7v9Xn9tz9z7u/lp9rjMVzOxCMysK/9+VZva8mR2T5LLx+0bp3z2pLnO6JXNiEfb5XSGYbTKzT8zsHjPrUprH3V9398PiFrseeM3dW7v73cA5QGdgf3c/N2UfKHH5qz34ht/6zrjP+JGZ/crM9ktXORuT+g7aZnaWmc0xs41mttbMppceV7PhBCFZqagBfsfd9wW6AKuAP1SR9yngTOBCYD/gCGAWcFIKyiV1YGY/Bu4EbiU68B4E/BE4qwar+Y677xv3d0X9l7TRmOTurYH2wHDgAGBWfBCs4GBgXoX3n7h7rKYbTmOt//bwGTsC3wO+CrxpZvukaftZycwOAcYD1xIdV7sT/VZ3Z7JcDZK719sfsAz4Ztz704l+pKXvHwH+X5j+JrAN6FaD9d8A/BvYBMwHhsfNuwR4A/gtsB5YCpwWN7878K+w7EvAPcCESrZzPFBcybz4z1AuH3BlKFcB0DyU5TOiE4H7gZaVLFf2vQE3A5OJduBNRAe9wXF5ewGvARvCvDPj5u0XllsDfAr8D5AT5uWG8qwFlgBjAQfykvje9wM2A+fW175RYd4lwBtx7w8P/6MvgIXAeXHz9gf+DmwE3gP+X4Vl7wKWh/mzgGPj5lX33f4UWBHmLQROqqS83wbeD9tYDtwcN++z8L1uDn9HJ1j+5or7Xvj/fAD8tuI+ArwClADbwzqfBHYCu8L70SHfpcACov3/n8DBcev38D9fBCwNaWcAc8K+9BbQv8L/6zpgLvAlMAloAexD9LvdHfcZD6zqdxKX1hpYCVwRl5awzES/l99WWH4K8OMwfSDwNNG+vhS4srLvl+gke174nK8BvSp8zhuJfrfrgYeBFvH/A6La9+pQ9mGE4xrR/vmzuHXlsOcYtY5oX2sf5hWG/8HFRPvIWuCmMO/UCv/PD0L698J3s4noN/vDJI9R5wBzKplX422F+deHz/858P3wWQ4J86o61nUAngvf/RfA64RjUkP4q9+VlT+QtwIeBcYn+lEAtwH/quH6zw07fg5wPrAF6BLmXRL+qT8gOphcHv5ZFua/Dfw+/LOOC//oeguAwM+B2UDH8P5OYCrRGX5rooP2rxKtn70D4HaiH1ku8CvgnTAvH1gM/AxoBpwYPsdhYf54ooNEa6If3CfsOTheBnwMdAtlepW4AEj0w32uih9NjCSCZTL7RoJ5lxCCGNEBdjnRDzIPGEh0sOgT5k8Mf62A3iFvfAC8iChI5hGdAf+HPQe0qr7bw8K6DgzvC4GvVLF/9CPaD/sT/eiHxS1X5YkFCQJgSL8FeLeSfeQ14PuVrYPowLyY6AQpj+jk5624+U50UtEeaBm+19XA0PBdXBz+R83j/l8ziX5v7YkOjpdV9/tI9DupkD6eqPZbZZmJfqPL2fP7bUcUeEt//7OAXxD9DnoQHbRPqfjdAIcSHSe+RfT7uT5ss1nc5/yIPb+LNyn/+46F7eQTHVvWAE8Q/cb6EO1PPUL+q4F32HMC/CfgyQr7xZ/D938EsIMQjCv+P0Pat4GvAAZ8A9gKDEziGNUjlOsO4ARg3+r2v2q2dSrR76gP0e/uMcoHwDup/Fj3K6KAmB/+ji39nzaEv/pdWbQzbSaK9jGiANQv0Y8i7AgT67i9OcBZYfoSYHHcvFbhn3QA0eW6GLBP3PwnKu4EcfOq2rniP8PxRDWG3xPVPvcL6Ub0o/tK3HJHs+fMu9z62TsAvhw3rzewLUwfG3bEnLj5T4Zlcol+UL3j5v2Q6L4RRLWIy+LmnUzyNcCRwH/qcd8o/ftB3P+uNACeD7xeYdk/Ab8Mn3EXIeCHeeVqgAm2ux44Ionv9hCigPBNIL+Gn+1O4I4wXVjd90rlAfAyYFEl+8hrVB0Anyec7IT3OUQHsYPDewdOjJt/H/B/Fba/EPhG3P/rorh5twP3V/f7SPQ7qZB+G/BSdWUm+g19BhwX5v0AeCVMDwU+q7DeG4GHK343RCemkytsYwVwfNznjP9dnA78O+5zbgNyw/vW4XscGpd/FntOfhYQd9WA6DbQLqLgXrpfFMTNnwmMqGqfqPAZ/wZclcz/gOhy82SigL09/D/2reW2HiIEtLjfiofX6o51txCdlB9Sk99Uuv5ScQ9wmLu3JToDugL4l5kdkCDfOqIdJGlmNirc2N1gZhuAvkRV7FL/KZ1w961hcl+is8b17r4lLu+nNdl2FdoCY4h2kC9DWkeiADwrrqwvhPRk/CdueivQIty3ORBY7u7x1/I/BboSfQ/NKP+5SudRumyFeclaB3Soh3tHw9y9bdzfnxPkORgYWvq9he9uJNGJTEeig0n854ifxsyuNbMFoXXkBqLLtwn3EeK+W3dfTHQGfzOw2swmmtmBiT6EmQ01s1fNbI2ZfUkUuDokyltDXYkuE9XGwcBdcd/ZF0QHp65xeZZXyH9the+5G9F+Uqrid7VvLcsWL/4zVlpmj46eE4ELQt4LgcfjljuwQtl/RnRvuqIDidvXw29nOZV/L59S/jtY5+4lYXpbeF0VN38be76Xg4Fn48q0gOjSdXy5kv5Ozew0M3vHzL4I6zudJPczd3/H3c9z945EJ87HATfVclsVjx3x09Ud635DVON+0cyWmNkNyZQ/XVL2GIS7l7j7M0Q7QKKWgi8DQ8ysIJn1mdnBRLXGK4havrUlunRhSSy+EmhX4eb7QclsNwnrie6lPGxmXw9pa4l+GH3iDvb7edQ4qC4+B7qZWfz/7SCiM9q1RGebByeYB9F30K3CvGS9TXQWOayG5a2N5USXxuMD5b7ufjnR2WyM6BJTqfhWxscS3cc7D2gX9pEvSW4fwd2fcPdjiL5DB35dSdYniC75dHP3/Ygu8ZRuw5P7mOWF/+l3iO6R1MZyovs28d9bS3d/Ky6PV8g/rkL+Vu7+ZBLbqu1n3Jeohl36Gasr85PAOeG3P5Tonl/pcksrLNfa3U9PsNnPiftNmJkR7TMr4vJU/F18XpvPF8p1WoVytXD3FdUuWeE7NbPmRJ/3t0DnsC9PI8l9udyK3d8DniGqMNRmWyup5DdHNcc6d9/k7te6ew+i/fvHZtZgGjmmLABa5Cyia/cLKs5395eJ7kk8a2aDzCzPzFqb2WVmdmmCVe5D9I9bE9b/Pfb8Q6vk7p8CRcD/WvToxTFE/4zqPkOLCn8Jdz6PHmEYGT7L0HCW+WfgDjPrFNbV1cxOSaa8VXiX6HLD9WaWb9Gzg98hupRcQnTJY1z4Hg8GfgyUNneeDFxpZgVm1o7onl9SQs32F8C9ZjbMzFqF7Z9mZrfX8TNV9BxwqJn9V9hGvpkdZWa9wmd8Brg5lOFwYFTcsq2JAuQaIM/MfgG0SWajFj1nd2I4GGwn+lGXVJK9NfCFu283syFEtZNSa4gaiPRIcrv5ZtaL6GB/ANHl9Nq4H7jRzPqE9e5nZlU9HvFn4LJQmzUz28fMvm1mrZPY1ipgf0vykQYza25mg4guq60namhSbZnd/X2i7/NB4J/uviHMmglsNLOfWvQ8ZK6Z9TWzoxJsfjLwbTM7yczyie4L7yBq9FNqbPhdtCeqSU5K5nMlcD/R7+/g8Hk6hmNgMlYBhXEnt82IrqKtAWJmdhrRbYtqmdkxZvaDuGPP4UQNgd6p5bYmA98zs15m1oroWACU1agrPdaZ2Rlmdkg4dm4k+k1V9rtKu1QEwL+b2WaiDzsOuNjd51WS9xyiM41JRGfqHwGDiWqH5bj7fOB3RLWRVUSNEN6sQbkuJDqL/ILoftL4avJ3JToIxv99pbLM7v4SUcONqeHH/lOiqv87ZrYxfKbDKls+Ge6+k2hHPo3ozOuPwCh3/zhk+RFRgFxCdE/yCaLr9xDtpP8kamk4myiQlDGzn5nZ81Vs+/dEAfV/iH4oy4lq438Ly480s8r+z6X+buWfA3w2wXY2Ef34RhCdif+HqCbWPGS5guiy5n+IbsY/SXRAI3y+54ka/3xKFMjKXSKtQnOi+1Nrw7o7ER0ME/lv4BYz20R0MJgcV/6tRPv9m+GS0FcrWcf54Xeygag2uQ4Y5O61qn24+7NE39PEsL99RLSfVJa/iOi+2j1EQWkx0b3YZLb1MdH3viR8xoSXiolO1DYR/ebGE90v+1rprYgky/wkUa3xibjtlxCd+A0gagG6lihI7hWQ3X0hUcOoP4R83yF6HGdnXLYngBeJfjdLiO4r18ZdRP/LF8PnfofomJOM0s4M1pnZ7PA7uJJo31pPdPyamuS6NhAdJz4M+9gLwLNE93FrvC13fx64m6jh3GKiYzDs+d1VdazrGd5vDsv90RvQM8+lLaxEGiUz+zVwgLtfnOmySONjZsuIGhftddItiYUrFh8RtRiu8XOoDUnWdoUmjZOZHW5m/cNluyHAaKKzWxFJETMbHm4ftSOqtf+9sQc/UACUxqc10eXbLUSXbH5H1MxaRFLnh0S3Pv5NdA/v8swWp37oEqiIiGQl1QBFRCQrKQCKiEhWUgAUEZGspAAoIiJZSQFQRESykgKgiIhkJQVAERHJSgqAIiKSlRQARUQkKykAiohIVlIAFBGRrKQAKCIiWUkBUEREspICoIiIZCUFQBERyUoKgCIikpUUAEVEJCspAIqISFZSABQRkaykACgiIllJAVBERLKSAqCIiGQlBUAREclKCoAiIpKVFABFRCQrKQCKiEhWUgAUEZGspAAoIiJZSQFQRESykgKgiIhkJQVAERHJSgqAIiKSlRQARUQkKykAiohIVlIAFBGRrKQAKCIiWUkBUEREspICoIiIZCUFQBERyUoKgCIikpUUAEVEJCspAIqISFZSABQRkayUl+kCpEqHDh28sLAw08UQEWlUZs2atdbdO2a6HOmQsgBoZg8BZwCr3b1vXPqPgCuAGPAPd78+pN8IjAZKgCvd/Z8hfRDwCNASmAZc5e5e3fYLCwspKiqq188kItLUmdmnmS5DuqTyEugjwKnxCWZ2AnAW0N/d+wC/Dem9gRFAn7DMH80sNyx2HzAG6Bn+yq1TRESkNlIWAN19BvBFheTLgdvcfUfIszqknwVMdPcd7r4UWAwMMbMuQBt3fzvU+sYDw1JVZhERyR7pbgRzKHCsmb1rZv8ys6NCeldgeVy+4pDWNUxXTBcREamTdDeCyQPaAV8FjgImm1kPwBLk9SrSEzKzMUSXSznooIPqXFgRadh27dpFcXEx27dvz3RRGp0WLVpQUFBAfn5+pouSMekOgMXAM+Fy5kwz2w10COnd4vIVAJ+H9IIE6Qm5+wPAAwCDBw+utqGMSDoVrZ7NtGXPs37HBto1b8vphacxuNPATBerUSsuLqZ169YUFhZiluh8WRJxd9atW0dxcTHdu3fPdHEyJt2XQP8GnAhgZocCzYC1wFRghJk1N7PuRI1dZrr7SmCTmX3Vor17FDAlzWUWqbOi1bOZvOgp1u/YAMD6HRuYvOgpilbPzmzBGrnt27ez//77K/jVkJmx//77Z33NOWUB0MyeBN4GDjOzYjMbDTwE9DCzj4CJwMUemQdMBuYDLwBj3b0krOpy4EGihjH/Bp5PVZlFUmXasufZtXtXubRdu3cxbZl257pS8KsdfW8pvATq7hdUMuuiSvKPA8YlSC8C+u69hEjjUVrzSzZdRFJPXaGJpEG75m1rlC6NxwsvvMBhhx3GIYccwm233ZYwzyWXXMJTTz2V5pJJdRQARdLg9MLTyM8p39ouPyef0wtPy1CJslNs9y62l2xje8lWtpdsI1bhsnRNlZSUMHbsWJ5//nnmz5/Pk08+yfz58+uptMltX2pPAVAkDQZ3Gsh5Pc8pq/G1a96W83qeo1agaRTbvYuY72LPk1ROzHfVKQjOnDmTQw45hB49etCsWTNGjBjBlCmJ2+nNmDGDr33ta/To0aOsNuju/OQnP6Fv377069ePSZMmAfDaa69xxhlnlC17xRVX8MgjjwBRN4+33HILxxxzDH/961+5++676d27N/3792fEiBG1/izZqMl2hi3S0AzuNFABL4NiHqs0PY/aPQu3YsUKunXb8wRXQUEB7777bsK8K1eu5I033uDjjz/mzDPP5JxzzuGZZ55hzpw5fPDBB6xdu5ajjjqK4447rtrttmjRgjfeeAOAAw88kKVLl9K8eXM2bNhQq8+RrVQDFJEsUdmjwbV/ZDhRv/yVta4cNmwYOTk59O7dm1WrVgHwxhtvcMEFF5Cbm0vnzp35xje+wXvvvVftds8///yy6f79+zNy5EgmTJhAXp7qNDWhACgiWaKyZv+1fxygoKCA5cv39OJYXFzMgQcemDBv8+bNy6ZLA2dlA9vk5eWxe/fusvcVn9fbZ599yqb/8Y9/MHbsWGbNmsWgQYOIxRLXdGVvCoAikhXyLHHtqLL0ZBx11FEsWrSIpUuXsnPnTiZOnMiZZ56Z9PLHHXcckyZNoqSkhDVr1jBjxgyGDBnCwQcfzPz589mxYwdffvkl06dPT7j87t27Wb58OSeccAK33347GzZsYPPmzbX+PNlG9WURyQp5Ofmwu/ReYNTVcJ7lRem1XWdeHvfccw+nnHIKJSUlXHrppfTp0yfp5YcPH87bb7/NEUccgZlx++23c8ABBwBw3nnn0b9/f3r27MmRRx6ZcPmSkhIuuugivvzyS9yda665hrZt29b682QbS2Js2UZp8ODBrgFxRZq2BQsW0KtXr0wXo9FK9P2Z2Sx3H5yhIqWVLoGKiEhWUgAUEZGspAAoIiJZSQFQRESyklqBiqTJq8vfZvyCZ1m7bR0dWu7PqF7DOaHb0ZkulkjWSuV4gA+Z2eow9l/FedeZmZtZh7i0G81ssZktNLNT4tIHmdmHYd7dpkGspBF6dfnb3PPBeNZsW4cDa7at454PxvPq8rczXTSRrJXKS6CPAKdWTDSzbsC3gM/i0noDI4A+YZk/mllumH0fMIZolPieidYp0tCNX/AsO0p2lkvbUbKT8QuezVCJpL4UFhbSr18/BgwYwODBiZ8e0HBIDVMqB8SdYWaFCWbdAVwPxHeZfhYw0d13AEvNbDEwxMyWAW3c/W0AMxsPDEOjwksjs3bbuhqlS2psjW1hc2wju72EHMtl37w2tMrbp/oFq/Hqq6/SoUOH6jPWs5KSEnJzc6vPKAmltRGMmZ0JrHD3DyrM6gosj3tfHNK6humK6SKNSoeW+9coXerf1tgWNu7awG6PxtDb7SVs3LWBrbEtadm+hkNqeNLWCMbMWgE3AScnmp0gzatIr2wbY4gul3LQQQfVopQiqTGq13Du+WB8ucugzXObMarX8AyWKrtsjm1k78OHszm2sU61QDPj5JNPxsz44Q9/yJgxYxLm03BIDU86a4BfAboDH4RLmwXAbDM7gKhm1y0ubwHweUgvSJCekLs/4O6D3X1wx44d67n4IrV3QrejueKIUXRsuT8GdGy5P1ccMUqtQNOotOaXbHqy3nzzTWbPns3zzz/Pvffey4wZMxLm03BIDU/avi13/xDoVPo+BMHB7r7WzKYCT5jZ74EDiRq7zHT3EjPbZGZfBd4FRgF/SFeZRerToM69OHi/9sR8F3mWT4cWnTNdpKySY7kJg12O1e0eWunwR506dWL48OHMnDkzYS0ulcMhzZgxg6lTp/J///d/zJs3T4EwSal8DOJJ4G3gMDMrNrPRleV193nAZGA+8AIw1r1sT70ceBBYDPwbNYCRRmjjzg2s2raCmO8CIOa7WLVtBRt3bshswbLIvnlt2PuuioX02tmyZQubNm0qm37xxRfp27dv0strOKTMSmUr0AuqmV9Y4f04YFyCfEVA8nuUSAO0dvsqvML9J8dZu30VbZq1zUyhskzpfb76bAW6atUqhg+P7uPGYjEuvPBCTj01+Se1NBxSZmk4JJE0+OTLvfqDKHPofjq/qy0Nh1Q3Gg5JRFIuzxIPulpZuoikngKgSBp0aNEZq3D/yTA1hBHJIDUVEkmD0vt8a7evKtcKVPf/RDJHAVAkTdo0a6uAJ9KA6BKoiIhkJQVAERHJSgqAImkS272L7SXb2F6yle0l24jt3pXpIkk9uOOOO+jTpw99+/blggsu2KvXFtBwSA2VAqBIGsR27wq9wJQ+d+vEfJeCYJpt3LmBJRsX8smXH7Fk48I698SzYsUK7r77boqKivjoo48oKSlh4sSJ9VPYJJSU1K0f02ynACiSBjGP1Shd6l+quqOLxWJs27aNWCzG1q1by/oGrUjDITU8agUqkhbOxFem8otHfk/xmpUUdOzCLZf8mBEnnpnpgmWNVHRH17VrV6677joOOuggWrZsycknn8zJJyca8U3DITVEqgGKpMHEV/7Of9/1Pyxf/TnuzvLVn/Pfd/0PE1/5e6aLljVKa37Jpidj/fr1TJkyhaVLl/L555+zZcsWJkyYkDCvhkNqeBQARdLgFw//lm07yjeO2LZjO794+LcZKlH2SUV3dC+//DLdu3enY8eO5Ofnc/bZZ/PWW28lzJvK4ZDGjh3LrFmzGDRoELGYLqsnSwFQJA2K1/6nRulS/1LRHd1BBx3EO++8w9atW3F3pk+fXqPOuTUcUmalrL5sZg8BZwCr3b1vSPsN8B1gJ9HYft9z9w1h3o3AaKAEuNLd/xnSBwGPAC2BacBV3lSHsJAmq2uHzhSv2TvYde2gvkDTJRXd0Q0dOpRzzjmHgQMHkpeXx5FHHsmYMWOSXl7DIWVWyoZDMrPjgM3A+LgAeDLwirvHzOzXAO7+UzPrDTwJDCEaEf5l4NAwIvxM4CrgHaIAeLe7VzsoroZDkobkvhce4tp7bi13GbRl8xb87oqfcfmpl2awZI2bhkOqGw2HlCLuPgP4okLai+5l7b7fAQrC9FnARHff4e5LiUZ/H2JmXYA27v52qPWNB4alqswiqXLxN8/nt2NvpKDjAZgZBR0P4Ldjb+Tib55f/cIikhKZbDJ0KTApTHclCoilikParjBdMV2kUWmVtw+XfGsE55xwer2NRi4idZORAGhmNwEx4PHSpATZvIr0ytY7BhgD0c1pkYakVd4+CngiDUjaA6CZXUzUOOakuMYsxUC3uGwFwOchvSBBekLu/gDwAET3AOux2CJ1tmzjEuZ+MZutsS20ytuH/u0HUtimR6aLJZK10voYhJmdCvwUONPdt8bNmgqMMLPmZtYd6AnMdPeVwCYz+6qZGTAKmJLOMovUh2Ubl/DemrfYGtsCwNbYFt5b8xbLNi7JcMlEslfKAqCZPQm8DRxmZsVmNhq4B2gNvGRmc8zsfgB3nwdMBuYDLwBj3b20l9fLgQeJGsb8G6i2BahIQzP3i9mUePmOi0u8hLlfzM5QiUQkla1AL3D3Lu6e7+4F7v4Xdz/E3bu5+4Dwd1lc/nHu/hV3Pyz+MQd3L3L3vmHeFXoGUBqj0ppfsunSeFx66aV06tSJvn377jXvD3/4A4cddhh9+vTh+uuv32t+xU6vJb3UcZxIGrTK2ydhsFOjmPRKxX3YSy65hCuuuIJRo0aVS3/11VeZMmUKc+fOpXnz5qxevbpO26mJWCymfkGToK7QRNKgf/uB5FpuubRcy6V/+4EZKlH2SdV92OOOO4727dvvlX7fffdxww03lPUB2qlTp4TLb968mXPOOYfDDz+ckSNHlvUPOn36dI488kj69evHpZdeyo4dO4BoOKS1a9cCUFRUxPHHHw/AzTffzJgxYzj55JMZNWoU8+bNY8iQIQwYMID+/fuzaNGiOn3OpkgBUCQNCtv04KiOXyur8bXK24ejOn5NrUDTKN33YT/55BNef/11hg4dWuUoD++//z533nkn8+fPZ8mSJbz55pts376dSy65hEmTJvHhhx8Si8W47777qt3mrFmzmDJlCk888QT3338/V111FXPmzKGoqIiCgoJql882qiOLpElhmx4KeBmU7vuwsViM9evX88477/Dee+9x3nnnsWTJEqIG7XsMGTKkLDgNGDCAZcuW0bp1a7p3786hhx4KwMUXX8y9997L1VdfXeU2zzzzTFq2bAnA0Ucfzbhx4yguLubss8+mZ8+e9f8hGznVAEUkK1R2vzVV92ELCgo4++yzMTOGDBlCTk5O2aXLePHDJOXm5hKLxSodJgnKD5VU1TBJF154IVOnTqVly5accsopvPLKK3X9SE2OAqCIZIV034cdNmxYWdD55JNP2LlzJx06dEhq2cMPP5xly5axePFiAB577DG+8Y1vANE9wFmzZgHw9NNPV7qOJUuW0KNHD6688krOPPNM5s6dW5eP0yQpAIpIVkjVfdgLLriAo48+moULF1JQUMBf/vIXIHo8YsmSJfTt25cRI0bw6KOP7nX5szItWrTg4Ycf5txzz6Vfv37k5ORw2WXRU2O//OUvueqqqzj22GPJzc2tdB2TJk2ib9++DBgwgI8//nivVqqSwuGQMk3DIYk0fRoOqW40HJKIiEgWUgAUEZGspAAoIiJZSQFQRESykgKgiIhkJQVAERHJSqkcD/AhM1ttZh/FpbU3s5fMbFF4bRc370YzW2xmC83slLj0QWb2YZh3tyX7II2ISIotX76cE044gV69etGnTx/uuuuusnk///nP6d+/PwMGDODkk0/m888/32t5DYeUWamsAT4CnFoh7QZgurv3BKaH95hZb2AE0Ccs80ezsi4b7gPGEI0S3zPBOkVEklK0eja3zBzHNa//hFtmjqNodd06ws7Ly+N3v/sdCxYs4J133uHee+9l/vz5APzkJz9h7ty5zJkzhzPOOINbbrmlPj5CUmKxWNq21ZilckDcGcAXFZLPAh4N048Cw+LSJ7r7DndfSjT6+xAz6wK0cfe3w0C44+OWERFJWtHq2Uxe9BTrd2wAYP2ODUxe9FSdgmCXLl0YODDqSq1169b06tWLFStWANCmTZuyfFu2bKm0FxgNh5Q56b4H2NndVwKE19IBsroCy+PyFYe0rmG6YrqISI1MW/Y8u3bvKpe2a/cupi17vl7Wv2zZMt5//32GDh1alnbTTTfRrVs3Hn/88UprgBoOKXMaSiOYRKdGXkV64pWYjTGzIjMrWrNmTb0VTkQav9KaX7LpNbF582a++93vcuedd5ar+Y0bN47ly5czcuRI7rnnnoTLlg6HlJOTUzYc0sKFC/caDmnGjBnVlqPicEi33norv/71r/n000/L0mWPdAfAVeGyJuF1dUgvBrrF5SsAPg/pBQnSE3L3B9x9sLsP7tixY70WXEQat3bN29YoPVm7du3iu9/9LiNHjuTss89OmOfCCy+sdOQGDYeUOekOgFOBi8P0xcCUuPQRZtbczLoTNXaZGS6TbjKzr4bWn6PilhERSdrphaeRn5NfLi0/J5/TC0+r9TrdndGjR9OrVy9+/OMfl5sXf89t6tSpHH744UmvV8MhpUcqH4N4EngbOMzMis1sNHAb8C0zWwR8K7zH3ecBk4H5wAvAWHcvCau6HHiQqGHMv4H6uWAvIlllcKeBnNfznLIaX7vmbTmv5zkM7lT78QDffPNNHnvsMV555RUGDBjAgAEDmDZtGgA33HADffv2pX///rz44ovlHpGojoZDSg8NhyQijZaGQ6obDYckIiKShRQARUQkKykAiohIVsrLdAEakumfvcFjH/+NddvWs3/LdvzX4cM46aBjMl0sERFJAQXAYPpnb3Dv3AnsLIl6ili7bT33zp0AoCAoItIE6RJo8NjHfysLfqV2luzisY//lpkCiYhISikABuu2ra9RuohIVcMhnX/++WXPBhYWFjJgwIC9ltdwSJmlS6DB/i3bsTZBsNu/ZbsEuUWkMXp1+duMX/Asa7eto0PL/RnVazgndDu61usrHQ5p4MCBbNq0iUGDBvGtb32L3r17M2nSpLJ81157Lfvtt199fISkxGIx8vJ0eK+OaoDBfx0+jGa55btJapabz38dPiwzBRKRevXq8re554PxrNm2DgfWbFvHPR+M59Xlb9d6nVUNh1TK3Zk8eTIXXHBBwnVoOKTMUQAMTjroGMb2v4gOLdthQIeW7Rjb/yI1gBFpIsYveJYdJTvLpe0o2cn4Bc/Wy/oTDYcE8Prrr9O5c2d69uyZcDkNh5Q5qiPHOemgYxTwRJqotdvW1Si9JiobDgngySefrLT2B3uGQwLKhkNq3br1XsMh3XvvvVx99dVVlqPicEjjxo2juLiYs88+u9IAnM1UAxSRrNCh5f41Sk9WVcMhxWIxnnnmGc4///xKl9dwSJmjACgiWWFUr+E0z21WLq15bjNG9Rpe63VWNRwSwMsvv8zhhx9e48uPGg4pPaoMgGZ2iJl9PUH6sWb2ldQVS0Skfp3Q7WiuOGIUHVvujwEdW+7PFUeMqlMr0KqGQwKYOHFilZc/K6PhkNKjyuGQzOw54GfuPrdC+mDgl+7+nVpt1Owa4PuAAx8C3wNaAZOAQmAZcJ67rw/5bwRGAyXAle7+z+q2oeGQRJo+DYdUNxoOqWqFFYMfgLsXEQWqGjOzrsCVwGB37wvkAiOAG4Dp7t4TmB7eY2a9w/w+wKnAH82s8tMeERGRJFQXAFtUMa9lHbabB7Q0szyimt/nwFnAo2H+o8CwMH0WMNHdd7j7UqKR4YfUYdsiIiLVBsD3zOwHFRPNbDQwqzYbdPcVwG+Bz4CVwJfu/iLQ2d1XhjwrgU5hka7A8rhVFIc0ERGRWqvuOcCrgWfNbCR7At5goBlQq6ZTZtaOqFbXHdgA/NXMLqpqkQRpCW9cmtkYYAzAQQcdVJviiYhIlqgyALr7KuBrZnYC0Dck/8Pd6/JAyTeBpe6+BsDMngG+Bqwysy7uvtLMugCrQ/5ioFvc8gVEl0wTlfcB4AGIGsHUoYwiItLEVRkAzax9mPwg/JVLd/cvarHNz4CvmlkrYBtwElAEbAEuBm4Lr1NC/qnAE2b2e+BAoCcwsxbbFRERKVPdPcBZRMGp9LV0uvR9jbn7u8BTwGyiRyByiGpttwHfMrNFwLfCe9x9HjAZmA+8AIx195LabFtEJBVKSko48sgjyw1tdPPNN9O1a9eEzweW0nBImVXdJdDuqdiou/8S+GWF5B1EtcFE+ccB41JRFhHJHo9Pf5abHrqNz9Z8zkEdD2TcpTcw8qTa9wRT6q677qJXr15s3LixXPo111zDddddV+f115SGQ0pOtV2hmVmemVmY7mZm55jZgJSXTESkHj0+/VnG3HE9n65egbvz6eoVjLnjeh6fXrfRIIqLi/nHP/7B97///Votr+GQMqe6rtB+QNQY5dMwPR04B5hkZj9NQ/lEROrFTQ/dxtYd28qlbd2xjZseuq1O67366qu5/fbbycnZ+3B6zz330L9/fy699FLWr997wG3QcEiZVF0N8GrgK8AxwJ3A19x9BHAkoI7lRKTR+GxNwsbjlaYn47nnnqNTp04MGjRor3mXX345//73v5kzZw5dunTh2muvTbiO0uGQcnJyyoZDWrhw4V7DIc2YMaPa8lQcDunWW2/l17/+NZ9++mlZuuxRXQDc6e7r3f0zYLG7rwVw963AzqoXFRFpOA7qeGCN0pPx5ptvMnXqVAoLCxkxYgSvvPIKF10UPdbcuXNncnNzycnJ4Qc/+AEzZyZuvK7hkDKnugDY0syONLNBQLMwPTC8r6qbNBGRBmXcpTfQqnn5WlCr5i0Zd+kNtV7nr371K4qLi1m2bBkTJ07kxBNPZMKECQCsXLmyLN+zzz5L3759K1vNXjQcUnpU10xoJfD7MP2fMN0lpP8nheUSEalXpa09U9EKNJHrr7+eOXPmYGYUFhbypz/9Kell44dDisViHHXUUeWGQxo9ejS33norQ4cOrXQdkyZNYsKECeTn53PAAQfwi1/8os6fqampcjikhAuYve/uR6aoPPVGwyGJNH0aDqluNBySiIhIFqpNAPxzvZdCREQkzWocAN39j6koiIhIbdT0No5E9L3pEqiINGItWrRg3bp1OpjXkLuzbt06WrTI7sb86ixORBqtgoICiouLWbNmTaaL0ui0aNEi63uHUQAUkUYrPz+f7t1T0me/ZAFdAhURkayUkQBoZm3N7Ckz+9jMFpjZ0WbW3sxeMrNF4bVdXP4bzWyxmS00s1MyUWYREWlaMlUDvAt4wd0PB44AFgA3ANPdvSfRqBM3AJhZb2AE0Ac4FfijmeWmolCPT3+WwpFDyTm5G4Ujh9Z5mBQREWm40h4AzawNcBzwFwB33+nuG4CzgEdDtkeBYWH6LGCiu+9w96XAYmBIfZcrVWOFiYhIw5SJGmAPYA3wsJm9b2YPmtk+QGd3XwkQXjuF/F2B5XHLF4e0epWqscJERKRhykQAzAMGAveFPkW3EC53VsISpCV86MfMxphZkZkV1bRZdCrGChMRkYYrEwGwGCh293fD+6eIAuIqM+sCEF5Xx+XvFrd8AZAwKrn7A+4+2N0Hd+zYsUaFSsVYYSIi0nClPQC6+3+A5WZ2WEg6CZgPTAUuDmkXA1PC9FRghJk1N7PuQE8g8ciSdZCKscJERKThytSD8D8CHjezZsAS4HtEwXiymY0GPgPOBXD3eWY2mShIxoCx7l5S3wVK91hhIiKSWTUeD7Cx0HiAIiI1p/EARUREmjgFQBERyUoKgCIikpUUAEVEJCspAIqISFZSABQRkaykACgiIllJAVBERLKSAqCIiGQlBUAREclKCoAiIpKVFABFRCQrZWo0iAZpa2wLm2Mb2e0l5Fgu++a1oVXePpkuloiIpIACYLA1toWNuzZQOtj8bi8J71EQFBFpgjJ2CdTMcs3sfTN7Lrxvb2Yvmdmi8NouLu+NZrbYzBaa2SmpKM/m2EZKg98eHtJFRKSpyeQ9wKuABXHvbwCmu3tPYHp4j5n1BkYAfYBTgT+aWW59F2Z3JWPsVpYuIiKNW0YCoJkVAN8GHoxLPgt4NEw/CgyLS5/o7jvcfSmwGBhS32XKsVx2luxg084v2bhzA5t2fsnOkh3k1H+sFRGRBiBTNcA7geuB3XFpnd19JUB47RTSuwLL4/IVh7R6ZW5sL9mGh8ugjrO9ZBvmVt+bEhGRBiDtAdDMzgBWu/usZBdJkFbxZl3puseYWZGZFa1Zs6ZG5dq068sapYuISOOWiRrg14EzzWwZMBE40cwmAKvMrAtAeF0d8hcD3eKWLwA+T7Rid3/A3Qe7++COHTvWqFAx31WjdBERadzSHgDd/UZ3L3D3QqLGLa+4+0XAVODikO1iYEqYngqMMLPmZtYd6AnMrO9y5Vl+jdJFRKRxa0jPAd4GTDaz0cBnwLkA7j7PzCYD84EYMNa9/ptmdmjRmVXbVpTdAwQwjA4tOtf3pkREpAEw94S30xq9wYMHe1FRUY2W2bhzA2u3ryLmu8izfDq06EybZm1TU0ARkQbIzGa5++BMlyMdGlINMOPaNGurgCcikiXUGbaIiGQlBUAREclKCoAiIpKVFABFRCQrKQCKiEhWUgAUEZGspAAoIiJZSQFQRESykgKgiIhkJQVAERHJSgqAIiKSlRQARUQkKykAiohIVkp7ADSzbmb2qpktMLN5ZnZVSG9vZi+Z2aLw2i5umRvNbLGZLTSzU9JdZhERaXoyUQOMAde6ey/gq8BYM+sN3ABMd/eewPTwnjBvBNAHOBX4o5nlZqDcIiLShKQ9ALr7SnefHaY3AQuArsBZwKMh26PAsDB9FjDR3Xe4+1JgMTAkrYUWEZEmJ6P3AM2sEDgSeBfo7O4rIQqSQKeQrSuwPG6x4pAmIiJSaxkLgGa2L/A0cLW7b6wqa4I0r2SdY8ysyMyK1qxZUx/FFBGRJiojAdDM8omC3+Pu/kxIXmVmXcL8LsDqkF4MdItbvAD4PNF63f0Bdx/s7oM7duyYmsKLiEiTkIlWoAb8BVjg7r+PmzUVuDhMXwxMiUsfYWbNzaw70BOYma7yiohI05SXgW1+Hfgv4EMzmxPSfgbcBkw2s9HAZ8C5AO4+z8wmA/OJWpCOdfeStJdaRESalLQHQHd/g8T39QBOqmSZccC4lBVKRESyjnqCERGRrKQAKCIiWUkBUEREspICoIiIZCUFQBERyUoKgCIikpUUAEVEJCspAIqISFbKRE8wDVbR6tlMW/Y863dsoF3ztpxeeBqDOw3MdLFERCQFFACDotWzmbzoKXbt3gXA+h0bmLzoKQAFQRGRJkiXQINpy54vC36ldu3exbRlz2eoRCIikkoKgMH6HRtqlC4iIo2bAmDQrnnbGqWLiEjjpgAYnF54Gvk5+eXS8nPyOb3wtAyVSEREUqnRBEAzO9XMFprZYjO7ob7XP7jTQE7sehzNc5sB0Dy3GSd2PU4NYEREmqhGEQDNLBe4FzgN6A1cYGa963MbyzYuYVPJWgZ0OpShXfoyoNOhbCpZy7KNS+pzMyIiKXHBg5fR6uxDsG8V0OrsQ7jgwcsyXaQGr1EEQGAIsNjdl7j7TmAicFZ9bmDuF7MpqTDQfImXMPeL2fW5GRGRenfBg5fx16ensW3TdgC2bdrOX5+epiBYjcYSALsCy+PeF4e0erM1tqVG6SIiDcWUaS9TEttdLq0ktpsp017OUIkah8YSAC1Bmu+VyWyMmRWZWdGaNWtqtIFWefvUKF1EpKEorfklmy6RxhIAi4Fuce8LgM8rZnL3B9x9sLsP7tixY4020L/9QHItt1xaruXSv70awYhIw9aydYsapUuksQTA94CeZtbdzJoBI4Cp9bmBwjY9OKrj18pqfK3y9uGojl+jsE2P+tyMiEi9O+v0b5KbV/5wnpuXw1mnfzNDJWocGkVfoO4eM7MrgH8CucBD7j6vvrdT2KaHAp6INDpPfv9+4DKmTHuZbZu207J1C846/ZshXSpj7nvdSmsSBg8e7EVFRZkuhohIo2Jms9x9cKbLkQ6N5RKoiIhIvVIAFBGRrKQAKCIiWUkBUEREslKTbQRjZmuAT2u5eAdgbT0WRySe9i9JpbruXwe7e80epG6kmmwArAszK8qWVlCSftq/JJW0fyVPl0BFRCQrKQCKiEhWUgBM7IFMF0CaNO1fkkrav5Kke4AiIpKVVAMUEZGs1CACoJkVmtk2M5sT3i8Lr8eb2XMp2uZbtVzuNTOr9xZW4Tt4LUwfa2bzzeyj+t6OiIhEGkQADP7t7gPStTF3/1q6tlXKrMKAg5Vw99eB01NcnKxW2UlXLdZzs5ldF6YfMbNzqsh7tZm1qs12qilDezN7ycwWhdd2If14M3skTJ9vZotTdUKZaRX/nyFtWXjNyhPpavItq+X6m9T+3pACYLz44dz3NbOnzOxjM3vczAzAzAaZ2b/MbJaZ/dPMuoT018zsDjObYWYLzOwoM3smfFn/r3SlZrY5vHYJeeeY2UdmdmzpfDP7nZnNNrPpZhb/YOi5ZjbTzD6Jy59rZr8xs/fMbK6Z/TCkH29mr5rZE8CHleUDSoAvUvR9SmJpPekCrgbq/YAA3ABMd/eewPTwvhx3nwR8PwXbbkjS/f9s0CfSDcDVNPD9vUEGQHc/Ku7tkURfZG+gB/B1M8sH/gCc4+6DgIeAcXHL7HT344D7gSnAWKAvcImZ7V9hcxcC/ww/nCOAOSF9H2C2uw8E/gX8Mm6ZPHcfEspVmj4a+DKU/SjgB2bWPcwbAtzk7r0ry+fuy9397OS/JalnZSddZna9mX1oZh+Y2W0h7Stm9kI44XrdzA6vycrN7ErgQODVcEI02szuiJv/AzP7fTiD/9jMHg0nSE+VnkVXdtIHnAU8GqYfBYaF6Z3Al7X4LpoKnUgn8d1k9f7u7hn/AwqBjxKkHw+8FPf+PuAiomC2kShYzQE+BF4MeV4Dvh6mT6yw/AxgQJjeHF6PAxYDN5fOC+klRIEOosA7J8H6OwOLw/RTwCdxZVoKnBw+w6tx602YL9nvRH8p3+dOA94CWoX37cPrdKBnmB4KvBKmbwauC9OPEJ2UVbbNZUCHML0P8G8gP7x/C+gXyuVx+9hDwHVAfsjTMaSfTzQwNMCGCttZX8n2jweey/R3n87/Z9zn/hIoIDrpfxs4pprv9DXg12H6KuBzoAvQHCgG9g/zSo8j1xKd5EI0aHfrMO3AyDD9C+CeuPX/LkyfDrwcpscA/xOmmwNFQPfwGbYA3avKV4vvLav398YwIvyOuOkSolHsDZjn7kdXs8zuCsvvDsuXcfcZZnYc8G3gMTP7jbuPT7DO+OdFStdZWh5CmX7k7v+MX8jMjifacakqnzQY3wQedvetAO7+hZntC3wN+GuoOEB00Kk1d99iZq8AZ5jZAqIDw4dmVggsd/c3Q9YJwJXAC0Qnfi+FMuQCK+tShiwz092LASy6T1gIbKDq73RqeP2Q6HizMiy/BOgGrIvL+x7wULg69Td3nxPSdwOTwvQE4Jm4ZUqnZ4XyQHTS3N/23FvbD+hJVLuZ6e5Lq8lXOj9ZWb2/N4YAmMhCoKOZHe3ub4ed7lB3n1fTFZnZwcAKd/+zme0DDATGE50pngNMJLpM+kY1q/oncLmZveLuu8zsUGBFsvncfUuCvJJ+RvmTHYj2hQ1e//eXHgR+BnwMPByXXnH7TtUnfavMrIu7rwyXiVbXczmbAp1IJ5bV+3uDvAdYHXffSRScfm1mHxBdSqztzejjgTlm9j7wXeCukL4F6GNms4gupd5SzXoeBOYDsy16fOFPJD7BSDafZMaLwKVx9yHau/tGYKmZnRvSzMyOqMW6NwGtS9+4+7tENYkLgSfj8h1kZqU//AuITr7KTvpCGfLNrE/IMxW4OExfTHTfW6pX1XdaI+FEerW7/xn4C9GJNOw5kYaanUjnh/UeGk7Ma5zPzLqa2fRqtpfV+3uDPvC6+2tE18pL318RNz2H6P5dxWWOr2L5+Hn7htdH2XNDteK6fg78vIr1ryVcunD33URnNz+rsJqKZagsnzQA7v6CmQ0AisxsJzCN6H81ErjPzP6H6P7EROCDGq7+AeB5M1vp7ieEtMlE957Xx+VbAFxsZn8CFgH3ufvOcLnrbjPbj+i3eycwD7gNmGxmo4HPgHNr+rmzUTXfaU0dD/zEzHYBm4FRIT3+RPpLontZVXmQ6Jgy26Jrf2vY08ijpvm6ALGqNpbt+3uD6ArNzLoR3fBcl4Jqd62Y2ebSIJmBbR8L/JHo+zg+E2Vo6sK9h+fcvW+Gy/EccIe7T091ucJltOvc/Yz6XnemNZT/Z0UZPo5cAXzm7lOrzZwmDW1/bxCXQD16BKBbQwl+sKeGmKFtv+7u/RT8UqoE2M/iHpxOJzNra2afANtKDwYp3t75RCdV66vL20hl9P/ZELn7PQ0l+DXU/b1B1ABFmhIze5ao6Xq8n6rlrzRFjXl/VwAUEZGs1CAugYqIiKSbAqCIiGQlBUBpssyss5k9YWZLLOpP8G0zG57pclUnNBj47zRt63iL+q+cZ2b/Ssc2RRoKBUBpksKzUX8DZrh7D486TR9B1B9kxbwN7XnYtkDKA6CZtSVqKXemu/dBzw9KllEAlKbqRKJRQe4vTXD3T939DwBmdomZ/dXM/g68aNEYY3+zqEf6d8ysf8hXNv5ZeP+RRT3YV9WL/W0WDWg818x+W1UhzayPRSMCzAn5exI96PuVkPYbM9vXopEEZlvUa/9Zccv/PJTjJTN70vaM1ZZMb/4XAs+4+2fh+1EXapJVGtqZr0h96QPMribP0UD/0AHwH4D33X2YmZ1I1B/sgGqWPwwY7e5vmtlDwH+H1+HA4e7uoZZVlcuAu9z9cTNrRtTp7w1A39LnYkMNdbi7bzSzDsA7ZjYVGETUfd+RRL/l2UQdK0PUC8dl7r7IzIYS1fROrLDtQ4F8iwZQbR3Kkaj/SpEmSQFQsoKZ3Us0BM5O3zPe5EvuXjp22jFEwQR3f8XM9g9dMFUlUS/2dwLbgQfN7B9AdSORvw3cZGYFRLWxRbanB/6y4gO3WtTZ8m6gK9FQXMcAU9x9W/iMfw+vyfbmn0cURE8CWgJvm9k77v5JNWUWaRJ0CVSaqnns6ZAYdx9LdKCPH5C0Yu/6FTlRX4rxv5MWFeaXy+/uMaIBkJ8m6pvxhaoK6e5PAGcC24B/htpnRSNDuQeFWuGqUI5EZYa43vzj/nolyFcMvODuW0K/tjOIBoUWyQoKgNJUvQK0MLPL49JaVZF/BlGgKe1DcG3oFX8ZIZCa2UDK93ixVy/2ofa1n7tPA64mXEY1s+Fm9quKGzWzHsASd7+bqJf7/lToRZ9orLfVYfisE4CDQ/obwHfMrEXY7rcBatCb/xTgWDPLC/cvhxJ1TCySFXQJVJqkcP9tGHCHmV1P1Fv+FuCnlSxyM/Cwmc0FtrJnuJWngVGhj8n3gPjLg3v1Yk8UrKaYWWkN7ZqQ9yvAxgTbPR+4yKJRBP4D3BLuSb5p0XBZzwO/Bv5uZkVEQ399HD7je+Fe4AfAp0Sjgn8Z1lttb/7uvsDMXgDmEl1afdDdP6rk+xFpctQVmkgtWA17sTezCcA17r6mnsuxr7tvDjW4GcAYd6+u8Y+IoBqgSFq4+0UpWvUDZtab6J7gowp+IslTDVBERLKSGsGIiEhWUgAUEZGspAAoIiJZSQFQRESykgKgiIhkJQVAERHJSk32OcAOHTp4YWFhposhItKozJo1a627d6w+Z+PXZANgYWEhRUVFmS6GiEijYmafZroM6aJLoCIikpUUAEVEJCspAIqISFZqsvcAE9m1axfFxcVs374900VpdFq0aEFBQQH5+fmZLoqISL3IqgBYXFxM69atKSwsxKyywbSlIndn3bp1FBcX07179+oXEJG0K1o9m2nLnmf9jg20a96W0wtPY3CngZkuVoOWVZdAt2/fzv7776/gV0Nmxv7776+as0gDVbR6NpMXPcX6HRsAWL9jA5MXPUXRao2OVZWsCoCAgl8t6XsTabimLXueXbt3lUvbtXsX05Y9n6ESNQ5ZFwBFRJqa0ppfsukSUQBMsxdeeIHDDjuMQw45hNtuuy1hnksuuYSnnnoqzSUTkcaqXfO2NUqXiAJgGpWUlDB27Fief/555s+fz5NPPsn8+fPTun0RaXpOLzyN/JzyLbTzc/I5vfC0DJWocVAArEJs9y62l2xje8lWtpdsI1bhGntNzZw5k0MOOYQePXrQrFkzRowYwZQpUxLmnTFjBl/72tfo0aNHWW3Q3fnJT35C37596devH5MmTQLgtdde44wzzihb9oorruCRRx4Boi7hbrnlFo455hj++te/cvfdd9O7d2/69+/PiBEj6vR5RKRhGNxpIOf1PKesxteueVvO63mOWoFWI6seg6iJ2O5dxDw+4Hn0fjfk5dTuWbgVK1bQrVu3svcFBQW8++67CfOuXLmSN954g48//pgzzzyTc845h2eeeYY5c+bwwQcfsHbtWo466iiOO+64arfbokUL3njjDQAOPPBAli5dSvPmzdmwYUOtPoeINDyDOw1UwKsh1QArEfNYjdKT4e57pVXWunLYsGHk5OTQu3dvVq1aBcAbb7zBBRdcQG5uLp07d+Yb3/gG7733XrXbPf/888um+/fvz8iRI5kwYQJ5eTr/EZHspQBYqb2DVdXp1SsoKGD58uVl74uLiznwwAMT5m3evPmeLYbAmSiAAuTl5bF79+6y9xWf19tnn33Kpv/xj38wduxYZs2axaBBg4jFah/QRUQaMwXASlX23Fvtn4c76qijWLRoEUuXLmXnzp1MnDiRM888M+nljzvuOCZNmkRJSQlr1qxhxowZDBkyhIMPPpj58+ezY8cOvvzyS6ZPn55w+d27d7N8+XJOOOEEbr/9djZs2MDmzZtr/XlERBozXQOrRJ7lVbgHuCe91uvMy+Oee+7hlFNOoaSkhEsvvZQ+ffokvfzw4cN5++23OeKIIzAzbr/9dg444AAAzjvvPPr370/Pnj058sgjEy5fUlLCRRddxJdffom7c80119C2bdtafx4RkcbMKrus1tgNHjzYKw6Iu2DBAnr16pX0OqKGMDGiy55GnuXVugFMU1DT709EGh8zm+XugzNdjnRQDbAKeTn55JG9AU9EpCnTPUAREclKCoAiIpKVFABFRCQrKQCKiEhWUgAUEZGspACYZoWFhfTr148BAwYweHDilsYaDklEJPX0GEQGvPrqq3To0CHt2y0pKSE3Nzft2xURaYhSVgM0s4fMbLWZfRSX1t7MXjKzReG1Xdy8G81ssZktNLNT4tIHmdmHYd7dVlnv0SmwNbaF1dtX8p9txazevpKtsS3p2rSGQxIRSbFU1gAfAe4Bxsel3QBMd/fbzOyG8P6nZtYbGAH0AQ4EXjazQ929BLgPGAO8A0wDTgWeT2G5gSj4bdy1gdLOr3d7SXgPrfL2qXS56pgZJ598MmbGD3/4Q8aMGZMwn4ZDEhFJrZTVAN19BvBFheSzgEfD9KPAsLj0ie6+w92XAouBIWbWBWjj7m971Gfb+LhlUmpzbCN7j/zgIb323nzzTWbPns3zzz/Pvffey4wZMxLm03BIIiKple5GMJ3dfSVAeO0U0rsCy+PyFYe0rmG6YnrK7faSGqUnq3T4o06dOjF8+HBmzpyZMJ+GQxIRSa2G0go00X09ryI98UrMxphZkZkVrVmzpk4FyrHEjUUqS0/Gli1b2LRpU9n0iy++SN++fZNeXsMhiYjUn3RfA1tlZl3cfWW4vLk6pBcD3eLyFQCfh/SCBOkJufsDwAMQjQZRl4Lum9em3D3AiLFvXptar3PVqlUMHz4cgFgsxoUXXsipp56a9PIaDklEpP6kdDgkMysEnnP3vuH9b4B1cY1g2rv79WbWB3gCGELUCGY60NPdS8zsPeBHwLtEjWD+4O7Tqtt2fQyHtDW2hc2xjez2EnIsl33z2tSpAUxjp+GQRJo+DYdUD8zsSeB4oIOZFQO/BG4DJpvZaOAz4FwAd59nZpOB+UAMGBtagAJcTtSitCVR68+UtwAt1Spvn6wOeCIiTVnKAqC7X1DJrJMqyT8OGJcgvQhI/kaZiIhIEhpKIxgREZG0UgAUEZGspAAoIiJZSQFQRESykgJgmt1xxx306dOHvn37csEFF+zVawtoOCQRkXRQAEyjFStWcPfdd1NUVMRHH31ESUkJEydOTNv2S0rq1o2biEhTogBYhY07N7Bk40I++fIjlmxcyMadG+q8zlgsxrZt24jFYmzdurWsb9CKNBySiEhqaTiASmzcuYFV21bgoSu0mO9i1bYVALRp1rZW6+zatSvXXXcdBx10EC1btuTkk0/m5JNPTphXwyGJiKSWaoCVWLt9VVnwK+U4a7evqvU6169fz5QpU1i6dCmff/45W7ZsYcKECQnzajikpqdo9WxumTmOa17/CbfMHEfR6tmZLpJIVlMArETMd9UoPRkvv/wy3bt3p2PHjuTn53P22Wfz1ltvJcyr4ZCalqLVs5m86CnW79gAwPodG5i86CkFQZEMUgCsRJ7l1yg9GQcddBDvvPMOW7duxd2ZPn16jTqX1nBIjde0Zc+za3f5k6ddu3cxbVnaurYVkQp0DawSHVp0LncPEMAwOrToXOt1Dh06lHPOOYeBAweSl5fHkUceyZgxY5JeXsMhNV6lNb9k00Uk9VI6HFIm1cdwSBt3bmDt9lXEfBd5lk+HFp1r3QCmKdBwSLV3y8xxCYNdu+Zt+cWQm9JfIJFKaDgkAaLWntkc8KT+nF54GpMXPVXuMmh+Tj6nF56WwVKJZDcFQJE0GNxpIBDdC1y/YwPtmrfl9MLTytJF6mrZxiXM/WI2W2NbaJW3D/3bD6SwTY9MF6tBUwAUSZPBnQYq4ElKLNu4hPfWvEVJGEd8a2wL762JWpgrCFZOrUBFRBq5uV/MLgt+pUq8hLlf6DGbqigAiog0cltjW2qULhEFQBGRRq5V3j41SpeIAmCaXXrppXTq1Im+ffvuNe8Pf/gDhx12GH369OH666/fa37FTq9FRAD6tx9IruWWS8u1XPq31z3nqqgRTJpdcsklXHHFFYwaNapc+quvvsqUKVOYO3cuzZs3Z/Xq1WkrUywWU7+gIo1YaUMXtQKtGR31qpCKZsXHHXccy5Yt2yv9vvvu44YbbijrA7RTp04Jl9+8eTPnnHMOH330EYMGDWLChAmYGdOnT+e6664jFotx1FFHcd9999G8eXMKCwspKiqiQ4cOFBUVcd111/Haa69x88038/nnn7Ns2TI6dOjATTfdxPe+9z127tzJ7t27efrpp+nZs2edPquIpE9hmx4KeDWkS6CVKG1WXHoTubRZ8bKNS1KyvU8++YTXX3+doUOHVjnKw/vvv8+dd97J/PnzWbJkCW+++Sbbt2/nkksuYdKkSXz44YfEYjHuu+++arc5a9YspkyZwhNPPMH999/PVVddxZw5cygqKqKgoKC+P6KISIOSkQBoZteY2Twz+8jMnjSzFmbW3sxeMrNF4bVdXP4bzWyxmS00s1PSUcZ0NyuOxWKsX7+ed955h9/85jecd955CUd/GDJkCAUFBeTk5DBgwACWLVvGwoUL6d69O4ceeigAF198MTNmzKh2m2eeeSYtW7YE4Oijj+bWW2/l17/+NZ9++mlZuog0DqkYwLupS3sANLOuwJXAYHfvC+QCI4AbgOnu3hOYHt5jZr3D/D7AqcAfzSrc7U2BdDcrLigo4Oyzz8bMGDJkCDk5Oaxdu3avfPHDJOXm5hKLxSodJgnKD5VU1TBJF154IVOnTqVly5accsopvPLKK3X9SCKSJqUDeJcO11Y6gLeCYNUydQk0D2hpZnlAK+Bz4Czg0TD/UWBYmD4LmOjuO9x9KbAYGJLqAqa7WfGwYcPKgs4nn3zCzp076dChQ1LLHn744SxbtozFixcD8Nhjj/GNb3wDgMLCQmbNmgXA008/Xek6lixZQo8ePbjyyis588wzmTt3bl0+joikUSoG8M4GaQ+A7r4C+C3wGbAS+NLdXwQ6u/vKkGclUNoKpCuwPG4VxSFtL2Y2xsyKzKxozZo1dSpnqpoVX3DBBRx99NEsXLiQgoIC/vKXvwDR4xFLliyhb9++jBgxgkcffRQzS2qdLVq04OGHH+bcc8+lX79+5OTkcNlllwHwy1/+kquuuopjjz2W3NzKK86TJk2ib9++DBgwgI8//nivVqoi0nClYgDvbJD24ZDCvb2ngfOBDcBfgaeAe9y9bVy+9e7ezszuBd529wkh/S/ANHevvDpD/QyHpM5ly9NwSCIN05KNCxMGuzzLp0ebw2q0Lg2HlFrfBJa6+xoAM3sG+Bqwysy6uPtKM+sClD4IVwx0i1u+gOiSacqpWbGINAapGMA7G2TiHuBnwFfNrJVF1/hOAhYAU4GLQ56LgSlheiowwsyam1l3oCcwM81lFhFpsNo0a0vnll3Js3wgqvl1btlV45lWI+01QHd/18yeAmYDMeB94AFgX2CymY0mCpLnhvzzzGwyMD/kH+te4fkEEZEspwG8ay4jPcG4+y+BX1ZI3kFUG0yUfxwwLtXlEhGR7KGeYEREJCspAIqISFZSAEyj5cuXc8IJJ9CrVy/69OnDXXfdVTbv5z//Of3792fAgAGcfPLJfP753g1dNRySiEj9UQBMo7y8PH73u9+xYMEC3nnnHe69917mz58PwE9+8hPmzp3LnDlzOOOMM7jlllvSVq5YLJa2bYmINBQKgFUoWj2bW2aO45rXf8ItM8dRtLpuHWF36dKFgQOjnmRat25Nr169WLFiBQBt2rQpy7dly5ZKe4EpHQ7p8MMPZ+TIkWX9gE6fPp0jjzySfv36cemll7Jjxw4g6gqttE/RoqIijj/+eABuvvlmxowZw8knn8yoUaOYN28eQ4YMYcCAAfTv359FixbV6bOKiDR0Gg+wEkWrZzN50VPs2h31rrB+xwYmL3oKgMGd6j7K8rJly3j//fcZOnRoWdpNN93E+PHj2W+//Xj11VcTLvf+++8zb948DjzwQL7+9a/z5ptvMnjwYC655BKmT5/OoYceyqhRo7jvvvu4+uqrqyzDrFmzeOONN2jZsiU/+tGPuOqqqxg5ciQ7d+6kpERPmohI06YaYCWmLXu+LPiV2rV7F9OWPV/ndW/evJnvfve73HnnneVqfuPGjWP58uWMHDmSe+65J+GyGg5JRKR+KABWYv2ODTVKT9auXbv47ne/y8iRIzn77LMT5rnwwgsrHblBwyGJiNQPBcBKtGvetkbpyXB3Ro8eTa9evfjxj39cbl78PbepU6dy+OGHJ71eDYckIlJzVd4DNLMfVzXf3X9fv8VpOE4vPK3cPUCA/Jx8Ti88rdbrfPPNN3nsscfo168fAwYMAODWW2/l9NNP54YbbmDhwoXk5ORw8MEHc//99ye93vjhkGKxGEcddVS54ZBGjx7NrbfeWu5+Y0WTJk1iwoQJ5Ofnc8ABB/CLX/yi1p9TRKQxqHI4JDMr7a7sMOAooo6pAb4DzHD376e2eLVXH8MhFa2ezbRlz7N+xwbaNW/L6YWn1UsDmMZKwyGJNH0aDilw9/8FMLMXgYHuvim8v5loHL8mbXCngVkd8EREmrJk7wEeBOyMe78TKKz30oiIiKRJss8BPgbMNLNnw/thwKMpKZFIE7Vs4xLmfjGbrbEttMrbh/7tB2rAZZEMSioAuvs4M3seOBZw4Hvu/n5KSybShCzbuIT31rxFSRjKcmtsC++teQtAQVAkQ2ryGEQJsDvuT0SSNPeL2WXBr1SJlzD3i7p1rycitZdUADSzq4DHgQ5AJ2CCmf0olQUTaUq2xrbUKF1EUi/ZGuBoYKi7/9LdfwF8FfhB6orVNFU1HNL555/PgAEDGDBgAIWFhWXPCcbTcEiNV6u8fWqULlJTW2NbWL19Jf/ZVszq7St1cpWEZBvBGNEl0FIlIU1qoHQ4pIEDB7Jp0yYGDRrEt771LXr37s2kSZPK8l177bXst99+aStXLBYjL0/9oqdS//YDy90DBMi1XPq312M2UndbY1vYuGsDURMN2O0l4b1OsqqSbA3wYeBdM7vZzP4XeAf4S+qK1TC8uvxtvvfi9Xxnymi+9+L1vLr87Tqtr6rhkEq5O5MnT+aCCy5IuA4Nh9Q4FbbpQd92R9AspxkAzXKa0bfdEWoAI/Vic2wjpcFvDw/pUplkW4H+3sxeA44JSU2+Feiry9/mng/Gs6MkevxxzbZ13PPBeABO6HZ0ndefaDgkgNdff53OnTvTs2fPhMtpOKTGaePODeTl5tB3/35laYaxcecG2jRrm7mCSZOw2xP/XitLl0hNW4E6WdIKdPyCZ8uCX6kdJTsZv+DZSpZIXmXDIQE8+eSTldb+QMMhNVZrt6/CK5yhO87a7asyVCJpSnIst0bpElEr0Eqs3bauRunJqmo4pFgsxjPPPMP5559f6fIaDqlxivmuGqWL1MS+eW3Yu1mGhXSpTEZagZpZWzN7ysw+NrMFZna0mbU3s5fMbFF4bReX/0YzW2xmC83slNputyY6tNy/RunJqGo4JICXX36Zww8/nIKCghqtV8MhNXx5ll+jdJGaaJW3D23y25bV+HIslzb5bdUAphrJBsD6bgV6F/CCux8OHAEsAG4Aprt7T2B6eI+Z9QZGAH2AU4E/mqW+Xj+q13Ca5zYrl9Y8txmjeg2v9TpLh0N65ZVXyh55mDZtWtn8iRMnVnn5szLxwyH169ePnJyccsMhXXXVVRx77LHk5lb+tU2aNIm+ffsyYMAAPv74Y0aNGlXzDyiV6tCiM1bhJ2MYHVp0zlCJpKlplbcPnVp04YCWBXRq0UXBLwlVDodUlikaF/BiIL4v0Efc/c4ab9CsDfAB0MPjNm5mC4Hj3X2lmXUBXnP3w8zsRgB3/1XI90/gZnevsklmfQyH9Orytxm/4FnWbltHh5b7M6rX8HppANNYaTikutm4cwNrt68i5rvIs3w6tOisBjDS4Gg4pApCK9B/AV8nqvnVpRVoD2AN8LCZHQHMAq4COrv7yrC9lWbWKeTvSvTYRanikLYXMxsDjAE46KCDalm8PU7odnRWBzypX22atVXAE2lAatIKdA7wFFEtcJ2Z1TbC5AEDgfvc/UhgC+FyZyUSXWpNWG119wfcfbC7D+7YsWMtiyciItkgqRpgaPH5S2AVe+7/OdC/FtssBord/d3w/imiALjKzLrEXQJdHZe/W9zyBcDntdiuiIhImWRrgFcBh7l7H3fv7+793L02wQ93/w+w3MwOC0knAfOBqUT3GQmvU8L0VGCEmTU3s+5AT2BmbbYtIiJSKtkOIJcDX9bjdn8EPG5mzYAlwPeIgvFkMxsNfAacC+Du88xsMlGQjAFj3dW9gYiI1E2VATC0/oQoSL1mZv8AdpTOd/ff12aj7j4HSNTK6KRK8o8DxtVmWyIi2eDx6c9y00O38dmazzmo44GMu/QGRp5U+8e2skF1l0Bbh7/PgJeAZnFprVNbtKarpKSEI488stzQRjfffDNdu3ZN+HxgKQ2HJCKJPD79WcbccT2frl6Bu/Pp6hWMueN6Hp9e964bm7Iqa4Du/r/pKkg2ueuuu+jVqxcbN5bvqf2aa67huuuuS3t5NBySSON200O3sXXHtnJpW3ds46aHblMtsApV1gDN7M7w+nczm1rxLy0lzKDHpz9L4cih5JzcjcKRQ+vlbKq4uJh//OMffP/736/V8hoOSUQq+mxN4obxlaVLpLrT/sfC629TXZCGpvSSQulZVeklBaBOZ1RXX301t99+O5s2bdpr3j333MP48eMZPHgwv/vd72jXrt1eeTQckohU1LH9/qxetzZhulSuyhqgu88Kr/9K9JeeImZGVZcUauu5556jU6dODBo0aK95l19+Of/+97+ZM2cOXbp04dprr024Dg2HJCIVHfbVnuTmlT+c5+blcNhXE48rKpHqLoF+aGZzE/x9aGZNeriAVFxSePPNN5k6dSqFhYWMGDGCV155hYsuugiAzp07k5ubS05ODj/4wQ+YOTPxo44aDklEKtqvexv6ndCHlq1bANCydQv6ndCH/bprOKSqVHcJNGubHB7U8UA+Xb0iYXpt/epXv+JXv/oVELXo/O1vf8uECRMAWLlyJV26dAHg2WefpW/fvkmvN344pEMOOSThcEinnXZa0sMhLVmyhLlz53LiiSfW9qOKSBp1aLk/fhh0PezAvdKlctVdAv209C8k9QzTq4EvUl66DBp36Q20al7+MmCr5i0Zd2lV3ZbW3vXXX0+/fv3o378/r776KnfccUfSy2o4JJHslorh27JBssMh/YBolIX27v4VM+sJ3O/uCR9cbwjqYzgkPVhanoZDqputsS1sjm1kt5eQY7nsm9dGY7ZJvamv4ds0HNLexgJDgHcB3H1R3HBFTdbIk4ZndcCT+rM1toWNuzZQOpDJbi8J71EQlHqh4dtqLtnOsHe4+87SN2aWRyVDEonI3jbHNrL3T8ZDuohkQrIB8F9m9jOgpZl9C/gr8PfUFSt1krnkK3vT91Y3uyvpv72ydBFJvWQD4A1Eo7h/CPwQmObuN6WsVCnSokUL1q1bp4N5Dbk769ato0WLFpkuSqOVY4kbIFWWLiKpl+w9wJvd/RfAnwHMLNfMHnf3kakrWv0rKCiguLiYNWvWZLoojU6LFi0oKCjIdDEarX3z2pS7Bxgx9s3Tc1oimZJsADzIzG5091+FMfz+CryfwnKlRH5+Pt27d890MSQLlTZ0UStQkYYj2QD4PaIBbG8ETgCed/fkH1QTEVrl7aOAJ9KAVDcg7sC4t3cBfwLeJGoUM9DdZ6eycCIiIqlSXQ3wdxXerwd6h3QH1FeWiIg0StUNiHtCugoiIiKSTtVdAr3I3SeY2Y8TzXf336emWCIiIqlV3SXQ0jv2rRPM08N0IiLSaFV3CfRP4fV/K84zs6tTVCaRJkmdq4s0LMn2BJNIwsuiyQoP079vZs+F9+3N7CUzWxRe28XlvdHMFpvZQjM7pS7bFcmEx6c/y5g7rufT1Stwdz5dvYIxd1zP49OfzXTRRLJWXQKg1XHbVwEL4t7fAEx3957A9PAeM+sNjAD6AKcCfzRT/1HSuNz00G1s3bGtXNrWHdu46aHbMlQiEalLAKz1PUAzKwC+DTwYl3wW8GiYfhQYFpc+0d13uPtSYDHR0EwijcZnqz+vUbqIpF51rUA3kTjQGdAyQXqy7gSup3zjms7uvhLA3VfGjTfYFXgnLl9xSBNpNFq1acmWjVsTpovUh9juXcQ8RnTINvIsj7yc/EwXq0Grsgbo7q3dvU2Cv9bunmw3auWY2RnAaneflewiiYpWybrHmFmRmRWpw2tpSHoO/Qq5eeV/brl5OfQc+pUMlUiakij47WLPodGJ+S5iu3dlslgNXq2CWB19HTjTzE4HWgBtzGwCsMrMuoTaXxdgdchfDHSLW74ASHjdyN0fAB4AGDx4sB7TkAZjwIB+ACx8ZxHbNm2nZesWHPbVnmXpInUR1fwSp+ehWmBl0h4A3f1G4EYAMzseuM7dLzKz3wAXA7eF1ylhkanAE2b2e+BAoCcwM83FFqmTUb2Gs3HnJroedmBZWvPcZozqpccgpD5Udr6vekBVMlEDrMxtwGQzGw18BpwL4O7zzGwyMB+IAWPdNYy2NC4ndDsagPELnmXttnV0aLk/o3oNL0sXqRuj8uYaUhlrqqOjDx482IuKijJdDBGRlNtzD7C8PMuvcUMYM5vl7oPrq2wNWV0egxARkQYgLyefPMtnT43PahX8sk1DugQqIiK1lJeTrwYvNaQaoIiIZCUFQBERyUoKgCIikpUUAEVEJCspAIqISFZSABQRkaykACgiIllJAVBERLKSAqCIiGQlBUAREclKCoAiIpKV1BeoSJpEPfbHiIatMfIsT50Vi2SQaoAiabBnuJrS4cecmO8itnvvIWxEJD0UAEXSIKr5JZ8uIqmnACiSFpUNPN00B6QWaQwUAEXSwmqYLiKppgAokgZ5lri9WWXpIpJ6CoAiaZCXk0+e5bOnxmfkWb5agYpkkE4/RdIkLyefPBTwRBoK1QBFRCQrKQCKiEhWSnsANLNuZvaqmS0ws3lmdlVIb29mL5nZovDaLm6ZG81ssZktNLNT0l1mERFpejJRA4wB17p7L+CrwFgz6w3cAEx3957A9PCeMG8E0Ac4FfijmeVmoNwiItKEpD0AuvtKd58dpjcBC4CuwFnAoyHbo8CwMH0WMNHdd7j7UmAxMCSthRYRkSYno/cAzawQOBJ4F+js7ishCpJAp5CtK7A8brHikJZofWPMrMjMitasWZOycouISOOXsQBoZvsCTwNXu/vGqrImSEvYf5S7P+Dug919cMeOHeujmCIi0kRlJACaWT5R8Hvc3Z8JyavMrEuY3wVYHdKLgW5xixcAn6errCIi0jRlohWoAX8BFrj77+NmTQUuDtMXA1Pi0keYWXMz6w70BGamq7wiItI0ZaInmK8D/wV8aGZzQtrPgNuAyWY2GvgMOBfA3eeZ2WRgPlEL0rHuXpL2UouISJOS9gDo7m9QeRf4J1WyzDhgXMoKJSIiWUc9wYiISFZSABQRkaykACgiIllJAVBERLKSAqCIiGQlBUAREclKCoAiIpKVFABFRCQrKQCKiEhWUgAUEZGspAAoIiJZSQFQRESykgKgiIhkJQVAERHJSgqAIiKSlRQARUQkKykAiohIVlIAFBGRrKQAKCIiWUkBUEREspICoIiIZCUFQBERyUoKgCIikpUaTQA0s1PNbKGZLTazG1KxjQsevIxWZx+CfauAVmcfwgUPXpaKzYiISAPQKAKgmeUC9wKnAb2BC8ysd31u44IHL+OvT09j26btAGzbtJ2/Pj1NQVBEpIlqFAEQGAIsdvcl7r4TmAicVZ8bmDLtZUpiu8ullcR2M2Xay/W5GRERaSAaSwDsCiyPe18c0soxszFmVmRmRWvWrKnRBkprfsmmi4hI49ZYAqAlSPO9EtwfcPfB7j64Y8eONdpAy9YtapQuIiKNW2MJgMVAt7j3BcDn9bmBs07/Jrl55b+O3Lwczjr9m/W5GRERaSAaSwB8D+hpZt3NrBkwAphanxt48vv3c+53Ty+r8bVs3YJzv3s6T37//vrcjIiINBB5mS5AMtw9ZmZXAP8EcoGH3H1efW/nye/fD9+v77WKiEhD1CgCIIC7TwOmZbocIiLSNDSWS6AiIiL1SgFQRESykgKgiIhkJXPf63G6JsHM1gCf1nLxDsDaeiyOSDztX5JKdd2/Dnb3mj1I3Ug12QBYF2ZW5O6DM10OaZq0f0kqaf9Kni6BiohIVlIAFBGRrKQAmNgDmS6ANGnavySVtH8lSfcARUQkK6kGKCIiWUkBUEREslKDCIBmVmhm28xsTni/LLweb2bPpWibb9VyudfMrN6bGIfv4LUwfayZzTezj+p7OyIiEmkQATD4t7sPSNfG3P1r6dpWKTPLTSafu78OnJ7i4mS1yk66arGem83sujD9iJmdU0Xeq82sVW22U00Z2pvZS2a2KLy2C+nHm9kjYfp8M1ucqhPKTKv4/wxpy8JrVp5IV5NvWS3X36T294YUAOOtiZve18yeMrOPzexxMzMAMxtkZv8ys1lm9k8z6xLSXzOzO8xshpktMLOjzOyZ8GX9v9KVmtnm8Nol5J1jZh+Z2bGl883sd2Y228ymm1l8zwjnmtlMM/skLn+umf3GzN4zs7lm9sOQfryZvWpmTwAfVpYPKAG+SNH3KYml9aQLuBqo9wMCcAMw3d17AtPD+3LcfRJNf7CvdP8/G/SJdANwNQ18f2+QAdDdj4p7eyTRF9kb6AF83czygT8A57j7IOAhYFzcMjvd/TjgfmAKMBboC1xiZvtX2NyFwD/DD+cIYE5I3weY7e4DgX8Bv4xbJs/dh4RylaaPBr4MZT8K+IGZdQ/zhgA3uXvvyvK5+3J3Pzv5b0nqWdlJl5ldb2YfmtkHZnZbSPuKmb0QTrheN7PDa7JyM7sSOBB4NZwQjTazO+Lm/8DMfh/O4D82s0fDCdJTpWfRlZ30AWcBj4bpR4FhYXon8GUtvoumQifSSXw3Wb2/u3vG/4BC4KME6ccDL8W9vw+4iCiYbSQKVnOAD4EXQ57XgK+H6RMrLD8DGBCmN4fX44DFwM2l80J6CVGggyjwzkmw/s7A4jD9FPBJXJmWAieHz/Bq3HoT5kv2O9Ffyve504C3gFbhffvwOh3oGaaHAq+E6ZuB68L0I0QnZZVtcxnQIUzvA/wbyA/v3wL6hXJ53D72EHAdkB/ydAzp5xMNDA2wocJ21ley/eOB5zL93afz/xn3ub8ECohO+t8GjqnmO30N+HWYvgr4HOgCNAeKgf3DvNLjyLVEJ7kQDdrdOkw7MDJM/wK4J279vwvTpwMvh+kxwP+E6eZAEdA9fIYtQPeq8tXie8vq/b0xDIi7I266hGgQXwPmufvR1Syzu8Lyu6kwCLC7zzCz44BvA4+Z2W/cfXyCdcY/MFm6ztLyEMr0I3f/Z/xCZnY80Y5LVfmkwfgm8LC7bwVw9y/MbF/ga8BfQ8UBooNOrbn7FjN7BTjDzBYQHRg+NLNCYLm7vxmyTgCuBF4gOvF7KZQhF1hZlzJkmZnuXgxg0X3CQmADVX+nU8Prh0THm5Vh+SVAN2BdXN73gIfC1am/ufuckL4bmBSmJwDPxC1TOj0rlAeik+b+tufe2n5AT6LazUx3X1pNvtL5ycrq/b0xBMBEFgIdzexod3877HSHuvu8mq7IzA4GVrj7n81sH2AgMJ7oTPEcYCLRZdI3qlnVP4HLzewVd99lZocCK5LN5+5bEuSV9DPKn+xAtC9s8Pq/v/Qg8DPgY+DhuPSK23eqPulbZWZd3H1luEy0up7L2RToRDqxrN7fG+Q9wOq4+06i4PRrM/uA6FJibW9GHw/MMbP3ge8Cd4X0LUAfM5tFdCn1lmrW8yAwH5ht0eMLfyLxCUay+SQzXgQujbsP0d7dNwJLzezckGZmdkQt1r0JaF36xt3fJapJXAg8GZfvIDMr/eFfQHTyVXbSF8qQb2Z9Qp6pwMVh+mKi+95Svaq+0xoJJ9Kr3f3PwF+ITqRhz4k01OxEOj+s99BwYl7jfGbW1cymV7O9rN7fG/SB191fI7pWXvr+irjpOUT37youc3wVy8fP2ze8PsqeG6oV1/Vz4OdVrH8t4dKFu+8mOrv5WYXVVCxDZfmkAXD3F8xsAFBkZjuBaUT/q5HAfWb2P0T3JyYCH9Rw9Q8Az5vZSnc/IaRNJrr3vD4u3wLgYjP7E7AIuM/dd4bLXXeb2X5Ev907gXnAbcBkMxsNfAacW9PPnY2q+U5r6njgJ2a2C9gMjArp8SfSXxLdy6rKg0THlNkWXftbw55GHjXN1wWIVbWxbN/fG0RfoGbWjeiG57oUVLtrxcw2lwbJDGz7WOCPRN/H8ZkoQ1MX7j085+59M1yO54A73H16qssVLqNd5+5n1Pe6M62h/D8ryvBx5ArgM3efWm3mNGlo+3uDuATq0SMA3RpK8IM9NcQMbft1d++n4JdSJcB+FvfgdDqZWVsz+wTYVnowSPH2zic6qVpfXd5GKqP/z4bI3e9pKMGvoe7vDaIGKNKUmNmzRE3X4/1ULX+lKWrM+7sCoIiIZKUGcQlUREQk3RQARUQkKykASpNlZp3N7AkzW2JRf4Jvm9nwTJerOqHBwH+naVvHW9R/5Twz+1c6tinSUCgASpMUno36GzDD3Xt41Gn6CKL+ICvmbWjPw7YFUh4AzawtUUu5M929D3p+ULKMAqA0VScSjQpyf2mCu3/q7n8AMLNLzOyvZvZ34EWLxhj7m0U90r9jZv1DvrLxz8L7jyzqwb6qXuxvs2hA47lm9tuqCmlmfSwaEWBOyN+T6EHfr4S035jZvhaNJDDbol77z4pb/uehHC+Z2ZO2Z6y2ZHrzvxB4xt0/C9+PulCTrNLQznxF6ksfYHY1eY4G+ocOgP8AvO/uw8zsRKL+YAdUs/xhwGh3f9PMHgL+O7wOBw53dw+1rKpcBtzl7o+bWTOiTn9vAPqWPhcbaqjD3X2jmXUA3jGzqcAgou77jiT6Lc8m6lgZol44LnP3RWY2lKimd2KFbR8K5Fs0gGrrUI5E/VeKNEkKgJIVzOxeoiFwdvqe8SZfcvfSsdOOIQomuPsrZrZ/6IKpKol6sb8T2A48aGb/AKobifxt4CYzKyCqjS2yPT3wlxUfuNWizpZ3A12JhuI6Bpji7tvCZ/x7eE22N/88oiB6EtASeNvM3nH3T6ops0iToEug0lTNY0+HxLj7WKIDffyApBV716/IifpSjP+dtKgwv1x+d48RDYD8NFHfjC9UVUh3fwI4E9gG/DPUPisaGco9KNQKV4VyJCozxPXmH/fXK0G+YuAFd98S+rWdQTQotEhWUACUpuoVoIWZXR6X1qqK/DOIAk1pH4JrQ6/4ywiB1MwGUr7Hi716sQ+1r/3cfRpwNeEyqpkNN7NfVdyomfUAlrj73US93PenQi/6RGO9rQ7DZ50AHBzS3wC+Y2Ytwna/DVCD3vynAMeaWV64fzmUqGNikaygS6DSJIX7b8OAO8zseqLe8rcAP61kkZuBh81sLrCVPcOtPA2MCn1MvgfEXx7cqxd7omA1xcxKa2jXhLxfATYm2O75wEUWjSLwH+CWcE/yTYuGy3oe+DXwdzMrIhr66+PwGd8L9wI/AD4lGhX8y7Deanvzd/cFZvYCMJfo0uqD7v5RJd+PSJOjrtBEasFq2Iu9mU0ArnH3NfVcjn3dfXOowc0Axrh7dY1/RATVAEXSwt0vStGqHzCz3kT3BB9V8BNJnmqAIiKSldQIRkREspICoIiIZCUFQBERyUoKgCIikpUUAEVEJCv9f9HEEXE1pVoNAAAAAElFTkSuQmCC\n"
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(2)\n",
    "out_index = wit_results.index.levels[0]\n",
    "cm = plt.cm.get_cmap('Greens')\n",
    "colors = cm(ages/max(ages))\n",
    "for i, key in enumerate(out_index):\n",
    "    data = wit_results.loc[key]\n",
    "    bic_data = data['bic'].to_numpy()\n",
    "    norm_bic = bic_data - np.min(bic_data)\n",
    "    like_data = data['likelihood'].to_numpy()\n",
    "    norm_like = like_data - np.min(like_data)\n",
    "    _ = ax[0].scatter(data['bic'].index,\n",
    "                     norm_bic,\n",
    "                     label=str(ages[i]) + ' hours',\n",
    "                     color=np.array(colors)[i])\n",
    "    _ = ax[1].scatter(data['likelihood'].index,\n",
    "                     norm_like,\n",
    "                     label=str(ages[i]) + ' hours',\n",
    "                     color=np.array(colors)[i])\n",
    "    _ = ax[0].set(xlabel='Groups, stage ' + str(i), ylabel='-BIC')\n",
    "    _ = ax[1].set(xlabel='Groups, stage ' + str(i), ylabel='Likelihood')\n",
    "_ = fig.suptitle('BIC and Likelihood: C. Elegans at Different Developmental Stages')\n",
    "_ = fig.set_size_inches(6, 8)\n",
    "_ = ax[0].legend(loc='best')\n",
    "_ = ax[1].legend(loc='best')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Once again, we find that cell type explains the topology the best. This is still not entirely surprising, since the\n",
    "connectivity of interneurons should be very different than that of sensory or motor neurons"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## C. Elegans Male and Hermaphrodite Conectomes (Worm Wiring) ##"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "# load chemical synapse connectome files\n",
    "path_h = '../json_connectomes/worm_wiring/connectome/Hermaphrodite/0.json'\n",
    "path_m = '../json_connectomes/worm_wiring/connectome/Male/0.json'\n",
    "ww_connectomes = []\n",
    "ww_connectomes.append(GraphIO.load(path_h)[0])\n",
    "ww_connectomes.append(GraphIO.load(path_m)[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "                                                                bic  \\\nC. Elegans Hermaphrodite ['hemisphere']               -45484.654884   \n                         ['cell_type0']               -37783.419558   \n                         ['hemisphere', 'cell_type0'] -41098.822390   \nC. Elegans Male          ['hemisphere']               -53425.404138   \n                         ['cell_type0']               -43612.375238   \n                         ['hemisphere', 'cell_type0'] -47806.496439   \n\n                                                         likelihood  n_params  \\\nC. Elegans Hermaphrodite ['hemisphere']               -22687.264567         9   \n                         ['cell_type0']               -18500.151558        64   \n                         ['hemisphere', 'cell_type0'] -17588.252151       484   \nC. Elegans Male          ['hemisphere']               -26655.512739         9   \n                         ['cell_type0']               -21399.507936        64   \n                         ['hemisphere', 'cell_type0'] -20541.786468       529   \n\n                                                            estimator  \nC. Elegans Hermaphrodite ['hemisphere']                SBMEstimator()  \n                         ['cell_type0']                SBMEstimator()  \n                         ['hemisphere', 'cell_type0']  SBMEstimator()  \nC. Elegans Male          ['hemisphere']                SBMEstimator()  \n                         ['cell_type0']                SBMEstimator()  \n                         ['hemisphere', 'cell_type0']  SBMEstimator()  ",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th></th>\n      <th>bic</th>\n      <th>likelihood</th>\n      <th>n_params</th>\n      <th>estimator</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">C. Elegans Hermaphrodite</th>\n      <th>['hemisphere']</th>\n      <td>-45484.654884</td>\n      <td>-22687.264567</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-37783.419558</td>\n      <td>-18500.151558</td>\n      <td>64</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-41098.822390</td>\n      <td>-17588.252151</td>\n      <td>484</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th rowspan=\"3\" valign=\"top\">C. Elegans Male</th>\n      <th>['hemisphere']</th>\n      <td>-53425.404138</td>\n      <td>-26655.512739</td>\n      <td>9</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['cell_type0']</th>\n      <td>-43612.375238</td>\n      <td>-21399.507936</td>\n      <td>64</td>\n      <td>SBMEstimator()</td>\n    </tr>\n    <tr>\n      <th>['hemisphere', 'cell_type0']</th>\n      <td>-47806.496439</td>\n      <td>-20541.786468</td>\n      <td>529</td>\n      <td>SBMEstimator()</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "groups = [[\"hemisphere\"], [\"cell_type0\"], [\"hemisphere\", \"cell_type0\"]]\n",
    "types = ['Hermaphrodite', 'Male']\n",
    "rows = {}\n",
    "for i in range(len(ww_connectomes)):\n",
    "    g_row = {}\n",
    "    for group in groups:\n",
    "        g_row[str(group)] = fit_apriori_sbm(ww_connectomes[i], group, plot_sbms=False)\n",
    "    rows['C. Elegans ' + types[i]] = g_row\n",
    "ww_results = pd.DataFrame.from_dict({(i,j): rows[i][j]\n",
    "                           for i in rows.keys()\n",
    "                           for j in rows[i].keys()},\n",
    "                       orient='index')\n",
    "ww_results.head(6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One last time, we'll plot -BIC and Likelihood, the lowest -BIC is again set to 0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "<Figure size 432x576 with 2 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcYAAAIaCAYAAAC6du3gAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABbkElEQVR4nO3de5xd873/8ddbhMQlRIRGQpM6oTJkIjM0KIIipa1QKspBf061Dq1QRR0qmjqH07REe6ooRYU0UinHrdII0SMqGYLELVEhNxERcY3cPr8/1nfGmsncb3sm834+Hvux1/6u9f2u79p77fVZ3++6KSIwMzOzzCaFroCZmVlb4sBoZmaW48BoZmaW48BoZmaW48BoZmaW48BoZmaW0yECo6RRku5oRL7fSbosDQ+VtLARZVTKJ2mOpKFNqVcj6tBXUkjatKXn1ZZIukTS72sZv4ukDyV1as16NVVH/T0LRdLpkv7eguVXbAfa8jrZ0t9Da82jPhodGCXNl/RJ+hFXSHpA0s658bdK+nnu82ZpBZgr6aOU/xZJfZu4DE1SW8CLiO9HxOjmnF9EFEXEY81ZZnOS9BNJD1ZJm1tD2ohWrttfJV2Y+9w7BYjq0j4XEf8ZEf9WU3kR8WZEbBUR61qovvtKelDSe5LelfS0pO/UM++tklan/1f567mWqGdLq2ljl7YBXylEndqqquukpMck1bgOtxW5nbVnqqRvn9bj+QWqWqM0tcX49YjYCugFLAV+Xcu0E4FvAN8GtgGKgTLgsCbWwZrXNOCA8j1WSZ8DOgODq6T9S5q23pqhhTMNODj3+SDg5WrS5kbEWy1cl1pJ2g94FHic7LvqAZwFfLUBxfx32kiWv4pboKrtTltuKSvTIXriarClpD1zn78NvF6oyjRWs/yAEbGKLPANqG582is8HDgmImZExNqIWBkR/xMRN9eQ52JJr0n6QNKLko7NjTtd0t8ljUmt1dclfTU3vp+kx1PeycD2jVmuqq3eKuN+mOrVR9LmqS5vSlqaumC71pCv6l7yZpJuT3WdI6k0N+0eaY/xvTTuG7lx26R8yyS9IenS8j+kpE6pPu9I+idwdAMWewZZIByUPh8ETAVeqZL2WkQslrSTpPtSi2iepO/m6jhK0kRJd0h6Hzg9Lc/PJT2ZWkH/K6mHpHGS3pc0o5ZehPKgXb7eHghcC5RWSZuWm395F1X5Hu0Zkt4EHlWVLslUt9GS/i/9Ho9Iqlh3JJ2avuvlki6ro8XzC+C2iLg6It6JTFlEfKvOX6CB0rpws6Qlkhal77d8J6aTpF+mdeF1SedUWebvSHopLe8/JX0vV+5QSQsl/UjS26n87+TGH5X+Ax+k+V7QxOX4f6kuK5T1Dnw+Ny4knS1pLjA3V7cLc3Ubnur0alofL8nl31fS9PRfWiLpN5I2q1L+D9N38I6kX6hKgFPN25vHJF0p6f+Aj4EvSNo/rcsr0/v+uen7qYbtU36dlHQl2fr8m/Rf+U2a5ouSJqdlfEVSjetUE3/fHsr+2+9LehrYtR4/4x+B03KfTwVur1KnGrft1dS/3svarCKiUS9gPvCVNLwFcBtwe278rcDP0/BVwOMNLP8EYCey4H0i8BHQK407HVgDfBfoRLYnvhhQGj8d+BWwOdlG/APgjhrmMxRYWMO4/DJUTAdcBjwD9EyfrwXuA7YDtgb+F/iv6sqv8r2NAlYBR6Xl+C/gqTSuMzAPuATYDDg0LcfuafztwL1pfn2BV4Ez0rjvk7Wkdk51mgoEsGkafzFwfy3f/VTgvDT8G+D/AVdWSbslDT8O/BboQhY4lwGH5ZZvDTA8/Y5dgcfScu1K1nPwYqr7V4BN03L9oYZ6bQ58AuydPs8GvgD8X5W0U3PzvyMN903fwe3Alqku5Wnl38tjwGvAbrm6XpXGDQA+BL6cfo8xadm+Uk09twDWAYc04f91K2ndq2Zc1Xr/BbghLdcOwNPA93LrwotAH6A78LcqeY9Ov4XIWt4fA4Nz6+5a4Gdk6+NRaXz3NH4JcGAa7l6er5r6ng78vY5tyPC0XuyR1oNLgSdz0wYwmWx97pqr209T3b5Ltu7dSfafKCL7b30h5S8BhqSy+wIvASOrlD81lb8L2Tr5b/Xc3jwGvJnmuSmwI7AC+Nf0+aT0uUdd26dqftvHyuuRPm8JLAC+k8oeDLwDFNXw3Tfl9x0PTEjz3BNYVN3vWKXefVP9OqXf8hWy//b8Bmzb/96YZW3OV1MD44fAe+nLXQzsVd0fG7gJGN+kisIsshZn+Zc3r8qGKIDPka3Ua4Etc+PvpHkC46K0Qv8d2CalK/2wu+by7Qe8Xl35bBgY/5YbNwD4JA0fCLwFbJIbf1fK0wn4FBiQG/c94LE0/Cjw/dy4I8j92erxXY8CJqXh54D+wLAqaaeRBd51wNa5vP8F3JorZ1qVsh8D/iP3+ZfAQ7nPXwdm1VK3x4BzyTZg5TsqV+XS1gOfz82/6gbnC7myytPyG6FLc+P/HXg4Df8UuKvKOrea6gNj71TuF5uwvt9KtmF/L/e6rWq9yTbCnwJdc3lPAqbm1oXv5cZ9pbZ1gSzInptbdz/JTwu8DQxJw2+m9a5bHctyOtl/8r0qr/V89l94iLRjlz5vQraRLv8tAzi0yv/2E6BT+rx1muZLuWnKgOE11GkkaX3OlT+sym8/JVf/arc3ufXmZ7nx/wo8XWV+01M5tW6fqDswngg8UaXsG4DL67le1ev3JdvGrCG3DgP/Sd2BcVOyna8jyf6X/0GVwFhN3llU3raXB8YmLWtTXk3tSh0eEduS7fmcAzyu7PhTVcvJjkPWW+q2mpW6Pt4j22PJd4lWHEOKiI/T4FZkeyIrIuKj3LRvNGTetdgWOJOsNbgypfUk+6OU5er6cEqvj/yxsI+BLqmbaydgQUSsz41/g2yjuz1Zq+WNasZRnrfKuIaYBnxZUneyVvFc4Elg/5S2Z5pmJ+DdiPighnpQpR7lluaGP6nm81Z11O0gsh2H8hM6/p5LWxARtS1vdfXJq/p7lNel0nea1rnlNZSxgmyj36B1vhpjImLb3Ou0aqb5PNne/pLc+ncDWctxg3pXGUbSVyU9lbqq3iNrNeT/Z8sjYm3uc/47+Waa/o3UNbhfLcvyVJVl2ZYssOaXY2xuGd4l2+msbV1aHp+dOPVJeq92XZK0m6T7Jb2lrFv/P9nwEEvV/8xOuc81bW+qy7sTG/7nyv8XTd0+fR74Uvn3lL6rk8kaBRtowu/bkyzINWY7cjtZgDsJ2OCs+3ps28s1aFmbU3MdY1wXEfeQtR6+XM0kfwP2ldSnPuWlYws3kQXbHulPNJvsj1KXJUB3SVvm0napz3zrYQXwNeAPkg5Iae+Q/QGLcn/6bSI7KakpFgM7VznOsQtZq/Udsr25z1czDrLvYOcq4xpiOlk355lk3ZRExPupTmcCiyPi9fR5O0lb11APyPYim9M0sgB4EPBESvs/4ICUVtcJQY2tzxKy7kgAlB1D7lHtDLIN53SywNHSFpC1GLfPrX/dIqIoja9Ub3LrhaTNgT+TdQvvmP5nD1K//xmRnS9wDFkQ/gtZt1tTluN7VYJn14h4Mj/LJpR/Pdnhhf4R0Y3sEEXV5az6n1ncgPLzdVtM5f9meXmLaPj2qeoyLyA7LJX/nraKiLOqZmzi77uMrGXbmO3In8m6cP9ZdSe1gdv2ei9rc2uWwKjMMWTHGV6qOj4i/kZ2fGCSpJJ0YHlrSd+X9P+qKXJLshViWSr/O2R7FXVKP8RM4Apll4h8max7rq5l6FLlVe3KE9mlFienZflSatHdBFwjaYdUVm9JR9anvrX4B1kX7YWSOiu79vHrZF3S68g2Qlem7/HzwPl8tnc2AfihshODupMdU6y3iPiE7Ds8n8+CD2Qts/NJwSciFpC1JP8rfWcDgTOAcY1Y3vp6kqzlfkp53SJiBdm6cgoNPFO2ASYCX1d2UsVmwBXUvoG5kOxkox9L6gEgqVjS+OasVEQsAR4Bfimpm6RNJO0q6eA0yQTg3LRObgtclMu+GVlvzzJgrbITSo6oz3zTf+tkSdtExBrgfbId48b6HfATSUWp/G0kndCE8qramqyOH0r6Itlxwqp+LKm7ssvOzgX+1Mh5PQjsJunbaVt3ItlhkvsbsX1aSnYcvdz9qex/TduFzpL2kbRHNXkb/fumbcw9wChJW0gaQOWTamrL+xHZORHVXWbSkG17Q5a1WTU1MP6vpA/JVrgrgdMiYk4N0x5PtsL8CVhJtpdQStaarCQiXiQ79jSdbMXYi9RyqadvA18i6465nCpnRVWjN1mrL/+q8QysiJhMdkD4PkklZBubecBTqZvmb8DuDahvdfNYTXZ5y1fJWoi/JTup5OU0yQ/IAuc/yQLWncAtadxNwF/JjgU+Q7aCV1B24ftDdVThcbKWQP76sydSWj74nER2fGExMIms/39yfZezoVJrrIzsDz+7jro153znkH3n48n2+j8gOx7zaQ3TP0m2cTgU+Kekd4Ebyf4D+Qu5a9sLv1CVr2N8p4bpTiXbCL5I1qsxkc+6cW8iC5zPA8+m+a8F1qUu8B+SBc8VZP+b++r6LnL+FZif1vnvk+2YNEpETAKuBsan8mbTsEtb6nIB2fJ9QPadVBf07iVbt2YBDwDVnjFfl4hYTtaz9COy7vYLga9FRPnv15Dt01jgeGVnw16XfrMjgBFk/7m3yL63zaupR1N/33PIulXfIjvm/Yf6ZoyImRHxWjXp9d62N2RZm1v5WVVm1gCStiI7gaR/6lZuF1Kr4XcRUbWrr0OTFGS/5bxC18UKryNfiGrWIJK+nrqVtiQ7bvMC2VnGbZakrsqu7dtUUm+yFsqkQtfLrC1zYDSrv2PIunQWk13CMiLafpeLyI6HriDrSn2J7NITM6uBu1LNzMxy3GI0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPL2bTQFWht22+/ffTt27fQ1TAza1fKysreiYieha5Ha+hwgbFv377MnDmz0NUwM2tXJL1R6Dq0lhbrSpV0i6S3Jc3OpW0nabKkuem9e27cTyTNk/SKpCNz6SWSXkjjrpOklL65pD+l9H9I6ttSy2JmZh1HSx5jvBUYViXtYmBKRPQHpqTPSBoAjACKUp7fSuqU8lwPnAn0T6/yMs8AVkTEvwDXAFe32JKYmVmH0WKBMSKmAe9WST4GuC0N3wYMz6WPj4hPI+J1YB6wr6ReQLeImB4RAdxeJU95WROBw8pbk2ZmZo3V2mel7hgRSwDS+w4pvTewIDfdwpTWOw1XTa+UJyLWAiuBHi1WczMz6xDayuUa1bX0opb02vJsWLh0pqSZkmYuW7askVU0azlvvb+UkisPYOn7bxe6KmYdXmsHxqWpe5T0Xr4VWAjsnJuuD7A4pfepJr1SHkmbAtuwYdctABFxY0SURkRpz54d4mxja2fGPDKWN99dwJjJYwtdFbMOr7UD433AaWn4NODeXPqIdKZpP7KTbJ5O3a0fSBqSjh+eWiVPeVnHA4+m45Bm7cpb7y/lrqcnsD6CO5+e4FajWYG15OUadwHTgd0lLZR0BnAVcLikucDh6TMRMQeYALwIPAycHRHrUlFnAb8nOyHnNeChlH4z0EPSPOB80hmuZu3NmEfGsj7t061fv96tRrMCU0drZJWWloYv8Le24q33l1Ly8wNYtfbTirQunbvwzH/8Hzt226GWnGatS1JZRJQWuh6toa2cfGPWIeVbi+XcajQrLAdGswJ6eM5kVq9bXSlt9brVPDT7kQLVyMw63L1SzdqS2ZfPKHQVzKwKtxjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyHBjNzMxyChIYJZ0naY6k2ZLuktRF0naSJkuam96756b/iaR5kl6RdGQuvUTSC2ncdZJUiOUxM7ONR6sHRkm9gR8CpRGxJ9AJGAFcDEyJiP7AlPQZSQPS+CJgGPBbSZ1ScdcDZwL902tYKy6KmZlthArVlbop0FXSpsAWwGLgGOC2NP42YHgaPgYYHxGfRsTrwDxgX0m9gG4RMT0iArg9l8fMzKxRWj0wRsQiYAzwJrAEWBkRjwA7RsSSNM0SYIeUpTewIFfEwpTWOw1XTd+ApDMlzZQ0c9myZc25OGZmtpEpRFdqd7JWYD9gJ2BLSafUlqWatKglfcPEiBsjojQiSnv27NnQKpuZWQdSiK7UrwCvR8SyiFgD3APsDyxN3aOk97fT9AuBnXP5+5B1vS5Mw1XTzczMGq0QgfFNYIikLdJZpIcBLwH3AaelaU4D7k3D9wEjJG0uqR/ZSTZPp+7WDyQNSeWcmstj1m5MLJtE8eghbH/+LhSPHsLEskmFrpJZh7Zpa88wIv4haSLwDLAWeBa4EdgKmCDpDLLgeUKafo6kCcCLafqzI2JdKu4s4FagK/BQepm1GxPLJjFywkV8suYTABauWMTICRcBcHzJsYWsmlmHpeyEzo6jtLQ0Zs6cWehqmAFQPHoIC1cs2iC9T/fePHfZUwWokVn1JJVFRGmh69EafOcbswJatKL6w+I1pZtZy3NgNCug3t13alC6mbU8B0azArrsqIvo2rlrpbSunbty2VEXFahGZtbqJ9+Y2WfKT7AZ/eDVLFqxmN7dd+Kyoy7yiTdmBeSTb8zMrE4++cbMzKyDcmA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPLKUhglLStpImSXpb0kqT9JG0nabKkuem9e276n0iaJ+kVSUfm0kskvZDGXSdJhVgeMzPbeBSqxTgWeDgivggUAy8BFwNTIqI/MCV9RtIAYARQBAwDfiupUyrneuBMoH96DWvNhTAzs41PqwdGSd2Ag4CbASJidUS8BxwD3JYmuw0YnoaPAcZHxKcR8TowD9hXUi+gW0RMj4gAbs/lMTMza5RCtBi/ACwD/iDpWUm/l7QlsGNELAFI7zuk6XsDC3L5F6a03mm4avoGJJ0paaakmcuWLWvepTEzs41KIQLjpsBg4PqI2Bv4iNRtWoPqjhtGLekbJkbcGBGlEVHas2fPhtbXzMw6kEIExoXAwoj4R/o8kSxQLk3do6T3t3PT75zL3wdYnNL7VJNuZmbWaK0eGCPiLWCBpN1T0mHAi8B9wGkp7TTg3jR8HzBC0uaS+pGdZPN06m79QNKQdDbqqbk8ZmZmjVKos1J/AIyT9DwwCPhP4CrgcElzgcPTZyJiDjCBLHg+DJwdEetSOWcBvyc7Iec14KFWXAYzs3bjrfeXUnLlASx9/+26J+7gNi3ETCNiFlBazajDapj+SuDKatJnAns2a+XMzDZCYx4Zy5vvLmDM5LH84psbbE4tx3e+MTPbyL31/lLuenoC6yO48+kJbjXWwYHRzGwjN+aRsayP7KT99evXM2by2ALXqG1zYDQz24iVtxZXr1sNwOp1q91qrIMDo5nZRizfWiznVmPtHBjNzDZiD8+ZXNFaLLd63Woemv1IgWrU9hXkrFQzM2sdsy+fUegqtDtuMZqZmeXUGhgl/YukA6pJP1DSri1XLTMzs8Koq8V4LfBBNemfpHFmZmYblboCY9+IeL5qYrrjTN8WqZGZmVkB1RUYu9QyrmtzVsTMzKwtqCswzpD03aqJks4AylqmSmZmZoVT1+UaI4FJkk7ms0BYCmwGHNuC9TIzMyuIWgNjRCwF9pd0CJ89xeKBiHi0xWtmZmZWALUGRknbpcHn0qtSekS823JVMzMza311daWWAQEovZOGSZ+/0EL1MjMzK4i6ulL7tVZFzMzM2oI6bwknaVNJSsM7Szpe0qAWr5mZmVkB1HVLuO8CbwNvpOEpwPHAnyRd1Ar1MzMza1X1uVxjV2Br4CXg8xHxjqQtgBnA1S1bPTMzs9ZVV1fq6ohYERFvAvMi4h2AiPgYWF17VjMzawsmlk2iePQQtj9/F4pHD2Fi2aRCV6lNq6vF2FXS3mQBdLM0rPSq7XZxZmbWBkwsm8TICRfxyZpPAFi4YhEjJ2RHwo4v8X1aqqOIqHmkNLWa5F7AEoCIOKSF6tViSktLY+bMmYWuhplZqygePYSFKxZtkN6ne2+eu+ypepcjqSwiSpuzbm1VXZdrbBD4JD3bHgOimVlHtGjF4galWz0u1zAzs/ard/edGpRujQuMNzV7LczMrEVcdtRFdO1c+SmBXTt35bKjfMVdTeo6+WYDEfHblqiImZk1v/ITbEY/eDWLViymd/eduOyoi3ziTS1qPfmmRWcsdQJmAosi4mvpxuR/AvoC84FvRcSKNO1PgDOAdcAPI+KvKb0EuJXsockPAudGHQvkk2/MzBquI518U8hjjOeS3TSg3MXAlIjoT3aHnYsBJA0ARgBFwDDgtymoAlwPnAn0T69hrVN1MzPbWBUkMErqAxwN/D6XfAxwWxq+DRieSx8fEZ9GxOvAPGBfSb2AbhExPbUSb8/lMTMza5RCtRivBS4E1ufSdoyI8usjlwA7pPTewILcdAtTWu80XDV9A5LOlDRT0sxly5Y1ywKYmdnGqdUDo6SvAW9HRFl9s1STFrWkb5gYcWNElEZEac+ePes5WzMz64gafFZqMzgA+Iako8huK9dN0h3AUkm9ImJJ6iZ9O02/ENg5l78PsDil96km3czMrNFavcUYET+JiD4R0ZfspJpHI+IU4D7gtDTZacC9afg+YISkzSX1IzvJ5unU3fqBpCHpeZGn5vKYmZk1SiFajDW5Cpgg6QzgTeAEgIiYI2kC8CKwFjg7ItalPGfx2eUaD6WXmZlZoxXsOsZC8XWMZmYN5+sYzczMOigHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzsxwHRjMzs5xWD4ySdpY0VdJLkuZIOjelbydpsqS56b17Ls9PJM2T9IqkI3PpJZJeSOOuk6TWXh4zM9u4FKLFuBb4UUTsAQwBzpY0ALgYmBIR/YEp6TNp3AigCBgG/FZSp1TW9cCZQP/0GtaaC2JmZhufVg+MEbEkIp5Jwx8ALwG9gWOA29JktwHD0/AxwPiI+DQiXgfmAftK6gV0i4jpERHA7bk8ZmZmjVLQY4yS+gJ7A/8AdoyIJZAFT2CHNFlvYEEu28KU1jsNV02vbj5nSpopaeayZcuadRnMzGzjUrDAKGkr4M/AyIh4v7ZJq0mLWtI3TIy4MSJKI6K0Z8+eDa+smZl1GAUJjJI6kwXFcRFxT0pemrpHSe9vp/SFwM657H2AxSm9TzXpZmZmjVaIs1IF3Ay8FBG/yo26DzgtDZ8G3JtLHyFpc0n9yE6yeTp1t34gaUgq89RcHjMzs0bZtADzPAD4V+AFSbNS2iXAVcAESWcAbwInAETEHEkTgBfJzmg9OyLWpXxnAbcCXYGH0svMzKzRlJ3Q2XGUlpbGzJkzC10NM7N2RVJZRJQWuh6twXe+MTMzy3FgNDMzyynEMUYzA9asWcPChQtZtWpVoatiVqFLly706dOHzp07F7oqBePAaFYgCxcuZOutt6Zv3774Nr/WFkQEy5cvZ+HChfTr16/Q1SkYd6WaFciqVavo0aOHg6K1GZLo0aNHh+/FcGA0KyAHRWtrvE46MJp1aG+99RYjRoxg1113ZcCAARx11FG8+uqrteYZOnQou+++O4MGDWLQoEEcf/zxAIwaNYoxY8a0RrXr5bHHHuNrX/tapbTTTz+diRMnFqhGlQ0dOpTmvHSsb9++vPPOOwDsv//+AMyfP58777yz2ebRUfgYo1kHFREce+yxnHbaaYwfPx6AWbNmsXTpUnbbbbda844bN47S0o3/kra1a9ey6aaF20w2dv5PPvkk8Flg/Pa3v93cVduoucVo1k5MLJtE8eghbH/+LhSPHsLEsklNKm/q1Kl07tyZ73//+xVpgwYN4sADD2xqVXnttdcYNmwYJSUlHHjggbz88ssV6UOGDGGfffbhpz/9KVtttRUAH374IYcddhiDBw9mr7324t57s7s7zp8/nz322IPvfve7FBUVccQRR/DJJ58AcN111zFgwAAGDhzIiBEjGlzHsrIyDj74YEpKSjjyyCNZsmQJkLXkLrnkEg4++GDGjh3L0KFDOe+88zjooIPYY489mDFjBscddxz9+/fn0ksvrShv+PDhlJSUUFRUxI033liRvtVWW/GjH/2IwYMHc9hhh5F/ws/dd9/Nvvvuy2677cYTTzwBwK233soJJ5zA17/+dY444gjeffddhg8fzsCBAxkyZAjPP/88AMuXL+eII45g77335nvf+x75m7WUf68XX3wxTzzxBIMGDeKaa65h3bp1/PjHP2afffZh4MCB3HDDDQ3+3jqEiOhQr5KSkjBrC1588cV6T3v3zHui94X9Y7vz+lS8el/YP+6eeU+j5z927NgYOXJkg/MdfPDBsdtuu0VxcXEUFxfHBRdcEBERl19+efziF7+IiIhDDz00Xn311YiIeOqpp+KQQw6JiIijjz467rzzzoiIuP7662PLLbeMiIg1a9bEypUrIyJi2bJlseuuu8b69evj9ddfj06dOsWzzz4bEREnnHBC/PGPf4yIiF69esWqVasiImLFihUb1HPq1KnRrVu3inoWFxdH9+7d4+67747Vq1fHfvvtF2+//XZERIwfPz6+853vVCzfWWedVWl5L7zwwoiIuPbaa6NXr16xePHiWLVqVfTu3TveeeediIhYvnx5RER8/PHHUVRUVJEOxB133BEREVdccUWcffbZFeWef/75ERHxwAMPxGGHHRYREX/4wx+id+/eFeWdc845MWrUqIiImDJlShQXF0dExA9+8IO44oorIiLi/vvvDyCWLVsWEVHxvU6dOjWOPvroimW54YYbYvTo0RERsWrVqigpKYl//vOfG3x31a2bwMxoA9vw1ni5K9WsHRj94NV8suaTSmmfrPmE0Q9ezfElx7Z6fWrrSv3www958sknOeGEEyrSPv30UwCmT5/OX/7yFwC+/e1vc8EFFwDZDvoll1zCtGnT2GSTTVi0aBFLly4FoF+/fgwaNAiAkpIS5s+fD8DAgQM5+eSTGT58OMOHD6+2LgceeCD3339/xefTTz8dgFdeeYXZs2dz+OGHA7Bu3Tp69epVMd2JJ55YqZxvfOMbAOy1114UFRVVTPuFL3yBBQsW0KNHD6677jomTcpa8QsWLGDu3Ln06NGDTTbZpKK8U045heOOO66i3PLh/HIBHH744Wy33XYA/P3vf+fPf/4zAIceeijLly9n5cqVTJs2jXvuyR5OdPTRR9O9e/dqv4O8Rx55hOeff77iOOvKlSuZO3duh740ozoOjGbtwKIV1T9Rrab0+igqKmqRE1HWr1/Ptttuy6xZs+qdZ9y4cSxbtoyysjI6d+5M3759Ky4Z2HzzzSum69SpU0VX6gMPPMC0adO47777GD16NHPmzKn38biIoKioiOnTp1c7fsstt6z0ubwOm2yySaX6bLLJJqxdu5bHHnuMv/3tb0yfPp0tttiCoUOH1njJQ/6sz/KyOnXqxNq1a6udf1RzP+vyMhp6BmlE8Otf/5ojjzyyQfk6Gh9jNGsHenffqUHp9XHooYfy6aefctNNN1WkzZgxg8cff7zRZQJ069aNfv36cffddwPZxvi5554DYMiQIRWtn/ITfiBrueywww507tyZqVOn8sYbb9Q6j/Xr17NgwQIOOeQQ/vu//5v33nuPDz/8sN513H333Vm2bFlFYFyzZg1z5sxp0HLmrVy5ku7du7PFFlvw8ssv89RTT1Wqa/kOyJ133smXv/zlBpV90EEHMW7cOCA703b77benW7duldIfeughVqxYsUHerbfemg8++KDi85FHHsn111/PmjVrAHj11Vf56KOPGrawHYADo1k7cNlRF9G1c9dKaV07d+Wyoy5qdJmSmDRpEpMnT2bXXXelqKiIUaNGsdNOWbAt776szsknn1xxucZXvvKVDcaPGzeOm2++meLiYoqKiipOprn22mv51a9+xb777suSJUvYZpttKsqbOXMmpaWljBs3ji9+8Yu11n3dunWccsop7LXXXuy9996cd955bLvttvVe9s0224yJEydy0UUXUVxczKBBgyrO5GyMYcOGsXbtWgYOHMhll13GkCFDKsZtueWWzJkzh5KSEh599FF++tOfNqjsUaNGMXPmTAYOHMjFF1/MbbfdBsDll1/OtGnTGDx4MI888gi77LLLBnkHDhzIpptuSnFxMddccw3/9m//xoABAxg8eDB77rkn3/ve9yq1VC3jx06ZFchLL73EHnvsUe/pJ5ZNYvSDV7NoxWJ6d9+Jy466qCDHF5vi448/pmvXrkhi/Pjx3HXXXRVBc2O11VZbNag12xZUt252pMdO+RijWTtxfMmx7S4QVlVWVsY555xDRLDttttyyy23FLpKZhtwYDSzVnPggQdWHG/sKNpba9F8jNHMzKwSB0YzM7McB0YzM7McB0YzM7McB0azDmxjf+yUJG6++eaKtGeffRZJddazrS2LtS6flWrWQUVs/I+d2muvvfjTn/7EGWecAWR32ykuLi5wraytc4vRrB156/2llFx5AEvff7vJZXWEx07tsssurFq1iqVLlxIRPPzww3z1q1+tGH/TTTexzz77UFxczDe/+U0+/vjjei+LbbwcGM3akTGPjOXNdxcwZvLYJpc1e/ZsSkpKGpU3f0u4H//4xxuMP/PMM/n1r39NWVkZY8aM4d///d8BOPfcczn33HOZMWNGxa3nALp06cKkSZN45plnmDp1Kj/60Y8qbp49d+5czj77bObMmcO2225bca/Vq666imeffZbnn3+e3/3udzXW9fjjj+fuu+/mySefZPDgwZVuAn7ccccxY8YMnnvuOfbYY49K3a51LYttvNyVatZOvPX+Uu56egLrI7jz6QlccPi57Nhth4LUpb08dgrgW9/6FieeeCIvv/wyJ510UqV7os6ePZtLL7204ibkVZ86Uduy2MbLLUazdmLMI2NZn1pR69evb3KrsaioiLKysuaoWiX5x06Vv1566aVa8+QfOzVr1ix23HHHGh87VX7T6wceeICzzz6bsrIySkpKarwZ9uc+9zk6d+7M5MmTOeywwyqNO/300/nNb37DCy+8wOWXX77Bo6IasyzW/rX7wChpmKRXJM2TdHGh62PWEspbi6vXrQZg9brV3Pn0hCYda+xIj5362c9+xtVXX02nTp0qpX/wwQf06tWLNWvWVDzCqb7LYhuvdh0YJXUC/gf4KjAAOEnSgOaez8SySex5xT70OH9n9rxiHyaWTWruWZjVKt9aLNfUVmNHeuzU/vvvX2136+jRo/nSl77E4YcfXuM8a1oW23i168dOSdoPGBURR6bPPwGIiP+qKU9DHzs1sWwSIydcxCdrPqlI69q5K9d+6+p2/6QDK6yGPHZqzyv2YcnKtzZI77XN55h9+YzmrlqL6YiPnWqP/Nip9q03sCD3eSHwpeacwegHr64UFAE+WfMJox90YLTW056CX2382ClrD9p7YFQ1aRs0gSWdCZwJVPuU69osWrG4QelmVrOO+Ngpa3/a9TFGshbizrnPfYANIlZE3BgRpRFR2rNnzwbN4HPb7NigdDMza9/ae2CcAfSX1E/SZsAI4L7mnEH/HXatNn23Hf6lOWdjHVR7PsZvGyevk+08MEbEWuAc4K/AS8CEiJjTnPOY+/Zr1aa/+va85pyNdUBdunRh+fLl3hBZmxERLF++nC5duhS6KgXVrs9KbYyGnpVq1lLWrFnDwoULN7io3KyQunTpQp8+fejcuXOldJ+VamYtrnPnzvTr16/Q1TCzKtp1V6qZmVlzc2A0MzPLcWA0MzPL6XAn30haBtR+h+KabQ+804zVMcvz+mUtrSnr2OcjomEXgrdTHS4wNoWkmR3lrCxrfV6/rKV5Hasfd6WamZnlODCamZnlODA2zI2FroBt1Lx+WUvzOlYPPsZoZmaW4xajmZlZTpsPjJL6SvpE0qz0eX56Hyrp/haa55ONzPeYpGY/4yt9B4+l4QMlvShpdnPPx8zM2kFgTF6LiEGtNbOI2L+15lVOUqf6TBcRTwBHtXB1OqyadsQaUc4oSRek4VslHV/LtCMlbdGY+dRRh+0kTZY0N713T+lDJd2ahk+UNK+ldjILrervmdLmp/cOuXNdx3TzG1n+RrW+t5fAmLcsN7yVpImSXpY0TpIAJJVIelxSmaS/SuqV0h+TdI2kaZJekrSPpHvSF/nz8kIlfZjee6VpZ0maLenA8vGSfinpGUlTJOUvej1B0tOSXs1N30nSLyTNkPS8pO+l9KGSpkq6E3ihpumAdcC7LfR92oZadUcMGAk0+4YCuBiYEhH9gSnpcyUR8Sfg31pg3m1Ja/+ebXrnug0YSRtf39tdYIyIfXIf9yb7kgcAXwAOkNQZ+DVwfESUALcAV+byrI6Ig4DfAfcCZwN7AqdL6lFldt8G/pr+VMXArJS+JfBMRAwGHgcuz+XZNCL2TfUqTz8DWJnqvg/wXUnlj1XYF/iPiBhQ03QRsSAijqv/t2TNqGJHTNKFkl6Q9Jykq1LarpIeTjthT0j6YkMKl/RDYCdgatpJOkPSNbnx35X0q7TH/7Kk29JO08Tyve6adgSBY4Db0vBtwPA0vBpY2YjvYmPhnet6fDcden2PiDb9AvoCs6tJHwpMzn2+HjiFLMi9TxbEZgEvAI+kaR4DDkjDh1bJPw0YlIY/TO8HAfOAUeXjUvo6sgAIWUCeVU35OwLz0vBE4NVcnV4HjkjLMDVXbrXT1fc78atF17evAk8CW6TP26X3KUD/NPwl4NE0PAq4IA3fSrajVtM85wPbp+EtgdeAzunzk8BeqV6RW79uAS4AOqdpeqb0E4Fb0vB7Veazoob5DwXuL/R335q/Z265VwJ9yBoJ04Ev1/GdPgZcnYbPBRYDvYDNgYVAjzSufBvyI7IdX4BOwNZpOICT0/BPgd/kyv9lGj4K+FsaPhO4NA1vDswE+qVl+AjoV9t0jfjeOvT63t6fx/hpbngd2fMlBcyJiP3qyLO+Sv71VHk+ZURMk3QQcDTwR0m/iIjbqykzf81LeZnl9SHV6QcR8dd8JklDyVZqapvO2oSvAH+IiI8BIuJdSVsB+wN3p4YGZBujRouIjyQ9CnxN0ktkG4wXJPUFFkTE/6VJ7wB+CDxMtjM4OdWhE7CkKXXoYJ6OiIUAyo5D9gXeo/bv9L70/gLZtmZJyv9PYGdgeW7aGcAtqSfrLxExK6WvB/6Uhu8A7snlKR8uS/WBbEd6oD47drcN0J+sNfR0RLxex3Tl4+urQ6/v7T0wVucVoKek/SJielohd4uIOQ0tSNLngUURcZOkLYHBwO1ke5fHA+PJulv/XkdRfwXOkvRoRKyRtBuwqL7TRcRH1UxrrUtU3gGCbD14L5r/+NXvgUuAl4E/5NKrzj+ofUdwqaReEbEkdTe93cz13Bh457p6HXp9b3fHGOsSEavJgtbVkp4j65Js7IHwocAsSc8C3wTGpvSPgCJJZWRdsj+ro5zfAy8Czyi7zOIGqt8pqe901voeAf5f7jjHdhHxPvC6pBNSmiQVN6LsD4Ctyz9ExD/IWh7fBu7KTbeLpPINwklkO2QVO4KpDp0lFaVp7gNOS8OnkR1Tt7rV9p02SNq5fjsibgJuJtu5hs92rqFhO9edU7m7pZ31Bk8nqbekKXXMr0Ov7+12oxsRj5H1x5d/Pic3PIvs+GDVPENryZ8ft1V6v43PDuZWLesy4LJayn+H1A0SEevJ9oguqVJM1TrUNJ0VWEQ8LGkQMFPSauBBst/pZOB6SZeSHf8YDzzXwOJvBB6StCQiDklpE8iOa6/ITfcScJqkG4C5wPURsTp1m10naRuy//S1wBzgKmCCpDOAN4ETGrrcHVEd32lDDQV+LGkN8CFwakrP71yvJDtWVpvfk21PnlHWh7iMz04uaeh0vYC1tc2so6/vbf6WcJJ2JjvYurwFmvCNIunD8uBZgHkfCPyW7PsYWog6bMzSsY37I2LPAtfjfuCaiJjS0vVK3XEXRMTXmrvsQmsrv2dVBd6GnAO8GRH31TlxK2lr63ub70qN7FKFndtKUITPWpQFmvcTEbGXg2KLWQdso9wF4a1J0raSXgU+Kd9ItPD8TiTb0VpR17TtVEF/z7YoIn7TVoJiW13f23yL0WxjIWkS2Sn2eRf5LGTbGLXn9d2B0czMLKfNd6WamZm1JgdGMzOzHAdGsxxJO0q6U9I/ld2LcbqkYwtdrzxl98eMdFp6edreKe2COvKOqmsas47OgdEsSdd9/QWYFhFfiOwm9CPI7qVZddpCXwP8ApWvfRtBw68nM7NqODCafeZQsqev/K48ISLeiIhfA0g6XdLdkv4XeETZ89/+ouzu/09JGpimq9QqU/ZUhb6q/YkBVyl7APXzksbUo65vAl1SC1fAMOCh3Dy/q+wJC89J+rOqef6dmvikBLONlQOj2WeKgGfqmGY/4LSIOBS4Ang2IgaS3RWkuntgVrU7cGPK8z7w75K2A44FilL6z2srIGci2d099k/1zt+3856I2CciisnuIHJGNflvJLuvZgnZkwt+W8/5mm3UHBjNaiDpf1KLa0YueXJElD/X7svAHwEi4lGgR7pNVW2qPjHgy2QBchXwe0nHAR/Xs4oTyALjSVS+xyTAnqkV+ALZbbwq3etTlZ+UMIvsvry9MDMHRrOcOXx2k2ci4mzgMCD/ENmqTzKoKsjuQ5n/b3WpMr7S9BGxluyB1X8mu6/lw/WpbES8BawBDid7Tl7ercA5EbEXWcu2S5XxFU9KyL32qM98zTZ2Doxmn3mU7LjdWbm0DY7N5Uwja42V33/xnfQEgvmkACtpMJXv/rHBEwNS622biHgQGAkMSnmPlfRfddT5p2R3E1lXJX1rYImypyycXDVTMz4pwWyjU+gz68zajIgIScOBayRdSPZkgo+Ai2rIMgr4g6Tnybo/yx9582fg1NRFOQN4NZdngycGkD1M9l5JXchaoeelaXcl62atrc5P1jDqMuAfwBtkZ7BuXc00zfGkBLONjm8JZ9ZKGvrEAEl3AOdFxLIWrZiZVeIWo1kbFRGnFLoOZh2RW4xmZmY5PvnGzMwsx4HRzMwsx4HRzMwsx4HRzMwsx4HRzMwsx4HRzMwsp8Ndx7j99ttH3759C10NM7N2pays7J2I6Fn3lO1fhwuMffv2ZebMmYWuhplZuyLpjULXobW4K9XMzCzHgdHMzCzHgdHMzCzHgdHMzCzHgdHMrAN46/2llFx5AEvff7vQVWnzHBjNzDqAMY+M5c13FzBm8thCV6XNc2A0M9vIvfX+Uu56egLrI7jz6QluNdahIIFR0raSJkp6WdJLkvaTtJ2kyZLmpvfuuel/ImmepFckHZlLL5H0Qhp3nSQVYnnMzNqyMY+MZX169u769evdaqxDoVqMY4GHI+KLQDHwEnAxMCUi+gNT0mckDQBGAEXAMOC3kjqlcq4HzgT6p9ew1lwIM7O2rry1uHrdagBWr1vtVmMdWj0wSuoGHATcDBARqyPiPeAY4LY02W3A8DR8DDA+Ij6NiNeBecC+knoB3SJiekQEcHsuj5mZUbm1WM6txtoVosX4BWAZ8AdJz0r6vaQtgR0jYglAet8hTd8bWJDLvzCl9U7DVdPNzCx5eM7kitZiudXrVvPQ7EcKVKO2rxD3St0UGAz8ICL+IWksqdu0BtUdN4xa0jcsQDqTrMuVXXbZpWG1NTNrx2ZfPqPQVWh3CtFiXAgsjIh/pM8TyQLl0tQ9Snp/Ozf9zrn8fYDFKb1PNekbiIgbI6I0Ikp79uwQN4c3M7NGavXAGBFvAQsk7Z6SDgNeBO4DTktppwH3puH7gBGSNpfUj+wkm6dTd+sHkoaks1FPzeUxMzNrlEI9duoHwDhJmwH/BL5DFqQnSDoDeBM4ASAi5kiaQBY81wJnR8S6VM5ZwK1AV+Ch9DIzM2s0RVR7WG6jVVpaGn4eo5lZw0gqi4jSQtejNfjON2ZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZmZjkOjGZtwFvvL6XkygNY+v7bdU9sZi3KgdGsDRjzyFjefHcBYyaPLXRVzDo8B0azAnvr/aXc9fQE1kdw59MT3Go0KzAHRrMCG/PIWNanB4avX7/erUazAnNgNCug8tbi6nWrAVi9brVbjWYF5sBoVkD51mI5txqtuU0sm0Tx6CFsf/4uFI8ewsSySYWuUpvmwGhWQA/PmVzRWiy3et1qHpr9SIFqZBubiWWTGDnhIhauWEQQLFyxiJETLnJwrIWiyt7qxq60tDRmzpxZ6GqYmbWK4tFDWLhi0Qbpfbr35rnLnqp3OZLKIqK0OevWVhWkxShpvqQXJM2SNDOlbSdpsqS56b17bvqfSJon6RVJR+bSS1I58yRdJ0mFWB4zs7Zq0YrFDUq3wnalHhIRg3J7IBcDUyKiPzAlfUbSAGAEUAQMA34rqVPKcz1wJtA/vYa1Yv3NzNq83t13alC6ta1jjMcAt6Xh24DhufTxEfFpRLwOzAP2ldQL6BYR0yPrD749l8fMzIDLjrqIrp27Vkrr2rkrlx11UYFq1PYVKjAG8IikMklnprQdI2IJQHrfIaX3Bhbk8i5Mab3TcNX0DUg6U9JMSTOXLVvWjIthZta2HV9yLNd+62r6dO+NEH269+bab13N8SXHFrpqbdamBZrvARGxWNIOwGRJL9cybXXHDaOW9A0TI24EboTs5JuGVtbMrD07vuRYB8IGKEiLMSIWp/e3gUnAvsDS1D1Kei+/wnkhsHMuex9gcUrvU026mZlZo7V6YJS0paSty4eBI4DZwH3AaWmy04B70/B9wAhJm0vqR3aSzdOpu/UDSUPS2ain5vKYmZk1SiG6UncEJqUrKzYF7oyIhyXNACZIOgN4EzgBICLmSJoAvAisBc6OiHWprLOAW4GuwEPpZWZm1mi+wN/MzOrkC/zNzMw6KAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGMzOzHAdGswKbWDaJ4tFD2P78XSgePYSJZZMKXSWzDm3TQlfArCObWDaJkRMu4pM1nwCwcMUiRk64CIDjS44tZNXMOiy3GM0KaPSDV1cExXKfrPmE0Q9eXaAamZkDo1kBLVqxuEHpZtbyHBjNCqh3950alG5mLa9ggVFSJ0nPSro/fd5O0mRJc9N799y0P5E0T9Irko7MpZdIeiGNu06SCrEsZo112VEX0bVz10ppXTt35bKjLipQjcyskC3Gc4GXcp8vBqZERH9gSvqMpAHACKAIGAb8VlKnlOd64Eygf3oNa52qmzWP40uO5dpvXU2f7r0Rok/33lz7rat94o1ZARXkrFRJfYCjgSuB81PyMcDQNHwb8BhwUUofHxGfAq9LmgfsK2k+0C0ipqcybweGAw+1ykKYNZPjS451IDRrQwrVYrwWuBBYn0vbMSKWAKT3HVJ6b2BBbrqFKa13Gq6avgFJZ0qaKWnmsmXLmmUBzMxs49TqgVHS14C3I6KsvlmqSYta0jdMjLgxIkojorRnz571nK2ZmXVEhehKPQD4hqSjgC5AN0l3AEsl9YqIJZJ6AW+n6RcCO+fy9wEWp/Q+1aSbmZk1Wqu3GCPiJxHRJyL6kp1U82hEnALcB5yWJjsNuDcN3weMkLS5pH5kJ9k8nbpbP5A0JJ2Nemouj5mZWaO0pVvCXQVMkHQG8CZwAkBEzJE0AXgRWAucHRHrUp6zgFuBrmQn3fjEGzMzaxJFVHtYbqNVWloaM2fOLHQ1zMzaFUllEVFa6Hq0hka3GCWdX9v4iPhVY8s2MzMrlKZ0pW6d3ncH9iE7FgjwdWBaUyplZmZWKI0OjBFxBYCkR4DBEfFB+jwKuLtZamdmZtbKmuOs1F2A1bnPq4G+zVCumZlZq2uOs1L/CDwtqfyx48PJbulmZmbW7jQ5MEbElZIeAg4ku/PMdyLi2SbXzMzMrACa6zrGdWT3PQ0q3//UzMysXWnyMUZJ5wLjgO3Jbvx9h6QfNLVcMzOzQmiOFuMZwJci4iMASVcD04FfN0PZZmZmrao5zkoVWVdquXVU/+QLMzOzNq85Wox/AP6RzkoV2YOFb26Gcs3MzFpdc5yV+itJjwFfTkk+K9XMzNqt5jwrNfBZqWZm1s75rFQzM7Mcn5VqZmaW47NSzczMcpr7rFTI7pXqs1LNzKxdaq6zUh8HDiBrKfqsVDMza7ea66zUWcCS8vIk7RIRbzZT2WZmZq2myYExnYF6ObCUz44vBjCwqWWbmZm1tuZoMZ4L7B4Ry5uhLDMzs4JqjrNSFwAr6zuxpC6Snpb0nKQ5kq5I6dtJmixpbnrvnsvzE0nzJL0i6chceomkF9K46yT5bFgzM2uSRrcYJZ2fBv8JPCbpAeDT8vER8asasn4KHBoRH0rqDPw9Pej4OGBKRFwl6WLgYuAiSQOAEUARsBPwN0m7RcQ64HrgTOAp4EFgGPBQY5fJzMysKS3GrdPrTWAysFkubeuaMkXmw/Sxc3oF2c3Hb0vpt5Fd9kFKHx8Rn0bE68A8YF9JvYBuETE9IgK4PZfHzMysURrdYoyIKxqbV1InoAz4F+B/IuIfknaMiCWp7CWSdkiT9yZrEZZbmNLWpOGq6dXN70yyliW77LJLY6ttZmYdQFO6Uq+NiJGS/pesxVdJRHyjprypG3SQpG2BSZL2rG1W1RVRS3p187sRuBGgtLS02mnMzMygaWel/jG9j2lsARHxXnpk1TBgqaReqbXYC3g7TbYQ2DmXrQ+wOKX3qSbdzMys0Rp9jDEiytL749W9asonqWdqKSKpK/AV4GXgPuC0NNlpwL1p+D5ghKTNJfUD+gNPp27XDyQNSWejnprLY2Zm1ihN6Up9geq7LkV2jk1NF/j3Am5Lxxk3ASZExP2SpgMTJJ1BdkLPCWQFzZE0AXgRWAucnbpiAc4CbgW6kp2N6jNSzcysSZSd0NmIjNLnaxsfEW80quAWVlpaGjNnzix0NczM2hVJZRFRWuh6tIamnJVaEfhSkOwfEX9L3aPNdQ9WMzOzVtXkO99I+i4wEbghJfUB/tLUcs3MzAqhOW4JdzbZI6feB4iIucAOteYwMzNro5ojMH4aEavLP0jalBquJzQzM2vrmiMwPi7pEqCrpMOBu4H/bYZyzczMWl1zBMaLgWXAC8D3gAcj4j+aoVwzM7NW1xxnj46KiJ8CN0F2H1RJ4yLi5GYo28zMrFU1R4txF0k/AZC0GXAPMLcZyjUzM2t1zREYvwPslYLj/cBjETGqGco1MzNrdU25Jdzg3MexZNcx/h/ZyTiDI+KZplbOzMystTXlGOMvq3xeAQxI6QEc2oSyzczMCqIpt4Q7pDkrYmZm1hY0pSv1lIi4Q9L51Y2PiF81vlpmZmaF0ZSu1C3T+9bVjPOdb8zMrF1qSlfqDen9iqrjJI1sQp3MzMwKpjku16hOtd2rZmZmbV1LBUa1ULlmZmYtqqUCo48xmplZu9SUs1I/oPoAKKBro2tkZmZWQE05+aa6s1HNzMzatZbqSq2RpJ0lTZX0kqQ5ks5N6dtJmixpbnrvnsvzE0nzJL0i6chceomkF9K46yT52KaZmTVJqwdGYC3wo4jYAxgCnC1pANlzHadERH9gSvpMGjcCKAKGAb+V1CmVdT1wJtA/vYa15oKYmdnGp9UDY0QsKb/BeER8ALwE9AaOAW5Lk90GDE/DxwDjI+LTiHgdmAfsK6kX0C0ipkdEALfn8piZmTVKIVqMFST1BfYG/gHsGBFLIAuewA5pst7Agly2hSmtdxquml7dfM6UNFPSzGXLljXrMpiZ2calYIFR0lbAn4GREfF+bZNWkxa1pG+YGHFjRJRGRGnPnj0bXlkzM+swChIYJXUmC4rjIuKelLw0dY+S3t9O6QuBnXPZ+wCLU3qfatLNzMwarRBnpQq4GXipyhM47gNOS8OnAffm0kdI2lxSP7KTbJ5O3a0fSBqSyjw1l8fMzKxRmvJ0jcY6APhX4AVJs1LaJcBVwARJZwBvAicARMQcSROAF8nOaD07ItalfGcBt5LdUOCh9DIzM2s0ZSd0dhylpaUxc+bMQlfDzKxdkVQWEaWFrkdrKOhZqWZmZm2NA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVmOA6OZmVlOQQKjpFskvS1pdi5tO0mTJc1N791z434iaZ6kVyQdmUsvkfRCGnedJLX2spiZ2calUC3GW4FhVdIuBqZERH9gSvqMpAHACKAo5fmtpE4pz/XAmUD/9KpappmZWYMUJDBGxDTg3SrJxwC3peHbgOG59PER8WlEvA7MA/aV1AvoFhHTIyKA23N5zMzMGqUtHWPcMSKWAKT3HVJ6b2BBbrqFKa13Gq6abmZm1mibFroC9VDdccOoJX3DAqQzybpc2WWXXZqvZmZNsGbNGhYuXMiqVasKXRWzCl26dKFPnz507ty50FUpmLYUGJdK6hURS1I36dspfSGwc266PsDilN6nmvQNRMSNwI0ApaWl1QZPs9a2cOFCtt56a/r27YvPG7O2ICJYvnw5CxcupF+/foWuTsG0pa7U+4DT0vBpwL259BGSNpfUj+wkm6dTd+sHkoaks1FPzeUxa/NWrVpFjx49HBStzZBEjx49OnwvRkFajJLuAoYC20taCFwOXAVMkHQG8CZwAkBEzJE0AXgRWAucHRHrUlFnkZ3h2hV4KL3M2g0HRWtrvE4WKDBGxEk1jDqshumvBK6sJn0msGczVs3MzDq4ttSVamat7K233mLEiBHsuuuuDBgwgKOOOopXX3211jxDhw5l9913Z9CgQQwaNIjjjz8egFGjRjFmzJjWqHa9PPbYY3zta1+rlHb66aczceLEAtWosqFDhzJz5sxmK69v37688847AOy///4AzJ8/nzvvvLPZ5tFRODCatRMTyyZRPHoI25+/C8WjhzCxbFKTyosIjj32WIYOHcprr73Giy++yH/+53+ydOnSOvOOGzeOWbNmMWvWrDYTaFrC2rVr2+X8n3zyScCBsbEcGM3agYllkxg54SIWrlhEECxcsYiREy5qUnCcOnUqnTt35vvf/35F2qBBgzjwwAObXN/XXnuNYcOGUVJSwoEHHsjLL79ckT5kyBD22WcffvrTn7LVVlsB8OGHH3LYYYcxePBg9tprL+69NzuPbv78+eyxxx5897vfpaioiCOOOIJPPvkEgOuuu44BAwYwcOBARowY0eA6lpWVcfDBB1NSUsKRRx7JkiVLgKwld8kll3DwwQczduxYhg4dynnnncdBBx3EHnvswYwZMzjuuOPo378/l156aUV5w4cPp6SkhKKiIm688caK9K222oof/ehHDB48mMMOO4xly5ZVjLv77rvZd9992W233XjiiScAuPXWWznhhBP4+te/zhFHHMG7777L8OHDGThwIEOGDOH5558HYPny5RxxxBHsvffefO973yO7z8ln8wS4+OKLeeKJJxg0aBDXXHMN69at48c//jH77LMPAwcO5IYbbmjw99YhRESHepWUlIRZW/Diiy/We9qBP/tSbHdenw1eA3/2pUbPf+zYsTFy5MgG5zv44INjt912i+Li4iguLo4LLrggIiIuv/zy+MUvfhEREYceemi8+uqrERHx1FNPxSGHHBIREUcffXTceeedERFx/fXXx5ZbbhkREWvWrImVK1dGRMSyZcti1113jfXr18frr78enTp1imeffTYiIk444YT44x//GBERvXr1ilWrVkVExIoVKzao59SpU6Nbt24V9SwuLo7u3bvH3XffHatXr4799tsv3n777YiIGD9+fHznO9+pWL6zzjqr0vJeeOGFERFx7bXXRq9evWLx4sWxatWq6N27d7zzzjsREbF8+fKIiPj444+jqKioIh2IO+64IyIirrjiijj77LMryj3//PMjIuKBBx6Iww47LCIi/vCHP0Tv3r0ryjvnnHNi1KhRERExZcqUKC4ujoiIH/zgB3HFFVdERMT9998fQCxbtiwiouJ7nTp1ahx99NEVy3LDDTfE6NGjIyJi1apVUVJSEv/85z83+O6qWzeBmdEGtuGt8WpL1zGaWQ0Wraj2Et0a01vauHHjKC0trXbchx9+yJNPPskJJ5xQkfbpp58CMH36dP7yl78A8O1vf5sLLrgAyHbQL7nkEqZNm8Ymm2zCokWLKrp0+/Xrx6BBgwAoKSlh/vz5AAwcOJCTTz6Z4cOHM3z48GrrcuCBB3L//fdXfD799NMBeOWVV5g9ezaHH344AOvWraNXr14V05144omVyvnGN74BwF577UVRUVHFtF/4whdYsGABPXr04LrrrmPSpKwFv2DBAubOnUuPHj3YZJNNKso75ZRTOO644yrKLR/OLxfA4YcfznbbbQfA3//+d/785z8DcOihh7J8+XJWrlzJtGnTuOeeewA4+uij6d694rkLNXrkkUd4/vnnK7q/V65cydy5czv0NYvVcWA0awd6d9+JhSsWVZveWEVFRS1yfHD9+vVsu+22zJo1q955xo0bx7JlyygrK6Nz58707du34lq6zTffvGK6Tp06VXSlPvDAA0ybNo377ruP0aNHM2fOHDbdtH6btIigqKiI6dOnVzt+yy23rPS5vA6bbLJJpfpssskmrF27lscee4y//e1vTJ8+nS222IKhQ4fWeC1g/nKI8rI6depU6Xhifv4RG96TpLyMhl5aERH8+te/5sgjj6x74g7MxxjN2oHLjrqIrp27Vkrr2rkrlx11UaPLPPTQQ/n000+56aabKtJmzJjB448/3ugyAbp160a/fv24++67gWxj/NxzzwEwZMiQitbP+PHjK/KsXLmSHXbYgc6dOzN16lTeeOONWuexfv16FixYwCGHHMJ///d/89577/Hhhx/Wu4677747y5YtqwiMa9asYc6cOQ1azryVK1fSvXt3tthiC15++WWeeuqpSnUt3wG58847+fKXv9ygsg866CDGjRsHZGfabr/99nTr1q1S+kMPPcSKFSs2yLv11lvzwQcfVHw+8sgjuf7661mzZg0Ar776Kh999FHDFrYDcIvRrB04vuRYAEY/eDWLViymd/eduOyoiyrSG0MSkyZNYuTIkVx11VV06dKFvn37cu211wLZiTg1tfpOPvlkunbNAvX222/P3/72t0rjx40bx1lnncXPf/5z1qxZw4gRIyguLubaa6/llFNO4Ze//CVHH30022yzTUV5X//61yktLWXQoEF88YtfrLXu69at45RTTmHlypVEBOeddx7bbrttvZd9s802Y+LEifzwhz9k5cqVrF27lpEjR1JUVFTvMvKGDRvG7373OwYOHMjuu+/OkCFDKsZtueWWzJkzh5KSErbZZhv+9Kc/NajsUaNG8Z3vfIeBAweyxRZbcNtt2UOILr/8ck466SQGDx7MwQcfXO19oAcOHMimm25KcXExp59+Oueeey7z589n8ODBRAQ9e/as6Nq2z6i6ZvrGrLS0NJrz2iGzxnrppZfYY489Cl2NVvXxxx/TtWtXJDF+/HjuuuuuijNQN1ZbbbVVg1qzbUF166aksoio/sDyRsYtRjNrNWVlZZxzzjlEBNtuuy233HJLoatktgEHRjNrNQceeGDF8caOor21Fs0n35iZmVXiwGhmZpbjwGhmZpbjwGhmZpbjwGjWgW3sj52SxM0331yR9uyzzyKpznq2tWWx1uXAaNaOvPX+UkquPICl77/d5LIiNv7HTu21116VLqgfP348xcXFBayRtQcOjGbtyJhHxvLmuwsYM3lsk8vqCI+d2mWXXVi1ahVLly4lInj44Yf56le/WjH+pptuYp999qG4uJhvfvObfPzxx/VeFtt4OTCatRNvvb+Uu56ewPoI7nx6QpNbjbNnz6akpKRReU8++eSKrtQf//jHG4w/88wz+fWvf01ZWRljxozh3//93wE499xzOffcc5kxYwY77fTZDdC7dOnCpEmTeOaZZ5g6dSo/+tGPKm6ePXfuXM4++2zmzJnDtttuW3Gv1auuuopnn32W559/nt/97nc11vX444/n7rvv5sknn2Tw4MGVbgJ+3HHHMWPGDJ577jn22GOPSt2udS2Lbbx8gb9ZOzHmkbGsT8Fi/fr1jJk8ll9888qC1KW9PHYK4Fvf+hYnnngiL7/8MieddFLF0+0h2zm49NJLK25CXvWpE7Uti2283GI0awfKW4ur160GYPW61U1uNRYVFVFWVtZcVayQf+xU+eull16qNU/+sVOzZs1ixx13rPGxU+WPZ3rggQc4++yzKSsro6SkpNJjm/I+97nP0blzZyZPnsxhhx1Wadzpp5/Ob37zG1544QUuv/zyDR4V1Zhlsfav3QdGScMkvSJpnqSLW2IeE8smsecV+9Dj/J3Z84p9mFg2qSVmY1ajfGuxXHmrsbE60mOnfvazn3H11VfTqVOnSukffPABvXr1Ys2aNRWPcKrvstjGq10HRkmdgP8BvgoMAE6SNKA55zGxbBIjJ1zEkpVvAbBk5VuMnHCRg6O1qofnTK5oLZZbvW41D81+pNFllj92avLkyey6664UFRUxatSoimN/5d2X1ckfY/zKV76ywfhx48Zx8803U1xcTFFRUcXJNNdeey2/+tWv2HfffVmyZEmlx07NnDmT0tJSxo0bV+/HTu21117svffedT52av/996+2u3X06NF86Utf4vDDD69xnjUti2282vVjpyTtB4yKiCPT558ARMR/1ZSnoY+dKh49pNonp/fp3pvnLnuqmhxm9ePHTnWMx061R37sVPvWG1iQ+7wQ+FLViSSdCZwJVPswz9osWrG4QelmVjM/dsrag/YeGFVN2gZN4Ii4EbgRshZjQ2bwuW12rOhGrZpuZg3TER87Ze1Puz7GSNZC3Dn3uQ/QrE25/jvsWm36bjv8S3POxjqo9nwowzZOXifbf2CcAfSX1E/SZsAI4L7mnMHct1+rNv3Vt+c152ysA+rSpQvLly/3hsjajIhg+fLldOnSpdBVKah23ZUaEWslnQP8FegE3BIRc5pzHrMvn9GcxZlV6NOnDwsXLmTZsmWFropZhS5dutCnT59CV6Og2nVgBIiIB4EHC10Ps4bq3Lkz/fr1K3Q1zKyK9t6VamZm1qwcGM3MzHIcGM3MzHLa9Z1vGkPSMqD2GzHWbHvgnWasjlme1y9raU1Zxz4fET2bszJtVYcLjE0haWZHuSWStT6vX9bSvI7Vj7tSzczMchwYzczMchwYG+bGQlfANmpev6yleR2rBx9jNDMzy3GL0czMLMeB0czMLKfNB0ZJfSV9ImlW+jw/vQ+VdH8LzfPJRuZ7TFKznwqdvoPH0vCBkl6UNLu552NmZu0gMCavRcSg1ppZROzfWvMqJ6lTfaaLiCeAo1q4Oh1WTTtijShnlKQL0vCtko6vZdqRkrZozHzqqMN2kiZLmpveu6f0oZJuTcMnSprXUjuZhVb190xp89N7h9y5rmO6+Y0sf6Na39tLYMzLP6NnK0kTJb0saZwkAUgqkfS4pDJJf5XUK6U/JukaSdMkvSRpH0n3pC/y5+WFSvowvfdK086SNFvSgeXjJf1S0jOSpkjK3w3iBElPS3o1N30nSb+QNEPS85K+l9KHSpoq6U7ghZqmA9YB77bQ92kbatUdMWAk0OwbCuBiYEpE9AempM+VRMSfgH9rgXm3Ja39e7bpnes2YCRtfH1vd4ExIvbJfdyb7EseAHwBOEBSZ+DXwPERUQLcAlyZy7M6Ig4CfgfcC5wN7AmcLqlHldl9G/hr+lMVA7NS+pbAMxExGHgcuDyXZ9OI2DfVqzz9DGBlqvs+wHcllT9vaF/gPyJiQE3TRcSCiDiu/t+SNaOKHTFJF0p6QdJzkq5KabtKejjthD0h6YsNKVzSD4GdgKlpJ+kMSdfkxn9X0q/SHv/Lkm5LO00Ty/e6a9oRBI4BbkvDtwHD0/BqYGUjvouNhXeu6/HddOj1PSLa9AvoC8yuJn0oMDn3+XrgFLIg9z5ZEJsFvAA8kqZ5DDggDR9aJf80YFAa/jC9HwTMA0aVj0vp68gCIGQBeVY15e8IzEvDE4FXc3V6HTgiLcPUXLnVTlff78SvFl3fvgo8CWyRPm+X3qcA/dPwl4BH0/Ao4II0fCvZjlpN85wPbJ+GtwReAzqnz08Ce6V6RW79ugW4AOicpumZ0k8ke2A3wHtV5rOihvkPBe4v9Hffmr9nbrlXAn3IGgnTgS/X8Z0+Blydhs8FFgO9gM2BhUCPNK58G/Ijsh1fyB6mvnUaDuDkNPxT4De58n+Zho8C/paGzwQuTcObAzOBfmkZPgL61TZdI763Dr2+t/cHFX+aG15H9uBlAXMiYr868qyvkn89VR7cHBHTJB0EHA38UdIvIuL2asrMXwxaXmZ5fUh1+kFE/DWfSdJQspWa2qazNuErwB8i4mOAiHhX0lbA/sDdqaEB2cao0SLiI0mPAl+T9BLZBuMFSX2BBRHxf2nSO4AfAg+T7QxOTnXoBCxpSh06mKcjYiGAsuOQfYH3qP07vS+9v0C2rVmS8v8T2BlYnpt2BnBL6sn6S0TMSunrgT+l4TuAe3J5yofLUn0g25EeqM+O3W0D9CdrDT0dEa/XMV35+Prq0Ot7ew+M1XkF6Clpv4iYnlbI3SJiTkMLkvR5YFFE3CRpS2AwcDvZ3uXxwHiy7ta/11HUX4GzJD0aEWsk7QYsqu90EfFRNdNa6xKVd4AgWw/ei+Y/fvV74BLgZeAPufSq8w9q3xFcKqlXRCxJ3U1vN3M9Nwbeua5eh17f290xxrpExGqyoHW1pOfIuiQbeyB8KDBL0rPAN4GxKf0joEhSGVmX7M/qKOf3wIvAM8ous7iB6ndK6judtb5HgP+XO86xXUS8D7wu6YSUJknFjSj7A2Dr8g8R8Q+ylse3gbty0+0iqXyDcBLZDlnFjmCqQ2dJRWma+4DT0vBpZMfUrW61facNknau346Im4CbyXau4bOda2jYznXnVO5uaWe9wdNJ6i1pSh3z69Dre7vd6EbEY2T98eWfz8kNzyI7Plg1z9Ba8ufHbZXeb+Ozg7lVy7oMuKyW8t8hdYNExHqyPaJLqhRTtQ41TWcFFhEPSxoEzJS0GniQ7Hc6Gbhe0qVkxz/GA881sPgbgYckLYmIQ1LaBLLj2ity070EnCbpBmAucH1ErE7dZtdJ2obsP30tMAe4Cpgg6QzgTeCEhi53R1THd9pQQ4EfS1oDfAicmtLzO9cryY6V1eb3ZNuTZ5T1IS7js5NLGjpdL2BtbTPr6Ot7m79XqqSdyQ62Lm+BJnyjSPqwPHgWYN4HAr8l+z6GFqIOG7N0bOP+iNizwPW4H7gmIqa0dL1Sd9wFEfG15i670NrK71lVgbch5wBvRsR9dU7cStra+t7mu1Iju1Rh57YSFOGzFmWB5v1EROzloNhi1gHbKHdBeGuStK2kV4FPyjcSLTy/E8l2tFbUNW07VdDfsy2KiN+0laDYVtf3Nt9iNNtYSJpEdop93kU+C9k2Ru15fXdgNDMzy2nzXalmZmatyYHRzMwsx4HRLEfSjpLulPRPZfdinC7p2ELXK0/Z/TEjnZZenrZ3Srugjryj6prGrKNzYDRL0nVffwGmRcQXIrsJ/Qiye2lWnbbQ1wC/QOVr30bQ8OvJzKwaDoxmnzmU7OkrvytPiIg3IuLXAJJOl3S3pP8FHlH2/Le/KLv7/1OSBqbpKrXKlD1Voa9qf2LAVcoeQP28pDH1qOubQJfUwhUwDHgoN8/vKnvCwnOS/qxqnn+nJj4pwWxj5cBo9pki4Jk6ptkPOC0iDgWuAJ6NiIFkdwWp7h6YVe0O3JjyvA/8u6TtgGOBopT+89oKyJlIdneP/VO98/ftvCci9omIYrI7iJxRTf4bye6rWUL25ILf1nO+Zhs1B0azGkj6n9TimpFLnhwR5c+1+zLwR4CIeBTokW5TVZuqTwz4MlmAXAX8XtJxwMf1rOIEssB4EpXvMQmwZ2oFvkB2G69K9/pU5SclzCK7L28vzMyB0SxnDp/d5JmIOBs4DMg/RLbqkwyqCrL7UOb/W12qjK80fUSsJXtg9Z/J7mv5cH0qGxFvAWuAw8mek5d3K3BOROxF1rLtUmV8xZMScq896jNfs42dA6PZZx4lO253Vi5tg2NzOdPIWmPl9198Jz2BYD4pwEoaTOW7f2zwxIDUetsmIh4ERgKDUt5jJf1XHXX+KdndRNZVSd8aWKLsKQsnV83UjE9KMNvoFPrMOrM2IyJC0nDgGkkXkj2Z4CPgohqyjAL+IOl5su7P8kfe/Bk4NXVRzgBezeXZ4IkBZA+TvVdSF7JW6Hlp2l3Jullrq/OTNYy6DPgH8AbZGaxbVzNNczwpwWyj41vCmbWShj4xQNIdwHkRsaxFK2ZmlbjFaNZGRcQpha6DWUfkFqOZmVmOT74xMzPLcWA0MzPLcWA0MzPLcWA0MzPLcWA0MzPL+f85R7jQrfz37wAAAABJRU5ErkJggg==\n"
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(2)\n",
    "out_index = ww_results.index.levels[0]\n",
    "colors = sns.color_palette(\"dark\")[2]\n",
    "for i, key in enumerate(out_index):\n",
    "    data = ww_results.loc[key]\n",
    "    bic_data = data['bic'].to_numpy()\n",
    "    norm_bic = bic_data - np.min(bic_data)\n",
    "    like_data = data['likelihood'].to_numpy()\n",
    "    norm_like = like_data - np.min(like_data)\n",
    "    if 'Male' in key:\n",
    "        shape = '^'\n",
    "    else:\n",
    "        shape = 'o'\n",
    "    _ = ax[0].scatter(data['bic'].index,\n",
    "                     norm_bic,\n",
    "                     label=key,\n",
    "                     color=np.array(colors),\n",
    "                     marker=shape,)\n",
    "    _ = ax[1].scatter(data['likelihood'].index,\n",
    "                     norm_like,\n",
    "                     label=key,\n",
    "                     color=np.array(colors),\n",
    "                     marker=shape,)\n",
    "    _ = ax[0].set(xlabel='Groups, ' + types[i], ylabel='-BIC')\n",
    "    _ = ax[1].set(xlabel='Groups, ' + types[i], ylabel='Likelihood')\n",
    "_ = ax[0].legend(loc='best')\n",
    "_ = ax[1].legend(loc='best')\n",
    "_ = fig.suptitle('BIC and Likelihood: Worm Wiring C. Elegans Hermaphrodite and Male')\n",
    "_ = fig.set_size_inches(6, 8)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For a third time, cell type groupings get the best BIC. This makes sense since this is the same species as the developmental data\n"
   ]
  }
 ],
 "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.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}