nn_analysis.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json\n",
    "import os\n",
    "import pickle\n",
    "from collections import defaultdict\n",
    "from glob import glob\n",
    "from pathlib import Path\n",
    "import matplotlib.ticker as mtick\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from matplotlib import rcParams\n",
    "from matplotlib.lines import Line2D\n",
    "from rdkit import Chem\n",
    "from rdkit.Chem import AllChem\n",
    "\n",
    "from plot_utils import *\n",
    "import pandas as pd \n",
    "from sklearn.ensemble import RandomForestClassifier\n",
    "from sklearn.model_selection import train_test_split\n",
    "import copy\n",
    "import re \n",
    "from sklearn.neighbors import NearestNeighbors\n",
    "from matplotlib.ticker import PercentFormatter\n",
    "from utils import *\n",
    "\n",
    "params = {\n",
    "    'font.size': 18,\n",
    "    'legend.fontsize': 18,\n",
    "    'text.usetex': False,\n",
    "}\n",
    "rcParams.update(params)\n",
    "\n",
    "from rdkit import RDLogger\n",
    "RDLogger.DisableLog('rdApp.*')\n",
    "\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def fit_clfs(df, is_aldh1=False, topliss_augmented=False):\n",
    "    \"\"\"\n",
    "    Args:\n",
    "        chid: which assay to use:\n",
    "        external_file:\n",
    "    Returns:\n",
    "        clfs: Dictionary of fitted classifiers\n",
    "        aucs: Dictionary of AUCs\n",
    "        balance: Two numbers showing the number of actives in split 1 / split 2\n",
    "        df1: data in split 1\n",
    "        df2: data in split 2\n",
    "    \"\"\"\n",
    "    # read data and calculate ecfp fingerprints\n",
    "\n",
    "    \n",
    "    if is_aldh1:\n",
    "        df = pd.read_csv(\"210415_LIT_PCBA_ALDH1_Best_AB_sets.csv\")\n",
    "        in_set_A = np.where(df['in_set_A']==1)[0]\n",
    "        in_set_B = np.where(df['in_set_B']==1)[0]\n",
    "        smiles = df[\"SMILES\"]\n",
    "        activity = df[\"Act\"]\n",
    "        index_test = []\n",
    "        for i in range(len(smiles)):\n",
    "            if i not in in_set_A and i not in in_set_B:\n",
    "                index_test.append(i)\n",
    "\n",
    "        in_test = np.array(index_test)\n",
    "\n",
    "        #in_test = np.where(df['in_set_A']==0 and df['in_set_B']==0)\n",
    "\n",
    "\n",
    "        X1 = np.array(ecfp(smiles[in_set_A]))\n",
    "        X2 = np.array(ecfp(smiles[in_set_B]))\n",
    "\n",
    "        y1 = np.array(activity[in_set_A])\n",
    "        y2 = np.array(activity[in_set_B])\n",
    "        smiles_test = np.array(smiles)[in_test]\n",
    "    else:\n",
    "        df['ecfp'] = ecfp(df.smiles)\n",
    "        df_train, df_test = train_test_split(df, test_size=0.1, stratify=df['label'], random_state=0)\n",
    "\n",
    "\n",
    "        df1, df2 = train_test_split(df_train, test_size=0.5, stratify=df_train['label'], random_state=0)\n",
    "        X1 = np.array(list(df1['ecfp']))\n",
    "        X2 = np.array(list(df2['ecfp']))\n",
    "\n",
    "        y1 = np.array(list(df1['label']))\n",
    "        y2 = np.array(list(df2['label']))\n",
    "\n",
    "        del df1['ecfp']\n",
    "        del df2['ecfp']\n",
    "        smiles_test = df_test.smiles\n",
    "    \n",
    "    \n",
    "    balance = (np.mean(y1), np.mean(y2))\n",
    "    clf_1 = RandomForestClassifier(\n",
    "        n_estimators=100, random_state=0, min_samples_leaf=1, max_depth=None)\n",
    "\n",
    "    clf_2 = RandomForestClassifier(\n",
    "        n_estimators=100, random_state=0)\n",
    "    clf_1.fit(X1, y1)\n",
    "    clf_2.fit(X2, y2)\n",
    "    test_scores = clf_1.predict_proba(ecfp(smiles_test))[:, 1]\n",
    "    \n",
    "    if topliss_augmented:\n",
    "        augmented_smiles_test = []\n",
    "        for i, s in enumerate(smiles_test):\n",
    "            if test_scores[i]>0.25:\n",
    "                augmented_smiles_test.extend(topliss_walk(s, n_combinations=2))\n",
    "            else:\n",
    "                augmented_smiles_test.extend([s])\n",
    "        augmented_smiles_test = np.array(augmented_smiles_test)\n",
    "        X_test = ecfp(augmented_smiles_test)\n",
    "    else:\n",
    "        augmented_smiles_test = smiles_test\n",
    "        X_test = ecfp(augmented_smiles_test)\n",
    "\n",
    "    scores_1 = clf_1.predict_proba(X_test)[:, 1]\n",
    "    scores_2 = clf_2.predict_proba(X_test)[:, 1]\n",
    "\n",
    "    return df1.smiles, df2.smiles, y1, y2, augmented_smiles_test, scores_1, scores_2\n",
    "\n",
    "def get_sims(smiles_1, smiles_2):\n",
    "    neigh = NearestNeighbors(10, 0.1, metric='jaccard', algorithm='brute')\n",
    "    neigh.fit(ecfp(smiles_2))\n",
    "    sims = []\n",
    "    for x in ecfp(smiles_1):\n",
    "        dist, idx = neigh.kneighbors(x.reshape(1, -1))\n",
    "        sims.extend(dist[0])\n",
    "    return sims"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Similarity of high scoring molecules in held-out test set with Split 1/2 molecules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns = ['AssayID', 'TargetName', 'TargetID']\n",
    "\n",
    "rows = [('CHEMBL3888429', 'JAK2', 'CHEMBL2971'),\n",
    "        ('CHEMBL1909203', 'EGFR', 'CHEMBL203'),\n",
    "        ('CHEMBL1909140', 'DRD2', 'CHEMBL217'),\n",
    "        ('CHEMBL1794345', 'Pla. Fal.', 'CHEMBL364'),\n",
    "        ('CHEMBL1741322', 'Cytochrome', 'CHEMBL3356'),\n",
    "        ('CHEMBL1963961', 'CCRF-CEM', 'CHEMBL382')]\n",
    "\n",
    "assay_info = pd.DataFrame(rows, columns=columns)\n",
    "\n",
    "chid2name = dict(zip(assay_info['AssayID'], assay_info['TargetName']))\n",
    "\n",
    "\n",
    "chids = ['CHEMBL3888429', 'CHEMBL1909203', 'CHEMBL1909140']\n",
    "optimizers = ['graph_ga', 'lstm_hc', 'mso']\n",
    "auc_title = False # whether or not to include auc in plot\n",
    "lw = 0.7 # line width\n",
    "nrows = 1 # for two optimizers\n",
    "ncols = len(chids) # for the number of tasks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABB4AAAK3CAYAAAAiZYwvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAABhuklEQVR4nO3deZhcZZn38e+NrGFpRomsAkJEFjOihiU6EEBx3Fryjq/rBAgyybyDLDKuIIqOCjgKUXDUCY5GQEZBQSYyo4RhEccFAurADBIkBGJQTAg0Iju53z/OaSg61Z3q6jpVXdXfz3XVVamzPOeuLvKk+PVznicyE0mSJEmSpCqs1+kCJEmSJElS7zJ4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kNYhIpZFREbEQcPs3ywiHiqPWbyOtnYujxt2HduI+PvymMci4i3ltkkR8eaImB8R/xsRD0fEIxFxW0R8ISJ2GMt7lKR2qelT1/XYss65m0TE30XE5RGxIiIejYg/lW1eVu57Xp3zZjdwve81UOfDEXF7RHw1Inav7qckSaNXp99aExED5fbLI+IjEbHjMOceVKfPWxMRD0TEzyPigxGx8SjOfTwiVkbEzRFxXkT89XDnl228PCI+GhHXRcTqiHgiIu6NiH+LiL9s1c9InROZw/7/jySKThzYCTg4M6+ps/8o4Gs1m6Zm5i3DtLUzcCdAZkad/R8F/gF4BHhLZv5Huf1vgHPLw+4AfgWsD+wLbAMMAH+ZmT8f9RuUpDaq6VN/CPx+hEOPycyHa847APgWsB3wJHATcDeQwPbANGBD4CGgv7a/jojZwNeBe4EfDHO9mzLz7HXUuTVFv/tc4FHgdZl57TrftCS1wTD91qbAtsDLgEnAGuCLwIcy89Gacw8Crgb+BHyn3Lw+8EJgOhDAYuCgzPzTkOvWO3c9oA/YDXhxef69wN9m5mVDzl8feKJ8OQBcD6wuz31Zuf0zmfnhUfw4NM6s3+kCpB4wu3y+h+IL8Wzg/aNtJCI+A3yQOl+aKTrjc4EvZOb/1JyzCfDPwOHAtyJit8x8Akka/86oF+bWExF/AfwnsAHwZeDUzFw55JhNgSOAk4HhRoH9OjNnj6XOiPgz4DLgAIr+15EPksabtfrXiNiI4vviZ4HjgV0i4rDMXDPk3FVD+8mIeBVFHzwNOA44Y5jrrnVuef6OwKnAu4FLI+IdmXnRkMNuAE4Hvl/7XTYiZlKEGR+KiCsz88rh3rTGN2+1kMYgInah+PL5J+CocvNfl8lto21ERHyRInS4H3jN0H8sMvMbmTm3NnQotz8C/C3wILAzRSItST2jHJr7LYrQ4fTMPGZo6ACQmX/KzC8DLwEqG/2VmfdT9NcALy7/HZCkcS0zH8vMrwIHUoysfRNwdIPn/hewoHz5hiaufXdmHg18imLkw1cj4rk1+5/MzH0z89Khv0DLzO/xzMjid4322ho/DB6ksTmSogO9JDOvAJZQ3PrwukZOjoj1gH8B3gOspLidY1RfmMvw4bby5fajOVeSusDhFH3bPRS/MRtRZg5k5u0V11R7O93WFV9LklomM28GPl++/PtRnDrY742lz/sEsBzYHJgzivN+WT77PbeLGTxITYqIoAge4JkU+Bvl8+wGzl8fuJBipMQK4MDM/FUTdTyHYrQDjHy/tCR1ozeVzxeNo1vJtqj58x86VoUkNefC8nn3iNiuwXMG+72m+7zMfBIYvMXi1aM4dUr57PfcLmbwIDXvYIoJfO6mmFAH4DyKSXv6a4eQDeO7wNspJps8IDN/3WQdRwGTKUZM/KTJNiRpvHp5+XxjR6t4tsGhxr8GlnayEElqwv8Cj5d/3rPBc15fPl8+xmv/YjTXjYjJPPMLvctGOFTjnJNLSs2bXT6fn+XyMJn524j4T+BQivvQvjjC+W8un9+dmXc2U0BEvAj4XPny5Mx8rJl2JKkDri4GjtV1bWYeVP55q/J5rXkdACLiVIpZ12v9OjPrTX42I4ZZzrjeSkN1rvV8ii/fnwX+CPzNYP8vSd0iM9dExGqK24PXWoJ4UDk6dxfgA8BfANcB54zx8qvK52GvW3P9oFiR6M+AHwOXjvHa6iCDB6kJEbE58Jby5TeG7P4GRfAwm5GDhx9RTPBzfkTMyMxR/dasXKv+3yiWKvrXcsIgSeoWIy2nOZoRYP3AK4Zsu5b6s66PtJzmcOoFJPcAMzLzN6NsS5LGi8GR70PD052GCWj/FTg8M5+q6Lr1fAZ4I0Xf/dcGvd3N4EFqztso1kL+SZ1JzC6hWGXiFRHxksy8Za2zC2+k+OL9Soovtgdm5l2NXLwMPv6DYhm3H9DAnBKSNM40upzmKorlMbeqtzMzpw3+OSLeBCwcoa1mltMcDEjWA15AsZLRdsCFZb/96Cjbk6SOKucH27J8uXrI7j9RLF8JsAnwUuDFwDspbnk7c4yXH+zLh153aI0fphhpcT/w2sy8e4zXVYcZPEjNmV0+7xQRP66zP2uOe3+9BjLzoYh4PXAFsB9F+DAjM5ePdOGImERxf90+FKMm/iozHx/pHEnqYr+gCB5eAXyzA9d/VkASEXtQzOuzD3Aao5sVXpLGg72ADcs/D/0F2aqhAW1EvBeYB5wREVdn5k1juPbLhrlu7fWOBU4HHgJen5n/PYbraZxwcklplCJiV4r73KBY1udVdR595f6/Lu+PqyszHwT+EriB4h7lq0aaXTgiNgK+R/Ebt58DbyqX05SkXvX98vltI/Wn7ZKZt/LMMnDHRsQunaxHkprwzvL5fzJznStFZObnKUb0rk8xx01Tyj78reXLK4c55t3A2cAjQP9ol5nX+GXwII3e7PL50syM4R7A7RST9rxupMYycwB4LXATxXJBV0XENkOPKzvriyjmj/gFRQL8x5a9K0kan86jmFNhe+DUDtcCQGYupBj1sAFwSofLkaSGRcRU4Pjy5Whum/gg8CRwSEQc2OTlT6UYwfYgsNbcZBHxDuBcihU3/k+Dt+OpSxg8SKMQEesBR5QvL1jH4YNDgmevq93MfIAiUPglxX10V5Wzp9de93yKlTD+h+Jet/tHUbokdaVyDoV3AU8Ap0TEF8vl1Z4lIjYA9m1jaR8tnw+PiKGrakjSuBIRG0XE0RS36U6iWJpy6ATpw8rMO4AF5ctRhcARsWNEnEsR1CbFim4PDDnmMIrvuk8Bb8vMH47mGhr/wslBpZFFxDJgJ2AGxf1wi4AHgG1GWr4yIqZQjHp4HNg2M1dHxM7AnVB/6bZypYqrgakU974dnJmrIuJ44AvlYT+gmN23nu9l5vdG9w4lqX1q+tSRVrUAOLv2PuKIOIhiVvVtKH7rdiMwONnYdhT3LG9JcU/wSZn5xZpzZ1MsyVa7TGejdR483G/dIuI/KEa1/Utm/k0j7UpSVYbpXydR9JsvBzYF1lDcynBS7eS4ZR97NXBXZu48TPs7Uny33RA4IDN/POTc2okp1wO2AF4E7AFEWdOczPz+kHafT9GfbwT8BvivYd7iqsysO3eaxr+O3yspdYFNyuc/AXPLP188UugAkJm/iYjrKX4D9y5GXlpz8Jz7IuLVFJ33S4BF5evn1hw20q0byyjmgJCk8e4v17H/exS3oAGQmdeUc+wcBbyJYqb1vSm+RK+kWOP9hxTLC99XQb31nELRJx8REZ/OzDvbdF1JGslg/5oUYexqimWG/ws4f10TmQ8nM++OiH8GjgM+DrxmyCGbAkeWf36C4paK31GMEv4P4JJhvj9PoggdoLjteMowJdzFMJO2a/xzxIM0gnLZygGKlHarNn6ZlSRJkqSe4BwP0siOpwgdfmnoIEmSJEmj560W0hARsRXweWB3inXj1wAnd7ImSZIkSepW3mohDVEzAeQAsBg4IzPrrjUsSZIkSRqZwYMkSZIkSaqMczxIkiRJkqTKOMdDB2y11Va58847d7oMSR124403rsrMyZ2uoxfZz0oaZF9bDftZSYMa6WcNHjpg5513ZvHixZ0uQ1KHRcRdna6hV9nPShpkX1sN+1lJgxrpZ73VQpIkSZIkVcbgQZIkSZIkVcbgQZIkSZIkVcbgoY0ioj8i5g8MDHS6FEmSJEmS2sLgoY0yc2Fmzu3r6+t0KZIkSdKo+Ys0Sc0weJAk9Qy/EEtStfxFmqRmGDxIknqGX4glSZLGn/U7XYA0ET366KOsXLmSRx99lCeffLLT5aiFNthgA57//OezxRZbdLoUSZIkaVwweJDabGBggHvvvZfJkyezzTbbsP766xMRnS5LLZCZPPLII6xYsQLA8EGSJEnCWy3aynuPBbBq1Sp22GEH/uzP/owNNtjA0KGHRASTJk1i++235w9/+EOny5EkSZLGBYOHNvLeYwE8/vjjbLLJJp0uQxXaZJNNeOKJJzpdhiRJkjQuGDxIHeAoh97m5ytJkiQ9w+BBkiRJUkO8dVhSMwweJEk9wy/EklQtbx2W1AyDB0ldacGCBUQE11xzzdPbrrnmGiKCBQsWdKwudZZfiCVJksYfl9OUxpl5i5Z0uoRnOfHQ3cZ0/tKlSznjjDP40Y9+xN13381GG23ENttsw7777svs2bM5+OCDW1RpfcuWLWPBggXMnDmTvffeu6FzHnroIc4880xuvPFGbrrpJlasWMGMGTOeFXJIkiRJaozBQxtFRD/QP2XKlE6XIrXF4sWLmTFjBhtssAFHHHEEe+21F4888gi33347V1xxBZtvvnlLg4cDDzyQRx55hA022ODpbcuWLeMTn/gEO++8c8PBw6pVq/j4xz/O1ltvzSte8QruvffeltUoSZIkTTQGD22UmQuBhdOmTZvT6VqkdvjEJz7Bww8/zC9/+Ute+tKXrrX/97//fUuvt95667HxxhuPuZ1tt92W5cuXs8MOOwCw2WabjblNSZIkaaJyjgdJlbn99tt53vOeVzd0ANhmm22e9ToimD17NldeeSX7778/kyZNYptttuGEE07goYceWuf1hs7xsGDBgqdHVBx11FFEBBHBQQcdNGI7G2200dOhgyRJkqSxMXiQVJldd92V++67j0suuaThc2666SZmzpzJ9OnT+dznPscBBxzA2WefzWGHHcaaNWtGdf0DDzyQk08+GYC5c+dy/vnnc/755/ORj3xkVO1IkiRJap7Bg6TKnHLKKWywwQa85S1vYbfdduPd7343X/7yl7n11luHPefmm2/mggsuYN68eRxzzDFcfPHFHH/88Vx11VVcdNFFo7r+LrvswqGHHgrA9OnTmTVrFrNmzXp6myRJGh2XLZbUDIMHSZWZPn06N954I0ceeSQDAwN8/etf55hjjmHPPffkwAMPZOnSpWud8+IXv5iZM2c+a9uHP/xhAC699NJ2lK0u5hdiSaqWyxZLaobBg6RKTZ06lQULFnDvvfeybNkyvvGNb3DAAQdw3XXXcdhhh/H4448/6/g99thjrTa23XZbttxyy7pBhVTLL8SSJEnjj8GDpLbZaaedOOKII7j22mt51atexS233ML111/f6bIkSZIkVcjgQVLbRQT77bcfACtWrHjWvnrzP/zud7/jgQceYJdddmnqWpIkSZI6x+Chjbz3WBPNokWLePLJJ9fa/sgjj3DFFVcAsOeeez5r32233cb3vve9Z237zGc+A7DW3A+N2GyzzQBYvXr1qM+VJEmSNHbrd7qAiSQzFwILp02bNqfTtUjtcOKJJ3Lffffx5je/malTpzJp0iSWL1/OhRdeyJIlSzjiiCOYOnXqs86ZOnUqs2bNYs6cObzoRS/i6quv5jvf+Q4zZszg7W9/+6hr2HPPPdl888350pe+xKRJk9hyyy15/vOfzyGHHDLieV/84hd54IEHAHjiiSe46667+NSnPgXAS1/6Uvr7+0ddiyRJkjQRGTxIqsxZZ53FZZddxo9//GO++93v8sADD9DX18ef//mf86EPfYjZs2evdc7LX/5yzjrrLD7ykY/wla98hS222IJjjz2W0047jfXWG/0grU022YRvfetbnHLKKbz3ve/lscceY8aMGesMHj73uc9x1113Pf162bJlfPSjHwXgyCOPNHiQJEmSGhSZ2ekaJpxp06bl4sWLO12GOuTWW2+tu3KDivkYjjzySBYsWNDpUsaskc85Im7MzGltKmlCsZ+VNMi+thr2s5IGNdLPOseDJEmSpIY4Z5mkZhg8SJIkSWpIZi7MzLl9fX2dLkVSFzF4kCRJkiRJlXFySUnjhnPOaKwioh/onzJlSqdLkSRJUskRD5KknuEQYEmSpPHH4KGNnIxHkiRJkjTRGDy0kb+JkyRJkiRNNAYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJklTP1acXD0ljYvAgqSstWLCAiOCaa655ets111xDRLBgwYKO1SVJUi9zeXhJzVi/0wVIGmK8peoHnzSm05cuXcoZZ5zBj370I+6++2422mgjttlmG/bdd19mz57NwQcf3KJC61u2bBkLFixg5syZ7L333g2dc9NNN3HhhRdy1VVXceeddwIwZcoUjjrqKObMmcMGG2xQYcWSJI1fmbkQWDht2rQ5na5FUvcweJBUmcWLFzNjxgw22GADjjjiCPbaay8eeeQRbr/9dq644go233zzlgYPBx54II888sizgoFly5bxiU98gp133rnh4OEf//EfufLKK5k5cyZz5szhqaee4vvf/z7vec97uOyyy/jBD35ARLSsbkmSJKmXGTxIqswnPvEJHn74YX75y1/y0pe+dK39v//971t6vfXWW4+NN954zO0cd9xxLFiw4FltHXvsscyaNYtvfvObXH755bzpTW8a83XUehHRD/RPmTKl06VIkiSp5BwPkipz++2387znPa9u6ACwzTbbPOt1RDB79myuvPJK9t9/fyZNmsQ222zDCSecwEMPPbTO6w2d42HBggVPj6g46qijiAgigoMOOmjEdl71qlfVDTDe/va3A3DLLbessxZ1RmYuzMy5fX19nS5FkiRJJYOHNnIyHk00u+66K/fddx+XXHJJw+fcdNNNzJw5k+nTp/O5z32OAw44gLPPPpvDDjuMNWvWjOr6Bx54ICeffDIAc+fO5fzzz+f888/nIx/5yKjaGfTb3/4WgK233rqp8yVJkqSJyOChjfxNnCaaU045hQ022IC3vOUt7Lbbbrz73e/my1/+Mrfeeuuw59x8881ccMEFzJs3j2OOOYaLL76Y448/nquuuoqLLrpoVNffZZddOPTQQwGYPn06s2bNYtasWU9vG42HHnqIz372s/T19XHYYYeN+nxJkiRpojJ4kFSZ6dOnc+ONN3LkkUcyMDDA17/+dY455hj23HNPDjzwQJYuXbrWOS9+8YuZOXPms7Z9+MMfBuDSSy9tR9lreeqpp5g1axZ33nknX/7yl3nuc5/bkTokSZKkbmTwIKlSU6dOZcGCBdx7770sW7aMb3zjGxxwwAFcd911HHbYYTz++OPPOn6PPfZYq41tt92WLbfcsm5QUbU1a9bw7ne/m8suu4xPf/rTvPOd72x7DZIkSVI3M3iQ1DY77bQTRxxxBNdeey2vetWruOWWW7j++us7Xdaw1qxZw9/8zd9w3nnnceqppz49X4QkSZKkxhk8SGq7iGC//fYDYMWKFc/aV2/+h9/97nc88MAD7LLLLk1dqxmDocPXv/51TjnlFD7+8Y831Y4kSZI00Rk8SKrMokWLePLJJ9fa/sgjj3DFFVcAsOeeez5r32233cb3vve9Z237zGc+A7DW3A+N2GyzzQBYvXp1w+dkJnPmzOHrX/86J598Mp/85CdHfV1JkiRJhfU7XYCk3nXiiSdy33338eY3v5mpU6cyadIkli9fzoUXXsiSJUs44ogjmDp16rPOmTp1KrNmzWLOnDm86EUv4uqrr+Y73/kOM2bM4O1vf/uoa9hzzz3ZfPPN+dKXvsSkSZPYcsstef7zn88hhxwy7Dkf+MAH+NrXvsZLX/pS9thjDy644IJn7d91112ZPn36qGuRJEmSJiKDB0mVOeuss7jsssv48Y9/zHe/+10eeOAB+vr6+PM//3M+9KEPMXv27LXOefnLX85ZZ53FRz7yEb7yla+wxRZbcOyxx3Laaaex3nqjH6S1ySab8K1vfYtTTjmF9773vTz22GPMmDFjxOBh8eLFAPzqV7/i8MMPX2v/kUceafAgSZIkNcjgQRpvDj6p0xW0zGtf+1pe+9rXjvq817zmNbzmNa8Z8ZjZs2evFVwcdNBBZOZax77hDW/gDW94Q8PXv+aaaxo+VpKkiSQi+oH+KVOmdLoUSV3EOR4kSZIkNSQzF2bm3L6+vk6XIqmLGDxIkiRJkqTKGDxIktSD5i1awrxFSzpdhiRJknM8SBo/6s3PII2G9x5LkiSNP454kCT1DO89liRJGn8MHtooIvojYv7AwECnS5EkSZIkqS0MHtrI38RpkLcU9DY/X0mSJOkZBg9Smz3nOc/hiSee6HQZqtCTTz7J+us7hY4kSZIEBg9S222++eY8+OCDnS5DFfrjH//Ixhtv3OkyJEmSpHHB4EFqs+c+97ncf//9rFq1iscff9xh+T0kM3n44YdZtWoVkydP7nQ5kiRJ0rjgWGCpzTbaaCN23HFHVq9ezbJly3jqqac6XZJaaKONNmLrrbd2xIMkSZJUMniQOmCjjTZi2223Zdttt+10KZIkSZJUKW+1kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJktSQiOiPiPkDAwOdLkVSFzF4kCRJktSQzFyYmXP7+vo6XUql5i1awrxFSzpdhtQzDB4kSZIkSVJlDB4kSZIkSVJlDB4kSeoVV59ePCRJksYRgwdJkiRJklQZgwdJUs9wtnVJkqTxx+BBktQzJsps65IkSd3E4EGSJEmSJFXG4EGSJEmSJFXG4KGNvPdYkiRJkjTRGDy0kfceS5IkSZImGoMHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSZJUmfU7XYAkSRqbeYuWAHCi/6pLkqRxyBEPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkiT1sqtPLx6SJEkdYvAgSZIkSZIqY/AgSZIkSZIqY/AgSZIkSZIqY/AgSZIkSZIqY/AgSZIkSZIqY/AgSZIkqSER0R8R8wcGBjpdiqQuYvAgSZIkqSGZuTAz5/b19XW6FEldxOBBGgfmLVrCvEVLOl2GJEmSJLWcwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSRPN1acXD0mSpDYweJAk9QyXeZMkSRp/DB4kST3DZd4kSZLGH4MHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSWqEk/NKTTF4kMYT/yGTJEmS1GMMHiRJkiRJUmUMHqROcaieJEmSpAmgo8FDRJwUERdHxNKIyIhYto7j94uIKyPijxHxYET8ICL2HuU1G2ojIl4RET+OiIci4taIeMcw7V0WEZePpgZJkiRJ3WPeoiXMW7Sk02VIXWv9Dl//NGA1cBOw5UgHRsT+wDXACuBj5eZjgesi4pWZefO6LtZoGxGxOfB94LfA+4GDgG9GxB2ZeUNNe28FDgH2WvdblSRJkiRp4ul08LBrZi4FiIhbgM1GOPZs4HHgwMxcUZ5zEXArcCbw2gau12gbrwS2AaZn5rKImA/sB8wEbijP27Js75TMvLvB96uJavCWioNP6mwdkiRJktRmHb3VYjB0WJeImALsA1w8GBiU568ALgZeExHbtLCNTcrn1eUxa4AHgE1rmvwscDdwTiPvQZIkSZKkiahbJpfcp3z+aZ19PwMCeEUL27gReAL4ZETsFBFHAi8FfgIQETOAI4E5ZSghSZIkSZLq6PStFo3arnxeUWff4LbtW9VGZi6PiOOBzwPHl/sWABdHxEbAfOBzmfnf66xckiRJkqQJrFtGPEwqnx+rs+/RIce0pI3M/AqwLTAd2DEzj8rMBD5K8XP7h4h4bkScFxG/jYhflJNN1hURcyNicUQsXrly5TpKlSRJkiSpN3TLiIeHy+eN6uzbeMgxLWsjM++nuA0DgIh4CfAB4HWZ+WhEXAo8D/grYF/g2xFxd2b+fOgFMnM+xUgJpk2bluuoVT1scCmmE7vlb58kSZIkjUG3/K/PPeVzvdspBrfVu4WiZW1ExHrAucAFmXl1RGwHvA44NDOvB66PiHcA7wbWCh4kSZIkSZqIuuVWixvK5+l19u0PJMWEkFW28R7ghcD7y9c7lM/La45ZDrxgHXVIkiRJkjRhdEXwkJm/ARYDby1HGgBQ/vmtwFWZ+fua7VtFxO4R0ddsG7Ui4gXAp4ETytsv4JkRFFNrDp1as12SJEmSpAmvo7daRMThwE7ly8nAhhFxSvn6rsw8v+bwE4Crgesi4pxy23EU4cn7hjR9LHAqcBTFahTNtFHrS8CPMvPbgxsy87cRcQ3whTK8eAWwF8XICEmSJEmSROfneDgamDFk2yfL52uBp4OHzPxJRBwEfKp8JPAT4K2Z+atGLtZMGxHxNuAgilBhqHcBXwb+AVgFHJ2Z1w53/YjoB/qnTJnSSLmSJEmSJHW9jgYPmXnQKI//KfDqBo77OPDxsbRRc/xFwEXD7PsdMHMUbS0EFk6bNm1Oo+doYnp65YtDd+twJZIkSZI0Nl0xx4MkSZIkSepOBg+SJEmSJKkyBg+SJEmSJKkyBg+SJEmSJKkyBg9tFBH9ETF/YGCg06VIkiRJktQWBg9tlJkLM3NuX19fp0uRJEmSVJWrTy8ekgCDB0mSJEmSVCGDB0mSJEmSVBmDB0mSJGkMIuKkiLg4IpZGREbEsnUcv19EXBkRf4yIByPiBxGx9yiv2VAbEfGKiPhxRDwUEbdGxDuGae+yiLh8NDVIUqPW73QBUi+bt2gJACf6N03SeDd4L/LBJ3W2Dqk7nQasBm4CthzpwIjYH7gGWAF8rNx8LHBdRLwyM29e18UabSMiNge+D/wWeD9wEPDNiLgjM2+oae+twCHAXut+q5I0ev7vkCRJkjQ2u2bmUoCIuAXYbIRjzwYeBw7MzBXlORcBtwJnAq9t4HqNtvFKYBtgemYui4j5wH7ATOCG8rwty/ZOycy7G3y/kjQq3mrRRi6nKUmS1HsGQ4d1iYgpwD7AxYOBQXn+CuBi4DURsU0L29ikfF5dHrMGeADYtKbJzwJ3A+c08h4kqRkGD23kcpoai3mLljx964YkSepK+5TPP62z72dAAK9oYRs3Ak8An4yInSLiSOClwE8AImIGcCQwpwwlJKkS3mohSZIktcd25fOKOvsGt23fqjYyc3lEHA98Hji+3LcAuDgiNgLmA5/LzP9eZ+WSNAaOeJAkSZLaY1L5/FidfY8OOaYlbWTmV4BtgenAjpl5VGYm8FGK/xf4h4h4bkScFxG/jYhflJNNriUi5kbE4ohYvHLlynWUKUnPcMSDJEmS1B4Pl88b1dm38ZBjWtZGZt5PcRsGABHxEuADwOsy89GIuBR4HvBXwL7AtyPi7sz8+ZB25lOMkmDatGm5jjol6WkGD5IkSVJ73FM+17udYnBbvVsoWtZGRKwHnAtckJlXR8R2wOuAQzPzeuD6iHgH8G7g58O1I0mj4a0WkiRJUnvcUD5Pr7NvfyApJoSsso33AC8E3l++3qF8Xl5zzHLgBeuoQ5IaZvAgjWdXn148huFKF1qXiDgpIi6OiKURkRGxbB3H7xcRV0bEHyPiwYj4QUTsPcprNtRGRLwiIn4cEQ9FxK3lb9jqtXdZRFw+mhokaTzKzN8Ai4G3liMNACj//Fbgqsz8fc32rSJi94joa7aNWhHxAuDTwAnl7RfwzAiKqTWHTq3ZLklj5q0WbRQR/UD/lClTOl2KpInjNIr1228CthzpwIjYH7iGYojux8rNxwLXRcQrM/PmdV2s0TYiYnPg+8BvKX7rdhDwzYi4IzNvqGnvrcAhwF7rfquS1BkRcTiwU/lyMrBhRJxSvr4rM8+vOfwE4GqKfvGccttxFL8QfN+Qpo8FTgWOoliNopk2an0J+FFmfntwQ2b+NiKuAb5QhhevoOhz3zPSe5ak0TB4aKPMXAgsnDZt2pxO1yJpwtg1M5cCRMQtwGYjHHs28DhwYGauKM+5CLgVOBN4bQPXa7SNVwLbANMzc1lEzAf2A2ZSDiOOiC3L9k7JzLsbfL+S1AlHAzOGbPtk+Xwt8HTwkJk/iYiDgE+VjwR+Arw1M3/VyMWaaSMi3kYR8tYLct8FfBn4B2AVcHRmXjtMO/4iTdKoGTxIrTZ4a8TBJ3W2DgkYDB3WJSKmAPsAXxsMDMrzV0TExcBREbHNcMN3m2hjk3L36vKYNRHxALBpTZOfBe4GzkGSxrHMPGiUx/8UeHUDx30c+PhY2qg5/iLgomH2/Y4i+G2kHX+RJmnUnONBkgRFYADw0zr7fgYExfDbVrVxI/AE8MmI2CkijgReSvEbOyJiBnAkMCcz1zT6JiRJkjT+OOJBkgQwOEFZvSXYBrfVW7qtqTYyc3lEHA98Hji+3LcAuDgiNqJYJ/5zmfnf66xcDRucjPZE//WXJElt5IgHSRLApPL5sTr7Hh1yTEvayMyvANtSLAm3Y2YelZkJfJTi36d/iIjnRsR5EfHbiPhFOdlkXRExNyIWR8TilStXrqNUSZIktYu/85AkATxcPm9UZ9/GQ45pWRvlcm4/G3wdES8BPgC8LjMfjYhLgecBfwXsC3w7Iu7OzJ8PvUBmzqcYKcG0adNyHbVKkqQmPT2C7tDdOlyJuoXBgyQJnlmvvd7tFIPb6t1C0bI2ImI94Fzggsy8ulzW7XXAoZl5PXB9RLwDeDewVvAgSZKk8clbLSRJUC5hSXHbw1D7UyzVdmPFbbwHeCHw/vL1DuXz8ppjlgMvWEcdkqSKRER/RMwfGBjodCmSuojBQxvZUUsarzLzN8Bi4K3lSAMAyj+/FbiqdinNiNgqInaPiL5m26gVES8APg2cUN5+Ac+MoJhac+jUmu2SpDbLzIWZObevr2/dB0tSyVst2sh1jyW1W0QcDuxUvpwMbBgRp5Sv78rM82sOPwG4GrguIs4ptx1HEVK/b0jTxwKnAkdRrEbRTBu1vgT8KDO/PbghM38bEdcAXyjDi1cAe1GMjJAkSVKXMHiQmtDRCXWuPr14Pvik9l9b3ehoYMaQbZ8sn68Fng4eMvMnEXEQ8KnykcBPgLdm5q8auVgzbUTE24CDKEKFod4FfBn4B2AVcHRmXttILZIkSRofDB6ksRgvIUC9OsZLbeqozDxolMf/FHh1A8d9HPj4WNqoOf4i4KJh9v0OmNloWxHRD/RPmTKl0VMkSZJUMed4kCT1DO89liRJGn8MHqRGXX36M6MIJEmSJEkN8VYLqcc8Pf9ERX+7n2n/u8UGb+WQulrVfYYkSZIjHiRJkiRJUmX8/YbUIhP5t4YdXeVDkiS1jZP4SmqGIx4ktY7zYEiS1NOcxFdSMwweJEmSJElSZQweJEk9IyL6I2L+wMBAp0uRJKknzFu05OnbaqVmGTy0kV+IJalaDgGWJEkafwwe2sgvxJIkSZKkicbgQZpInPxRkiRJw/C2ClXF4EGSJEmSJFXG4EGSJEmSJFXG4EGSJElSQ5wsXVIzDB4kSVJ9zgkjaQgnS5fUjPU7XYCk6g1OEnSif+PV4yKiH+ifMmVKp0uRJGlCefr75qG7dbgSjUeOeJAk9Qx/EydJkjT+GDxIkiRJkqTKGDxIkiRJUtWuPt25czRhGTxIkiRJkqTKGDxIGta8RUueniioaab7kiRJ0oTmHPeSKuFKGpIkaUIY/AXLwSd1tg5pHHPEgyRJkiRJqoy/i2wj15fvTv7mXpIkSRq9p79HH7pbhytRpznioY1cX16SJEnSeNGS+bykBhg8SJIkSWpIRPRHxPyBgYFOlyKpixg8SJJ6hl+IJalajuCV1AyDB0lSz/ALsSRJ0vhj8CBp3a4+/ZmlosbA+wglSZKkicfgQRJgKCBJkiSpGgYPkiRJkiSpMgYP0kTXglsoJEmSJGk4Bg+SJEmSJKky63e6AEmdMTifw4md6AUGR1kcfFIHLi5JkiSpnRzxIEmSJEmSKmPwIEmSJEmSKuOtFpKezdsgpAmvo7diSdJE4vcuTRCOeJAkSZIkSZUxeJAk9YyI6I+I+QMDA50uRZJ6kv2spGYYPEiSekZmLszMuX19fZ0uRZJ6kv2spGYYPEiSJEmSpMoYPEiSJEnSeHH16c9MOin1CIMHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSZJUGYOHNnL5IUmSJEnSRGPw0EYuP9RlnNRHkiRJksbM4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGS1DMioj8i5g8MDHS6FElSl5i3aAnzFi3pdBlSTzN4kCT1jMxcmJlz+/r6Ol2KJEmSSgYPkiRJkiSpMgYPkiR1o6tPLx6S1Ebe0iapGQYPkjrGeyolSeou3tImqRkGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZKknuAcYuOTwYMkSZIkuVqQVBmDB0mSJEmSVBmDB0mSJEmSVBmDB0mSJEmSVBmDB0mSJEmSVBmDB0md52ROkiRJUs8yeJAkSZIkSZUxeJAkSZIkSZUxeJAkSZIkSZUxeJAkSZIkSZUxeJAk9YyI6I+I+QMDA50uRZIkSSWDB0lSz8jMhZk5t6+vr9OlSJKkOuYtWsK8RUt6/pp6NoMHSZIkSZJUGYMHSZIkSZJUma4JHiLipIi4OCKWRkRGxLIm2zkiIn4REY9ExL0R8dWImFznuHdFxK8j4qGIuC4iXlbnmC0iYnlEfLCZWiRJkiRJ6nVdEzwApwGHAHcA9zfTQEScCHwDGABOAP4ZeAdwTURsWnPcfsAFwC+B9wMbA5dHxOZDmjwdWAmc1Uw9kiRJkiT1uvU7XcAo7JqZSwEi4hZgs9GcHBFbAZ8CbgBenZlPldtvAP6NIog4rTx8JrAMeGdmZkT8EFgK7A8sKs97JTAHmJ6ZT47pnUmSJEmS1KO6ZsTDYOgwBjOBScA5g6FD2e5CilBhVs2xmwAPZGaWr1eXz5sCRMSGwLnA2Zl54xjrkiRJE4Qzq0uSJqKuCR5aYJ/y+ad19v0M2D0iNqs5Zu+IODwidqIYKfE4MBgyfJginPhYhfVKkiRJktT1uulWi7HarnxeUWffCiDKY5YAFwGvB84r9z8GnJCZyyNid+Bk4M2Z+XCjF4+IucBcgB133LGpNyBJkiRJKgyOIDvx0N06XInWZSKNeJhUPj9WZ9+jtcdkYTawIzAd2C4zvxIRAcwHLs7MKyJiakRcGRG/i4irImLqcBfPzPmZOS0zp02evNYiGpIkSdK4FxH9ETF/YGCg06VI6iITKXgYHJ2wUZ19Gw85BoDMXJ6ZP8vMwTke5gJ7ACeWK1xcCdwFvJFinogr66x8IUmSJPWEzFyYmXP7+vo6XYqkLjKRgod7yuft6+zbHsiaY9YSEdsCnwH+PjNXAf3AlsBxmXkTcHz5+k2tK1mSJEmSpO42kYKHG8rn6XX27Q/clpkPjXD+F4HrM/P88vUOwOrBeR7K59XAC1pUryRJkiRJXa8ng4eI2DEido+IDWo2XwY8AhwbEc+pObYf2AX45gjtzQReB/xtzeZ7gMkRsXV5zNbAZEYYNSFJkiRJ0kTTNataRMThwE7ly8nAhhFxSvn6rpqRCFCsRjEDeCGwDCAzV0bER4HPUczF8K8Ut1i8D/g18PlhrrsFxWiHUzPzzppdlwMPApdGxPnALGCg3D7ce+gH+qdMmdLgu5YkSZIkqbt104iHo4FPlo/nU8ynMPj66EYayMwzgaOA5wJnA39HsXTmjBFuszgDuBeYN6St+ymW3Fwf+CywIfDGcvtw13cyHkmSJEnShNI1Ix4y86BWHJuZC4AFo2jrmBH2/RzYt9G2JEmSJEmaaLppxIMkSZIkSeoyBg+SJEmSJKkyBg+SJEmSJKkyBg9tFBH9ETF/YGCg06VI49a8RUuYt2hJp8uQJEmS1CIGD23kqhaSJEmSpInG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGS1DNcPUiSJGn8MXhoI78QS1K1XD1IkiRp/DF4aCO/EHeHeYuWMG/Rkk6XIUmaAPw3R5I0EYw6eIiIv6iiEEmSJEmS1HuaGfHwo4j434h4X0RMbnlFkiRJkiSpZzQTPHyofP4s8NuI+E5EvC4iooV1SZrorj69eEiSJEnqaqMOHjLzs5m5J3AA8E3gL4HLgbsi4hMRsXNrS5QkSZIkSd2q6cklM/O/MvPdwLbA3wIrgI8Cv4mIKyLibRGxQYvqlCRJkiRJXWjMq1pk5kOZ+VXgr4ALyjZfA3yL4laMD0TEc8Z6HUmSJEmS1H3GFDxExHoR8eaIuAy4C5gF/Bg4Ang78GvgDOALYy20F0REf0TMHxgY6HQpkiRJkiS1RVPBQ0S8KCLOAH4LXAq8EjgH2DMzD8zMCzLz4sycAfwz8M6WVdzFMnNhZs7t6+vrdCmSJEmSJLXF+qM9ISKuowgaArgWeB/w3cx8fJhTrgP+X9MVSpIkSdI4M2/REgBOHPX/UUkTTzN/TV4MnAXMz8zbGzj+SuDgJq4jSZIkSZK6XDPBw/aZ+USjB2fmSoqREZIkSZIkaYJpZo6HHSKif7id5QSKOzdfkiRJkiRJ6hXNjHj4NPACYOEw+98H3E2xsoUkSZJG4en7xg/drcOVSJLUGs2MePgL4Icj7L8COLC5ciRJkiRJUi9pJnh4PvD7Efb/Adi6uXIkSZIkSVIvaSZ4eADYdYT9U4A/NlWNJEmSJEnqKc0ED9cBcyJim6E7ym1/A/x4rIX1onLizfkDAwOdLkWSJEmSpLZoJnj4NLAZ8IuIeH9EvKZ8vB/4RbnvtFYW2Ssyc2Fmzu3r6+t0KZIkSeqAiDgpIi6OiKURkRGxrMl2joiIX0TEIxFxb0R8NSIm1znuXRHx64h4KCKui4iX1Tlmi4hYHhEfbKYWSVqXUQcPmflL4P9SrIjxjxQTTf6w/PNzgLdm5uIW1ihJkiT1itOAQ4A7gPubaSAiTgS+AQwAJwD/DLwDuCYiNq05bj/gAuCXwPuBjYHLI2LzIU2eDqwEzmqmHklal2aW0yQzvx8ROwJ/Cbyo3LwEuCIzH2lVcZIkSVKP2TUzlwJExC0Uo4UbFhFbAZ8CbgBenZlPldtvAP6NIogYHH08E1gGvDMzMyJ+CCwF9gcWlee9EpgDTM/MJ8f0ziRpGE0FDwBlwPC91pUiSZLUO+YtWgLAiYfu1uFKNJ4Mhg5jMBOYBJwzGDqU7S6MiKXALJ4JHjYBHsjMLF+vLp83BYiIDYFzgbMz88Yx1iVJw2pmjgdJkiRJnbFP+fzTOvt+BuweEZvVHLN3RBweETtRjJR4HBgMGT5MEU58rMJ6Jam5EQ8R8Q7gOIrbLJ5X55DMzKZHU0iSJEmqa7vyeUWdfSuAKI9ZAlwEvB44r9z/GHBCZi6PiN2Bk4E3Z+bD1ZYsaaIbdTgQER8AzgDuo0hV72t1UZIkSZLqmlQ+P1Zn36O1x5S3WMyOiI8C2wNLMnN1RAQwH7g4M6+IiKnAPGAv4FaKcOLmoY1HxFxgLsCOO+7Ywrckqdc1MyrhPcDPKSazcSJJSZIkqX0GRydsBAz9Lr7xkGMAyMzlwPKaTXOBPYC/Kle4uBL4PvBB4BjgyoiYkpl/HNLOfIrAgmnTpiWS1KBm5njYBrjA0EGSJElqu3vK5+3r7NseyJpj1hIR2wKfAf4+M1cB/cCWwHGZeRNwfPn6Ta0rWdJE10zw8BuKzkiSJElSe91QPk+vs29/4LbMfGiE878IXJ+Z55evdwBWD87zUD6vBl7Qonolqang4Uzg6JrZciVJkiS1WETsGBG7R8QGNZsvo7jF4tiIeE7Nsf3ALsA3R2hvJvA64G9rNt8DTI6IrctjtgYmM8KoCUkarWbmeHgK+APw64j4GnBnue1ZMvO8odsmuvIfhP4pU6Z0uhRJkiR1QEQcDuxUvpwMbBgRp5Sv76oZiQDFahQzgBcCywAyc2U5WeTnKOZi+FeKWyzeB/wa+Pww192CYrTDqZl5Z82uy4EHgUsj4nxgFjBQbpeklmgmeFhQ8+dThjkmeWbZHpUycyGwcNq0aXM6XYskSZI64miKMKHWJ8vna4HzWYfMPDMi7gNOBM6mCA4uAj48wm0WZwD3UqxeUdvW/RHxeuAc4LMUq1q8MTPvr9eIv0iT1IxmgoeDW16FJEmSNAFk5kGtODYzF/DsXwiuq61jRtj3c2DfBtvxF2mSRm3UwUNmXltFIZIkSZIkqfc0M7nk0yJio4jYPiI2bFVBkiRJkiSpdzQVPETEyyPiKuCPwN3AX5Tbnx8R/xkRr2lhjZIkSZIkqUuNOniIiL2B64BdGTKBZGb+AdgEOLIVxUmSJEmSpO7WzIiHf6BY13cv4MNADNn/nzQ4OY0kSZIkSeptzQQPBwDnlkv1ZJ39dwPbjakqSZIkSeNORPRHxPyBgYFOlyKpizQTPGwMjNTTbNFkLZKkDouIkyLi4ohYGhEZEcuabOeIiPhFRDwSEfdGxFcjYnKd494VEb+OiIci4rqIeFmdY7aIiOUR8cFmapEktU5mLszMuX19fZ0uRVIXaSZ4uAN4xQj7DwH+t7lypA65+vTiIek0in78DuD+ZhqIiBOBb1CE1CcA/wy8A7gmIjatOW4/4ALgl8D7KYLtyyNi8yFNng6sBM5qph5JkiR11vpNnHMh8NGIuAj4RbktASLifcDrKL5oSpK6z66ZuRQgIm4BNhvNyRGxFfAp4Abg1Zn5VLn9BuDfKP59OK08fCawDHhnZmZE/BBYCuwPLCrPeyUwB5iemU+O6Z1JFZq3aAkAJx66W4crkSRp/GlmxMPngJ8BPwR+RBE6zIuIFcA/UnxZ/FLLKpQktc1g6DAGM4FJwDmDoUPZ7kKKUGFWzbGbAA9k5uB8QavL500BImJD4Fzg7My8cYx1SZIkqUNGHTxk5uPAoRTDYh8BHgV2A1YBHwTelJlrWlmkJKlr7FM+/7TOvp8Bu0fEZjXH7B0Rh0fEThQjJR4HBkOGD1OEEx+rsF5JkiRVrJlbLSiHu84rH5JUKYcwd5XBVY1W1Nm3gmIJ5u2AJcBFwOuB88r9jwEnZObyiNgdOBl4c2Y+XG3JatjgXDgHn9TZOiRJUldp5lYLNcnlhyRNAJPK58fq7Hu09pgszAZ2BKYD22XmVyIigPnAxZl5RURMjYgrI+J3EXFVREytd+GImBsRiyNi8cqVK1v6piRJUmfNW7Tk6V9GqfuMesRDRBzRyHGZed66j5pYynucF06bNm1Op2uRpIoMjk7YiOJ2vFobDzkGgMxcDiyv2TQX2AP4q3KFiyuB71PczncMcGVETMnMPw5pZz5FYMG0adMSSZIkjQvN3GqxgGJCyRiyfeiXPIMHSZp47imftwd+M2Tf9hT/VtzDMCJiW+AzwHGZuSoi3gVsWb5+OCKOBw4H3gT8a4trlyStQ0T0A/1TpkzpdCnqEd5SOzE0c6vFwRRrvB9c8zgU+DvgZuDnwGtaVaAkqavcUD5Pr7Nvf+C2zHxohPO/CFyfmeeXr3cAVg/O81A+rwZe0KJ6JUmjkJkLM3NuX19fp0uR1EVGPeIhM68dZtd/RsQ3gOuBlwNXj6UwSdL4FhE7UszXcEdmPlFuvgw4Gzg2Ii4cXFKz/A3ZLsBHR2hvJvA64CU1m+8BJkfE1pl5b0RsDUxmhFETkiRJGl9aOrlkZj4GXEBxD64kqcuUS1ueEhGnUPwPft/g64g4fMjh5wG3UtxCAUBmrqQIF/almIthbkR8guK2iF8Dnx/multQjHY4NTPvrNl1OfAgcGlE/B1wCTBQbpckSWqKE1W2V1PLaa7DY9R8CZUkdZWjgRlDtn2yfL4WOJ91yMwzI+I+4ESK0Q8PUiyd+eERbrM4A7iXIcs0Z+b9EfF64BzgsxRBxxsz8/7G3o4kSZI6raXBQzkp2P8D7lzXsZKk8SczD2rFsZm5gGIy4kbbGnakXGb+nGIExTo56ZkkSd3PCSd7TzPLaV41zK7nArsDGwJHjqUoSZKa4bLFkiRJ408zIx52Ye2lM5NilvFLgC9m5k/GWpgkSZIkSep+zaxqsXMFdUiSJEmSpB7U0lUtJEmSJPWuiOiPiPkDAwOdLkVSFzF4kCRJktSQzFyYmXP7+vo6XYqkLtLM5JJrWHuOh3XJzKxi6U5JkiRJkjSONRMGnAe8HHgJcBvFmuoAewK7ATcDN7WkOkmSJEmS1NWaCR6+CbwFmJmZ/1a7IyJmAucD78vMK8deniRJjYuIfqB/ypQpnS5FkiRJpWbmePgk8M9DQweAzPweMB/41BjrkiRp1Lz3WJIkafxpJnj4c+COEfb/BpjaXDmSJEmSJKmXNBM83A+8doT9rwNcX0eSJEmSJDUVPFwIHBYR/xIRe0TEc8rHHhHxNeBNFPNASJIkSZKkCa6ZySVPAaYARwGzgTXl9vWAABaWx0iSJEmSpAlu1MFDZj4G/J+IeC0wE3hhuWspcFlmXtG68nqLs61LkiRJkiaaZkY8AFAGDIYMo5CZC4GF06ZNm9PpWiRJkiRJaodm5nh4WkRMiYhXRYTrlkmq3tWnFw9pGBHRHxHzBwac41iSqmA/K6kZTQUPEfGmiLgDuA34EfCKcvvzI+I3EfF/W1ijJEkNycyFmTm3r888XJKqYD8rqRmjDh4i4iDgUmA18AmKCSUByMw/AHcA72hNeZIkSZIkqZs1M+LhY8CvgP2Af6qz/6fAy8dSlCRJkiRJ6g3NBA/7AN/MzDXD7P8tsE3zJUmSJEmSpF7RTPCwHvDYCPu3Ah5vrhxJkiRJktRLmgkebgUOGGH/myhuxZAkSZIkSRNcM8HDvwD/NyKOrjk/I2JSRJwNTAfmt6pASZIkSZLUvdYf7QmZ+eWIeBVwLnAmkMC/As8DngN8PTO/2dIqJUmSJElSVxp18ACQmbMi4rvALGB3iiU1fw6cl5nfbWF9kiRJkiSpi40qeIiITYC3Ardl5qXApZVUJUmSJEmSesJo53h4jOIWi5dVUIskSWMSEf0RMX9gYKDTpUiSJKk0quAhM9cAy4EtqilHkqTmZebCzJzb19fX6VKklpq3aAnzFi3pdBmSJDWlmVUtvgEcHhEbtboYSZIkSeOXI8skNaOZySV/AvwV8MuI+BJwO/Dw0IMy80djrE2SJEnSOJKZC4GF06ZNm9PpWiR1j2aCh0U1f/4CxXKataLc9pxmi5IkSZIkSb2hoeAhIvYFfpOZq4Gjqi1JkiRJkiT1ikZHPPwUOBy4MDO/ERGbAfOBT2Xm/1ZWnSRJkiRpzAYnqD3x0N06XIkmokYnl4whrzcC3g5s09pyJEmSJElSL2lmVYtBQ8MISZIkSZKkZxlL8CBJkiRJkjQigwdJkiRJklSZ0Syn+YaIGJzTYRLFkplvjYi96xybmTlvrMVJkiRJkqTuNprg4V3lo9bfDnNsAgYPkqS2ioh+oH/KlCmdLkWSJEmlRoOHgyutQpKkFsjMhcDCadOmzel0LZIkSSo0FDxk5rVVFyJJkiRJknqPk0tKkiRJkqTKGDxIkiRJkqTKGDxIkiRJkqTKdFXwEBHrRcSJEfHriHg0IpZHxJkRsWmD558UERdHxNKIyIhYNsKx7yqv81BEXBcRL6tzzBZlDR8cw9uSJEmSJKlndVXwQLFE51nA/wLHARcDxwMLI6KR93IacAhwB3D/cAdFxH7ABcAvgfcDGwOXR8TmQw49HVhZ1iRJkiRJkoZodDnNjouIvSjChksy8y012+8EzgbeAVy4jmZ2zcyl5Xm3AJsNc9xMYBnwzszMiPghsBTYH1hUnv9KYA4wPTOfbPJtSZIkSV0jIvqB/ilTpnS6FEldpJtGPLwTCODzQ7afCzwMzFpXA4OhQwM2AR7IzCxfry6fNwWIiA3L656dmTc22KYkSZLU1TJzYWbO7evr63QpkrpINwUP+wBrgOtrN2bmoxS3ROzTwmv9FNg7Ig6PiJ2ATwGPA4Mhw4cpwomPtfCakiRJkiT1nK651QLYDliVmY/V2bcCeGVEbJiZj7fgWhcBrwfOK18/BpyQmcsjYnfgZODNmflwC64lSZIkSVLP6qYRD5MoAoB6Hq05ZsyyMBvYEZgObJeZX4mIAOYDF2fmFRExNSKujIjfRcRVETF1uDYjYm5ELI6IxStXrmxFmZIkSZIkjXvdFDw8DGw0zL6Na45pmcxcnpk/y8zBOR7mAnsAJ5YrXFwJ3AW8kWLyySvrrHwx2Nb8zJyWmdMmT57cyjIlSZIkSRq3uil4uAfYKiLqhQ/bU9yG0YrbLOqKiG2BzwB/n5mrgH5gS+C4zLyJYlnPLYE3VVWDJEmSJEndppuChxso6t23dmNEbAzsDSyu+PpfBK7PzPPL1zsAqwfneSifVwMvqLgOSdIwIqI/IuYPDAx0uhRJkiSVuil4+DaQwHuHbJ9DMbfDNwc3RMSu5SSQLRERM4HXAX9bs/keYHJEbF0eszUwudwuqWLzFi1h3qIlcPXpxUPCZd4kSZLGo65Z1SIzb46IfwKOjYhLgH+nmG/heOBa4MKaw/8T2AmI2jYi4vByOxQhwYYRcUr5+q6a0Qy152xBMdrh1My8s2bX5cCDwKURcT4wCxgot9cVEf1A/5QpUxp705IkSZIkdbmuCR5K7wWWUUzy+EZgFXAO8LHMXNPA+UcDM4Zs+2T5fC2wVvAAnAHcC8yr3ZiZ90fE68vrfxa4FXhjZt4/3MUzcyGwcNq0aXMaqFWSJEmSpK7XVcFDZj4FnFk+Rjpu52G2H9TENY8ZYd/PGTLnhCRJkiRJekY3zfEgSZLUcU/PMSNJkhpi8CBJkiRJkipj8CBJkiRJkipj8NBGri8vSeo5LmkrSZLWweChjVxfXpIkSZI00Rg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8tJHLaUqSJEmSJhqDhzZyOU1JkiRJ0kRj8CBJkiRJkipj8CBJkiSpId46LKkZBg+SpJ7hF2JJqpa3DktqhsGDJKln+IVYkiRp/DF4kCRJkiRJlTF4kCRJkiRJlTF4aCPvPZYkSZIkTTQGD23kvceSJEmSpInG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFVm/U4XIEmSxrd5i5YAcKLfGiRJUhP8CiFJkiRJLWBQK9XnrRZtFBH9ETF/YGCg06VIkiRJktQWBg9tlJkLM3NuX19fp0uRJEmSJKktDB4kSZIkSVJlDB4kSZIkSVJlDB4kSZIkSVJlnG9VUs94Zibp7xYbDj6pg9VIkiRJAkc8SJIkSZKkChk8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8tFFE9EfE/IGBgU6XIkmSJElSWxg8tFFmLszMuX19fZ0uRZJ6kgGvJElqxrxFS5i3aEmny+hZBg+SpJ5hwCtJkjT+GDxIkiRJkqTKGDxIkiRJkqTKGDxIkiRJkqTKGDxIkiRJkqTKGDxIkiRJkqTKGDxIkiRJkqTKGDxIkiRJkjTEvEVLmLdoSafL6AkGD5IkSZIkqTIGD5IkSVKbRMR6EXFiRPw6Ih6NiOURcWZEbNrg+SdFxMURsTQiMiKWjXDsu8rrPBQR10XEy+ocs0VZwwfH8LYkaUQGD5IkSVL7zAPOAv4XOA64GDgeWBgRjXw3Pw04BLgDuH+4gyJiP+AC4JfA+4GNgcsjYvMhh54OrCxrkqRKrN/pAiRJkqSJICL2oggbLsnMt9RsvxM4G3gHcOE6mtk1M5eW590CbDbMcTOBZcA7MzMj4ofAUmB/YFF5/iuBOcD0zHyyybclSevkiAdJkiSpPd4JBPD5IdvPBR4GZq2rgcHQoQGbAA9kZpavV5fPmwJExIbldc/OzBsbbFOSmmLw0EYR0R8R8wcGBjpdiiRJktpvH2ANcH3txsx8lOKWiH1aeK2fAntHxOERsRPwKeBxYDBk+DBFOPGxFl5TkuoyeGijzFyYmXP7+vo6XYokSZLabztgVWY+VmffCmCrciRCK1wEnFc+llHcUnFCZi6PiN2Bk4H/l5kPN9pgRMyNiMURsXjlypUtKlPSRGDwIEmSJLXHJKBe6ADwaM0xY5aF2cCOwHRgu8z8SkQEMB+4ODOviIipEXFlRPwuIq6KiKkjtDk/M6dl5rTJkye3okxJE4TBgyRJktQeDwMbDbNv45pjWiYzl2fmzzJzcI6HucAewInlChdXAncBb6SYfPLKOitfSNKYGDxIkiRJ7XEPxe0U9cKH7Sluw3i8qotHxLbAZ4C/z8xVQD+wJXBcZt5EsaznlsCbqqpB0sRk8CBJkiS1xw0U37/3rd0YERsDewOLK77+F4HrM/P88vUOwOrBeR7K59XACyquQ9IEY/AgSZIktce3gQTeO2T7HIq5Hb45uCEidi0ngWyJiJgJvA7425rN9wCTI2Lr8pitgcnldklqmfU7XYAkSZI0EWTmzRHxT8CxEXEJ8O8U8y0cD1wLXFhz+H8COwFR20ZEHF5uhyIk2DAiTilf31UzmqH2nC0oRjucmpl31uy6HHgQuDQizgdmAQPldklqGYMHSZIkqX3eS7G85VyKCR1XAecAH8vMNQ2cfzQwY8i2T5bP1wJrBQ/AGcC9wLzajZl5f0S8vrz+Z4FbgTdm5v3DXTwi+oH+KVOmNFCqJBUMHiRJkqQ2ycyngDPLx0jH7TzM9oOauOYxI+z7OUPmnFhHWwuBhdOmTZsz2jokTVzO8SBJkiRJkipj8CBJkiRJkipj8CBJkiRJkipj8CBJkiRJkipj8CBpQpm3aAnzFi3pdBmSJEnShGHwIEmSJKkhEdEfEfMHBgY6XYqkLmLwIEmSJKkhmbkwM+f29fV1uhRJXcTgQZIkSZIkVcbgQZIkSZIkVcbgQZIkSZIkVcbgQZIkSZIkVcbgoY2cBViSJLVKveWBXTJYkjQeGTy0kbMAS5IkSZImGoMHSZIkSZJUGYMHSZIkSQ3x1mFJzTB4kCRJktQQbx2W1AyDB0nS0yJivYg4MSJ+HRGPRsTyiDgzIjZt8PyTIuLiiFgaERkRy0Y49l3ldR6KiOsi4mV1jtmirOGDY3hbkiRJ6iCDB0lSrXnAWcD/AscBFwPHAwsjopF/M04DDgHuAO4f7qCI2A+4APgl8H5gY+DyiNh8yKGnAyvLmiRJktSF1u90AZKk8SEi9qIIGy7JzLfUbL8TOBt4B3DhOprZNTOXlufdAmw2zHEzgWXAOzMzI+KHwFJgf2BRef4rgTnA9Mx8ssm3JUmSpA5zxIMkadA7gQA+P2T7ucDDwKx1NTAYOjRgE+CBzMzy9eryeVOAiNiwvO7ZmXljg21KkiRpHDJ4kCQN2gdYA1xfuzEzH6W4JWKfFl7rp8DeEXF4ROwEfAp4HBgMGT5MEU58rIXXlCRJUgd4q4UkadB2wKrMfKzOvhXAKyNiw8x8vAXXugh4PXBe+fox4ITMXB4RuwMnA2/OzIdbcC1JkiR1kCMeJEmDJlEEAPU8WnPMmGVhNrAjMB3YLjO/EhEBzAcuzswrImJqRFwZEb+LiKsiYupwbUbE3IhYHBGLV65c2YoyJUmS1AIGD5KkQQ8DGw2zb+OaY1omM5dn5s8yc3COh7nAHsCJ5QoXVwJ3AW+kmHzyyjorXwy2NT8zp2XmtMmTJ7eyTElSKSL6I2L+wMBAp0uR1EUMHiRJg+4BtoqIeuHD9hS3YbTiNou6ImJb4DPA32fmKqAf2BI4LjNvoljWc0vgTVXVIEkaWWYuzMy5fX19nS5FUhcxeJAkDbqB4t+FfWs3RsTGwN7A4oqv/0Xg+sw8v3y9A7B6cJ6H8nk18IKK65AkSVILGTxIkgZ9G0jgvUO2z6GY2+GbgxsiYtdyEsiWiIiZwOuAv63ZfA8wOSK2Lo/ZGphcbpckSVKXcFULSRIAmXlzRPwTcGxEXAL8O8V8C8cD1wIX1hz+n8BOQNS2ERGHl9uhCAk2jIhTytd31YxmqD1nC4rRDqdm5p01uy4HHgQujYjzgVnAQLldkiRJXcLgQZJU673AMopJHt8IrALOAT6WmWsaOP9oYMaQbZ8sn68F1goegDOAe4F5tRsz8/6IeH15/c8CtwJvzMz7G3kjkiRJGh8MHjShzVu0BIAT/ZsgAZCZTwFnlo+Rjtt5mO0HNXHNY0bY93OGzDkxkojoB/qnTJky2jIkSZJUEed4kCT1jIkw2/q8RUueDk0lSZK6gcGDJEmSJEmqjMGDJEmSJEmqjMGDJEmSJEmqjMGDJEmSpIZERH9EzB8YGOh0KZK6iMGDJEmSpIZMhEl8JbWewYMkSZIkSaqMwYMkqWc4BFiSJGn8MXiQJPUMhwBLkiSNPwYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkiRJkiSpMgYPkqSeERH9ETF/YGCg06VIkiSpZPAgSeoZmbkwM+f29fV1uhRJkiSVDB4kSZIkSVJlDB7aJCLmRsTiiFi8cuXKTpcjTRjzFi1h3qIlnS5DkqSe4C1tkpph8NAmmTk/M6dl5rTJkyd3uhxJkiRp1LylTVIzDB4kSZIkSVJlDB4kSZIkSVJlDB4kSZIkSVJlDB7ayMl4JEmSJEkTjcFDGzkZjyRJkiRpojF4kCRJkiRJlTF4kCT1DG9pkyRJGn8MHiRJPcNb2tRq8xYtYd6iJZ0uQ5KkrmbwIEmSJEmSKmPwIEmSJEmSKmPw0EbeeyxJkiRJmmgMHtrIe48lSZIkSRONwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSarM+p0uQJLa4urTyz+8paNlSJLUzSKiH+ifMmVKp0uR1EUc8dBGLqcpSZKkbuYqbZKaYfDQRnbUklQtA15JkqTxx+BBktQzDHglSZLGH4MHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSZJUGYMHSZIkSZJUGYOHNnKZN0mSJEnSRGPw0EYu8yZJkiRJmmgMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJPSMi+iNi/sDAQKdLkSRJUsngoY38QixJ1crMhZk5t6+vr9OlSJIkqWTw0EZ+IZYkSZIkTTQGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIa4vLwkpph8CBJkiSpIS4PL6kZBg+SJEmSJKkyBg+SJEmSJKkyBg+SJEmSJKkyBg9t5GQ8kiRJkqSJxuChjZyMR5IkSZI00Rg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJJ6hqsHSZIkjT8GD5KknuHqQe0zb9ES5i1a0uky1IBWflZ+7pKkZhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8SJIkSZKkyhg8tFFE9EfE/IGBgU6XIkmSJElSWxg8tFFmLszMuX19fZ0uRZIkSZKktjB4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCT1jIjoj4j5AwMDnS5FkiRJJYMHSVLPyMyFmTm3r6+v06VIkiSpZPAgSZIkSZIqY/AgSZIkSZIqY/DQRt57LEmSJEmaaAwe2sh7jyVJktStImJuRCyOiMUrV67sdDmSuojBgyRJkqR1ysz5mTktM6dNnjy50+VI6iIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIkqTIGD5IkSZIaEhH9ETF/YGCg06VI6iIGD5IkSZIakpkLM3NuX19fp0uR1EUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiRJUmUMHiRJkiQ1xOU0pcbMW7SEeYuWdLqMccPgQZIkSVJDXE5TUjMMHiRJkiRJUmUMHiRpDBxGJ0mSJI3M4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJEmSJFXG4EGSJElSQyKiPyLmDwwMdLoUSV3E4EGSJElSQzJzYWbO7evr63QpkrqIwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZKkhkREf0TMHxgY6HQpkrqIwYMkSZKkhmTmwsyc29fX1+lSJHURgwdJkiRJklQZgwdJkiRJklQZgwdJkiRJklQZgwdJkiRJklQZgwdJkiRJklQZgwdJkiRJklQZgwdJkiRJklSZyMxO1zBhREQ/0A+8Hbi9wdO2AlZVVpQ6yc+2N43mc90pMydXWcxEEhFzgbnlyxcDtzV4qn8Xe5efbW8a7edqX1uBiFgJ3DWKU/z72Jv8XHtXS7/TGjyMcxGxODOndboOtZ6fbW/yc+0+fma9y8+2N/m5dic/t97k59q7Wv3ZequFJEmSJEmqjMGDJEmSJEmqjMHD+De/0wWoMn62vcnPtfv4mfUuP9ve5OfanfzcepOfa+9q6WfrHA+SJEmSJKkyjniQJEmSJEmVMXiQJEmSJEmVMXhos4hYLyJOjIhfR8SjEbE8Is6MiE3bcb6q04LPNod5PFR17RpeRJwUERdHxNLy81jWZDtHRMQvIuKRiLg3Ir4aEa4rXwH72d5lP9ub7Ge7k31tb7Kf7U3joZ91joc2i4gvAMcDlwL/AewBHAdcB7wmM9dUeb6q04LPNstjh07k8kRmfrv1FasR5eeyGrgJeAXwYGbuPMo2TgTOAq4FLgR2AP4euAvYNzP/1MqaJzr72d5lP9ub7Ge7k31tb7Kf7U3jop/NTB9tegB7AWuA7w7ZfhyQwLuqPN/H+P1sy2MTWNDp9+Jjrc9ll5o/3wIsG+X5WwF/Aq4HnlOzvb/8zE/u9HvspYf9bO8+7Gd792E/230P+9refNjP9u5jPPSz3mrRXu8EAvj8kO3nAg8Dsyo+X9Vp2WcTERtGxGatK01jkZlLx9jETGAScE5mPlXT7kJgKf69bTX72d5lP9uj7Ge7kn1tb7Kf7VHjoZ81eGivfShSxOtrN2bmo8Avy/1Vnq/qtOqz+b8UHfsfI+IPEXFORPS1slC13eBn/9M6+34G7O4/zC1lP9u77Gc1HPvZ9rOv7U32sxrOmPtZg4f22g5YlZmP1dm3AtgqIjas8HxVpxWfzfXAxyk66yOBq4Bjgev8wtTVtiufV9TZt4LiNwvb1dmn5tjP9i77WQ3Hfrb97Gt7k/2shjPmfnb9VlekEU0C6v1FBni05pjHKzpf1RnzZ5OZ+w3ZdF5E/DfwaeCE8lndZ1L5XO+/j0eHHKOxs5/tXfazGo79bPvZ1/Ym+1kNZ8z9rCMe2uthYKNh9m1cc0xV56s6VX02n6Xo3N/YTFEaFwY/93r/ffj3tvXsZ3uX/ayGYz/bfva1vcl+VsMZcz9r8NBe91AMUar3gW1PMbRppGR3rOerOpV8Npn5xGDbY6xPnXNP+bx9nX3bU8wEfE+dfWqO/Wzvsp/VcOxn28++tjfZz2o4Y+5nDR7a6waKn/m+tRsjYmNgb2BxxeerOpV8NuX5OwD3jrE+dc4N5fP0Ovv2B27LzIfaWE+vs5/tXfazGo79bPvZ1/Ym+1kNZ8z9rMFDe32bIg1675Dtcyjuifnm4IaI2DUidm/2fLXdmD7biHjeMO1+kmIuloUtq1SViYgdI2L3iNigZvNlwCPAsRHxnJpj+4Fd8O9tq9nP9i77WdnPjh/2tb3JflaV9bORmS0vVsOLiHMoZna9FPh3YA/geOC/gEMyc0153DJgp8yMZs5X+43ls42IeRRp4dXA3cBmwBuAg4GfAwdn5iNtezN6WkQcDuxUvjwO2BA4s3x9V2aeX3PsNcAM4IWZuaxm+/uAzwHXAP9KMSTtfcByYB9/E9da9rO9y362N9nPdif72t5kP9ubxkU/m5k+2vgAnlN+QLdRzAq6AjgL2GzIccuKj6e5831012cLHAb8sDznUeBPFOslnwxs3On3NpEfZeeawzyuGebYneu0Mxv4Vfn5/gH4GvD8Tr+/XnzYz/buw362Nx/2s935sK/tzYf9bG8+xkM/64gHSZIkSZJUGed4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF4kCRJkiRJlTF40LgXEQsiIjtdx3gQEcsi4poK2j0oIjIiZtds27nc9vFWX69sf3bZ/kFVtC+pcfazz7CflVQF+9ln2M9OTAYPWqfyL1Ojj507XW8VImLLiPi4nUq1ImLv8ue8c6drkdrJftZ+tl3sZzVR2c/az7aL/Wx963e6AHWFw4e8PgCYC8wHrhuyb2UF158D/L8K2h2NLYFTyz9f07kyeDFQRVr+I2AT4IkK2h7O+cC3gMdrtu1N8XO+BljWxlqkTrOftZ+tgv2s9Az7WfvZKtjPNsjgQeuUmRfUvo6I9Sk66p8O3VfR9Z+gvR3IuJWZj1XU7hrg0SraHioiNs/MP2bmU8BT7bimNN7Zz44f9rNSb7KfHT/sZycmb7VQS0TEehHxkYj4UUT8PiIej4i7I+LLEfG8Icc+fa9VRLwpIm6IiEcj4ncR8dnyH4La49e6J25wW0Q8r/zzqoj4Y0R8LyK2KY+ZGxG3lm3/OiIOq1P3+hHxoYj43/K4+yLi0oiYWnPMQcCd5ctTa4bhLRtNO+v4+W1c/jxui4iHI+KBiLg5Ij475Li17okb3BYRL42IKyPioYj4Q0ScWda1cUR8LiJWlLX9KCL2GNLGWvfEjVDrMRFxRdne4+XndkG94WRlmwsi4tUR8eOIeAhYWO571j1xUdx79/Xy1Ktrfs4LIuL/lH+eM0xN/xMRv4mIWFf9Ureyn7WftZ+VqmU/az9rP1sdRzyoVTYEPgB8F7gM+BOwD3A08BcR8YrMfHzIOW8AjgG+AnwNOAx4P3A/cFqD1/0B8FvgY8AU4Hjg0oi4hCLF/heK5PN44DsRsVtm3llz/jeBtwGLgC8D2wDvAX4aEQdk5i+AW4ETgXnApcAl5bkPjbKdkfwT8G7gPOAsir+bLwIOafDnsEN57W8D3wFeC/w98CSwF8WwszOArSh+xt+LiD3KZHi03g/8DDgbWA28BPgb4JCImJqZ9w05fhrwFuBc4BsjtHsJsC3F53Yaxc8d4A7gBuD3FD+jc2tPioj9gT2Bj2Smkzapl9nP2s/az0rVsp+1n7WfrUpm+vAxqgcwm+K+rNk12wLYpM6xR5fHvq1m287ltj8BOw9p4xbgd0PaWFD8p7r2NuCfhmw/q9x+N7BFzfY/L7efXrPt0HLbt4Go2f5Sig7uujo1f7zOe2y4nRF+pquBf2/guGXANXW2JfDWIdtvBNZQ/MNZW9fx5fF/WbPtoDqfad33DGxap65Xl8d+cMj2LB+vGeG/o4NG2laz77Ry355Dtp9b/py36/TfDR8+WvWwn13rPdrP2s/68NHSh/3sWu/RftZ+ttKHt1qoJbLwCEBEPCeKWXO3Aq4qD9mvzmnfy8xltW0AVwPbRMRmDV7680NeD04OdF5mPljT9n8DD1KkroP+T/n86fLag8f+imL41F9ExOQGamhFOwPAXhHxkgauV8+KzLx4yLYfU/zjd05tXTzzM3oRTcjMP8HTwxH7ys/5VxTvod7n/KvMvLKZaw1xLkVHffTghojYFHg78B+ZeU8LriGNW/azY27Hfnbd7Gc1odnPjrkd+9l1m7D9rMGDWiYi3hYRPwceoRhethJYWu7+szqnLK2zbXBY0/Pq7KtnaBv3l8931jn2/iHtvpAiQb21zrH/U3PMurSinfdS/Ixujog7IuKrEXFYRDT6d3S491tv3+D2Rn/GzxIRh5T35f0JeIDic14J9FH/c17SzHWGymJI4ZXA4RGxQbn5bcDmwFdbcQ1pvLOftZ/FflaqlP2s/Sz2s5UweFBLRMRfUQzNAjgB6KcYsvW6clu9/9ZGmgG2oUlVsphJtp7hto/LyVoy8zKKoWCHU6Tqrwa+B1wTERs20MRIP8uW/SwiYh/gCop7/j5McR/jayk+6/uo/zk/PNrrjGA+MBl4c/n6aIp75S5v4TWkccl+dmzsZxtmP6sJy352bOxnGzYh+1knl1SrHE4x6c3Bmfn0X8yI2L1zJa3TUoqOZQ/gv4fs27N8HkxXk+GNpp1hZeZq4ALggnI22zOAD1J0hkOHnXXKu4DnAK/PmkmNyiFi9dLh0Rrp5wzF/X1/AI6OiFuAVwGfycwnW3Btabyzn7WftZ+VqmU/az9rP1sRRzyoVZ6i+Ev29H9TZWdzSscqWrfvlc8n1S5bU96X9mbgx5m5stw8OOPvc8fYzloG7yGs3Vbew/aLEa7ZKYNp89B0+WRa05+M9HMmizWwFwB/CZxabv6XFlxX6gb2s/az9rNStexn7WftZyviiAe1yncolpi5KiLOAzYAZgKTOlnUSDJzUURcBLwD+LOI+D7PLBs0uGTR4LH3RcRvgHdExB3AvcCfMnPhaNoZxubA7yLi3yg65z9Q3EP3dxT3ry1s2Zseu0splmL694iYDzxOMSztz4FVLWj/Bor7Cz8SEX9Gcd/dnZn585pjzqVY6uqdwLWZeXsLrit1A/tZ+1n7Wala9rP2s/azFXHEg1oiM79FsV7tZsDnKIZU3UaR5I1nf01xb9dOwJkU6zBfC0zPtdcq/mvgdoplcP4VOKfJdoZ6mGI24xdSdEBfphjq92/AfuNpdtvM/C+Kf5D/BHwS+DjF5Eszym1jbf9uirWNN6H4OfwrxT9Ytcf8hmK2aJgA6bA0yH7Wfhb7WalS9rP2s9jPViaevSqJJI1/EfHvwHSKtY4f6XQ9ktRr7GclqVoTrZ91xIOkrhIRUyh+83DBROikJand7GclqVoTsZ91xIOkrhAR+1HMtHx8+bxHZi7raFGS1EPsZyWpWhO5n3XEg6Ru8XfA14AtgL+eKJ20JLWR/awkVWvC9rOOeJAkSZIkSZVxxIMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSaqMwYMkSZIkSarM/wcSJKj9cnS++wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 1080x720 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axes = plt.subplots(nrows=nrows, ncols=ncols, figsize=(5*ncols, 10*nrows), sharex=True)\n",
    "for c_idx, chid in enumerate(chids):\n",
    "    ax = axes[c_idx]\n",
    "    assay_file ='assays/processed/' + str(chid) + '.csv'\n",
    "    df = pd.read_csv(assay_file)    \n",
    "    smiles_1, smiles_2, y1, y2, smiles_test, scores_1, scores_2 = fit_clfs(df, topliss_augmented=True)\n",
    "    top_percentile = np.percentile(scores_1, 95)\n",
    "    similarities_1 = get_sims(np.array(smiles_test)[np.where(scores_1>top_percentile)[0]], np.array(smiles_1))\n",
    "    similarities_2 = get_sims(np.array(smiles_test)[np.where(scores_1>top_percentile)[0]], np.array(smiles_2))        \n",
    "    alpha=0.5\n",
    "    bins = np.linspace(0,1,50)\n",
    "    ax.hist(x= ([1-x for x in similarities_1], [1-x for x in similarities_2]), bins=bins, alpha=alpha, label=['Split 1', 'Split 2'], \n",
    "            weights=[np.ones(len(similarities_1)) / len(similarities_1), np.ones(len(similarities_2)) / len(similarities_2)], log=True)\n",
    "    ax.set_yscale('log')\n",
    "    \n",
    "    \n",
    "    ax.set_title(chid2name[chid])\n",
    "    low, high = ax.set_ylim(0, 1)\n",
    "    ax.set_yticks([0.0001, 0.001, 0.01, 0.1])\n",
    "    ax.set_yticklabels([\"\", \"0.1%\", \"1.0%\", \"10.0%\"])\n",
    "    \n",
    "    ax.set_title(chid2name[chid])\n",
    "\n",
    "\n",
    "axes[0].legend()\n",
    "\n",
    "axes[0].set_ylabel('Frequency')\n",
    "    \n",
    "for i in range(ncols):\n",
    "    axes[i].set_xlabel('Tanimoto similarity')\n",
    "\n",
    "xtext = -0.4\n",
    "plt.tight_layout()\n",
    "plt.savefig('figures/similarities_top_5.pdf', bbox_inches='tight')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.7"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}