5-Assess in-situ cell type calling in smFISH simulation.py.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Assess in-situ cell type calling in smFISH simulation\n",
    "\n",
    "In this step we assess the cell type calling accuracy that we can achieve in with cell expression vector identified in the smFISH simulation. This estimated accuracy can be interpret as the accuracy that would be obtained in the experiment while taking into account difficulties span by the limited number of markers and the RNA-cell assignement.  \n",
    "\n",
    "This accuracy depend of the previous choice of the marker list and the RNA assignement method that were used. A high cell type accuracy on smFISH simulation would confirm that our marker genes and RNA assignement method are relevant to detect thetargeted cell type\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### RNA assignement error estimation\n",
    "\n",
    "Assignement of RNA to cell is error prone in tissue due to the complexity of the cell shape, the possible high density of the nuclei and the lack of cytoplasm staining. Error in the RNA assignememnt to cell lead to error in the estimated RNA profile of the spatilly resolve cell and hence may lead to incorrect cell type calling.  \n",
    "    In this notebook we start by estimation the RNA assignememnt error using the mean jaccard index (also called IoU at the molecule level). In our context the jaccard index is\n",
    "$$\n",
    "  \\frac{GT\\cap P}{GT\\cup P}\\\n",
    "$$\n",
    "Where GT is the ground truth set of RNA of the cell and P is the predicted set of RNA of the cell"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scanpy as sc\n",
    "from comseg.utils.preprocessing import run_sctransform, sctransform_from_parameters\n",
    "from  sklearn.utils.random import  sample_without_replacement\n",
    "from pathlib import Path\n",
    "from sklearn import metrics\n",
    "from sklearn.metrics import confusion_matrix, ConfusionMatrixDisplay\n",
    "from matplotlib import pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "mean jaccard index is 0.8194918686079485\n"
     ]
    }
   ],
   "source": [
    "dict_result_jc = np.load(\"/home/tom/Bureau/phd/pipeline_review/test_set/dict_result_jc.npy\", allow_pickle=True).item()\n",
    "\n",
    "jaccard = []\n",
    "for image in dict_result_jc:\n",
    "    unique_cell = np.unique(dict_result_jc[image]['gr_cell'])\n",
    "    for cell in unique_cell:\n",
    "        df = dict_result_jc[image]\n",
    "        index_gr = df[df['gr_cell'] == cell].index\n",
    "        index_pred = df[df['pred_cell'] == cell].index\n",
    "        intersection = len(set(index_gr).intersection(set(index_pred)))\n",
    "        union = len(set(index_gr).union(set(index_pred)))\n",
    "        if union != 0:\n",
    "            jaccard.append(intersection / union)\n",
    "\n",
    "print(f'mean jaccard index is {np.mean(jaccard)}')\n",
    "      \n",
    "      "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The jaccard index is not equal to 1, hence estimated RNA profiles have some errors. We will try to classify theses estimated RNA profile to see if theses errors are small enought to not be detrimental for the cell type calling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### cell type classification method\n",
    "To classify each cell from its expression vector, we compute the cosine similarity between its RNA profile and the cell type centroids defined by scRNA-seq data "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We first compute the cell type centroid from the scRNA-seq data used for smFISH simulation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "optimisation can ba done here compute sctransform and centroind only once\n",
      "Writing temporary files...\n",
      "Running scTransform via Rscript...\n",
      "\n",
      "Attaching package: arrow\n",
      "\n",
      "The following objects are masked from package:feather:\n",
      "\n",
      "    read_feather, write_feather\n",
      "\n",
      "The following object is masked from package:utils:\n",
      "\n",
      "    timestamp\n",
      "\n",
      "Calculating cell attributes from input UMI matrix: log_umi\n",
      "Variance stabilizing transformation of count matrix of size 34 by 3996\n",
      "Model formula is y ~ log_umi\n",
      "Get Negative Binomial regression parameters per gene\n",
      "Using 34 genes, 3996 cells\n",
      "  |======================================================================| 100%\n",
      "Found 2 outliers - those will be ignored in fitting/regularization step\n",
      "\n",
      "Second step: Get residuals using fitted parameters for 34 genes\n",
      "  |======================================================================| 100%\n",
      "Calculating gene attributes\n",
      "Wall clock passed: Time difference of 0.6580091 secs\n",
      "Reading output files...\n",
      "Clipping residuals...\n"
     ]
    }
   ],
   "source": [
    "selected_gene = ['Atp6v0d2', 'Abcg1',# AM\n",
    "             'Rtkn2',  'Igfbp2', #AT1\n",
    "             'Sftpc','Cxcl15', #AT2,\n",
    "            'Cd79a', #B_cells\n",
    "             'Ms4a2', 'Fcer1a', #Basophils\n",
    "             'Ccdc153', #Ciliated\n",
    "             'Scgb3a2', 'Scgb1a1',#Club\n",
    "             'Cst3',#DC\n",
    "             'Cdh5', 'Clec14a',  #EC\n",
    "             'Inmt', 'Pcolce2', # Fibroblasts\n",
    "             'C1qc', 'C1qa', 'C1qb', # 'C3ar1', #IM\n",
    "             'Upk3b',# Mesotheliocytes\n",
    "             'Ifitm6','Plac8',# Monocytes\n",
    "            'Ms4a4b', 'Ccl5', 'Hcst', # NK_T_cells\n",
    "             'Gzma', 'Ncr1',# NK_cells\n",
    "             'S100a9',# Neutrophils\n",
    "             'Mmrn1',#Platelets\n",
    "           'Acta2','Myh11', # SMC\n",
    "             'Cd3g', 'Cd3d' #T_cells\n",
    "             ]\n",
    "\n",
    "path_to_rnaseq = \"../test_set/test_anndata.h5ad\"\n",
    "anndata = None\n",
    "normalize = True\n",
    "sample_size = 3999\n",
    "cell_type_column = \"cell_ID\"\n",
    "random_state= 1\n",
    "debug_path = None\n",
    "\n",
    "#todo implement the function without ssam package\n",
    "print('optimisation can ba done here compute sctransform and centroind only once')\n",
    "if anndata is None:\n",
    "    single_cells = sc.read(path_to_rnaseq)\n",
    "    random_index = sample_without_replacement(len(single_cells),\n",
    "                                              sample_size,\n",
    "                                   random_state = random_state)\n",
    "    anndata = single_cells[random_index, selected_gene].copy()\n",
    "    anndata.X = anndata.X.toarray()\n",
    "if type(anndata.X) != type(np.array([])):\n",
    "    anndata.X = anndata.X.toarray()\n",
    "#selected_genes_idx = [list(scrna_genes_qc).index(g) for g in ds.genes]\n",
    "scrna_unique_clusters =np.unique(list(anndata.obs[cell_type_column]))\n",
    "scrna_centroids = []\n",
    "## sc transform normalize the anndata\n",
    "# first remove the zero row (cell with no expression) #https://github.com/satijalab/sctransform/issues/86\n",
    "zero_rows_index  = list(np.where(np.sum(anndata[:, selected_gene].X, axis=1) == 0)[0])\n",
    "count_matrix = np.delete(anndata[:, selected_gene].X, zero_rows_index, axis=0)\n",
    "if normalize:\n",
    "    count_matrix, param_sctransform = run_sctransform(count_matrix, debug_path = debug_path)\n",
    "else:\n",
    "    param_sctransform  = None\n",
    "#print(count_matrix.shape)\n",
    "label_cell_type = np.delete(np.array(anndata.obs[cell_type_column]), zero_rows_index)\n",
    "\n",
    "for cl in scrna_unique_clusters:\n",
    "    #print(cl)\n",
    "    mean = np.mean(count_matrix[label_cell_type == cl], axis=0)\n",
    "    #print(mean)\n",
    "    scrna_centroids.append(mean)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/home/tom/Bureau/phd/pipeline_review/test_set/anndata/07_CtrlNI_Pdgfra-Cy3_Serpine1-Cy5_004.h5ad\n",
      "/home/tom/Bureau/phd/pipeline_review/test_set/anndata/07_CtrlNI_Pdgfra-Cy3_Serpine1-Cy5_006.h5ad\n",
      "cell type accuracy0.9402985074626866\n"
     ]
    }
   ],
   "source": [
    "y_true_cell_type  = []\n",
    "y_pred_cell_type = []\n",
    "for anndata_path in Path(\"/home/tom/Bureau/phd/pipeline_review/test_set/anndata/\").glob('*h5ad'):\n",
    "    print(anndata_path)\n",
    "    anndata = sc.read(anndata_path)\n",
    "    anndata.X = anndata.X.toarray()\n",
    "    for i in range(len(anndata)):\n",
    "        y_true_cell_type.append(anndata.obs['cell_type'][i])\n",
    "        expression_vector = anndata.X[i]\n",
    "        expression_vector = expression_vector.reshape(1, len(expression_vector))\n",
    "        norm_expression_vectors = sctransform_from_parameters(\n",
    "                    np.array(param_sctransform),\n",
    "                    expression_vector)\n",
    "        \n",
    "        correlation_array = metrics.pairwise.cosine_similarity(norm_expression_vectors,\n",
    "                                                                   np.array(scrna_centroids))[0]\n",
    "        index_cluster_max = np.argmax(correlation_array)\n",
    "        pred_rna_seq = scrna_unique_clusters[index_cluster_max]\n",
    "        \n",
    "        y_pred_cell_type.append(pred_rna_seq)\n",
    "\n",
    "acc = metrics.accuracy_score(y_true=y_true_cell_type,\n",
    "                             y_pred=y_pred_cell_type)\n",
    "print(f\"cell type accuracy{acc}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The obtain cell type accuracy is 0.94"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<sklearn.metrics._plot.confusion_matrix.ConfusionMatrixDisplay at 0x7f7e9078b880>"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAn0AAAJjCAYAAAB5r6cmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8WgzjOAAAACXBIWXMAAAsTAAALEwEAmpwYAAB8AklEQVR4nO3deZxWZf3/8dd7YAARB2QJkUCxDDMy1EnDrDArs41Wy2yxRbI0MzNbtLJFssUWl/JLlEuGlJlpWUpSKFJioLglWD9FTUAFA9yAYebz++OcwZthNuY+5z5z3/N+Ph73w/s+577f5zPHYeaa61zXdRQRmJmZmVltqyu6ADMzMzPLnxt9ZmZmZn2AG31mZmZmfYAbfWZmZmZ9gBt9ZmZmZn2AG31mZmZmfUD/ogswMzMz602OOGznWPtEc0WOteTOTddHxBsrcSw3+szMzMxKrH2imVuvH1+RY/Ub8++RFTkQbvSZmZmZbSOAFlqKLiNzHtNnZmZm1ge4p8/MzMxsG0FzuKfPzMzMzKqQe/rMzMzMSiRj+qLoMjLnnj4zMzOzPsCNPjMzM7M+wJd3zczMzNrwki1mZmZmVpXc02dmZmZWIgiawxM5zMzMzKwKuafPzMzMrA0v2WJmZmZmVck9fWZmZmYlAmh2T5+ZmZmZVSP39JmZmZm14TF9ZmZmZlaV3NNnZmZmViLA6/SZmZmZWXVyT5+ZmZlZG7V351339JmZmZn1Ce7pMzMzMysRhNfpMzMzM7Pq5EafmZmZWR/gy7tmZmZmpQKaa+/qrnv6zMzMzPoC9/SZmZmZlQi8ZIuZmZmZVSn39JmZmZltQzSjoovInHv6zMzMzPoA9/SZmZmZlQigxbN3zczMzKwauafPzMzMrA2P6TMzMzOzquSePjMzM7MSgXv6zMzMzKxKuafPzMzMrI2WcE+fmZmZmVUh9/SZmZmZlfCYPjMzMzOrWm70mZmZmfUBvrxrZmZmViIQzTXYL1Z7X5GZmZmZbcc9fWZmZmZteMkWMzMzM6tK7ukzMzMzK+ElW8zMzMysarmnrwoN6D84dqofmkt2bNyUS66ZmVlPbeRpNsemCna9ieaovX4xN/qq0E71Q5nywo/lkt18z/Jccs3MzHpqUcwruoSa4EafmZmZWYkAWmpwBFztfUVmZmZmth339JmZmZm14dm7ZmZmZlaV3NNnZmZmViKiNmfv1t5XZADU1zfzw/Nu4PwLr+enP7uOYz50d+bHaJy6gVkLlnHRwns56sRHqyY773zXXky+ay8m37UXk+/arSfc6KsASW+XFJL2SV/vmb7+Vsl7RkpqknR+FsdsaqrjS59/DScefwQnHv8GGhtXM/HFa7OIBqCuLjhhxiOcccwEjps6kcOmrWP83ht7fXbe+a69mHzXXky+ay8m37VXRguqyKOS3OirjKOBm9P/tnoAeHPJ6/cA92R3SLFxYz0A/fu30K9/SzIHPSMT93+GlSsGsPqhgWxpqmP+1cOYcsT6Xp+dd75rLybftReT79qLyXft1lNu9OVM0hDgUOBjwPtKdj0D3CupMX39XuA3WR67rq6F8y6cy+wrruH220azfNmIzLJH7NbE4ysHbH29ZlU9I8c09frsvPNdezH5rr2YfNdeTL5rt57yRI78TQOui4j7JK2VdCDQep11DvA+SY8CzcBKYPf2QiRNB6YDDKpv6NaBW1rq+PTxb2DnnTdzxpkL2WPP9Ty4Ip/bt5mZmdWKAJprsF+s9r6i3udoksYd6X9LL/FeB7yepAfw152FRMTMiGiMiMYB/QbvUAFPPz2AO+94Hgc2rtqhz3Vm7ep6Ru2+eevrkWOaWLOqvtdn553v2ovJd+3F5Lv2YvJdu/WUG305kjQceC0wS9IK4PPAUZCM3IyIzcAS4HPAb7M8dsPQjey8c/IPa8CALex/wKP89+Hu9RB2x/Klgxk7YTOjx22if30LU6et45a52fQi5pmdd75rLybftReT79qLyXftlZAs2VKJRyX58m6+3g38MiI+0bpB0o3AuJL3nAPcGBFPSNnN4hk+fCOfO+1W6uoCKVhw0zhuXdTuleMeaWkWF5w+lhmz76euH8ydM5wH7xvU67PzznftxeS79mLyXXsx+a7dekoRGU7ptG1I+hvwnYi4rmTbScCRwLiImNTm/ccCjRFxYme5Q3caE1Ne+LEcKobme5bnkmtmZtZTi2IeG+KJiq1v8sKXDo5zrn5RRY719hfcsSQiGrt+Z/nc05ejiDisnW3nAud28P6LgYvzrcrMzMz6Ijf6zMzMzNpojsounFwJnshhZmZm1ge4p8/MzMysRCCv02dmZmZm1ck9fWZmZmZttFR4Db1KcKOvCsXGTbktrdKvIbsFnNvTvGFDrvlmZmbWPjf6zMzMzEr43rtmZmZmVrXc02dmZmZWIpDX6TMzMzOz6uSePjMzM7M2WmqwX6z2viIzMzMz244bfWZmZmZ9gBt9Naxx6gZmLVjGRQvv5agTH800++Sz7mP2wlv4yTVLMs1tlWfteee79mLyXXsx+a69mHzXnq8IaI66ijwqyY2+nEh6u6SQtI+kRZKWSnpI0uPp86WS9pR0lqSHJT2V5fHr6oITZjzCGcdM4LipEzls2jrG770xs/wbrhrNV46blFleqbxrzzPftReT79qLyXftxeS7duspN/ryczRwM3B0RBwcEZOBrwK/jojJ6WMF8AfgoKwPPnH/Z1i5YgCrHxrIlqY65l89jClHrM8s/+7FQ3lyfT7zgPKuPc98115MvmsvJt+1F5Pv2itBtFToUUlu9OVA0hDgUOBjwPs6e29E3BIRq7KuYcRuTTy+csDW12tW1TNyTFPWh8lF3rXnme/ai8l37cXku/Zi8l279ZSXbMnHNOC6iLhP0lpJB0ZEPoPfzMzMLFMBFR9vVwm19xX1DkcDc9Lnc9LXZZE0XdJiSYub2NTl+9eurmfU7pu3vh45pok1q+rLLaMi8q49z3zXXky+ay8m37UXk+/arafc6MuYpOHAa4FZklYAnweOklTWhfuImBkRjRHRWM/ALt+/fOlgxk7YzOhxm+hf38LUaeu4Ze7QckqomLxrzzPftReT79qLyXftxeS79spopq4ij0ry5d3svRv4ZUR8onWDpBuBV1WyiJZmccHpY5kx+37q+sHcOcN58L5BmeWfds4y9nv5Ohp23cKl8xdx2Xl7MPfK3TLJzrv2PPNdezH5rr2YfNdeTL5rt55SRBRdQ02R9DfgOxFxXcm2k4AXA4uAxog4sWTfd4H3A7sDK4FZEXFmZ8do0PA4WIfnUD30a2jIJbdV84YNueabmVntWRTz2BBPVGyq6/hJDXHabxsrcqxPv/hvSyKiw4NJGgdcCowmGW44MyJ+nF5Z/DWwJ7ACOCoi/tfZsdzTl7GIOKydbeeWvLy4zb7TgNNyLsvMzMyq0xbgcxFxm6RdgCWS/gIcC8yLiLMlfRH4IvCFzoLc6DMzMzNro9Lj7TqSLuu2Kn3+pKR7gbEkK4VMTd92CTCfLhp9veMrMjMzM7NOSdoT2J9kuNjoknV+V5Nc/u2Ue/rMzMzMSgTQUrl1+kZKWlzyemZEzGz7pvTGD1cCJ0fEhtJFQSIiJHU5ScONPjMzM7PirOlsIgeApHqSBt+vIuJ36eZHJY2JiFWSxgCPdXUgX941MzMz24ZortCjy0qSLr2fA/dGxA9Kdl0DfDh9/mHg6q6y3NNnZmZm1nu9EvggcJekpem2LwNnA7+R9DHgQeCoroLc6LNt5L2OXtPrDsw1v/4G3+LYzMxqR0TcDB12Ce7Qor1u9JmZmZmVqPBEjoqpva/IzMzMzLbjnj4zMzOzNrozyaLauKfPzMzMrA9wT5+ZmZlZiQh5TJ+ZmZmZVSf39JmZmZm10eyePqsmjVM3MGvBMi5aeC9HnfhoVeXPPuc3zDrrKmZ+8/f89OtdLjK+w/KsvZrPu2svJt+1F5Pv2ovJz7t265gbfTmR9HZJIWkfSYskLZX0kKTH0+dLJe0r6VpJyyTdI+nsrI5fVxecMOMRzjhmAsdNnchh09Yxfu+NWcXnng9wyrePZPpX3s4nvzYt09w8a6/m8+7ai8l37cXku/Zi8ivxuyMLAbSgijwqyY2+/BwN3AwcHREHR8Rk4KvAryNicvp6BfD9iNgH2B94paQjszj4xP2fYeWKAax+aCBbmuqYf/UwphyxPovoiuTnKc/aq/m8u/Zi8l17MfmuvZj8av7dUQvc6MuBpCHAocDHgPd19L6IeCYi/pY+3wzcBjw/ixpG7NbE4ysHbH29ZlU9I8c0ZRFdkfwAvnfa9Vz49at589RlmeVCvrVX83l37cXku/Zi8l17Mfl5154d0Rx1FXlUkidy5GMacF1E3CdpraQDI6LTm8JKGga8FfhxB/unA9MBBjE443J7n898682s+d/ODNvlWb73het4eNUw7ly+W9FlmZmZVS339OXjaGBO+nxO+rpDkvoDlwPnRsT97b0nImZGRGNENNYzsMsC1q6uZ9Tum7e+HjmmiTWr6rtZftfyzl/zv50BWPfkTty8ZA/22evxzLLzrL2az7trLybftReT79qLyc+79qwk995VRR6V5EZfxiQNB14LzJK0Avg8cJSkzv7PzgT+HRE/yqqO5UsHM3bCZkaP20T/+hamTlvHLXOHZhWfa/6gAU3sNKhp6/PGSSt54L+7ZpIN+dZezefdtReT79qLyXftxeTnXbt1zpd3s/du4JcR8YnWDZJuBF7V3pslfQsYCnw8yyJamsUFp49lxuz7qesHc+cM58H7BlVF/q5Dn+Ubn5kHQL+6YN4/9uKfd2Uy1BHIt/ZqPu+uvZh8115MvmsvJj/v2rPUXIP9YoqIomuoKZL+BnwnIq4r2XYS8GJgEdAYESem258PPAwsAzalbz8/ImZ1dowGDY+DdXge5eeu6XUH5ppff0OnQyfNzKwKLYp5bIgnKnYtdMxLdo0PX16Z37PfedmVSyKisRLHck9fxiLisHa2nVvy8uKS7f+FCi/SY2ZmZn2SG31mZmZmJYLKT7KohNq7YG1mZmZm23FPn5mZmVkbLTXYL1Z7X5GZmZmZbcc9fWZmZmYlIqDZY/rMzMzMrBq5p88qKu919O678KDcsl90/K25ZZuZWe/i2btmZmZmVpXc02dmZmZWIlmnr/b6xWrvKzIzMzOz7binz8zMzKyN5hq8S6p7+szMzMz6APf0mZmZmZUIPHvXqkzj1A3MWrCMixbey1EnPlpV+Vlnj770fvb6/G3s8Y27tm4bcc1/2eObdzH+W3cz9sfL6Lduc9nHAZ/3ovJdezH5rr2YfNduPeFGX42qqwtOmPEIZxwzgeOmTuSwaesYv/fGqsjPI3vDlJE88umJ22z73+vH8OBXXspDZ0zi6ZcOY8S1j5R1DPB5LyrftReT79qLyXftlZDM3q3Eo5Lc6GtDUrOkpZLukHSbpEMyzD5T0qnp84slvTur7LYm7v8MK1cMYPVDA9nSVMf8q4cx5Yj1VZGfR/azezfQPHjb0QwtO/Xb+lybW0Dld+X7vBeT79qLyXftxeS7duspN/q292xETI6IlwFfAr5ddEE9MWK3Jh5fOWDr6zWr6hk5pqkq8vOufZtj/f5hJnxpKQ23rmXtW8eWn+fzXki+ay8m37UXk+/arafc6OtcA/C/zt4g6QuS7kp7Bs9Ot71A0nWSlkhaIGmfLjLOlvQvSXdK+n6G9VsX1r59HA98ezIbDhrBsPkeW2JmZokWVJFHJXn27vZ2krQUGASMAV7b0RslHQlMAw6OiGckDU93zQSOj4h/SzoY+ElHOZJGAO8A9omIkDSsg/dNB6YDDGJwl1/E2tX1jNr9uYkJI8c0sWZVfZef66488/OuvT1PHjSCseffx9q3Pr+sHJ/3YvJdezH5rr2YfNduPeWevu21Xt7dB3gjcKnU4WCv1wEXRcQzABHxhKQhwCHAFWnj8f9IGo8dWQ9sBH4u6Z3AM+29KSJmRkRjRDTWM7DLL2L50sGMnbCZ0eM20b++hanT1nHL3KFdfq678szPu/ZW9Y8+N3h4yB3/Y/PoQWVn+rwXk+/ai8l37cXku/b8RUBzqCKPSnJPXyci4h+SRgKjgMe6+bE6YF1ETO7mMbZIOgg4HHg3cCKd9C52V0uzuOD0scyYfT91/WDunOE8eF/5jZpK5OeRvdus/zD4vifp99QWJnzxdta+9fnsfPc6Bjy6EQRNwwfy2Pv37JW1VyrftReT79qLyXftxeTnXbt1ThFRdA29iqSnImJI+nwf4GZgdEQ0t/PeNwJfBV7Xenk37e37O/DDiLgi7SXcLyLukHQm8FREfF/SxcAfgeuAwRHxmKShwP0RMaKzGhs0PA7W4Rl+1bXjvgsPyi37Rcffmlu2mZl1bFHMY0M8UbFusZEvHhlvvmRaRY516cG/WBIRjZU4lnv6ttc6pg9AwIfba/ABRMR1kiYDiyVtBv4EfBk4BvippDOAemAOcEcHx9sFuFrSoPR4p2T1hZiZmZm1cqOvjYjo1/W7tnn/2cDZbbY9QDIesO17zyx5fmzJrvy6p8zMzGyHBPJt2MzMzMysOrmnrxskvRT4ZZvNmyLi4CLqMTMzs3xVeg29SnCjrxsi4i5gctF1mJmZmfWUG31mZmZmJQI8ps/MzMzMqpN7+qymeC09MzPLQkvUXr9Y7X1FZmZmZrYd9/SZmZmZlQqv02dmZmZmVcqNPjMzM7M+wJd3zczMzEoEtbk4s3v6zMzMzPoA9/SZmZmZteGJHFZVGqduYNaCZVy08F6OOvHRqsp37cXku/Zi8l17MfmuvZj8vGu3jtV0o09Ss6Slku6QdJukQyp03Islvbud7btL+m36fKqkP+ZVQ11dcMKMRzjjmAkcN3Uih01bx/i9N1ZFvmsvJt+1F5Pv2ovJd+3F5Odde1Zab8NWiUcl1XSjD3g2IiZHxMuALwHfLrKYiFgZEds1BvMwcf9nWLliAKsfGsiWpjrmXz2MKUesr4p8115MvmsvJt+1F5Pv2ovJz7t261ytN/pKNQD/A5A0RNK8tPfvLknT0u07S7o27Rm8W9J70+2HS7o9fe8vJA1Mt6+Q9N10+62SXlhyvFdL+ruk+1t7/STtKenutoVJek3aI7k0Pc4u5X6xI3Zr4vGVA7a+XrOqnpFjmsqNrUi+ay8m37UXk+/ai8l37cXk5117lmqxp6/WJ3LsJGkpMAgYA7w23b4ReEdEbJA0ErhF0jXAG4GVEfFmAElDJQ0CLgYOj4j7JF0KfBL4UZq1PiJeKulD6ba3pNvHAIcC+wDXAL/tpM5TgRMiYqGkIWl925A0HZgOMIjBO3oezMzMrI+r9Z6+1su7+5A06C6VJEDADEl3AjcAY4HRwF3A6yV9R9KrImI9MBF4ICLuSzMvAV5dcozLS/47pWT77yOiJSL+lWZ3ZiHwA0knAcMiYkvbN0TEzIhojIjGegZ2+YWvXV3PqN03b309ckwTa1bVd/m57soz37UXk+/ai8l37cXku/Zi8vOuPStBZXr5PKYvJxHxD2AkMAo4Jv3vgRExGXgUGJQ27A4gafx9S9JXuxPdwfNNJc87/b8aEWcDHwd2AhZK2qcbx+3U8qWDGTthM6PHbaJ/fQtTp63jlrlDy42tSL5rLybftReT79qLyXftxeTnXbt1rtYv726VNqT6AWuBocBjEdEk6TBgj/Q9uwNPRMRlktaRNMS+C+wp6YUR8R/gg8CNJdHvBc5O//uPHtb2goi4C7hL0stJLgkv60lWq5ZmccHpY5kx+37q+sHcOcN58L5B5URWLN+1F5Pv2ovJd+3F5Lv2YvLzrj1LtXhHDkVE1++qUpKaSXrtIOlt+3JEXJuO4/sDMARYDLwCOJLkUu73gBagCfhkRCyWdDjwfZJG8j/T7ZskrQB+nX52E3B0RPxH0sXAHyOidXmWpyJiiKQ90+2TJE0FTo2It0g6DzgsPe49wLERUdpTuI0GDY+DdXgm58jMzKy3WxTz2BBPVKwVNnSf0XHIzPdV5FjXvebcJRHRWIlj1XRPX0T062D7GrYdf9dqBXB9O++fB+zfwWG+FxFfaPP+Y9u8HpL+dwUwKX0+H5ifPv90R1+DmZmZVVj4jhxmZmZmVqVquqcvbxGxZ9E1mJmZWbZa78hRa9zTZ2ZmZtYHuNFnZmZm1gf48q6ZmZlZG768a2ZmZmZVyT19Zr1Ev4aG3LKbN2zILdvMrNa03oat1rinz8zMzKwPcE+fmZmZWRvhnj4zMzMzq0bu6TMzMzNrowX39JmZmZlZFXJPn5mZmVmJCK/TZ1WmceoGZi1YxkUL7+WoEx+tqnzX3r6Tz7qP2Qtv4SfXLMk0t5XPezH5rr2YfNdeTH7etVvH3OgDJO0maY6k/ydpiaQ/SXq1pN+m+6dK+mP6/G2SvthF3lRJh/SgjhWSRvbsq9hWXV1wwoxHOOOYCRw3dSKHTVvH+L03ZhGde75r79gNV43mK8dNyiyvlM97MfmuvZh8115Mft61ZylCFXlUUp9v9EkScBUwPyJeEBEHAl8CIiLe3fb9EXFNRJzdRexUYIcbfVmauP8zrFwxgNUPDWRLUx3zrx7GlCPWV0W+a+/Y3YuH8uT6fEZl+LwXk+/ai8l37cXk5127da7PN/qAw4CmiLiwdUNE3AE8LOnutm+WdKyk89Pnb5W0SNLtkm6QNFrSnsDxwGclLZX0KkmjJF0p6Z/p45Xp50dImivpHkmzILupQiN2a+LxlQO2vl6zqp6RY5qyis8137UXw+e9mHzXXky+ay8mv3p+RiZ35KjEo5Lc6INJQE8HSN0MvCIi9gfmAKdFxArgQuCHETE5IhYAP05fvxx4FzAr/fzXgJsj4iUkvY3je/5lmJmZmXXMs3fL83zg15LGAAOABzp43+uAfZMryQA0SBoCvBp4J0BEXCvpfx0dSNJ0YDrAIAZ3Wdja1fWM2n3z1tcjxzSxZlV9l5/rrjzzXXsxfN6LyXftxeS79mLyq+lnpO/IUZvuAQ7s4WfPA86PiJcCnwAGdfC+OpIewcnpY2xEPLUjB4qImRHRGBGN9Qzs8v3Llw5m7ITNjB63if71LUydto5b5g7dkUMWlu/ai+HzXky+ay8m37UXk1/NPyNrgXv64K/ADEnTI2ImgKT9gO58Fw4FHkmff7hk+5NAQ8nrucCnge+l+ZMjYilwE/B+4FuSjgR2LePr2EZLs7jg9LHMmH0/df1g7pzhPHhfR23S3pXv2jt22jnL2O/l62jYdQuXzl/EZeftwdwrd8sk2+e9mHzXXky+ay8mP+/arXOKiKJrKJyk3YEfkfT4bQRWACcDV0XEJElTgVMj4i2SjgUaI+JESdOAHwL/I2k8vjwipkp6EfBboIWksXcvcAHwYpKG9k0RcbykEcDlwFjg78AbgAMjYk1n9TZoeByswzP7+q136NfQ0PWbeqh5w4bcss3M8rYo5rEhnqjY9dadXzQmXnLuRypyrH8e+e0lEdFYiWO5pw+IiJXAUe3smpTunw/MT59fDFycPr8auLqdvPuA/dpsfm8771tL0tAzMzMzy5UbfWZmZmalIrkVW63xRA4zMzOzPsA9fWZmZmZttGR3v4Rewz19ZmZmZn2Ae/rMzMzMSgRenNnMzMzMqpR7+sx6iTzX0stzDUDwOoBmVmtEi3v6zMzMzKwauafPzMzMrA2v02dmZmZmVck9fWZmZmZtePaumZmZmVUl9/SZmZmZlYhwT5+ZmZmZVSk3+mpY49QNzFqwjIsW3stRJz5aVfmuvZj8k8+6j9kLb+En1yzJNLeVz3vls/POd+3F5Lt26wk3+jIgaTdJcyT9P0lLJP1J0osk3d2Nzz6VR011dcEJMx7hjGMmcNzUiRw2bR3j995YFfmuvbj8G64azVeOm5RZXimf98pn553v2ovJd+2V0RKqyKOS3OgrkyQBVwHzI+IFEXEg8CVgdJF1Tdz/GVauGMDqhwaypamO+VcPY8oR66si37UXl3/34qE8uT6fob4+75XPzjvftReT79qtp9zoK99hQFNEXNi6ISLuAB5ufS3pWEnnl7z+o6SpJa9/KOkeSfMkjcqiqBG7NfH4ygFbX69ZVc/IMU1ZROee79qLy8+Tz3vls/POd+3F5Lv2ykgmc+T/6A5Jv5D0WOkVRElnSnpE0tL08aauctzoK98koJwBUDsDiyPiJcCNwNfae5Ok6ZIWS1rcxKYyDmdmZmZV5mLgje1s/2FETE4ff+oqxEu2FK8F+HX6/DLgd+29KSJmAjMBGjS8y78N1q6uZ9Tum7e+HjmmiTWr6ssuthL5rr24/Dz5vFc+O+98115MvmuvjN60ZEtE3CRpz3Jz3NNXvnuAA7t4zxa2PdeDOnlvJnf7W750MGMnbGb0uE30r29h6rR13DJ3aBbRuee79uLy8+TzXvnsvPNdezH5rt1KnCjpzvTy765dvdk9feX7KzBD0vS0Nw5J+wGl38UrgE9JqgPGAgeV7KsD3g3MAd4P3JxFUS3N4oLTxzJj9v3U9YO5c4bz4H2dtTV7T75rLy7/tHOWsd/L19Gw6xYunb+Iy87bg7lX7pZJts975bPzznftxeS79vwFqmRP30hJi0tez2xtT3Thp8A3STqLvgmcA3y0sw8oujuK0DokaXfgRyQ9fhtJGnknA1dFxKR0hu9l6f57gV2BMyNifrpky0zgDcBjwHsj4vHOjteg4XGwDs/ni7Ga1K+hIdf85g0bcs03s75tUcxjQzxRsVbYoBeOjT2/+4mKHGv5u762JCIau3pfenn3jxGx3bpane0r5Z6+DETESuCodnZNSvcHcEwHnx2SY2lmZmbWA729S0zSmIhYlb58B9Dl2sBu9JmZmZn1YpIuB6aSXAr+L8lKH1MlTSZpn64AuuyadKPPzMzMrFT0utm7R7ez+ec7muPZu2ZmZmZ9gHv6zMzMzNrq7YP6esA9fWZmZmZ9gHv6zMzMzNroTWP6suJGn1kf4HX0bEd5bUez2uPLu2ZmZmZ9gHv6zMzMzNqoxRuWuafPzMzMrA9wT5+ZmZlZiaA2J3K4p8/MzMysD3BPn5mZmVmpANzTZ9WkceoGZi1YxkUL7+WoEx+tqnzXXky+ay8mv5prP/ms+5i98BZ+cs2STHNb+bxXPjvv/Lxrt4650ZczSc2Slkq6R9Idkj4nqa5k/0GSbpK0XNLtkmZJGlzucevqghNmPMIZx0zguKkTOWzaOsbvvbHc2Irku/Zi8l17MfnVXDvADVeN5ivHTcosr5TPe+Wz887Pu/YsRVTmUUlu9OXv2YiYHBEvAV4PHAl8DUDSaOAK4AsRMTEi9geuA3Yp96AT93+GlSsGsPqhgWxpqmP+1cOYcsT6cmMrku/ai8l37cXkV3PtAHcvHsqT6/MZKeTzXvnsvPPzrt0650ZfBUXEY8B04ERJAk4ALomIf5S857cRUXZ/94jdmnh85YCtr9esqmfkmKZyYyuS79qLyXftxeRXc+1583mvfHbe+VX1/RgVelSQG30VFhH3A/2A5wGTgHwGwpiZmZmV8OzdKiFpOkkvIYPoesjf2tX1jNp989bXI8c0sWZVfWb15Jnv2ovJd+3F5Fdz7Xnzea98dt751fP9KK/TZ+WTtBfQDDwG3AMc2J3PRcTMiGiMiMZ6Bnb5/uVLBzN2wmZGj9tE//oWpk5bxy1zh5ZVe6XyXXsx+a69mPxqrj1vPu+Vz847v5q/H2uBe/oqSNIo4ELg/IgISecDt0q6NiIWpe95J7Cw3HF9Lc3igtPHMmP2/dT1g7lzhvPgfYPK/yIqkO/ai8l37cXkV3PtAKeds4z9Xr6Ohl23cOn8RVx23h7MvXK3TLJ93iufnXd+3rVnqgbvvauoxTsK9yKSmoG7gHpgC/BL4AcR0ZLunwJ8l2SMXwtwE/DZiHimo8wGDY+DdXjepZtZH9avoSHX/OYNG3LNt9qyKOaxIZ6o2PXWgXs9P3b/xgkVOdaKD355SUQ0VuJY7unLWUT062L/P4BXVagcMzMz60r43rtmZmZmVqXc6DMzMzPrA3x518zMzKytGpzy4J4+MzMzsz7APX1mZmZm2/FEDjMzMzOrQu7pM7Pe7RX75Zt/y5355lcpr6NnfZ7H9JmZmZlZNXJPn5mZmVlb7ukzMzMzs2rknj4zMzOzUgH4NmxmZmZmVo3c02dmZmbWRnhMn5mZmZlVIzf6aljj1A3MWrCMixbey1EnPlpV+a69mPxqrh2grq6FC77/R77x5b9mnu3zXvnsvPNdezH5edeemajQo4Lc6KsQSc2SlpY8vphur5d0tqR/S7pN0j8kHVnu8erqghNmPMIZx0zguKkTOWzaOsbvvbH8L6QC+a69mPxqrr3V29+8jIf/OzTTTPB5LyI773zXXkx+JX4OWMfc6KucZyNicsnj7HT7N4ExwKSIOAB4O7BLuQebuP8zrFwxgNUPDWRLUx3zrx7GlCPWlxtbkXzXXkx+NdcOMHLE0xx04CP8+YYXZpbZyue98tl557v2YvLzrj1Toco8KsiNvgJJGgwcB3w6IjYBRMSjEfGbcrNH7NbE4ysHbH29ZlU9I8c0lRtbkXzXXkx+NdcOcPxHFzPr0gOIHH6I+rxXPjvvfNdeTH7etVvnPHu3cnaStLTk9beBe4GHIqLLm1xKmg5MBxjE4FwKNKtWBx/4X9atH8R/7h/Bfi9ZXXQ5ZlYDVIOzdzts9Ek6j06GGEbESblUVLuejYjJpRskdftO8hExE5gJ0KDhXX4rrl1dz6jdN299PXJME2tW1Xe72CLzXXsx+dVc+777PMYrXv5fXn7AIwyob2bw4CZO+8zNfPfHh2aS7/Ne+ey88117Mfl5126d6+zy7mJgSScPK99/gPGSGrIOXr50MGMnbGb0uE30r29h6rR13DI3uwHueea79mLyq7n2i351AB847l18+Ph38u0fvIo77totswYf+LwXkZ13vmsvJj/v2q1zHfb0RcQlpa8lDY6IZ/Ivqe+IiGck/Rz4saRPRMRmSaOAqRFxRTnZLc3igtPHMmP2/dT1g7lzhvPgfYOyKTznfNdeTH411543n/fKZ+ed79qLya+anwMFLKdSCYoulpyWNAX4OTAkIsZLehnwiYj4VCUKrBWSmoG7SjZdFxFflDQA+BbwTmAj8DTw1Yi4vqOsBg2Pg3V4rvWa9Rqv6PYoiJ655c58882sbItiHhviiYpNdR24x/NjzOmfqcixHvzEaUsiorESx+rORI4fAUcA1wBExB2SXp1nUbUoIvp1sH0zcFr6MDMzs8JVfjmVSujWki0R8XCbTc051GJmZmZmOelOT9/Dkg4BQlI98BmSpUbMzMzMalMNjunrTk/f8cAJwFhgJTA5fW1mZmZmVaLLnr6IWAMcU4FazMzMzHqHvtjTJ2kvSX+Q9LikxyRdLWmvShRnZmZmZtnozuXd2cBvgDHA7sAVwOV5FmVmZmZWqKjQo4K6M5FjcET8suT1ZZI+n1dBZmbb8Dp6ZmaZ6Ozeu8PTp3+W9EVgDkmb9L3AnypQm5mZmVnlBTW5Tl9nPX1LSL7s1q/6EyX7AvhSXkWZmZmZWbY6u/fuhEoWYmZmZtZbqAZn73ZnTB+SJgH7AlvvihwRl+ZVlJmZmZllq8tGn6SvAVNJGn1/Ao4Ebgbc6DMzM7PaVIM9fd1ZsuXdwOHA6oj4CPAyYGiuVZmZmZlZprrT6Hs2IlqALZIagMeAcfmWZVlonLqBWQuWcdHCeznqxEerKt+1F5Pv2ovJd+3F5Lv2YvLzrt061p1G32JJw4CfkczovQ34R55F7ShJzZKWljz2lPT3dN9USX8sI3tPSXd3sG++pMYeZL5d0r49rak76uqCE2Y8whnHTOC4qRM5bNo6xu+9sSryXXsx+a69mHzXXky+ay8mP+/arXNdNvoi4lMRsS4iLgReD3w4vczbmzwbEZNLHisi4pAdCZDUL6/i2vF2kjGSuZm4/zOsXDGA1Q8NZEtTHfOvHsaUI9ZXRb5rLybftReT79qLyXftxeTnXbt1rsNGn6QD2j6A4UD/9HmvJumpkpcNkq6VtFzShZLqWt8j6RxJdwBTJJ0i6e70cXLJ5/tL+pWkeyX9VtLgdo73U0mLJd0j6esl28+W9C9Jd0r6vqRDgLcB30t7JV8g6aSS98zJ4usfsVsTj68csPX1mlX1jBzTlEV07vmuvZh8115MvmsvJt+1F5Ofd+1ZUlTmUUmdzd49p5N9Abw241rKsZOkpenzByLiHW32H0TSs/YgcB3wTuC3wM7Aooj4nKQDgY8AB5MsSL1I0o3A/4CJwMciYqGkXwCfAr7f5hinR8QTaY/hPEn7AY8A7wD2iYiQNCwi1km6BvhjRPwWIL3jyYSI2JReSjczMzPLVGeLMx9WyULK9GxETO5k/60RcT+ApMuBQ0kafc3Alel7DgWuioin0/f9DngVcA3wcEQsTN93GXAS2zf6jpI0neScjiFpZP4L2Aj8PB1X2NHYwjuBX0n6PfD79t6QZk8HGMR2HY3bWbu6nlG7b976euSYJtasqu/yc92VZ75rLybftReT79qLyXftxeTnXXumavA2bN2ZyFEL2nagtr7eGBHNZXweAEkTgFOBwyNiP+BaYFBEbCHpZfwt8BaSXsb2vBm4ADgA+Kek7RrjETEzIhojorGegV0WvHzpYMZO2MzocZvoX9/C1GnruGVudivt5Jnv2ovJd+3F5Lv2YvJdezH5eddunevWHTlqwEFpw+xB4L3AzHbeswC4WNLZJJd33wF8MN03XtKUiPgH8H6SxalLNQBPA+sljSZZwHq+pCHA4Ij4k6SFwP3p+58EdgFIxxeOi4i/SboZeB8wBFhXzhfc0iwuOH0sM2bfT10/mDtnOA/eN6jrD/aCfNdeTL5rLybftReT79qLyc+79swENbk4syKq/6uS9FREDGlvm6SpwDdIGlovBP4GfCoiWtp+TtIpwEfTl7Mi4keS9iTpoVsMHEhyyfaDEfGMpPnAqRGxWNLFwCHAw8B6ksvC1wNXk9y+TsD3I+ISSa8kWQJnE0kj7+ckC14LuCwizu7s623Q8DhYh+/4iTIzM6tCi2IeG+KJil1vHThuXIz93GcrcqwHPvu5JRGxw8u/9UR3bsMm4Bhgr4j4hqTxwG4RcWvu1XVT2wZf6baImA+8ujufi4gfAD9os20FsE8Hn59a8vzYDso7qJ3PLWTbJVsO7eCzZmZmVoTq7xPbTnfG9P0EmAIcnb5+kmT8mZmZmZlVie6M6Ts4Ig6QdDtARPxP0oCuPmRmZmZWrSq9hl4ldKenryldey4AJI0CWnKtyszMzMwy1Z1G37nAVcDzJJ1FMnN1Rq5VmZmZmRUpKvSooC4v70bEryQtAQ4nmV369oi4N/fKzMzMzCwz3Zm9Ox54BvhD6baIeCjPwszMzMwKU4Nj+rozkeNaki9dJOvNTQCWAy/JsS4zs6rXf9zzc83f8vB/c803s9rSncu7Ly19LekA4FO5VWRmZmZmmdvh27BFxG2SDs6jGDMzM7OiKWpzyZbujOk7peRlHXAAsDK3iszMzMwsc93p6dul5PkWkjF+V+ZTjpmZmVkvEBW71W/FdNroSxdl3iUiTq1QPWZmZmaWgw4bfZL6R8QWSa+sZEFmZmZmhavBMX2d3ZHj1vS/SyVdI+mDkt7Z+qhEcVaexqkbmLVgGRctvJejTny0qvJdezH5rr3y+SOf9yzf/sk/+OmcG/nJ5Tfytvc+kFl2K5/3ymfnne/arSe6cxu2QcBa4LXAW4C3pv+1HSTpqfS/e0oKSd8q2TdSUpOk87M4Vl1dcMKMRzjjmAkcN3Uih01bx/i9N2YRnXu+ay8m37UXk9/cLGb9eF8++b7X8LmPvZK3vPtBxk14MpNs8HkvIjvvfNdeGa0zePN+VFJnjb7npTN37wbuSv97T/rfuytQW617AHhzyev3kJzfTEzc/xlWrhjA6ocGsqWpjvlXD2PKEeuzis8137UXk+/ai8n/39pB/L/lQwF49pn+PLxiCCNGZfdL0Oe98tl557t266nOGn39gCHpY5eS560PK88zwL2SGtPX7wV+k1X4iN2aeHzlgK2v16yqZ+SYpqzic8137cXku/bi8ls9b8wz7PWi9Sy/Z1hmmT7vlc/OO9+1V0hU6FFBnc3eXRUR36hYJX3THOB9kh4FmknWP9y92JLMrAiDdtrC6Wcv4Wc/3Jdnn64vuhwzq0GdNfpqb4Ga3uc64JvAo8CvO3ujpOnAdIBBDO4yeO3qekbtvnnr65FjmlizKrtfJHnmu/Zi8l17cfn9+rXw5bOX8LfrxvL3+WMyywWf9yKy88537RVQo3fk6Ozy7uEVq6KPiojNwBLgc8Bvu3jvzIhojIjGegZ2mb186WDGTtjM6HGb6F/fwtRp67hl7tBsCs8537UXk+/ai8oPPnPGnTy8Ygi/v3yvjDKf4/Ne+ey881279VSHPX0R8UQlC+nDzgFujIgnpOw6V1uaxQWnj2XG7Pup6wdz5wznwfsGVUW+ay8m37UXk7/vy/7H4W96hAf+vQvn/XIBAJf8dCKL//68TPJ93iufnXe+a6+QGuzpU0QNflW9lKSnImKIpD2BP0bEpDb7jwUaI+LEznIaNDwOljtizXq7/uOen2v+lof/m2u+WW+xKOaxIZ6o2LCzQWPHxfhPnVKRY/37jFOWRERj1+8sX3fuvWsZiYgh6X9XAJPa2X8xcHFFizIzM7Pt1WCfWHcWZzYzMzOzKudGn5mZmVkf4Mu7ZmZmZm30tSVbzMzMzKxGuNFnZmZm1ge40WdmZmbWB3hMn5mZmVlbNTimz40+M7OcePFkM+tN3OgzMzMzKxWevWtmZmZmVco9fWZmZmZtuafPzMzMzKqRe/rMzMzM2nJPn5mZmZlVkqRfSHpM0t0l24ZL+oukf6f/3bWrHDf6aljj1A3MWrCMixbey1EnPlpV+a69mHzXXky+ay8m37UXk5937VkQyezdSjy66WLgjW22fRGYFxF7A/PS153KrdEnKSRdVvK6v6THJf0xw2MMk/SpktdTdzRf0sWS3p0+nyVp3wzre3uWeTuiri44YcYjnHHMBI6bOpHDpq1j/N4bqyLftReT79qLyXftxeS79mLy8669VkXETcATbTZPAy5Jn18CvL2rnDx7+p4GJknaKX39euCRjI8xDPhUV2/qroj4eET8K6s8kv8BhTT6Ju7/DCtXDGD1QwPZ0lTH/KuHMeWI9VWR79qLyXftxeS79mLyXXsx+XnXnqmo0KPnRkfEqvT5amB0Vx/I+/Lun4A3p8+PBi5v3SFp5/Qa9a2Sbpc0Ld3+knTbUkl3Sto73X6KpLvTx8lpzNnAC9L3fi/dNkTSbyUtk/QrSUo/f6CkGyUtkXS9pDFti5U0X1Jj+vxoSXelx/tOyXveKOk2SXdImiepLr2ePirdXyfpP5JeA7wN+F5a3wvSx3VpDQsk7ZN+5j3pce6QdFMWJ37Ebk08vnLA1tdrVtUzckxTFtG557v2YvJdezH5rr2YfNdeTH7etVepkZIWlzym72hARHSrCZn37N05wFfTS677Ab8AXpXuOx34a0R8VNIw4FZJNwDHAz+OiF9JGgD0k3Qg8BHgYJJL7Ysk3Uhy/XpSREyG5PIusD/wEmAlsBB4paRFwHnAtIh4XNJ7gbOAj7ZXtKTdge8ABwL/A+ZKenua9zPg1RHxgKThEdGSXsY+BvgR8Drgjoi4UdI1wB8j4rdp7jzg+Ij4t6SDgZ8ArwW+ChwREY+k58LMzMz6hjUR0diDzz0qaUxErEo7sh7r6gO5Nvoi4k5Je5L08v2pze43AG+TdGr6ehAwHvgHcLqk5wO/SxtIhwJXRcTTAJJ+R9J4vKadw94aEf9N37cU2BNYB0wC/pJ2/PUDVrXz2VYvB+ZHxONpzq+AVwPNwE0R8UD69bVeX/8FcDVJo++jwEVtAyUNAQ4BrkhrABiY/nchcLGk3wC/a6+gtOU/HWAQgzspPbF2dT2jdt+89fXIMU2sWVXf5ee6K898115MvmsvJt+1F5Pv2ovJz7v2zFTHbdiuAT5MctXzwyTtkE5VYvbuNcD3Kbm0mxLwroiYnD7GR8S9ETGb5LLos8CfJL12B4+3qeR5M0nDVsA9Jcd6aUS8oWdfzvYi4mGSFvdrgYOAP7fztjpgXUkNkyPixennjwfOAMYBSySNaOcYMyOiMSIa67e2FTu2fOlgxk7YzOhxm+hf38LUaeu4Ze7QMr7KyuW79mLyXXsx+a69mHzXXkx+3rXXKkmXk3SKTZT0X0kfI2nsvV7Sv0muMp7dVU4lFmf+BUlj56708mur64FPS/p0RISk/SPidkl7AfdHxLmSxpNcFr6JpCfsbJIG3DuADwJPArt0o4blwChJUyLiH5LqgRdFxD0dvP9W4FxJI0ku7x5Ncnn4FuAnkiaUXN5t7e2bBVwG/DIimtNtW+uLiA2SHpD0noi4Ih1ruF9E3CHpBRGxiOSy9ZEkjb+13fi6OtTSLC44fSwzZt9PXT+YO2c4D943qJzIiuW79mLyXXsx+a69mHzXXkx+3rVnqhf19EXE0R3sOnxHcpSM/cuepKciYkibbVOBUyPiLems3h+RXPKsAx5It3+RpEHXRDIb5f0R8YSkU3huDN6siPhRmjmbpGH4Z+Da1vx03/nA4oi4WNJk4FxgKElj90cR8TNJF5OOu5M0P/38YklHA18maWReGxFfSDOPBGakNT8WEa9Pt9eTNNQOiohl6bZXkowB3AS8G2gBfgqMAeqBORHxjfRy9d7pseYBJ0cn/2MaNDwO1g79fzYzM6tai2IeG+IJdf3ObOw0ZlxM+OgpFTnWvTNOWdLDMX07LLdGX1+Tzvr9YUS8qss3l8mNPjMz60sKafR9pEKNvm9XrtHne+9mIO2d/CTJDF4zMzOzXseNvgxExNl0YwClmZmZVYcqmL27w3zvXTMzM7M+wD19ZmZmZm25p8/MzMzMqpF7+szMzMxKdetOttXHjT4zM9tOv4aGXPObN2zINd/MtudGn5mZmVkbnr1rZmZmZlXJPX1mZmZmbbmnz8zMzMyqkRt9ZmZmZn2AL++amZmZteGJHFZVGqduYNaCZVy08F6OOvHRqsp37cXku/Zi8qu59pPPuo/ZC2/hJ9csyTS3lc975bPzzs+7dutYTTT6JIWky0pe95f0uKQ/FlzXsZJ2L+LYdXXBCTMe4YxjJnDc1IkcNm0d4/feWBX5rr2YfNdeTH411w5ww1Wj+cpxkzLLK+XzXvnsvPPzrj1TUaFHBdVEow94Gpgkaaf09euBRwqsp9WxQCGNvon7P8PKFQNY/dBAtjTVMf/qYUw5Yn1V5Lv2YvJdezH51Vw7wN2Lh/Lk+nxGCvm8Vz477/y8a7fO1UqjD+BPwJvT50cDl7fukDRc0u8l3SnpFkn7pdvPlPQLSfMl3S/ppJLPnCLp7vRxcsn2D6U5d0j6paRdJD0gqT7d35C+fg/QCPxK0lJJO0k6UNKNkpZIul7SmPQzJ0n6V5o7J4uTMWK3Jh5fOWDr6zWr6hk5pimL6NzzXXsx+a69mPxqrj1vPu+Vz847v2q+HyvVy+eevh6bA7xP0iBgP2BRyb6vA7dHxH7Al4FLS/btAxwBHAR8TVK9pAOBjwAHA68AjpO0v6SXAGcAr42IlwGfiYgngfk81+B8H/C7iLgCWAwcExGTgS3AecC7I+JA4BfAWelnvgjsn9Z3fFYnxMzMzKxVzczejYg7Je1J0sv3pza7DwXelb7vr5JGSGq9seS1EbEJ2CTpMWB0+v6rIuJpAEm/A15F0ia/IiLWpFlPpBmzgNOA35M0Fo9rp8SJwCTgL5IA+gGr0n13kvQI/j7N2I6k6cB0gEEM7vJ8rF1dz6jdN299PXJME2tW1Xf5ue7KM9+1F5Pv2ovJr+ba8+bzXvnsvPOr5ftR6aPW1FJPH8A1wPcpubTbDZtKnjfTg4ZwRCwE9pQ0FegXEXe38zYB90TE5PTx0oh4Q7rvzcAFwAHAPyVtV0NEzIyIxohorGdglzUtXzqYsRM2M3rcJvrXtzB12jpumTt0R7+0QvJdezH5rr2Y/GquPW8+75XPzju/mr8fa0HN9PSlfgGsi4i70gZYqwXAMcA30+1rImJD2uPWngXAxZLOJmmsvQP4ILAZuErSDyJiraThJb19lwKzgW+W5DwJ7JI+Xw6MkjQlIv6RjgF8EXAvMC4i/ibpZpLLw0OAdT08BwC0NIsLTh/LjNn3U9cP5s4ZzoP3DSonsmL5rr2YfNdeTH411w5w2jnL2O/l62jYdQuXzl/EZeftwdwrd8sk2+e98tl55+dde6ZqcJ0+RVT/VyXpqYgY0mbbVODUiHiLpOEkDcK9gGeA6enl4DOBpyLi++ln7gbeEhErJJ0CfDSNmxURP0rf82Hg8yS9grdHxLHp9t2AB4AxEbEu3fYuYAbwLDCF5BLvucBQkgb3j4CLgb+l2wRcFhFnd/b1Nmh4HKzDd/AsmZl1X7+Ghq7fVIbmDRtyzbfasijmsSGeqNgV18Gjx8UL339KRY51149OWRIRjZU4Vk00+noDSe8GpkXEB/M+lht9ZpY3N/qsNymi0bf30ZVp9N3548o1+mrt8m4hJJ0HHAm8qehazMzMzNrjRl8GIuLTRddgZmZmGarBC6G1NnvXzMzMzNrhnj4zMzOzttzTZ2ZmZmbVyI0+MzMzsz7Al3fNzMzMSgWoBi/vutFnZmbbyXsdvTzXAfQagGbtc6PPzMzMrK0a7OnzmD4zMzOzPsA9fWZmZmZt1OKYPvf0mZmZmfUB7ukzMzMza8s9fWZmZmZWjdzoq2GNUzcwa8EyLlp4L0ed+GhV5bv2YvJdezH5rr19J591H7MX3sJPrlmSaW4rn/di8vOuPSuKyjwqyY2+GlVXF5ww4xHOOGYCx02dyGHT1jF+741Vke/ai8l37cXku/aO3XDVaL5y3KTM8kr5vBeTn3ft1rmaafRJCknnlLw+VdKZ6fMzJZ2aPh8k6S+t+9pkjJC0NH2slvRIyesBGdR4rKTz29aUh4n7P8PKFQNY/dBAtjTVMf/qYUw5Yn1V5Lv2YvJdezH5rr1jdy8eypPr8xl67vNeTH7etWcmKviooJpp9AGbgHdKGtnRG9KG25XAkog4s+3+iFgbEZMjYjJwIfDD1tcRsTmnunMxYrcmHl/5XDt1zap6Ro5pqop8115MvmsvJt+1F8PnvZj8av6eqQW11OjbAswEPtvB/v7Ar4F/R8QXszqopA9JulPSHZJ+mW4bJelKSf9MH6/sIuMkSf9Kc+ZkVZuZmZn1UA329NXaki0XAHdK+m47+04D/hIRJ2d1MEkvAc4ADomINZKGp7t+TNJLeLOk8cD1wIs7ifoiMCEiNkka1sGxpgPTAQYxuMva1q6uZ9Tuz3VOjhzTxJpV9V1/Ud2UZ75rLybftReT79qL4fNeTH41f8/Uglrq6SMiNgCXAie1s/tm4BBJL8rwkK8FroiINenxn0i3vw44X9JS4BqgQdKQTnLuBH4l6QMkPZbbiYiZEdEYEY31DOyysOVLBzN2wmZGj9tE//oWpk5bxy1zh3b/Kysw37UXk+/ai8l37cXweS8mv1q+Z0Rtzt6ttZ4+gB8BtwEXtdl+E3AJ8GdJh0bEqhxrqANeERHbTEmS1NH73wy8GngrcLqkl0ZEu42/7mppFhecPpYZs++nrh/MnTOcB+8bVE5kxfJdezH5rr2YfNfesdPOWcZ+L19Hw65buHT+Ii47bw/mXrlbJtk+78Xk5127dU4RtbHktKSnImJI+vy7wPuAX0TEmelM3aci4vuSPgkcD7wmItZ1krf1M5285yXAVcCUiFgraXhEPCFpNnB7RHwvfd/kiFgq6VigMSJObM0HfgCMj4gVkuqBB4F9O6utQcPjYB3ezTNjZtb79GtoyC27ecOG3LKtGItiHhviiQ57TrK286hxsc87TqnIsW772SlLIqKxEseqqcu7Jc4B2p3FGxE/JWmoXSOprD8vIuIe4CzgRkl3kDTgILm83JhOzPgXSSOzI/2AyyTdBdwOnNtZg8/MzMwqoAYnctRMT19f4p4+M6t27umzHVFIT9/bK9TTN6tyPX21OKbPzMzMrCyqwU6xPtvokzQCmNfOrsMjYm1P32tmZmbWG/XZRl/aWJuc9XvNzMysyhUw3q4SanUih5mZmZmV6LM9fWZmZmYdqfTCyZXgRp+ZmdWUfi+ZmGt+8z3Lc803y4sbfWZmZmZt1WBPn8f0mZmZmfUB7ukzMzMza6MWx/S5p8/MzMysD3BPn5mZmVlb7ukzMzMzs2rkRl8Na5y6gVkLlnHRwns56sRHqyrftReT79qLyXft7Tv5rPuYvfAWfnLNkkxzAerrm/nheTdw/oXX89OfXccxH7o782NU63nPOz/v2jMRyZi+SjwqyY2+GlVXF5ww4xHOOGYCx02dyGHT1jF+741Vke/ai8l37cXku/aO3XDVaL5y3KTM8ko1NdXxpc+/hhOPP4ITj38DjY2rmfji7G6lXs3nvZprt871yUafpJB0TsnrUyWdmT4/U9Kp6fNBkv7Sui+D4x4r6fy2x8nDxP2fYeWKAax+aCBbmuqYf/UwphyxviryXXsx+a69mHzX3rG7Fw/lyfV5DT0XGzfWA9C/fwv9+rdkOoarms97NdduneuTjT5gE/BOSSM7eoOkAcCVwJKIOLNShWVlxG5NPL5ywNbXa1bVM3JMU1Xku/Zi8l17MfmuvTh1dS2cd+FcZl9xDbffNprly0Zkll3N572aa89UVOhRQX210bcFmAl8toP9/YFfA/+OiC92FiTpQ5LulHSHpF+m20ZJulLSP9PHK7vIOEnSv9KcOT34eszMbAe1tNTx6ePfwIeOfgsvmvgEe+zpHierbX15yZYLgDslfbedfacBf4mIkzsLkPQS4AzgkIhYI2l4uuvHwA8j4mZJ44HrgRd3EvVFYEJEbJI0rINjTQemAwxicGdlAbB2dT2jdt+89fXIMU2sWVXf5ee6K898115MvmsvJt+1F+/ppwdw5x3P48DGVTy4YmgmmdV83qu59qwIL85cUyJiA3ApcFI7u28GDpH0oi5iXgtcERFr0swn0u2vA86XtBS4BmiQNKSTnDuBX0n6AEkvZHv1zoyIxohorGdgF2XB8qWDGTthM6PHbaJ/fQtTp63jlrnZ/DDLO9+1F5Pv2ovJd+3FaBi6kZ13ThofAwZsYf8DHuW/Dzdkll/N572aa7fO9eWePoAfAbcBF7XZfhNwCfBnSYdGxKodzK0DXhER20xJktTR+98MvBp4K3C6pJdGRLuNv+5qaRYXnD6WGbPvp64fzJ0znAfvG1ROZMXyXXsx+a69mHzX3rHTzlnGfi9fR8OuW7h0/iIuO28P5l65WybZw4dv5HOn3UpdXSAFC24ax62Lds8kG6r7vFdz7ZmK2uvqU9TgF9UVSU9FxJD0+XeB9wG/iIgz05m6T0XE9yV9EjgeeE1ErGsn5yXAVcCUiFgraXhEPCFpNnB7RHwvfd/kiFgq6VigMSJObD0O8ANgfESskFQPPAjs297xWjVoeByswzM6G2ZmldevIbtete2MG5NfNtB8z/Jc8217i2IeG+KJDntOsjZkxLiY9MaTK3KsRbNPXRIRjZU4Vp+9vFviHKDdWbwR8VOSRt01krb7UyQi7gHOAm6UdAdJAw6SS8aN6cSMf5E0HDvSD7hM0l3A7cC5nTX4zMzMLH+1uDhzn7y829rLlz5/FJ6bGdF2eZb09Tbb2uy/hORScOm2NcB723nvxcDF7Rzn0O7WbmZmZtYTfbLRZ2ZmZtahAtbQqwQ3+rpB0ghgXju7Do+I7O7bY2ZmZpYTN/q6IW3YTS66DjMzM6sMtRRdQfY8kcPMzMysD3BPn5mZmVlbNTimzz19ZmZmZn2Ae/rMzKzimjdsyC/8nhyzgX4vmZhbthd+7j18710zMzMzq0pu9JmZmZn1Ab68a2ZmZlYqgKi967vu6TMzMzPrA9zTZ2ZmZtaGJ3KYmZmZWVVyo6+GNU7dwKwFy7ho4b0cdeKjVZXv2ovJd+3F5Lv2YvLzzK6vb+aH593A+Rdez09/dh3HfOjuTPN93isgKvSooKpt9EkKSeeUvD5V0pk9zBom6VMZ1ranpHb/hUv6hqTXpc/nS2rM6ril6uqCE2Y8whnHTOC4qRM5bNo6xu+9sSryXXsx+a69mHzXXkx+3rU3NdXxpc+/hhOPP4ITj38DjY2rmfjitZlk+7xbT1Vtow/YBLxT0sgMsoYB7Tb6JGU67jEivhoRN2SZ2Z6J+z/DyhUDWP3QQLY01TH/6mFMOWJ9VeS79mLyXXsx+a69mPy8awexcWM9AP37t9Cvf0tmvTo+7/kTyZi+SjwqqZobfVuAmcBn2+6QNErSlZL+mT5emW4/U9KpJe+7W9KewNnACyQtlfQ9SVMlLZB0DfAvSYMkXSTpLkm3Szos/fyxkq5Oe+z+LelrJWX0k/QzSfdImitpp/QzF0t6d5t6+6Xb706Psd3XtKNG7NbE4ysHbH29ZlU9I8c0lRtbkXzXXky+ay8m37UXk5937QB1dS2cd+FcZl9xDbffNprly0Zkkuvzbj1V7bN3LwDulPTdNtt/DPwwIm6WNB64HnhxJzlfBCZFxGQASVOBA9JtD0j6HBAR8VJJ+wBzJb0o/exBwCTgGeCfkq4F1gB7A0dHxHGSfgO8C7isg+NPBsZGxKT0+MO6+fWbmVkv1dJSx6ePfwM777yZM85cyB57rufBFUOLLsu6I8Lr9PU2EbEBuBQ4qc2u1wHnS1oKXAM0SBqyg/G3RsQD6fNDSRtsEbEMeBBobfT9JSLWRsSzwO/S9wI8EBFL0+dLgD07Odb9wF6SzpP0RmC7G0dKmi5psaTFTWzqsvi1q+sZtfvmra9Hjmlizar6Lj/XXXnmu/Zi8l17MfmuvZj8vGsv9fTTA7jzjudxYOOqTPJ83q2nqrrRl/oR8DFg55JtdcArImJy+hgbEU+RXBIu/ZoHdZL7dDeP3/ZPgdbXpS2zZjrpVY2I/wEvA+YDxwOz2nnPzIhojIjGegZ2WdTypYMZO2Ezo8dton99C1OnreOWudn9hZlnvmsvJt+1F5Pv2ovJz7v2hqEb2XnnpHEzYMAW9j/gUf77cEMm2T7vlVGLY/qq/fIuEfFEevn0Y8Av0s1zgU8D3wOQNDntdVsBvCXddgAwIX3/k8AunRxmAXAM8Nf0su54YDnJJeDXSxoOPAu8Hfjojn4N6WSUzRFxpaTldHwZuNtamsUFp49lxuz7qesHc+cM58H7Omvj9p58115MvmsvJt+1F5Ofd+3Dh2/kc6fdSl1dIAULbhrHrYt2zyTb5916SlGl16wlPRURQ9Lno4EHgO9GxJlpI+oCknF8/YGbIuL4dDLF1cBYYBEwBTgyIlZImg3sB/wZuBY4NSJaG4iDgJ8CjSS9hadExN8kHUvS0BsKPB+4LCK+nk4O+WPJGL1TgSFpbRen+34raT5wKtAEXMRzvZBfiog/d/S1N2h4HKzDyzuBZmbWI/1eMjG37OZ7lueWXc0WxTw2xBOq1PF2Gfb82P/Vn6nIsRb84bQlEZHL8m1tVW1PX2uDL33+KDC45PUa4L3tfOZZ4A0d5L2/zab5Jfs2Ah/poJT/RsTb22StIJnc0fr6+yXPjy15PrXkYwd0kG9mZmZWtqpt9JmZmZnlpRbvvetGXxki4mLg4oLLMDMzM+uSG31mZmZmpQJoqb2uvlpYssXMzMzMuuBGn5mZmVkf4Mu7ZmZmZm3V3tVdN/rMzMx2RJ5r6fVryOauHR1p3rDdXT6tD3Gjz8zMzKyNWlyyxWP6zMzMzPoA9/SZmZmZtdWLblMraQXwJNAMbOnpbdvc6DMzMzPr/Q5LbzPbY270mZmZmbXhMX1mZmZmlqWRkhaXPKa3854A5kpa0sH+bnGjr4Y1Tt3ArAXLuGjhvRx14qNVle/ai8l37cXku/Zi8qu59pPPuo/ZC2/hJ9csyTS3VTWf90xEBR+wJiIaSx4z26no0Ig4ADgSOEHSq3vyZfW5Rp+kZklLJd0t6QpJg9PtT3XxuWGSPtXNY2SW1VN1dcEJMx7hjGMmcNzUiRw2bR3j995YFfmuvZh8115MvmsvJr+aawe44arRfOW4SZnllarm816rIuKR9L+PAVcBB/Ukp881+oBnI2JyREwCNgPHd/Nzw4CsGmpZZrVr4v7PsHLFAFY/NJAtTXXMv3oYU45YXxX5rr2YfNdeTL5rLya/mmsHuHvxUJ5cn8+w/Go+71kRoIiKPLqsRdpZ0i6tz4E3AHf35Ovqi42+UguAF5ZukDRE0jxJt0m6S9K0dNfZwAvSXsLvpe/9vKR/SrpT0tfbO0AH79kmS9IYSTeV9EC+qtwvbMRuTTy+csDW12tW1TNyTFO5sRXJd+3F5Lv2YvJdezH51Vx73qr5vNeo0cDNku4AbgWujYjrehLUZ2fvSupPcm287YnbCLwjIjZIGgncIuka4IvApIiYnH7+DcDeJF2sAq6R9OqIuKnkGO2+p52szwHXR8RZkvoBg3P6ss3MzKw7WoouIBER9wMvyyKrLzb6dpK0NH2+APh5m/0CZqSNsxZgLEkru603pI/b09dDSBp4N3XjPQ+1yfon8AtJ9cDvI2Jpm/2ks3WmAwzqRptw7ep6Ru2+eevrkWOaWLOqvsvPdVee+a69mHzXXky+ay8mv5prz1s1n3frXF+8vNs6pm9yRHw6Ija32X8MMAo4MO2JexQY1E6OgG+XZL0wItprQHb1HtLewVcDjwAXS/pQO++Z2Tqzp56BXX6Ry5cOZuyEzYwet4n+9S1MnbaOW+YO7fJz3ZVnvmsvJt+1F5Pv2ovJr+ba81bN5z1LvWVMX5b6Yk9fV4YCj0VEk6TDgD3S7U8Cu5S873rgm5J+FRFPSRoLNKUzazp9T9ssSXsA/42In0kaCBwAXFrOF9HSLC44fSwzZt9PXT+YO2c4D97XXtu19+W79mLyXXsx+a69mPxqrh3gtHOWsd/L19Gw6xYunb+Iy87bg7lX7pZJdjWfd+ucohfdW64SJD0VEUM62p6O4/sDyaXYxcArgCMjYoWk2cB+wJ8j4vOSPgN8PI14CvhARPy/0mN08p6tWSSzcD5P0iB8CvhQRDzQ0dfQoOFxsA4v91SYmVkv06+hIdf85g0bcs3Py6KYx4Z4QpU6XkPD8+PljSdU5Fh//duXl/T0Xro7qs/19LXX4Cvdnt7XbkoH73l/m9c/Bn7c2TE6ec/722y6pKvazczMrAKeWzi5pvTFMX1mZmZmfU6f6+kzMzMz61xADQ5/c0+fmZmZWR/gnj4zMzOzNlR7HX3u6TMzMzPrC9zTZ2ZmZtaWx/SZmZmZWTVyT5+ZmVkvUa2LJ9ecALUUXUT23NNnZmZm1ge4p8/MzMysLY/pMzMzM7Nq5J4+MzMzs7Zqr6PPPX1mZmZmfYF7+szMzMzakMf0WTVpnLqBWQuWcdHCeznqxEerKt+1F5Pv2ovJd+3F5Lv2YvLzrt065kZfRiSdLukeSXdKWirpYEnzJT0kSSXv+72kp0pev0jSnyT9W9Jtkn4jaXS59dTVBSfMeIQzjpnAcVMncti0dYzfe2O5sRXJd+3F5Lv2YvJdezH5rr2Y/Lxrt8650ZcBSVOAtwAHRMR+wOuAh9Pd64BXpu8bBowp+dwg4FrgpxGxd0QcAPwEGFVuTRP3f4aVKwaw+qGBbGmqY/7Vw5hyxPpyYyuS79qLyXftxeS79mLyXXsx+XnXnqmIyjwqyI2+bIwB1kTEJoCIWBMRK9N9c4D3pc/fCfyu5HPvB/4REX9o3RAR8yPi7nILGrFbE4+vHLD19ZpV9Ywc01RubEXyXXsx+a69mHzXXky+ay8mP+/arXNu9GVjLjBO0n2SfiLpNSX75gGvltSPpPH365J9k4AlFazTzMzMuhJAS4UeFeRGXwYi4ingQGA68Djwa0nHprubgZtJGnw7RcSKnhxD0nRJiyUtbmJTl+9fu7qeUbtv3vp65Jgm1qyq78mhK57v2ovJd+3F5Lv2YvJdezH5eddunXOjLyMR0Zxemv0acCLwrpLdc4Bzgd+0+dg9JI3F7uTPjIjGiGisZ2CX71++dDBjJ2xm9LhN9K9vYeq0ddwyd2j3vphuyDPftReT79qLyXftxeS79mLy8649KyJQVOZRSV6nLwOSJgItEfHvdNNk4EGSy7cAC4BvA5e3+ehs4EuS3hwR16ZZrwaeKHdcX0uzuOD0scyYfT91/WDunOE8eN+gciIrlu/ai8l37cXku/Zi8l17Mfl5126dU9Tg4oOVJulA4DxgGLAF+A/Jpd7fAqdGxOI2738qIoakz/cBfgS8AGgC7gQ+ExEdLl7UoOFxsA7P/gsxMzPrhRbFPDbEE+r6ndkYuvPu8YoXT6/IseYu+fqSiGisxLHc05eBiFgCHNLOrqkdvH9IyfNlwBvzqczMzMws4UafmZmZWVs1eCXUEznMzMzM+gD39JmZmZmVal2nr8a4p8/MzMysD3BPn5mZmVkblV5DrxLc02dmZmbWB7inz8zMakq/hoaiS+ix5g0bii6hx65fuTS37IOOeCa37A65p8/MzMzMqpEbfWZmZmZ9gC/vmpmZmW0jfHnXzMzMzKqTe/rMzMzMSgXu6TMzMzOz6uRGXw1rnLqBWQuWcdHCeznqxEerKt+1F5Pv2ovJd+3F5J981n3MXngLP7lmSaa5eWdDdZ33xx6p5/PvfgHHvWYfjps6katmjdxm/28vHMURu09m/dp+ZR0ncy0VelSQG301qq4uOGHGI5xxzASOmzqRw6atY/zeG6si37UXk+/ai8l37cXl33DVaL5y3KTM8iqVXW3nvV//YPpXV/KzG5fx4z/+mz9cPJIH7xsIJA3C227cheeN3ZxV+dYJN/rakDRC0tL0sVrSIyWvB5SZfayk89PnZ0o6NZuqtzdx/2dYuWIAqx8ayJamOuZfPYwpR6yvinzXXky+ay8m37UXl3/34qE8uT6foe15ZlfbeR8xegt77/csAIOHtDDuhZtYs6oegP87cywfO2MlUialZ0oRFXlUkht9bUTE2oiYHBGTgQuBH7a+joiq+VNkxG5NPL7yuTbqmlX1jBzTVBX5rr2YfNdeTL5rLy6/WlXzeV/98AD+3907sc8Bz/D36xoYuVsTL3hJdr2U1jk3+sok6UOS7pR0h6RfpttGSbpS0j/Txyu7yDhJ0r/SnDmVqdzMzKxynn26jm9+fE+O/8Yj9OsXzDlvNB/6/Kqiy+pYRGUeFeQlW8og6SXAGcAhEbFG0vB0149JeghvljQeuB54cSdRXwQmRMQmScM6ONZ0YDrAIAZ3Wdva1fWM2v25jsmRY5q2dqdnIc98115MvmsvJt+1F5dfrarxvG9pgm9+fE9e+87/ceib1vPAvYNY/dAAPvm6fQB4fFU9JxwxkXP/dB/Dn7elrGNZx9zTV57XAldExBqAiHgi3f464HxJS4FrgAZJQzrJuRP4laQPAO1+t0fEzIhojIjGegZ2WdjypYMZO2Ezo8dton99C1OnreOWuUO7/5UVmO/ai8l37cXku/bi8qtVtZ33CPjB58Yzbu9NvOsTjwMw4cUb+c1d93Dprf/i0lv/xagxTVxw/fLe0+ALoCUq86gg9/Tlow54RURsM1BBHY9UfTPwauCtwOmSXhoRZX3ntzSLC04fy4zZ91PXD+bOGc6D9w0qJ7Ji+a69mHzXXky+ay8u/7RzlrHfy9fRsOsWLp2/iMvO24O5V+7W67Or7bzfc+vOzPvtcCa8+Fk++bqJAHzkSys56PAnsyrZuklRgytOZ0XSmcBTEfH9Dva/BLgKmBIRayUNj4gnJM0Gbo+I76XvmxwRSyUdCzRGxImt2cAPgPERsUJSPfAgsG9ErOuorgYNj4N1eHZfqJlZDenX0FB0CT3WvGFD0SX02PUrl+aWfdARD7P4jo0Vm+M7dNBuccj4D1fkWNf9+7tLIqKxEsdyT18ZIuIeSWcBN0pqBm4HjgVOAi6QdCfJOb4JOL6DmH7AZZKGAgLO7azBZ2ZmZtYTbvR1IiLO7MZ7LgEuabNtDfDedt57MXBxO9mH9rxKMzMzy1wNXgn1RA4zMzOzPsA9fd0gaQQwr51dh0fE2krXY2ZmZraj3OjrhrRhN7noOszMzKxCfHnXzMzMzKqRe/rMzMzMSrUuzlxj3NNnZmZm1ge4p68KPcn/1twQv31wBz4yEliTUzl5Zued79qLyXftxeT3ndrX55zfe7Lzzt+h7H5jcs3fY4fTyxIQLZU9ZAW40VeFImLUjrxf0uK8VvvOMzvvfNdeTL5rLybftReT79qLy7ftudFnZmZm1pZn75qZmZlZNXJPX98ws0qz88537cXku/Zi8l17Mfmuvbj8nqvR2buKGuy+NDMzM+upoQNGxyG7HV2RY1338I+XVGpso3v6zMzMzNqqwU4xj+kzMzMz6wPc02dmZmYVJakOGBIRG4qupUPu6TOzzkjaOf1hhqQXSXqbpPoyM4d39sim8nxIGiRpu3UlJY2SNCiD/CMkvbud7e+W9Ppy89vJrZPUkHVuJUjaVdJ+Oebnem7yzK+2/6+SXilp5/T5ByT9QFJZixen/yb3bWf7vu39G+7hMWZLakhrvxv4l6TPZ5Ft3eNGXw2RdG5njxyP+9WMco6Q9DFJe7bZ/tEs8kvyvpv+4KmXNE/S45I+kFH8TcAgSWOBucAHgYvLzFwCLE7/+zhwH/Dv9PmSMrMBkDRC0nmSbpO0RNKPJY3IIPpc4FXtbD8U+GEG+V8Fbmxn+3zgGxnk5/aLqhL/XiXNT2sfDtwG/EzSD7LITvNz/SWeZ37O2Z9JsyXp5+m/qzdkkZ36KfCMpJcBnwP+H3BpmZnnkdwho60RwI/LzG61b9qz93bgz8AEkp+RvVAkPX2VeFSQG3215XiSX6Yrea6RUPrIy8fLDZA0AzgdeCkwT9KnS3afWG5+G29If/C8BVgBvBDI6heVIuIZ4J3ATyLiPcBLygmMiAkRsRdwA/DWiBgZESNI6p9bdsWJOcBjwLuAd5M0KH+dQe6BEfG7thsj4irg1RnkD4yIx9vJXwPsnEE+5PeLqhL/Xoemtb8TuDQiDgZel1E25P9LPM/8PLM/mma/Adg1zT07o2yALZEsvTENOD8iLgB2KTPzhRFxU9uNEbEAyKqHuD698vF24JqIaCJZHMUqxGP6assY4D3Ae4EtJL+0fxsR68oNltTRuAsBO5WbD7wV2D8itkg6E5gtaa+I+Gx6jCy1ft+/GbgiItZLmR1CkqYAxwAfS7f1yyj7FRFxXOuLiPizpO9mlD0mIr5Z8vpbkt6bQe7gTvZl8Udng6T+EbGldGP6iyWL70vY9hfV+RHRlNH3S27/Xkv0lzQGOIrkj6qstXdusvwlnmd+ntmt3yBvAn4ZEfcowx8ywJOSvgR8AHi1kiElZQ0jofNGY7nZrf6P5A/tO4Cb0kvSvXdMXw1yT18NiYi1EXFhRBwGfAQYRnLJIou/XtcBe0dEQ5vHLsCqDPK3/uJOf+m9leQX+hXAgAzyS/1R0jLgQJJexVHAxoyyTwa+BFyV/qDfC/hbRtkrJZ0hac/0cTpJL1EW5kp6Xzq2qU7SUcD1GeQ+JumgthslvZykN7FcvyO5ZLm1V0/SEODCdF8WWn9R7cxzv6jWlxua87/XVt8g+f/4/yLin+n3478zzG/v3GT5SzzP/Dyzl0iaS9Lou17SLkBLRtmQ/KGwCfhYRKwGng98r8zM/0h6U9uNko4E7i8zG4CIODcixkbEmyLxIHBYFtmZC6ClpTKPCvLizDVI0gHA0cDrSS4TnRMR/yoz81sk3fG3trPvOxHxhTLz/wh8LyJubLP9W8CXIyLTP1DSMU7rI6JZ0mCgIf3h2WulNX+N5LJokIwf/EZEPFFG5pNplkh++bX+BKoDnoqIsga3pw2+35CMa2y9ZNkIfAh4X0QsKjO/P/AtkiEGD6abxwM/B76SXj4qi6QJEfFAyWuRXArLpPGUx7/XIrXX81ot+Vllpz1vk4H7I2JdOj52bETcWW52mr/dz9xyfw5L2hu4Fvg72/5bnQK8JSLuKyP7lM72R0Rm40yzMrT+eXHIyPdU5FjXrf5JxRZndqOvhkj6Bskly3tJxmhdl9UPR0njI+KhLLI6yN8HeDAinm1n39iIeCSDY7yzs/3tjT3bgew/0MnYlIh4W0+z2znWzhHxdFZ5eZM0GvgUMCnddA/J5bTHMjzGTiRjMwH+0973URnZt0XEAW22LYmIA8vMze3fa8kxXkQy6H90RExSMnv3bRHxrYzyRwMzgN0j4kglsz+nRMTPy8zNrZGQc/YBne2PiNt6mt3mOO19T94ZEWWNvZM0EHg/2/5bnR0RZV0JkfS1zvZHxNfLyc/D0PrnxSEjtlsYIBfXPfpT35HDeuQM4AHgZeljRjqMREBLRLysjOzfA53+QCvT7LY/xFpl0eBLvbWTfUF5lwO/X8Znu0XSIcAsYAgwPp2594mI+FQG2e1OqmhvYPeOiohHga+ll9Fpb+JFT0k6LSK+GxHPStonIq4o2TcjIr5cRvY+JJNwhrb5g6EBKHu5GTr/9xrl/gJP/YxkktL/kYTeKWk2Se9oFi4GLuK58YL3kYxNLKvRR/mTEorKPqeTfQG8tpxwSZ8k+QNqL0mlvYa7kPTQlSUiNpH8/8xUb2zU9VVu9NWWCe1sEzCOZJxZObKeTFHpfCLiIzlmt7dsSNZ+CBwBXJMe846OGms9UDp7eRBwEMklnnJ/SYnkkvQJpBNaJDUD50VEFkuqvA9onczyJeCKkn1vBHrc6AMmksyQHsa2fzA8CRzX3gd2UHv/XrM2OCJubTOHIMvexJER8Zt0UgHpRKzmckPzbCTknJ33+LTZJDONvw18sWT7k+UM84Bthnpst4vkj5AeD/VQF0sQRcRJPc3OVQ1eCXWjr4akg2IBkLQ/STf9e0h6E64sM35sZ/9wM/hHm3d+3pd1fhMRR0m6i3Z+cGbUa0NEPNzmF3jZv2DT3G16QSWNA36UQfRngVcCB7WOi0snE/xU0mcjoty1+tTB8/Ze75CIuBq4WtKUiPhHOVkd5D/Y9bvKtkbSC0i/J5UsZJ3FxKtWT6fj1VrzX0EGk1zybCTknP3aiPhrR0NJyhlCkn5+Pcn5PVrSoSST6y6SNLLt2NMemAfsRnLF49cZf3/muWSY7QA3+mpIOn7n6PSxhuQyizL66/NZOv6Hm8WfQ3nnQ76XdT6T/vctOR7j4fQSbyhZauIzJOPB8vBf4MUZ5HwQeH0k6+YBEBH3K1kMey7lL9AcHTxv73VPvUPSPSTfo9eRrFn22Yi4rJzQPHtWSpwAzAT2kfQIyR+Ax2SQ2+oUkp7nF0haCIwi+UOzXHk2EvLMfg3wV9ofSlLuEJKt0jFyjSS90ReRrHBwGckfWD0SEW+XNJRkTceZSu6Y82tgTrm9iBFxSelrSYMjWc+0Fwtoqb2ePk/kqCGSWoAFJNP4/5Nuuz+ShX3Lzd5u4HC6/VUkszBP6M35tUDSSJKV8V9H0jCYC5xU7g/kNPs8nmuAtM48XBERZd2pRNLdETFpR/ftQH4z8DTPrRfZ+otEwKCIKHt9MUlLI2KypHeQNOpPAW4qc4xsRbT2/ihZ0qYuIp7MoEeoNH8gSW/zRJJzvjw9zqYs8kuOk1sjoToaINuStBTYH7gtIvZPt5U9kaMkv45k6MS5wIysZtcqWcP05yT33M10XHLWhtaPikOGvasix7puzf95Iof1yDtJ/qH+TdJ1JDMCsxort7n1SQ6XjiuRv1WeMxrTyzrfAZ5Hcu6z7LWZGBHb9NJIeiWwMIPsxSXPtwCXR0QWuZt7uK9bIiKrha8709pwzGMx77xdCRzQZrb3b0nWqMzCP9I/1u5p3SDpNjKa9FXaSCD7yUt5Zg8kubvNnpT8ns1oHCvA5ogIpYtJq2SdynKkVxKOJrl14s3AOyK5I0dWfkR+45KzFRBR2TX0KsGNvhoSEb8Hfp/+AJhGslDw8yT9lGSx4HJu2fXh9JJCHpeOK5FfKs8Zjd8luVVaHpddz2P7X6btbdthbS+/ZOhlav9uLiKbGbCV8Acli3k/C3xS2S7mnYu8Zx5L2g0YC+yU/pHW2gpuoPO7sOyoH5FfIyHP7KtJxt4tIVlEOWu/kfR/wDBJxwEfJfm51mOSHgT+R9JZMJ10wk/rMjRZLTeT17hk6x43+mpQ+lf9bJJbme1K0mP2Bcq7T+u9JJeO31Jy6fiz5dZawfxSec5ofDTrBl/aI3EIMKrNZJQGyrzFW0cTT1qVe7moQj1xuYqILyq53V3rYt5Pk/xR1ZvlPfP4COBYkjtBlF7620B5M6a3k2cjIcfs50fEGzPK2k5EfF/S60nO90TgqxHxlzJjHyD5WXBE+mj7c6GsmfypSo5LLl8Njulzo6/GRcT/SAZyzywzKs9Lx5XIL5X5jMaS3pTFkn5Nsq7h1r/wy5y1N4DkElR/tp2MsgEod/XQdwKjgYfbbB8H9Oo7lFTY7sDr0sHtrS4tqpiulMw8flXGl+da8y8BLpH0rojIdPhFG3k2EvLM/rukl0bEXRnlbSdt5JXb0Ct1GvBwRKwCkPRhkkvUK4AzMzrG8STjkscCj5B0RPT58dqV5IkctkNKLh0fTfKX36WUf+m4YvnpMfYiaQQfQnI54wHgAxGxoozMzhY0jYj4aE+zS46xR9bLfCi5/d2X2v5ykvRSkgHcnS1o3Sekww6mAvsCfwKOBG6OiMos118GSf8GlpLM8PxzZPwDX9L9JOMGf5HHkIYOJi99JiLW9sbskp7z/sDeJPes3URGC27nOeM7HYv5uoh4Ir3MPQf4NMmkrhdXw/d7lob2HxVTdqlMh/71637u27BZ71dy6fi9EXF4FeZvndGYdXYe0rFkp5GM1dra4xQRPb7sIumfEfHyDvbdFREv7Wl2rUh/kb8MuD0iXqbk1mOXRcTrCy6tS0quXb6OZMzXy0nvgxxl3Ee1Tf4uJD30HyGZ9f0LkiU+2hvHWfMk7dHZ/qz/aMuSpDtaZ6RLugB4PCLOTF8vjYjJGRzjEpKG9br09a4k95ou+4/irNVqoy/Tm9hb3xIR/4uImXk0yPLMlzRD0rCIeDpdwmJXSVndi3QvSX+Q9LikxyRdLSmrOy/8ClhGcieHr5NcdvlnmZnDOtm3U5nZteLZSKbxbZHUADxGcvm714vEXyLiaJKxfB8GbpV0YzpWtNz8JyPiZxFxCMm44a8BqyRdIumFXXy8S2nOsJLXu0r6Rbm5eWVHxIOtD2AEyVWLtwEjcuilP1TSR9LnIzP4OdNPUuuQr8NJ1htsldVQsP1aG3ywdfjR/hllZysCWloq86ggN/qsLzqynR88b8ooezZJb8oYknFgV5BcJsnCiEhuZN8UETemfx2XO7h6cTr7bxuSPo5X0W+1OG0c/IzknNwGZH6HjjxIGiHpM5IWA6eSXK4bCXyO5Hu13Px+kt4m6SqS2bDnAHsBfyC5FF6uPBsJuWVL+ipwCUnDbyRwkaQzsshO879G0shuvb1m6+LM5bgcuFHS1SQz1Rekx3ohGdxlJVWX9u6RZg/Hcwsqyifb+qJ+kgZGuoCspJ2AgRllD46IX5a8vkzS5zt8945pSv+7StKbgZXA8DIzTwauknQMzzXyGkl+ibyjzOyaULJu24XpJKOGiLizs8/0Iv8Afgm8PSL+W7J9saQLM8j/N/A34HsR8feS7b9VNsuf1EnaNW2QZd1IyDP7GOBlEbExzT6bZGxlJlcUSP5t7k/yBwgRsTK91N5jEXGWpHkkf7DOLRn/WUfyx0IWzgH+Ian1HtnvAc7KKNu6wY0+64t+BcwrmXzxEZK/yrPwZ0lfJOndC+C9wJ/SXyhEeXfP+JaS2yR9jmR9vgaSe9v2WEQ8Chwi6TCg9e4Y10bEXzv5WJ+gdH2yjvZltW5ZziZ2NHkjIr6TQf5+EfFUB/ll3y+bfBsJeWavJBl327qe40CS2apZyWVx5oi4pZ1tmYz/TLMuTXudW69QvDMi/tW6v7QR3ivU4JwHT+SwPknSG0kGuAP8JSKuzyi3s9tbRWRwSzyrDEl/62R3lDOBplIk/QV4T5uB83Mi4oiM8nMfmC9pX55rJPw1y0ZCXtmSfk8yceYvJH/8vR64leSe1mU3iCWdSjI7+PXAt0km6lweEeeWk1s0dXA7ziIM7TcypuxcmcULrn/yYt+GzSxPEXEdcF17+yT9IyJ6NMg9IrKatFFaz2kR8V1te3/c0mNm0aNibUQ+d4OptFFtx61Jel6G+duNi1Nyh47MpA2xf3Wwex5l3JEmx+yr0ker+T3MaVfkszhzb9Cr7m8YFZ5kUQlu9Jltr8e3qVKyyOsngdbxTPOB/4uIpg4/1LXW9c8Wd/ouy4WSBZk/BRxK0uheAFzYOl6rl2uWND4iHoKtS4pkeXknz3Fx3ZFnI6HH2RFxiaQBwIvSTcvL/BmwDUnfiYgvULI4c8m2auZLjzlzo89se+X84PkpUA/8JH39wXTbx3tcTMQf0v/mdX9c69ylJLcvOy99/X6SyRHvKayi7jsduFnSjSSNmFeR3Fc1K6Xj4kRyh5hKDszPs5HQ42xJU0nGCa8gOS/jJH04Im7KpLLksm7bBt6R7WzrFST1j4isbnVZIVGTY/rc6DPL1stbFzhN/VXSHeUESvoDnd8f923l5FuXJkXEviWv/yapo0uCvUpEXJdOSHlFuunkiFiTYX7pwPygzcD8Puwc4A0RsRxA0otIlkQ5sJxQSZ8k6XXeS1LpDPJdgIXlZOfsVrp3qbxXXd6tRW70mW2vnB88zZJeEBH/D7be8q3cm7h/v8zPW3luk/SK1pmNkg6mui61H8Jzww0A/phxfj3P/ZupzyJwB3qGMm0kSBoQEZszyK5vbfBBMgM2HfpRrtnAn0kmb3yxZPuTZa4MkLfunstcFvrvkQBa3NNnVlOU3H9zbZtlLT5YRuTnSXqC7if5QbcHyZIwPRYRN5bzeesZPXcf1Xrg75IeSneNJ7kzSq+Xrg/3cpJligA+I+mQiPhyRvmfIbnTx5Uk3++XSZoZEed1/skudbdnaIcbCZK+GhHfaGf7UOBqkvss9yi7xBJJs3huweRjyOAPhYhYT7JQ8tEA6aScQcAQSUNax272QqMkndLRzoj4Qfrf3txwrQlu9FmfIekVwNnAE8A3ScZljSQZjP6hdEYvEXF3T48REfMk7U0yow6SAdybyqz7NxFxVEkjpO0xy7qJu3XoLUUXkIE3AZMjuY1c6xIrtwOZNPqAjwEHR8TTaf53SBaELrfR162eoR42Eg6VdFZEnL71YMn9lK8HfldmdqvjgROA1pn1C3hunG/ZJL0V+AHJXX8eI/nj8l6S+3L3Rv2AIVTb5dvw7F2zanY+yS+7oST3lTwyIm6RtA/JeJt2l3DZEZLeA1wXEXemt106QNK3ylzI9zxJr2T7Rsg4YHUZudaJ0nulSnoZySQIgAURUdY4zQobRvKHDiTf+1kS2w5faCabX+zd6hnqobeR3DHkBxFxSvpH2p+B70dE2XcpkdQPuCMi9iFpmOXhWyTjNG+IiP3TxdU/kNOxsrCqvd5Vqzzfe9f6kv4RMTcirgBWt47RiogsL9V9JSKelHQoyeWhn5PM3i3HF4ANUXIz97RBsh74YZnZ1oX0EuavgOelj8skZXVbqrx9G7hd0sVpL98Ssp1dexGwSNKZks4EbiH5ni9Xa8/QLh08eixdaucdwJ6SLgduAD6fRYMvzW8Glksan0VeB5oiYi3JVYq6iPgbye0Te6vq6uEjuaQSLVGRRyW5p8/6ktK++mfb7MvqX15rr8ebgZ9FxLWSyr3f5uiIuKvtxoi4S9KeZWZb1/K6hJm7iLhc0nyScX0AX4iIzHqHI+IH6XIwr0w3fSQibs8gOreeoZIexEXAaSSXXie0bi+zF7HVrsA9km4Fnm7dmOFM+3WShgA3Ab+S9FjpcXqh3jNBo49zo8/6kpdJ2kDyV+dO6XPS1z1ekLmNRyT9H8k6Wt+RNJDye9SHdbJvpzKzrWt5XcLMjaR9ImJZyf2D/5v+d3dJu5c53KCtpcAq0t8npYtBl6Fb57eHt0or7Sk8t51tWfhKxnltTSO5r+9nSSaJDAV67eXTqpygEeExfWbVLCL6VeAwRwFvJBkftE7SGJIZveVYLOm4iPhZ6UZJHye5XGf5ar2E2XpbrbeTzSXMPH2OZFbtOe3sC56732xZ0svcXwMe5bnGcADlTi7qbs/QDt8qLSK+vuPl7LA3tb07RtpDnMlM/NZe55QXbbduU9TgitNmRStZSgGAcno+0pmFVwGbea6R1wgMAN6R5eU625akOpIB8xtJbsMGyUSOLC5hVj1J/yG59L22oOPfHhE7dK9fSV/tZHdExDfLLAtJt0XEAW223VnuTHtJT5I0qlt7Qlt/gYuk9oZy8u05DRoeB9e9viLHuqHlN0sioiJjMt3TZ5YhSW8j6V1pXUqhdU23Hi+lEBGPAoekM/QmpZuvjYi/llmudSEiWiRdkDYssrwkmitJ7+xsf0T8rrP9O+BhkglFRelJr0V7Y992Jhm7OYJkOaceKbljxgvauWPG33ua2yoisr4MbVVC0huBH5NMcpoVEWf3JMeNPrNsfZOcllJIZ+j9LYss2yHzJL0L+F1Uz6WRt3ayLyhZj65M9wPzJV0LbF2PMqPJELmIiK2XvCXtAnyGZAH1ObR/OXxH5HrHDEmDSNYAfCFwJ/CL6runre2odBmgC0jGiv8X+Keka3pyy0M3+syy1RQRayVtXUpB0o+KLsrK8gngFJJb7G1Mt/XqS2kRUdZdYHbAQ+ljQPrIXRa3SpM0nOT/6TEkY+IO6MGEkO203jFD0hfa7MrqjhmXAE0kM47fRHIF4TNlZlpHes9EjoOA/0TE/QCS5pBM5nGjz6xg1baUgnWhmi+ppeNBZwC7R8SRkvYFpkREJhNRWidFpN/zRMRTWeTmeas0Sd8D3gnMBF6aVc1tXMtzY+8GAROA5ZR/x4x9I+KlAJJ+TnK7Oqt9Y0mGUrT6L3BwT4I8kcMsQ5J2JlkDsI7nllL4VVED3S0b6Ri5Q0l+kS+IiN8XW1H3SPozyezj0yPiZZL6A7e3NhwyyJ9EcjvD4emmNcCHIuKeMnPnAv/s6FZp5azhJ6mF5FL0FrYdE5jbZIh06ZxPRcTHy8zZZoJIexNGLBuSriO5TWclDCKZLNZqZkTMLKnl3cAbW79/JH2QZALViTt6IPf0mWWoZCmFlnSc09oqGgdm7ZD0E5IxVJenm46X9PqIOKHAsrprZET8RtKXACJii6Tmrj60A2YCp6TjTZE0FfgZcEiZubndKi0iKn4nqoi4TVKPembaaF1rFLZdb9SzdzMWEW8suoYSj5DcdrPV89NtO8yNPrMMSHoFcDbJPU6/SdL7MZLkNkkfioiy7+trhXkt8OLWxnt6O7OyerIq6GlJI0h7tNLv0yxn2+7c2uADiIj5aW93WSJio6R3AL9Ob5V2CHByRFzVxUd7hTb3Da4jWUtwZbm5FVpr1HqffwJ7S5pA0th7H/D+ngS50WeWjfOBL5Nczv0rcGRE3CJpH5IeIjf6qtd/SJbeeTB9PS7dVg1OAa4hWUJkITAKeHeG+fdL+grJHzmQzFS/v9zQCt0qLU+l40C3kIzxu7KgWqzKpT30J5IMb+hHMmu7R394ekyfWQYkLY2IyenzeyPixSX7dngBWSuepD+Q9JANJbl37a3p64OBWyNianHVdV86jm8iySXA5RHRlGH2rsDXKVm4Gjiz3Jmwkr7W2f4K3VWjbJIGR8QzRddh1so9fWbZKJ3b/2ybff7Lqjp9v+gCMnIQsCfJz/sDJBERl2YRnDbuTsoiq01uVTTqOiJpCsmt+oYA4yW9DPhERHyq2Mqsr3NPn1kG0sHxT5MOrgZa/7oXMCgi6ouqzfouSb8EXgAsJbk3LiQD/stqqEm6prP9EfG2MvNzv1VaniQtIrmMfk1rL7+kuyNiUuefNMuXe/rMMuAB1rVH0s0RcWjJ/U637qJ6Zko2kqztlvVf91NI1g27nGTcXY8WSe5EbrdKq5SIeFja5rRkOWvarEfc6DMza98xUN2LMwN3A7sBqzLO3Y3kllBHk8wivBa4vNz1+VrlfKu0SnhY0iFASKonqf/egmsy8+VdM7P2lC58K+nKiHhX0TV1V8kklF2AySSTUErvjVvW5dc2xxpI0vj7HvD1iDg/o9y2t0r7cRa3SqsESSOBHwOvI+kFnQt8xou0W9Hc02dm1r7Sa3N7FVZFz+Q+CSVt7L2ZpMG3J3AukMk6ehW6VVpuImINaU+xWW/inj4zs3a06emryttdSfpORHyhq209yL0UmAT8CZgTEXeXk9dOfsVvlZaFap+AYrXPjT4zs3Z0MSO71zY8SrXXWJV0Z0TsV2ZuC89NtqiaRlneJH2unc1bJ6BExJAKl2S2DTf6zMxqjKRPAp8iuSz9/0p27QIsjIgPFFJYH1IyAeVjwG+AcyLisWKrsr7OjT4zsxojaSiwK/Bt4Islu56MiCeKqapvqOYJKFb73OgzM6th6d0gXpW+XBARdxRZTy1rMwHlgmqbgGK1z40+M7MaJekkYDrwu3TTO4CZEXFecVXVrmqdgGJ9hxt9ZmY1StKdwJSIeDp9vTPwj3IncphZdaorugAzM8uN2Pb2X81kf8s0M6sSXpzZzKx2XQQsknQVSWNvGvDzYksys6L48q6ZWQ2TdABwaPpyQUTcXmQ9ZlYcN/rMzGqUpBcA/42ITZIOA14KXBoR64qtzMyK4DF9Zma160qgWdILgQuBccDsYksys6K40WdmVrtaImILydpx50fE54ExBddkZgVxo8/MrHY1SToa+BDwx3RbfYH1mFmB3OgzM6tdHwGmAGdFxAOSJgC/LLgmMyuIJ3KYmdUwSTsB4yNiedG1mFmx3NNnZlajJL0VWApcl76eLOmaQosys8K40WdmVrvOBA4C1gFExFJgr+LKMbMiudFnZla7miJifZttLYVUYmaF823YzMxq1z2S3g/0k7Q3cBLw94JrMrOCuKfPzKx2fRp4CbAJuBzYAJxcZEFmVhzP3jUzMzPrA3x518ysxnQ1Qzci3lapWsys93Cjz8ys9kwBHia5pLsIULHlmFlv4Mu7ZmY1RlI/4PXA0cB+wLXA5RFxT6GFmVmhPJHDzKzGRERzRFwXER8GXgH8B5gv6cSCSzOzAvnyrplZDZI0EHgzSW/fnsC5wFVF1mRmxfLlXTOzGiPpUmAS8CdgTkTcXXBJZtYLuNFnZlZjJLUAT6cvS3/IC4iIaKh8VWZWNDf6zMzMzPoAT+QwMzMz6wPc6DMzMzPrA9zoM7NeS1KzpKWS7pZ0haTBZWRdLOnd6fNZkvbt5L1TJR3Sg2OskDSyu9vbvOepHTzWmZJO3dEazazvcqPPzHqzZyNickRMAjYDx5fulNSjZaci4uMR8a9O3jIV2OFGn5lZb+ZGn5lViwXAC9NeuAXp/WX/JamfpO9J+qekOyV9AkCJ8yUtl3QD8LzWIEnzJTWmz98o6TZJd0iaJ2lPksblZ9NexldJGiXpyvQY/5T0yvSzIyTNlXSPpFl043Znkn4vaUn6melt9v0w3T5P0qh02wskXZd+ZoGkfTI5m2bW53hxZjPr9dIevSOB69JNBwCTIuKBtOG0PiJeni5IvFDSXGB/YCKwLzAa+Bfwiza5o4CfAa9Os4ZHxBOSLgSeiojvp++bDfwwIm6WNB64Hngx8DXg5oj4hqQ3Ax/rxpfz0fQYOwH/lHRlRKwFdgYWR8RnJX01zT4RmAkcHxH/lnQw8BPgtT04jWbWx7nRZ2a92U6SlqbPFwA/J7nsemtEPJBufwOwX+t4PWAosDfwapL7zTYDKyX9tZ38VwA3tWZFxBMd1PE6YF9pa0deg6Qh6THemX72Wkn/68bXdJKkd6TPx6W1rgVagF+n2y8Dfpce4xDgipJjD+zGMczMtuNGn5n1Zs9GxOTSDWnj5+nSTcCnI+L6Nu97U4Z11AGviIiN7dTSbZKmkjQgp0TEM5LmA4M6eHukx13X9hyYmfWEx/SZWbW7HvikpHoASS+StDNwE/DedMzfGOCwdj57C/BqSRPSzw5Ptz8J7FLyvrnAp1tfSJqcPr0JeH+67Uhg1y5qHQr8L23w7UPS09iqDmjtrXw/yWXjDcADkt6THkOSXtbFMczM2uVGn5lVu1kk4/Vuk3Q38H8kVzGuAv6d7rsU+EfbD0bE48B0kkupd/Dc5dU/AO9oncgBnAQ0phNF/sVzs4i/TtJovIfkMu9DXdR6HdBf0r3A2SSNzlZPAwelX8NrgW+k248BPpbWdw8wrRvnxMxsO74Nm5mZmVkf4J4+MzMzsz7AjT4zMzOzPsCNPjMzM7M+wI0+MzMzsz7AjT4zMzOzPsCNPjMzM7M+wI0+MzMzsz7AjT4zMzOzPuD/A164ZCtwK5fXAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 720x720 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "cm = confusion_matrix(y_true_cell_type, y_pred_cell_type, labels=scrna_unique_clusters)\n",
    "disp = ConfusionMatrixDisplay(confusion_matrix=cm,\n",
    "                              display_labels=scrna_unique_clusters)\n",
    "fig, ax = plt.subplots(figsize=(10,10))\n",
    "disp.plot(ax=ax, xticks_rotation = \"vertical\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "the confusion matrix above give information about which cell type is miss-classified. As example, some T-cell are miss-identify as NK-cell or Neutriphilis.  \n",
    "\n",
    "For the purpose of this tutorail we only use a very small number of cells but those notebooks can be run on larger scale."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}