{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Predicting T cell Epitopes in cockroach antigens\n",
    "\n",
    "Bla g allergens are major targets of IgE responses associated with cockroach allergies. Epitopes in cockroach antigens have been studied by the group of Alessandro Sette.\n",
    "This examples shows how we can use epitopepredict to compare their experimental results to the predictions. Twenty-five distinct Bla g regions are recognized by T cell responses from allergic individuals Ag\n",
    "\n",
    "This workflow computes sequence regions with high bovine T cell epitope potential from selected antigens. This is done by finding promiscuous binders (present in multiple alleles) using NetMHCIIpan class II binding prediction methods. Cutoffs are allele-specific. NetMHCIIpan (Nielsen et al., 2008) is an artificial neural network algorithm trained on binding data for multiple MHC-II alleles. MHC-II alleles used are \n",
    "\n",
    "## References\n",
    "\n",
    "* C. Oseroff et al., “Analysis of T Cell Responses to the Major Allergens from German Cockroach: Epitope Specificity and Relationship to IgE Production,” J. Immunol., vol. 189, no. 2, pp. 679–688, 2012. http://www.jimmunol.org/content/189/2/679\n",
    "*  Oseroff C., J. Sidney, et al., \"Molecular determinants of T cell epitope recognition to the common Timothy grass allergen.\" J. Immunol. 2010. http://www.jimmunol.org/content/185/2/943"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os, math, time, pickle, subprocess\n",
    "from importlib import reload\n",
    "from collections import OrderedDict\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "pd.set_option('display.width', 130)\n",
    "import epitopepredict as ep\n",
    "from epitopepredict import base, sequtils, plotting, peptutils, analysis\n",
    "from IPython.display import display, HTML, Image\n",
    "%matplotlib inline\n",
    "import matplotlib as mpl\n",
    "import pylab as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  Alleles\n",
    "Each peptide was predicted for the capacity to bind to a panel of 20 HLA class II alleles (DPA1*0103/DPB1*0201, DPA1*0201/DPB1*0101, DPA1*0201/DPB1*0501, DPA1*0301/DPB1*0402, DQA1*0101/DQB1*0501, DQA1*0301/DQB1*0302, DQA1*0401/DQB1*0402, DQA1*0501/DQB1*0301, DRB1*0101, DRB1*0301, DRB1*0401, DRB1*0405, DRB1*0701, DRB1*0802, DRB1*1101, DRB1*1302, DRB1*1501, DRB3*0101, DRB4*0101, DRB5*0101) using the consensus prediction described by Wang et al. (35). Peptides with predicted binding scores in the top 20% for a given allele were considered potential binders, and the number of HLA molecules each peptide was predicted to bind was enumerated. All peptides predicted to bind seven or more HLA molecules were selected for synthesis and further study."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "HLA-DPA10103-DPB10201\n",
      "HLA-DPA10201-DPB10101\n",
      "HLA-DPA10201-DPB10501\n",
      "HLA-DPA10301-DPB10402\n",
      "HLA-DQA0101-DQB10501\n",
      "HLA-DQA10301-DQB10302\n",
      "HLA-DQA10401-DQB10402\n",
      "HLA-DQA10501-DQB10301\n",
      "HLA-DRB1*0101\n",
      "HLA-DRB1*0301\n",
      "HLA-DRB1*0401\n",
      "HLA-DRB1*0405\n",
      "HLA-DRB1*0701\n",
      "HLA-DRB1*0802\n",
      "HLA-DRB1*1101\n",
      "HLA-DRB1*1302\n",
      "HLA-DRB1*1501\n",
      "HLA-DRB3*0101\n",
      "HLA-DRB4*0101\n",
      "HLA-DRB5*0101\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>locus_tag</th>\n",
       "      <th>translation</th>\n",
       "      <th>description</th>\n",
       "      <th>type</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>NAIEFLNNIHDLLGIPHIPVTARKHHRRGVGITGLIDDIIAILPVD...</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>MIGLKLVTVLFAVATITHAAELQRVPLYKLVHVFINTQYAGITKIG...</td>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Bla-g-4</td>\n",
       "      <td>AVLALCATDTLANEDCFRHESLVPNLDYERFRGSWIIAAGTSEALT...</td>\n",
       "      <td>Bla-g-4</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>MAPSYKLTYCPVKALGEPIRFLLSYGEKDFEDYRFQEGDWPNLKPS...</td>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>Bla-g-6.0101</td>\n",
       "      <td>MDELPPEQIQLLKKAFDAFDREKKGCISTEMVGTILEMLGHRLDDD...</td>\n",
       "      <td>Bla-g-6.0101</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>Bla-g-6.0201</td>\n",
       "      <td>MDEIPAEQVVLLKKAFDAFDREKKGCISTEMVGTILEMLGTRLDQD...</td>\n",
       "      <td>Bla-g-6.0201</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>Bla-g-6.0301</td>\n",
       "      <td>MADEQLQLPPEQISVLRKAFDAFDREKSGSISTNMVEEILRLMGQP...</td>\n",
       "      <td>Bla-g-6.0301</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>Bla-g-7</td>\n",
       "      <td>MDAIKKKMQAMKLEKDNAMDRALLCEQQARDANIRAEKAEEEARSL...</td>\n",
       "      <td>Bla-g-7</td>\n",
       "      <td>CDS</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      locus_tag                                        translation   description type\n",
       "0       Bla-g-1  NAIEFLNNIHDLLGIPHIPVTARKHHRRGVGITGLIDDIIAILPVD...       Bla-g-1  CDS\n",
       "1       Bla-g-2  MIGLKLVTVLFAVATITHAAELQRVPLYKLVHVFINTQYAGITKIG...       Bla-g-2  CDS\n",
       "2       Bla-g-4  AVLALCATDTLANEDCFRHESLVPNLDYERFRGSWIIAAGTSEALT...       Bla-g-4  CDS\n",
       "3       Bla-g-5  MAPSYKLTYCPVKALGEPIRFLLSYGEKDFEDYRFQEGDWPNLKPS...       Bla-g-5  CDS\n",
       "4  Bla-g-6.0101  MDELPPEQIQLLKKAFDAFDREKKGCISTEMVGTILEMLGHRLDDD...  Bla-g-6.0101  CDS\n",
       "5  Bla-g-6.0201  MDEIPAEQVVLLKKAFDAFDREKKGCISTEMVGTILEMLGTRLDQD...  Bla-g-6.0201  CDS\n",
       "6  Bla-g-6.0301  MADEQLQLPPEQISVLRKAFDAFDREKSGSISTNMVEEILRLMGQP...  Bla-g-6.0301  CDS\n",
       "7       Bla-g-7  MDAIKKKMQAMKLEKDNAMDRALLCEQQARDANIRAEKAEEEARSL...       Bla-g-7  CDS"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "m2_alleles = ep.get_preset_alleles('mhc2_supertypes')\n",
    "my_alleles = ['HLA-DPA10103-DPB10201', 'HLA-DPA10201-DPB10101', 'HLA-DPA10201-DPB10501', 'HLA-DPA10301-DPB10402', \n",
    "            'HLA-DQA0101-DQB10501', 'HLA-DQA10301-DQB10302','HLA-DQA10401-DQB10402', 'HLA-DQA10501-DQB10301', \n",
    "            'HLA-DRB1*0101', 'HLA-DRB1*0301', 'HLA-DRB1*0401', 'HLA-DRB1*0405', 'HLA-DRB1*0701', 'HLA-DRB1*0802', 'HLA-DRB1*1101', \n",
    "            'HLA-DRB1*1302','HLA-DRB1*1501', 'HLA-DRB3*0101', 'HLA-DRB4*0101', 'HLA-DRB5*0101']\n",
    "for i in my_alleles:\n",
    "    print(i)\n",
    "prots = ep.fasta_to_dataframe('cockroach_allergens.fa')\n",
    "prots"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "took 66.336 seconds\n",
      "predictions done for 8 sequences in 20 alleles\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>1-log50k(aff)</th>\n",
       "      <th>Affinity</th>\n",
       "      <th>Core</th>\n",
       "      <th>HLA</th>\n",
       "      <th>Identity</th>\n",
       "      <th>Pos</th>\n",
       "      <th>Rank</th>\n",
       "      <th>allele</th>\n",
       "      <th>core</th>\n",
       "      <th>name</th>\n",
       "      <th>peptide</th>\n",
       "      <th>pos</th>\n",
       "      <th>rank</th>\n",
       "      <th>score</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>175</th>\n",
       "      <td>0.590</td>\n",
       "      <td>84.63</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>HLA-DPA1*0103-DPB1*0201</td>\n",
       "      <td>EFKNFLNFL</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>EFKNFLNFLQT</td>\n",
       "      <td>174</td>\n",
       "      <td>1.0</td>\n",
       "      <td>84.63</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>363</th>\n",
       "      <td>0.590</td>\n",
       "      <td>84.63</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>HLA-DPA1*0103-DPB1*0201</td>\n",
       "      <td>EFKNFLNFL</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>EFKNFLNFLQT</td>\n",
       "      <td>362</td>\n",
       "      <td>1.0</td>\n",
       "      <td>84.63</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>174</th>\n",
       "      <td>0.554</td>\n",
       "      <td>124.33</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>HLA-DPA1*0103-DPB1*0201</td>\n",
       "      <td>EFKNFLNFL</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>PEFKNFLNFLQ</td>\n",
       "      <td>173</td>\n",
       "      <td>3.0</td>\n",
       "      <td>124.33</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>362</th>\n",
       "      <td>0.554</td>\n",
       "      <td>124.33</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>HLA-DPA1*0103-DPB1*0201</td>\n",
       "      <td>EFKNFLNFL</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>PEFKNFLNFLQ</td>\n",
       "      <td>361</td>\n",
       "      <td>3.0</td>\n",
       "      <td>124.33</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>423</th>\n",
       "      <td>0.542</td>\n",
       "      <td>142.12</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>NaN</td>\n",
       "      <td>HLA-DPA1*0103-DPB1*0201</td>\n",
       "      <td>LYALFQEKL</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>ELYALFQEKLE</td>\n",
       "      <td>422</td>\n",
       "      <td>5.0</td>\n",
       "      <td>142.12</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    1-log50k(aff) Affinity Core  HLA Identity  Pos Rank                   allele       core     name      peptide  pos  rank  \\\n",
       "175         0.590    84.63  NaN  NaN      NaN  NaN  NaN  HLA-DPA1*0103-DPB1*0201  EFKNFLNFL  Bla-g-1  EFKNFLNFLQT  174   1.0   \n",
       "363         0.590    84.63  NaN  NaN      NaN  NaN  NaN  HLA-DPA1*0103-DPB1*0201  EFKNFLNFL  Bla-g-1  EFKNFLNFLQT  362   1.0   \n",
       "174         0.554   124.33  NaN  NaN      NaN  NaN  NaN  HLA-DPA1*0103-DPB1*0201  EFKNFLNFL  Bla-g-1  PEFKNFLNFLQ  173   3.0   \n",
       "362         0.554   124.33  NaN  NaN      NaN  NaN  NaN  HLA-DPA1*0103-DPB1*0201  EFKNFLNFL  Bla-g-1  PEFKNFLNFLQ  361   3.0   \n",
       "423         0.542   142.12  NaN  NaN      NaN  NaN  NaN  HLA-DPA1*0103-DPB1*0201  LYALFQEKL  Bla-g-1  ELYALFQEKLE  422   5.0   \n",
       "\n",
       "      score  \n",
       "175   84.63  \n",
       "363   84.63  \n",
       "174  124.33  \n",
       "362  124.33  \n",
       "423  142.12  "
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "reload(base)\n",
    "P = base.get_predictor('netmhciipan')\n",
    "b = P.predict_sequences(prots, alleles=my_alleles, threads=8)#, show_cmd=True)\n",
    "b.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The promiscuity versus ranking approach\n",
    "\n",
    "There are two broadly different approaches to computational selection of potentially dominant epitopes. The promiscuity method finds high scoring binders above the percentile (allele-specific) cutoff in greater than n alleles. The ranked approach calculates the median rank over all alleles and sorts the binders in ascending order. One may then take the top percentage of these binders as suited. We use the former method here as it lends itself better to the detection of clusters of binders.  Though they do give similar results as shown in the venn diagram."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAFoCAYAAACborbBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy86wFpkAAAACXBIWXMAAAsTAAALEwEAmpwYAABSxUlEQVR4nO3dd5hcZ3X48e+ZsrNN2pVWvVvFsiX3buNKsYMLmGYMDjamBRKSEAgQSljrR0+ABJIQCEloDmDjQMDGgHvDBRdZtmRZlmSt+kpabdGWmdkp5/fHe9carbbvzNw7M+fzPPPs7syde8/M7Nxz3y6qijHGGJNPIb8DMMYYU34suRhjjMk7Sy7GGGPyzpKLMcaYvLPkYowxJu8suRhjjMk7Sy6mZInIp0XkP/2OYzJE5CYRuXmExzeIyMUT3PcPROQLE41tokTk3SLyyCSe/x0R+ft8xuTtd8T32uRXxO8AjJkoVf2S3zEUmqqu9juGYlPVD/odg5k8K7mYohERu5gpAX5+TiIS9uvYJr8suZhJEZEWEfmUiLwgIh0i8n0RqfYeu1hEdonIJ0WkFfi+iMRE5J9FZI93+2cRiQ3a/hMisl9E9orI1SJyuYi8JCLtIvLpnGO/Us0hItUicrOIHBSRThF5UkRme49N9+La48X4f979R1XfiIiKyHLv9wYR+ZGIHBCR7SLyWREJDT629/cS77mRnH2/LCLdIrJNRK4b4W2sFpFbvG2fEZGTB72/r8055q1eTN1eldkZOdue6j2/W0RuAaoHvbYrReRZ7/15VEROGnScT4rIc0CviES8v3d7+9skIq8Z5n+gSUR+LSKHROSPwLJBjx8nInd7n98mEbkm57EfiMi/i8idItILXJJbnSciG0XkypztI97ncZr39znea+kUkXWSU4UoIseIyINe/HcDM0b4DEyeWXIx+XAdcBnupHIs8Nmcx+YA04HFwAeAzwDnAKcAJwNnDbF9NTAf+BzwPeBPgdOBC4C/F5FjhojhBqABWAg0AR8E4t5jPwZqgdXALOCfxvi6/sXb51LgIuB64MbRniQidcC3gNer6hTgPODZEZ7yRuDnuPfpJ8D/iUh0mG3fAPwMaAR+Dfyrd8wq4P9wr3W6t7+35MR0KvDfwJ/h3p/vAr8eSOyedwBXePteBnwYONN7DZcBLcPE9G9AApgLvMe75b4Xd3uvaxZwLfBtEVmV8/x3Al8EpgCD22p+6sU14DKgTVWfEZH5wG+AL3iv+W+B/xWRmd62PwGexiWVz+P+R0yRWHIx+fCvqrpTVdtxJ4nck0EWaFbVpKrGcYno/6nqflU9AKwB3pWzfQr4oqqmcCfRGcA3VbVbVTcAL+CS0mAp3ElzuapmVPVpVT0kInOB1wMfVNUOVU2p6oOjvSBx1TPXAp/yjt0CfH1QrCPJAieISI2q7vViH87Tqnqb95q/gUuu5wyz7SOqeqeqZnCJZOC9OAeIAv/svcbbgCdznvcB4Luq+oT3/vwQSA46zre8zzEOZIAYsEpEoqraoqpbBwfjvU9vAT6nqr2quh74Yc4mVwItqvp9VU2r6lrgf4G35WzzK1X9g6pmVTUx6BA/Ad4gIrXe3+/EJRxwFx13eu9HVlXvBp4CLheRRcCZwN97/3sPAbcP856aArDkYvJhZ87v24F5OX8fGHTCmOdtM9z2B70TJxwueezLeTwO1A8Rw4+B3wM/86q//sG7+l8ItKtqx5hfjTMDd7IeHOv80Z6oqr3A23Glp70i8hsROW6Ep7zy/qlqFtjFke9Jrtac3/twVWoRb/vdeuRMtLmxLwY+5lUfdYpIJ+69yT1ObhxbgI8ANwH7ReRnIjJUTDNxHYMG/w/kHvfsQce9DldCPeq4g3lxbASu8hLMG3AJZ2Dfbxu07/NxJah5QIf3WQwVlykwSy4mHxbm/L4I2JPz9+Bpt/fgTgrDbT8h3tX6GlVdhauGuhJXjbUTmC4ijUM8rRdXXQaAiOSe8NpwpaHBse4e6rkcebJEVX+vqq/DnehexFXvDeeV989r01nA+N+TvcB8EZFB8Q7YiSsRNubcalX1pznbHPFZqepPVPV83HugwFeHOO4BIM3R/wO5x31w0HHrVfVDwx13CANVY28EXvASzsC+fzxo33Wq+hXv/ZjmVcsNFZcpMEsuJh/+QkQWiMh0XJvKLSNs+1PgsyIyU0Rm4NpVJj32QEQuEZETvWqaQ7jEkFXVvcBvcfX800QkKiIXek9bB6wWkVPEdUK4aWB/XunpVuCLIjJFRBYDH82J9VngQhFZJCINwKdyYpktIm/0TmxJoAdXTTac00XkzV4J5CPecx4f51vwGO4k/1fea3wzrj1rwPeAD4rI2eLUicgVIjJlqJ2JyEoRebXXJpPAlRiPeg3e+/QL4CYRqfXaUnLbNu4AjhWRd3lxRUXkTBE5fhyv7WfApcCHOFxqAfdZXCUil4lIWFynjotFZIGqbsdVka0RkSoROR+4ahzHNJNkycXkw0+Au4CXga24BtbhfAH3pX8OeB54ZpTtx2oOcBsusWwEHsRVlYFrJ0nhShD7cSdwVPUl4P8B9wCbObox+S9xJZSXvcd+gmsUx6vfv8V7HU/jTqIDQrhEtAdox3UGyL1SH+xXuGq0Di/WN3vtL2Omqv3Am4F3e8d8O+6kP/D4U8D7cR0AOoAt3rbDiQFfwZXgWnGN8Z8aZtsP46oqW4EfAN/POW43LjFci3s/WnEloNhRexn+te3FJc/zyLlwUdWduNLMp3ElqJ3Axzl8XnsncDbu/WgGfjTWY5rJE1sszEyGiLQA71PVe/yOxRgTHFZyMcYYk3eWXIwxxuSdVYsZY4zJOyu5GGOMyTtLLqYoZBJTx0/gWAWdal5y5h8rV15X5Ge9ebn+SkRqROR2EekSkZ/7EE/ePlOxqfeLwpKLKQpVXa2qD+R7vzLJtUPGsP8HROR9hdp/sUjOBJhj9AngflWdoqrfAt4KzAaaVPVtIz91cgr9mZrisORiTAFJ6S4zsBjYMOjvl1Q17VM8psRYcjHjMrhKSI6cHn2GiNzhzfPULiIPy+Ep6sczdfxpIrLWe+zn4qajP6pKxBvl/R3gXBHp8eaWGjBN3Jxe3SLyhIgsy3nesFPAD9r/F3EzMf+rt/9/zXn4tSKy2Xut/zYw7Yp31f0HEfknETmIG7keE5GvicgOEdknbqXFmpzjDDsV/hAxrc6JfZ94SxAMrjbyRqrv8n7/MW7qk9u91/EJ7/43eO99p1dCO967/z7gkpzX/VPcTApv9/5+7xBx3eR9Vjd77/nzInKsuOUY9ovIThG5NGf7BhH5L3HLKuwWkS+IG2U/0c/0PHHLLHR5P8/LeewYGWbqfRlhqQYzSapqN7uN+YabB2p5zt8/AL7g/f5l3Ikh6t0u4HCPxBbgtd7vN+GmFLkcCHvPe9x7rAo3weBfe/t4M9A/cIwh4nk3bqZgBsV0EDf9SQT4H+Bn3mN1uJHcN3qPnYobhb5qmP0/gBskOvg9uAM3Nf0i3OjwP8mJJ40b3R8BanBT/P8aNy38FNzsvF/2tj8VN2vA2d57cYP3XsWGiGUKbs6sj+FmTp4CnD34c/D+vhjYlfP3K++/9/exuNkHXue9z5/AjdqvGup1e5/ZzSP8Xwx8ppd5r/tHwDbcdEBR3OwA23K2/yVu2v863Oj/PwJ/NsHPdDqHZzeI4OYh68BV4YEb3f8N3KwAFwLdA68FtwTB7bh54sK4pR2m+v09K4eblVxMPqVwEzUuVjeR5MPqfYOHMNLU8RHc9O8pVf0F7sQzXr9U1T+qq8b5H9z6MTC2KeDH4iuq2qmqO4D7c/YPsEdV/8U7dgI33f3fqGq7uulQvoSbDgXGNhX+gCuBVlX9uqom1C0F8MQ44x7wduA3qnq3uqlmvoZLhOeN/LQRPaxuws40bj2Zmbj3aWD5hCUi0uiVDC4HPqJumv79uAR87bB7dob7TK8ANqvqj73P9Ke4qX6uktGn3h9yqYZJvAfGU6r1wSaY/hF3BXuXV0v0H+pmqB3KeKaOH3ZK9hEM3v/ANP2vTAGf83iEw/OQTXb/cGS8M3FXxU/L4QmLBXeVPBDPDSLylznPqWLoKfcX4uZuy4cjlj5Q1ayI7GQMSwqMYPDSCG169PIJ9d6xo7jlCAa2DzH65zzcez54GQc4vDzCcFPvD8zi/GPv95+Jmzn7ZuAzOs653czRrORixquPYaaa966kP6aqS3HrbnxUhlkadwRDTR2/cLiNGX269sHGMgX8ZPY/+DltuBPr6pzjNajqwIlxLFPh58a+dJhjjrgEwBCv44ilD7z3eyGHlxQopJ240tmMnNc8VVVXDxPraAYv4wCHl0cYcep9HX6pBjNJllzMeD0LvNNrfP0T3Iy/wCsN08u9E1UXbjXDkaaaH8pj3vM+LG699Ddy5NTxg+0DFohb5ncsxjsF/D6GP6GPSt3iX98D/klEZgGIyHwRuczbZDxT4d8BzBWRj3idBKaIyNneY8/iVmCcLm5dmo+M8jpuBa4QkdeIW1TtY7gT/qMTfa1jpW6W47uAr4vIVBEJicgyERn4XxrvZ3on7jN9p/c/83ZgFXCHjjL1vgyzVEN+Xmlls+RixuuvcV/OTtyKgv+X89gK3PT1Pbgk8W1VvX88O9fDU8e/1zvGn+JOqslhnnIfrstsq4i0jWH/450C/pvAW0WkQ0S+NfZXcoRP4hrLHxeRQ7j3aKUXz5inwvdifx3u/W/FLRNwiffwj3Hr07TgTtyD19T5Mm4dnU4R+VtV3YR7b/8FV7q6CrjKe/+L4Xpc9d8LuNd9G669Dsb/mR7ElTg+hmv0/wRwpaoOPHekqfdHWqrBTILNLWYCT0SeAL6jqt8fdWNjTCBYycUEjohcJCJzvCqOG4CTgN/5HZcxZuyst5gJopW4NoE63CqQb/Xq6Y0xJcKqxYwxxuSdVYsZY4zJO0suxhhj8s6SizHGmLyz5GKMMSbvLLkYY4zJO0suxhhj8s6SizHGmLyz5GKMMSbvLLkYY4zJO0suxhhj8s6SizHGmLyz5GKMMSbvLLkYY4zJO0suxhhj8s6SizHGmLyz5GKMMSbvLLkYY0wBiMgDIvI+v+PwiyUXY4wJKBG5RkQeFZE+EXlgiMdPEZGnvcefFpFTih/l0Cy5GGPMEEQk4tNxZ+f82Q78M/CVIbarAn4F3AxMA34I/Mq733eWXIwxxiMiLSLySRF5DugVkc+KyFYR6RaRF0TkTTnbvltEHhGRr4lIh4hsE5HXD7PfuSLynIh8fJjHa0XkXSJyH3D/wP2qeo+q3grsGeJpFwMR4J9VNamq3wIEePWE34A8mlRyEZGLRWTXGLf9gYh8YYLHuUlEbp7Ic4NKRC4QkU05f7eIyGv9jMkYA8A7gCuARmATcAHQAKwBbhaRuTnbnu1tMwP4B+C/RERydyYixwAPAv+qqv846LFzReR7wG7geuC/gNPHGOdq4DlV1Zz7nvPu952VXIpERFRElg/8raoPq+rKIh5/iRdDT87t73Me/7iItInIBhE5Mef+V4nI/xUrTmMC4FuqulNV46r6c1Xdo6pZVb0F2AyclbPtdlX9nqpmcNVSc4Hcaq1VuJJIs6r+x8CdXlvKi8APgG3Aiar6OlX9H1WNjzHOeqBr0H1dwJRxvNaC8aVOMahEJKKqab/jyBcRiQJTVLU95+7Gwa/RuxJ7L7AUd/X0ZeBKr87568C1RQrZmCDYOfCLiFwPfBRY4t1VjyulDGgd+EVV+7xCS33O49cBW4DbBh1jATAfuBNYl7ufcegBpg66byrQPYF95d2oJRcROU1E1np1jj8XkVuGq94SkeO97ned3hXwGwZtMkNE7vb29aCILM557jdFZKeIHPJ6PVwwlhcwUDUnIp/2rrxbROS6nMdjXp3oDhHZJyLfEZGaQc/9pIi0At8XkbC3r4F61qdFZKG3/XFe/O0isklErsk5zg9E5N9E5Dfe854QkWXeYw95m63zSgxvH6lKUURCIvJ3XgwHReRWEZk+lvfDe/4JIvJ1YBfwujE8ZRGwVlUPAffgkgzAR4Bfq2rLWI9tTBlQAO/89D3gw0CTqjYC63HtGmN1E9AG/EREwq8cQPUbuORyL/AZYJeI/JOInDqOfW8AThpUDXeSd7/vRkwu4nod/BJXdJsO/BR40zDbRoHbgbuAWcBfAv8jIrlVP9cBn8dl/meB/8l57EngFO84PwF+LiLVY3wdc7x9zgduAP4j57hfAY719r3c2+Zzg547HVgMfAB3lfIO4HLcVcB7gD4RqQPu9mKbhbua/7aIrMrZ17W4etlpuKuVLwKo6oXe4yerar1XvB7JXwJXAxcB84AO4N9GeoKITBORPxeRJ3GfQRZ49RDH2u4l1O+LyMAV2BbgRBFpBF4LbPAS6rXA10aJ1ZhyVYdLNAcARORG4IRx7iMFvM3b149E5JVzrqoeUtX/UNXzcN/1BHC7iNw7sI13sVuNq2UKiUi1d64FeADIAH/lXUR/2Lv/vnHGWBiqOuwNuBDX0CQ59z0CfMH7/WJgl/f7BbiiXShn258CN3m//wD4Wc5j9bg3ZuEwx+7AnYzBZf+bh9nuYiAN1OXcdyvw97grjF5gWc5j5wLbcp7bD1TnPL4JeOMQx3k78PCg+76Lq0sdeH3/mfPY5cCLOX8rsHxQ3Lty/m4BXuv9vhF4Tc5jc3H/pJEh4poK/Azo9F735UB4iO3qgTNw/6SzccX03+c8/g7gGeC3uET7C+A13ut+ENflccFI/y92s1up33K/h97fX8R1B24DvuF9F97nPfZu4JFBz3/le447+Q9sW42rFfhB7jlyiOOHgHNz/n63t8/c2w9yHj8VeBqIe9/fU/1+Dwduo7W5zAN2q/cqPDtH2HanqmZz7tuOKykc9VxV7RGR9oHnicjf4ur953lv4FSOrNscSYeq9g467jxgJlALPJ1TchQgnLPtAVVN5Py9ENg6xDEWA2eLSGfOfRHgxzl/59ab9nFk3et4LAZ+KSK572UGlxR2D9o2iruaaseVBtera1w8gqr2AE95f+7zrnL2isgUVe1W1Z/iLgYQkSuAJLAWVx+8GngDrhRj7S9lRNZICPd/OtW7VeH+pyLeLff33L8zuAuzlHcb6vd+XLtAlzYf8R0LLFVdMujvz+CqrYba9ge4ZJF7n+T8fnHO7wlcrcBox88Cj410jEHbr2XsvcuKarTksheYLyKSk2CGO/nuARaKSCgnwSwCXsrZZuHALyJSj6uO2uO1r3wCd6W8QVWzItLB2Os2p4lIXU6CWYSrG23DZfTVqjr4pDxAB/29E1jmPX/w/Q+q6ljaMCZrJ/AeVf3DaBuq6kHgBBE5E7gReEZcH/0fAbd5SWXIp3o/j6ga9dqjvgS8HliBu2A45FW3fXpCr8b4StZIBNeVduoQt3rG14Yw0RiSuJ5MR920WVOFPr4pvtGSy2O4K5QPi8i/4/p+n4Ur7g32BO5q/RNeY/KrgKuAM3O2uVxEzgf+iGt7eVxVd4rr+prG1W1GROTvOLoXxGjWiMincf3Or8RVV2XF9SH/JxH5sKruF5H5wAmq+vth9vOfwOdF5AW8tghcaeEO4Csi8i5cNRS4dpweVd04hvj24RrKt4xh2+8AXxSRG1R1u4jMBM5T1V8N9wRVfRJ4UkT+Btde827gWyJyjar+TkTOxlWdbca1CX0LeEBVB3dl/Cyu2L1HRBRYKW7E8CXAy2OI3fhM1kgNrip1DofbFP0edhDDtVXOGvyArJE+YD/uAnUv0K7NOviiz5SYEZOLqvaLyJtxJ9wv4+rj78BVmQy17VXAt4FP4Q0KUtUXczb7CdCMa/d4BvhT7/7fA7/DlXJ6gX9i+Oq3obTi2mj24BLcB3OO+0lcA/7jXgP2buDfvWMO5Ru4L8JduGq5F4E3qepBEbnUe/wbuC/rOlwHgLG4CfihVzL4AO7LNJxv4q4m7xKRed62t+DaPUakqklv21u858a8h5biSiSzgEO4zgnvyH2uiBwHXAqc4+1rr4h8Bdf7ZD+u/aV8uLrSsHfLUKLd0GWNNHI4kcxh/BdmfqvFdfVd4v2dlDXSiks0e4CDlmxKj+g4PzMReQL4jqp+vzAhjY+IXIxr7F/gcyjGT65HTX3ObYr3sw7XThAe4jb4aj6LaycYuCUH/d2LK/11Ad0c2b5YNLJGwrgq5qW48RJj7VVZqvpxJf89wDZt1kM+x2PGYNTkIiIX4XpQteG6En8HWKqqewsf3ugsuVQYV/IbqF5p4nC7QbEHBGdxg9W6OJxwuoCDuNJjXnkN7wMJZTGu4b1SHcRV0b6szUdV65qAGMsXciWui2sd7gN9a1ASiylzboaAmbhEMvBzoj3w8i2EayRvwHUgOcz1gtz7ym3s03kcuRuXUBbgEsoSKjuh5GrybmfKGmnDtSNu0eaJvc+mMMZdLWZMQYnMwp1IF+Iaogvek6kIOnHtgq7B+shu80eRNTITOB44hsNtZmZkWdyMFJuBFm0+uju+KS5LLsZfbkqM+biEsgjXuFvuOnG1ANtwXckHSilLcWOWjupRZcYljuuEskGb819FacbGkospPtf4vgiXUBZQwROoJqrpuONSdn32eBZuqqHR73jKTBrX2/M5bR52vJcpEEsupjhct98FwHG4Bmm/x134qreO1EsnkdqzkGoNu/dif4T+h6aSun0asfZI5SbcAsjiSorPavMRM4SbArLkYgpLpBaXUI4jOI3xvknUkN5wKv17F1JDaOj2pAzoc7XEb5tO+Lk6a3PJs124JDPUyo4mjyy5mMJwo/pPwDVKV3QpBaC/isyLJ5PYeQw1AyWVsdgXof/X08j8vpFYMmTvYx7tBx7VZh1pMLOZBEsuJr9ElgIn47oOV7x0mOxLJxJvWUF1NnLEhKnjkhQyD04l8eMZxDqtyiyfXgL+qM3a53cg5caSi8kPkUW4eeSa/A4lKFqW0/fiyVSlq/KXDFKQvb+BuCWZvErhZgB/Tpv9mXWhHFlyMZPj5i87kyPXDa9o3VPpX3su2UPTCzctSwqy9zSQuHkGsUOTKBGZIxwCHtNm3e53IOXAkouZGDfY8SzcujkGVwW28RQS25cP31ifb/1C9i5XkqnuC1uSyZNduPaYTr8DKWWWXMz4iEzHJZVFo21aSfYuIP78mUT6q4mOvnX+9YZI/3AG/b+dVhGDUIshi1t87xmrKpsYSy5mbNxI+tOBk7DeX69Ih8k+ew6J1kXBOKnvqCLxjbnI1mrrwpwn+4H7bCbm8bPkYkbn2lUuwE3SaDwdTSSfOh9J1gZrQsks6ENTiH93NrEeqyrLhxTwB23Wl0bd0rzCkosZnkgVbuGw4/wOJUgU9MWT6dt6HLXFaluZiD4h85+zSN7dGIxSVRnYCjyszdrvdyClwJKLGZobr3IelTGR5Jj11ZJ68iIy3Y2ls0DX03XEvzaXKivF5EUPrpqs1e9Ags6SizmSm67lfA4vOWs8B+aQeOp8oplo6Z2ku0Okvz6X9NP1pZMUA0xx42KetuWXh2fJxRwmMhd4DVZaOcrmVfRtOrF4XYwL5d6p9P3bbKpTNpVMPuwB7rJqsqFZcjGOyMm4wZB20smRCZFdex6J1oXlk3DbIvTftAC2x4LVEaFEdQK/1Wbt9juQoLHkUulco/1FuAkmTY5EDenHLyHT01B+3Xr7hew/zyH58FRq/I6lDMSB32mzHvA7kCCx5FLJ3IDI12FdjI/S1Uj/Y68hlM95wYLo1430fm82dX7HUQbSwL02dcxhllwqlchy4EIqeBXI4RycSeKPFxPNVMicXS9UE1+zgCqbPmbSFDc32Xq/AwkCSy6VSORc4ES/wwiiffOIP3U+sfGsuVIODobp//uFsNPaYfLheeDxSu9JZsmlkoiEcO0rK/wOJYh2HkPfurNKv0fYRCWFTPMCUhtqrbtyHrTgqskyfgfiF0sulcLNDfY6bMLJIW09jr6Np1CDVGZiGZCC7FfnkXxiijX058F24O5KnfjSkkslcD3CLgPm+h1KEG09jr6Np5ZPV+PJyoD+yxzi9zbYe5IHW3Ej+ivuRFtR9coVSaQGuBJLLENqWW6JZbAwyF+1UvOmdmzp38lbhus4U3Gs5FLOROqBK7CuxkPavZj42nOprvSqsJHcNp3eH860rsp5sF6b9VG/gygmK7mUK5EG4I1YYhlS63zia8+xxDKat7ZT964D9PodRxk4QdbIWX4HUUyWXMqRm3zycrArzqEcmE3i6VcRq9ReYeN1TTt1VkWWF6fIGjnV7yCKxZJLuXGN968HpvgdShB1T6X/yQuJVto4lsl6zwFqX9dpCSYPzpQ1coLfQRSDfcHKietufCnQ5HcoQZSMkXn81Ui2Qkbe59tf7KPmgkPE/Y6jDJwna2SJ30EUmiWXciEiwCXAPL9DCaKsoE9cTCpZQ9TvWEpVGOSje4md3kPC71jKwMWyRsq6PdSSS/k4D1jqdxBBtfZc4oem28jzyYpA6NN7qDomga1hMjlVwKWyRsr2YseSSzkQORVY7XcYQbV5FX17F9tYlnypUkKf34U0pkn7HUuJm4abjqksWXIpdSIrcIt8mSEcmENi00k2lUm+NWSIfnEn6bBiA+UmZ6mskZP8DqIQLLmUMpFpwAV+hxFUyWrST7+KiI1lKYxF/VR/fI818OfB2bJGyq6t1JJLqRKJ4CaitPVYhqCgf7yQdLkv9uW3V/VQ+0YbAzNZArxW1ki934HkkyWX0nUB0Oh3EEG16ST6upqsAb8YbjxAzfF9JP2Oo8RV4xJM2XSTt+RSikRWYmuyDKujieSW460Bv1jCIJ/eQ6g6S0VOLZ9Hs4DT/Q4iXyy5lBq37v2r/A4jqNJhsk+dj9jULsXVmCH6sb02/iUPTpI1UhaDoC25lBLXzvJarJ1lWOvPJJGstaV6/XBOD7WXdFn7yySFgItkjZT8xZEll9Ji7SwjaJ9BctcS63bspz/fR2xGysa/TNIMoOS7J9t6LqVCZDFuNUkzhKygD1xJf189Mb9jGauP/I4bdh7ixFiY7p+8hTUAf3En7+9KMAegP0NNVZj4zW/m8/5GOj6bYyQ+usQ6U0xSGrhNm/WQ34FMlFWvlAJXHWbtLCPYfAJ9ffWltcTAxUt4tL6K+/97LTcO3Pdvl/O9gd//7h7eWhMpvXEkK5JUv72N3ltmlNbnETAR3AqWd/gdyERZtVhpOB0oqz7w+dRbR2rL8aVXHXb1cWyeUTv0QlxZha0dnPH6FTxZ7Ljy4Zp2aqx6bNLmyRo5zu8gJsqSS9C53mEn+h1GkK19FelyW5/l15tYUR3h0Fnz2e93LBNRpYT+qtUmt8yDc2SNlGS3+rL6QpapC7DPaVi7FxHvbCq9UstoHt7BmatmlmapZcCpfdSe2VN61XoBU4Wb8bzk2EkryESOA2b7HUZQKejGU8tv4a9kmlBLJ6e9YWVpJxeAv2glHLXBlZO1VNbILL+DGC9LLkElUg2c7XcYQdZyLPFEGY5p+cVGjp8ao3X1TDr9jmWymjJU/WmblV7y4Cy/AxgvSy7BdQ6UTrfaYkuHyW46sbQTy4fv5H1feIhP9vQz+20/56vfesL1CHx8N2eeOIs/+h1fvlzVQc3MFCm/4yhx82SNLPQ7iPGwcS5BJNIEvMXvMIJs48n0bl1lXV1LxWP1xL80v/zaxorsIPALbS6Nk7aVXILJFv8aQTJGZtuxNkivlJzdQ/X8pPUem6Qm4Bi/gxgrSy5BIzILWOR3GEG26SQS2Uj5NeSXsxDI+w+Q8TuOMnCa3wGMlSWX4Cm5hrtiSkXJ7LKpRUrS6b3UHJOw0sskTZc1stTvIMbCkkuQiMwBym6503x6+TgrtZSy9++30kselETpxZJLsJzqdwBBlhW0ZYX1oCtlJ8apWRG3VSsnabqskSV+BzEaSy5BITIDKKmuhsW2YzmJVMwmWy1172qz0kserPY7gNFYcgmOkijq+mnL8VYdVg5O6rNxL3kwX9bIVL+DGIkllyAQmQIs8TuMINu7gHiirrQHTRonDPLWdmvYz4Pj/Q5gJJZcgqFkp9Uulm3H+h2ByaeLu6iO2Zxjk3WsrJHAnsMDG1jFEBFgpd9hBFmymnT7TOt+XE5qlfClXST8jqPE1RDgGg9rHPXfIqAk12solpblJAnZVC/l5o3tRG6f5ncU49BHhG/zcbJEUMLM52n+lNv5d95LF4sRMjTSwru5mVjROi0cD7xcpGONi80t5jeRy4DFhdj1QzDtXXBjD0wFeB089DO47wR4/37cOu0JqKmG+H6Cu077PW+gv9zaW/ozsUw8VZ9JZ6PZjEbIZCNkNKxZDR/1hYyG+iUWiYeqwvFQLJwIRcP9ZdOx4bMLSKyrK5FSaRboJkYDSZKE+Rc+zkXcQg91XMR6AL7D+5jDZt7Mg0WM7GfarIeKeLwxsZKLn0RqKWD34xhkb4LbboQd2yF2Cnz257BxPYfXab8I3lpPcKdEb59BMlFXWmNbshrS7uS0VE9/Y6YvNUXj6TqNp6eEEulaSaZrwv2Z6rASDsNEe79lNRJKZ6KhZDYa7s9GQwmtqzqUbYi1SUN1W3hKVUc0EkqXRJX36zvJriuVMmkIaPDG6KQIo97nd4mXWABmsI0eil0eOx54osjHHJUlF3+tpIDtXmdD19nQBbAYkrNg70vQCOwFyADPwBk/hW8UKobJ2rYy2I2+qqKHktNT7YnZmY74bO1Kzgj39jdElVABS1ohSWerIulsFXFvlfqD8flHRBULx1N1VZ2ZqbH2TEPsoDTE2iJTY+1REZXCxTV+p/cSCyuaEQIV17DSCN/gsySYySIe4Ey2vfJYkjDbOYfzuaXIUa2UNfKkNmugviuWXPxVtIb8e6BpDyx8J4e/DN+EFXVw6EqCuU57VtD984JVHZbORrOtPYuTB/oWaGdiRqi3v6FKCQcqRhBJZmqjyXhttD1+eDahsKQy02r298+q3amz61ui9VWHoj4GCUC1Ej63m8QjU0ukaiyC8gk+z0Fq+BEfYj3zOIE9APyQd9LEZs5lS5GjqgbmAzuLfNwRWXLxi8g8vLaQQtsJsevggx+EW4/hcA+dW+HM8wnuUrptc0hmIv6fdHr7p6T2dC9LtfYsCXUmZ8QgVJLrkmQ0Gm7rm1/T1jefF9rOoSocT02v2ZueXbeTWXU7otWRuC/ng9d1oY8EejjgEJqIM4tNbGA1J7CHm7mSJPW8l5t9imgRllyMpyjrMnRD+CL44IXwxD/C2oH7+yD0HJz2j/CFYsQxEXsW+1Mlpip6MD4nubt7eWZ/78KqRLo+Cvh+lZ9v/ZmaaGvP0mhrj5tkty7amVwwdXN6UcOLsWImmhP7iFVnySZCAR8a0Uo9UTI0EaeXKPtYxSn8jl9xPvtYzQf5BmH86iG1CPiDT8cekvUW84vIO4AphTxEBjgdbpwCvQ/DrbmPrYHV/wWv3wFfK2QMk3HXm0j1VxfvpN6dbOzf2nFSek/30lhGq8qmR9b4ZXVa9f7EooYXdf6UrdXhUKbgJ/1/mU38rsaAr1T5HPP5PTeihFCEBTzFdfyGm/h3Yhwk7DX2z+cZruM3PkT4c23WDh+OOyRLLn4QaQSuKfRhvgnLPwIfnwG7BXdF9SH45RpYfxq8+wR4+UfwUKHjmIjO6SQfuazwvcQy2XB2d/fyxLbOE0KHkk2+V8EFTUjSmdl1O5KLG14IzazbU7D3Z0MN8b9bFPDkEnx/1GZ91u8gBlhy8YPIScA5focRZBtOpW/bcYUbXGqllPGrCsdTxzSu71867fmafHd17hey16xASqbXWDC1arP+2u8gBlhy8YPIldiiYCO67yqSffX5L7ns7V4Sf6n9dLFSysSFpT+zuOHFxIqmtdVV4WTeEnNJDagMJgV+pM0aiPVyrEG/2ESq8EbHm6H1V5HJd2LZ17Mo/kLb2aGe/mlW9TJJGa0Kv9x5Ul1L16rs/Clb+lY2PV1VE+2d9LnknB4yJTOgMpgENyi72F2hh2TJpfjmYxOGjqhtNv2Qn/r3/b0LEi8cOIfu/umWVPIsq5HQzkPH1e48dKzOrW/pO27Gk5H6qq4Jj/k5tdfOR3mwCEsuFWuR3wEE3YG5k++C3NY3N7HhwLkcSs6wapaCC8nenqW1e3uO0Xn1L/edMOsPsVgkMe7qsrkpqqamyRyK2KJwk7BQ1ohos//tHZZcis+WMh5Fx4yJ/18eSk7rX7fvokxnYpaVVIpOZE/Pstp9vYsyxzY907ds2nM145luJgRyZi/JextslvBJiAENQKfPcVhyKSqRemx6/RGlw2R7pox/ypdMNpx9oe3seEvnqloKOq+XGU1Go+GNbWfXtnSu6j959oPZ8XRhPr0X7m0oZHQVoQlLLhVnht8BBN3B2fQTGl+Pof29CxLPtl4USmbqrDk4QOLpKVWP776SmbU7+06e/dCYGv2XJqxKLA9mAFv9DsKSS3FZchlF2+yxL7KUTFdn1u27MLmvd4mVBgPsQN/C2nu3XZtdNv253uOanqodqapsToqqkpolOZgCcZ6x5FJcgfjQg+xQ49h60rV0Ht/3woGzYxmtssRSApRwaEv7qXX7ehYnzpr/u3BttGfIaX3CIEsTJDfXlNYaPgETiPOMdYktrkB86EHWM3XkapFkujrz8I6r48/vv6DWRtaXnu7+6dX3t1wT2nloRd9w2xwfD/YaPiUgJmukoPMWjoUll2Jxq07aVfYIMiGyyZrhJ6o82DcncX/LNWo9wUpbViPhZ1svqX1y96XxdDZ6VCI5LuHbzMLlpMnvACy5FI+VWkbR3UiKYeraNx08re/RXVfGUtlqq8otE629S2ru23ZNuiM+64jpSpZZo34++H6+seRSPL5/2EHXNf3o6pBUpirz6M4r4i8dPKMWQtbIW2aSmbqqR3a+oerFtjN6B+6blSIqaqWXSfL9fGPJpXh8/7CDrmvakcmlK9HUf3/L2zIH4/OtGqyshWRz+2l1j+68Ip7ORrIRCM1Kjb3XoBmSVYtVEN8b2IKuZ+rh/8eWzuP7Ht5xdSSZqbMBkRXiYHx+zYPb35JKpGvTC/pJ+x1PiauVNeJrSd+SS/HYAL9RxGtdXfv6/ef2Pr//glolbP+fFaYv1RB7oOWt1HfNTvkdS4kT8jT560TZl7cYREJg61SMJhETeWrPa/u2dZ5oibiCpbLVkczLfxai58y437GUOF97p1pyKQ7rgjyKeCiafnjvVem9PUvtvTLU92eF3Z+ppuP1w46HMaOy5FIB7Ep8BClqM/dWfS7ZkZhjDfcGgPp0PARhYf+f19L29t7Rn2GG4Ot5x5JLcVhyGUaSqZmH+af0geh8G9tgXlGTSR4ez3TwT+s4+DYrwYyflVwqgCWXIfRTn3mEb2T6mBdLRBM25Yd5RSyTOvLc1HZ9rSWYcbPkUgEsuQySpjr7KP+QjjO7CiBRlfA7JBMgITQkOmg1xbbrazn4FkswY2fJpQJYI3WONLHso3y1v4eFr8x8mwqnbES2OUI0mz66NNv27lra32QJZmwsuVSAYSdjrDRZwvoE/y95iKVHdM3OitWKmSNFNDP0BceB99TScYUlmNHZOJcKYHNiAYroH2mOd7DqqH/6TMhm+zBHimTTw5dm97+/hp7TrS51ZL6e3y25FIe9z8BzfLivjVOHLKpnQ1lLwOYI0ewwJRcAwsKeT0VJLu4vXkRmPOykVxwV/z5v44q+nVw6bMcGqxYzgw3Z5pJLY2F2fkFIN9o8ZEOzucUqQEVflbdxYmID7x+x/teSixksqmPIGZnGKDu/lCEbs3+ggLHkUhwV+z73MTP1FJ+JQHjEBJsN2bnBHKkqM0KbS67+hTF2f8baX45mJZcKUJEllzSx7ON8KZumbtTVI63kYial79Ra2t5h08QEiCWX4qjI93ktf5voY05s9C1BRSsyAZvhpULjXNH64Ntr6Ts+OfqGphgq8qTng4p7n3fwur59nDPmQVzhTNgGUZoj9Ici47zgCAt7Ph0iU2v92h2rFqsAFfXP3svs1Hr+bEwllgGR7DivUk3ZS4dGbqcbUqYxyp6/s9JLAFhyKY6KaWzMEtYn+VwmS2xcsxxH01EruZgj9IciEzs/9Z1aS/sbbQS/zyy5FEfFXEm9wHv6elg07lU3reRiBkuPu1osx4Ebqm2Apb/nHUsuxVERJZd2jk+2cOWEJsuLZqLWoG+OkJbwJM5P0RB7PlHpXRB9Lb1ZcimOsk8uWcK6lo8JhCaUJKJpm9vTHCk1mZILQP+i6gqfQTnu58EtuRRH2VeLbebavoG1WSYikp3kicSUlSySVZHJ/0+0XRcjPb1Sp4ex5FIByrrk0sfM1BbeMqnpvaNpqxYzhyXD0fz0sNRYmL0fTeVlX6XHqsUqQFmXXNbyt2klOqn/pdr+WvtfNK/ojVTnr/t+38k1dL/K16t4n1jJpQKUbcllNxcNuT7LeNUma627mHnFoWhdfhvjW/8iSra60hr4LblUgLJsVMwS1g28Ly9JIazhUDQdrdS6cTPIoWhtfsc9ZadEaHtHpZVeLLmUPdU+yrBqbCtvjvfTmLduXtWpaksuBoCuqvr8t8F1XlFNemolzZZhbS4VotPvAPIpTXV2C28d1xQvo6lJ1lRatYUZRle0Lv/nJo2Fabuh7C7yRmAllwrR4XcA+bSJ6+IZasc1xcto6pJ1NgWMAeBQ1ejLNExI16urSc2shN5j3do8ltXWCscaUYunbJJLkqmZ7Vw+6Ub8weoTBagKCbhEOhH5+Yaff1zRiKqGZ9bNfPryFZffrqrcufnOqw/0HThdkOyCqQsefM3S19znd7zFkEGyfZHqAp2bIiH2vz/J/C+V+6jdg34HYMmleMomuWzkPcksVROa5mUkdckCVIUEXCwcS7/p+Dd9o76qPpnKpMK3bLjl4xv2b1h/MH5wbiKdmHb9ydd/LiQhbetrm+J3rMXSF6lOAxMekDuqnrPdvGOx7YU7hv98Ty4V92X2UVkklyRTM7u5cNwTU47F1PjUcr+aPIqIUF9VnwRIZVPhrGbDIsL2zu0XnTn/zDtCElKAGbUzuv2NtHgOxqYWuDonJBx4d7l3HvE9uVjJpVhUexHpp5BXZEWwlbcklGhdIfYdS8fCsVQslYwmKyrJZLIZ+Z/n/+ezyXRy5twpcx9YNXPVtsd2PjbzhQMvnPFgy4OnRkKRnvMWnvezxY2L9/sdazG01kwv/EF6T60hNSNNtK1cz4G+JxcruRRXp98BTEaGaHYHl+W1h9hgU/sKfdUaPOFQWK8/+frPX7P6mk8eSh5asrV967ysZiNhCaffdfK7vnTMtGMefmTHIzf4HWex7K2dUYQTflhof1u59hzr12b1vaRryaW4fL+amIwd/EkiTYF68Xim90yv2O7IDdUN8WnV0za93PHy6qpwVcfxM49/BuDs+Wevjafj8/2OrxgySLYt1lCckmvXq6vJxsrx/y0Q5xlLLsXV6ncAk7GVqwt+RTm9Z3q5VlMM6WDfwfquRFcNQDwVj7bH21c1Vje2NtU2Pbu1fetKgOf2PXdsLByriCqxzqopqbzMhjwWWh2m8/JynJopEMmlor7IAbDX7wAmqpWz4wlm5b378WCNfY1RFEWoiG7J7fH2hsd2PXYjEFJVmVU366kz55/5fFeia8tvt/z2vT989oevDYfCyXMWnPMjv2Mthv01jWmgoFWvR2h/Y4Tpvyza4Yqkze8AwJJLcan2INID1Psdynht4w1FOU4kGwnVJmuTfdV9xTvB+GhF04rdK5pWfGHw/Q3VDfFrT7j2X/2IyU97apqKW5uSaaqi+9w4Ux4r+IVTEQWi5GLVYsW3x+8AxivJ1MxBVhek+/FQGvsaK2n+J5Njb+2M4vcU7Lyi6IcsoATQ7ncQYMnFDyWXXHZwWRLCRaummnloZrEOZQKkLxxLFW5k/kgHXh0jU1suFzR7tFkDMY2SJZfi2+V3AOO1k9fmdQ6x0czumh1DCcQXxBTPjrrZ/f4cORLi0GvKpVtyYM4vllyKzU2/H4hi61h0cUx/H/OK2v4RS8fCUxJTfDrRGL9snjq/qBcxR+h6TbmcCy25VLidfgcwVi1c6csMsrM7Z1fcYMpK1i/hTGtNk3+dOJLHxEjNKPX/uU5t1h6/gxhgycUfO/wOYKz2cbYv09XM65hXUVPAVLo9tTOSRRvfMqSQ0HVZqVeNBeq8YsnFH61Ar99BjKaTFf39FGm09CAN8YaqqlRVJay7YYCtU+b7P67p0IWlPjSjxe8Aclly8YOqAlv9DmM0u7nQ15P7jO4ZllwqQAbJbq+f7f+4ptS8Uq4aSwD7/A4ilyUX/2z2O4DR7OMsX6/k5nXM8/9q1hTc/pppyXQoEoxzUc+5pdqRZEdQuiAPCMYHWolUDxLgXmNJGtPF7iU22KyuWbFQNlQu4w/MMF6unxeck2L3uaV6QdPidwCDWXLx1xa/AxjOXvy/ggtrODSnc06pN7KaEaQllHmpYWHRZn8YVWJlDA0HJ9mNTZyANeaDJRe/BbZqbB9n+R0CAMfsP6bUG1nNCFrq5yRTQakSA9CqEPHjSu2CZpM2a+CWDgjOh1qJVHsJ6HQwXSwPRFfg6b3Tq2qSNb6XokxhrJu2PHgXDz1nlVpV7It+BzAUSy7+C1zVWILp6X4aA5FcABYeXFiqPXjMCDqq6pMHqxuCt+x338n+zRQwfnu0WQ/5HcRQLLn4bysQqCvzg5wYqC7ASw4ssbnGytD6xmOCWULoXxCYC6sxCGSpBSy5+E81Bbzgdxi52jgpUPW3sXQs3NTTVI4rBlaslIQzm6cuCE5Dfi6NhUnOD9QF3zASwDa/gxiOJZdgWA8E5iqug+MCVw++ZP+SUu0iaoawbcrc4IxtGUriuFKoit2szRqY88Zgwf1wK4mbKTkQPceyhLWH+YGrB5/bOTcWS8UCVV1nJu7ZIDbk54qvKoVq2MBWiYEllyBZB/63K3SzKFXMhcHGShBZ1rrMkksZ2F3TFO+MTQncBcwREsuC3qi/T5u1w+8gRmLJJShUuwjAKNtulgS2OmDJgSXV0XQ0sPGZsXl85uqgn7hLoVF/nd8BjMaSS7A863cAh1jie+lpOGENh47Zf0ypDXAzOXbVzowHsvvxYBoL0z83qCXlNm3WFr+DGI0llyBRPYDPgyq7WRS4KrFcy/YtqwlnwoFtxDTDU9BHS6HUMqB/QVBLyU/5HcBYWHIJnrV+HryXeYH+8keykdDitsXWLbkE7aydFfy2llypOUEsxe/XZg3cPGJDseQSNKq78XEZ5Dgzg17XzIq9K6pD2VCgxuKYkSnoY7NWB7uH2GD9AZqt+bCSKLWAJZegegwo+skzwfS0Eg38/0RVpiq8oH1B3O84zNhtr5sd76qqL51SC0BqdtC+C63arLv8DmKsgvbmGQDVTmBDsQ+bYFrJtGUcv+v4amt7KQ1Z0Mdnrgp8ifgoqZlBa38smVILWHIJsqdx0zsUTZLpJVPVVJWpCi/bt8zaXkrAxobFfYeq6ksvuWSmBan9cY82ayBnUB+OJZegUu0HnizmIRNMD2Id87BWtK6otVH7wRYPV6WemLmqxu84JiQzJRKghcNKqtQCllyC7kXgYLEOVmrJJaQhWbVrVVC7ixrgsZmr04GeQ2xEISFTF4TS/FZt1la/gxivEv3QK4SqAo8W63BJphfrUHmzoH1BTWNvozXuB9D+6sb4lqkLSrPUMiBb43dySVLEc0A+WXIJOtW9FGlBsSSNQWvAHJNTWk4J23ovwZJBsvfMPb20uh4PJVvr9//VE9qsJXnxZMmlNPwB6C30QTLECn2IV3yEj9zwFt7ytXfyzubBj32ZL7/uKq76bgst9WPZ15TElKpFbYv68h+lmahnmo5N9ERrS68Rf7BsrZ8llz3arIGe+XgkpX9lUQlUk4g8AFxRyMNkizjG5WIufrSe+vv/m/++Mff+9ayftpWtq2qoaR/P/lbvWl2zr3FfKhlNlv4JrcR1VNUnn52+YmzVYb/7yA0c2nki4Vg3b/nJmlfuf/QfL2HvMxeDKI1Lnuc1X/rfAoU7Mv9KLhngYZ+OnRdWcikVbuT+84U8RJZI0b5IV3P15hnMOKo09m2+fc11XPe/jHP5gUg2Ejr95dMzVj3mrwySvXvuGaIiY6tiXXLxo5zxoW8dcd/G/11J24uncNX3Ps9bfnITp733rkLEOibZWr+OvFabtcuvg+eDJZfS8kdgXFf046E+r+PyQ354cj31nZdwyYRGITf1NFUvObCkJOuny8Wjs05IjGv+sOOu3kztoIuMl++7iGWX/ZaqetcTcNqy7rwGOR4ZX0ouHQRghvTJsuRSSlQzwH0UaEnkLBHfkksHHVX3cu/r/4a/+fVk9rN61+qa2mStTcvvg5a6OfGNjUsmf6mf7JrNgQ0r+OX1f8evbvxbtvx+cR7CmyBfWg4e0mb1u5fapFlyKTWq7RRocKX62AT3HM/N7KV3xl/z139/Ddd8KUFi2sf5+Gc2s3nqePYT0pCcsfUMseqx4uqO1PTfN/fU/PQI0WyIVF8dV//gK6x++22s/a8/w69zrfQX+4jrtVn3FfughWAN+qVI9TlEFgLz/Q4lXy7iot0XcdHfDvx9Ddd86R/4hy8tYUnPePfVEG+oWtG6om/z3M2+VZhXkgyS/d38s8jbYMmqug4WnLMWCcHyP2nh2e9n6dxez7Rjxv2/MGmSKmZpfh/weBGPV1BWcild9wF5/bKF6C/a1f6H+fD7vsAXPtlDz+y38bavfotvvSqf+1+5Z2XNlPgUm3usCB6ddUKiIzY1fzMezzj+WVrXrQRg1xOz0GyExsXFTywAUrTZhRLAPeVQHTZA3CBwU5JEpgNvBPLS/fYPfDXeQYnOAzWEeDSefmD1A6TDaSuhF0hL3ey+u+afNfES4p0ffh89rceS6a8nEutm0QW/5tT3Ps69f3cDfQcXIqE0x119G6veuimPYY/dvC/FmfJYob8TCvym1CamHI0ll1Insgi4DJh08f0J1vQd4LSyqkpqq29LPHbsYzFk8u+POVJXtDb5v4svipbu3GFjMP+mBPVPVxf4KH/UZn22wMcouvL9p6gUqjtwi4tNWpjy68U7o2dG9apdq8rvhfmsLxxL/Xrhq8JlnVigGG0u28sxsYAll/Kguh54YbK7iRR3+ZiiWbZ/We289nk2PUye9Es4c/vC8zQeqS7/6sZQQa9LDgH3F/IAfrLkUj4eBSa1BGqEeNnWkZ7acqo18OdBBsn+dsE56ZJbsniiogcLtWBYGrhbm7XofZ2LxZJLuVDNAvfgRvdOSIyOsm2XCGlIznnpnEhVqsoWF5sgBb137unJfTXTizfDqa+ySrijUMnlIW3Woq3V5AdLLuXErV75W1xxe9xq2Ve2yQWgOl0dOWvLWdlQNlQ23T2L6bGZq+MtU+aWTW/CUYX6MogW4jvxiDZrUZbR8JMll3Kj2gPcDox70rsa9pf9/8O0vmmxs7ac1S9ZsQQzDuumLetdP21pWfUkHFX4UCGmWXpSm3XS7aOloOxPJhVJtReXYDrH87Ra9hWqCiBQZnbPrD7j5TOSNkXM2DzXuLT3iZmr6vyOo+jCnfm+AFmnzbo2z/sMLEsu5Uq1D7iDcSSYGB1hyFTECXdO15yaU1tOTViCGdmTTSt7H5+1uvISC0CkI5//Gxu1WZ/I4/4Cz5JLOXMJ5nbG2MgvqFTRky5sUMGxoH1BzYk7TrQxMENQ0Ednru5b23RsZSYWgOiBfCWXrcAjedpXybDkUu5U47gSzJjWganmYEGm8w+qJW1Lao/fdbyNgcmRBX1w9snximtjGSz2cj4a83cA92tz5U2FYsmlEhxOMG2jbVrPzopr6F6+b3ntsXuOPWpVzEqURbL3zD0j8VLDospOLADVWyY7SHQPbixLxX2nwJJL5VBNAL8GXh5pswbKvofkkFbuXVl34vYT+yq5DSYtoezv5p/VX1HdjYeVzlK1ZzITwr4M/FabtaJqAnKV//QN5jDVNHAPIqcBZwy1SSObK/aCY0nbktpYOhZ/+pinYxrSinof+sKx1G8WnKMdsamFnqSxNETb0khmorMQPK/Nmpf5/kpZRX2BjEf1GeBu3BQUR2hga2VM6zGMuZ1za8576bxUOBOumCvO1upp8VuXXBLK65ospa5q10Q6tijwmCUWx5JLpVLdBvwfgxYci5AIxeiomB5jQ5neOz12wYsXZCphqpjnG4/pvX3hq6r7w9GKGOM0ZrGXx1s9mgHu1WZ9vhDhlCJLLpVMtR34BdCae3cdu8r+pDqaKYkpVRduvFBqkjVlObFgWkKZe+ecFn9s1gl1KlLW0/5MSPWW8Zwbk7jFvkZsz6w0llwqnWvovwPYMHDXdDZWZO+WwWpSNZGLNl4UbupuKquxMD2R6v5fLLows3XqfGu4H1JWqV031irCbuBX2qyto25ZYSy5GDejsuofgN8B8Rk8a1UknmgmGj7vpfNqjt1zbG859CTbUTur77bFF4c7Y1OsfWU40dZ+wn1j+Q604hJLZ4EjKknWW8wcproDkdsaeekCIb1IKfNVBsdh5d6VdU3dTYmnlj0VSUVSJfe9SYai6Ydmn5TaNmWejV8ZTe360TpzKPAM8EwlDo4cK1F7b8wQ3ivrL72cHQtjZK0UkyMRTaSfWP5E+lDtoZLpsttSN7vvwTmnxJLhKvssx2LeV+JM+cNwVYY9wH1WDTY6Sy5mSCKcNp3EKZ/gmdRq2q1uPoei+vzC5/u2z9oe6Hm3kqFo+sHZJ6dsUOR4ZJXl12UJ9wyViLfhFvlKFjuqUmTJxQxJhFnA1QCvZUf8/WyoqiVjV745Whta4+sWr4v0R/snM5K7IFrq5sQfmHNKlXUxHqfonn6W/tng9qg0bvzKRj9CKlWWXMyQRBDgBqAKYCrJzJ/zfPJcWmtCYF1XPelQOrtu8brEnul7AtGW0R2p6f/DrBMzO+pnW2llIhru6mXOv+SWSNtx41cmvHx4pbLkYoYlwmuAZbn3LaS7/31syJxGm528chyYciCx9pi14WQ06UspJhGKpp+ccVz/iw2La2zcyiQs+GyCunXVQBZYj1s5smJma8gnSy5mWCIsAS4d6rFldCY/wIbsKjosyXjSoXR2w8IN8R1NO2qR4pTu0hLKPD9taWLt9BU16ZD17psU6cuw4toQovuBh7VZx7RMhRmaJRczLBHCwPXAsFfjJ9KWeB8vsJTS6T1VaB11Hcmnj3la4rF4wcaSZEG3TF0Qf3zGqlgiErN2lXyoXdvJws89p836ot+hlANLLmZEIrwaWD7admeyL/4+NoTm0RcrQliBp6i2zGxJbJq3KZrPcTFZ0B11s+OPz1wVPVRVH7iOBCVJUbbPj7Nl0QO6v2mX3+GUC0suZkQiLAYuG+v2F7E7/m42hmeQsBHguKqyzXM2x1+e/XJ1NjTxMUP9Es5saliUWDd9eawvUl1ygzgDq60xzh9PCtM1NQvcrIpNfZQnllzMiEQI4arGxpwsBNXT2Z+8khY9mbZYhMpaG2UoyUgys3H+xsTOpp3jao/pjtT0Pz9taXpjw+LqTChc8e9j3rRPTfDMajjQNFCd+6IqD/kaU5mx5GJGJcLFwLETeW49/ZnL2JG8lB1hqzKD3qre1PpF61P7p+6vGSnJ7KueFl87fQXWpTjPOqckWLsKWmcObiO8Q5U9vsRUpiy5mFGJsAC4fLL7WUFn/5VsS59Ha6y6wgdk9sR6UpvmbUrtnba3WsWV7Hoi1f1bp8xPbWxYVGXtKXnWMTXJcyuz7BkyWXepckvRYypzllzMmIhwLTA1H/uKksleyJ7E5WyX5XRWV/KgzH2RbOLWBYf6160KRfY1TA3EQMyykUXZMzvOhhVh2htHKjU/poot8pVnllzMmIhwAnBevvc7k77UBexJncYBWUFnRUwx00Z1/xPMTt/Pgsgmprm2rFAmy/IdCVa+HKG+cF2YK0IqnOHlhQk2Lo8RH7XzQxrXkF+Wi8L5yZKLGRMRqoDrGGHMSz4sozN5BvvTp3IgvJyuaDnMytxFVeplpqY2Mk0fZW7VdqaO/B42dSQ5dluGBftiRMo/2eaForRNS7BlsbJjbjXZMXd+2KjKwwWNrUJZcjFjJsJ5wAlFOx6qx9GROoP9qVNoCx1DVywa8J5n/YSyO6nv38S0zHqmhzbQFG1ngl2HQ5ksi/YmWNEiNHVWF2vUf0npqeln24IUW5bESMQm8j7fpoqNxC8ASy5mzERoAN7u1/HDZPVYOvuX0ZVdQrcuoEfm0BueRjLqR7tNGskepDq9hYb0BppYz/RIC1OjSgHm9qpJpFm0J8XCvdDUGSMU7CRbUD21SXbPSrNtQYSOEdtSRrNXldvzFpc5giUXMy4ivB5Y6HccucJkdR69qTn0ZebQp3Po1VnEZQaJUCPJUBiVEEqYrAhIGCXk3ed+HpmYkoQy3VRluqjKdhDLtlOtbVTTRo20URNqozp0kOpwHz5NZx9JZ1nQmmTRHmX2wfKvOsuidDQm2DlH2T4vSl9tvqpm71Hl5TztywxiycWMS766JQdNmKxGyGoWIUUpDVZUZXpXP3P3p5ndFqKpq6rkk00W5dCUfg5MT9M6I0zrzCrSeZ+U8xBwq43ILxxLLmbcRLgamOV3HGYYjV39zNufpqkTGrvD1PZVEQpwe02iKkV3XZq2aVlaZ0Y4MC1KpuAzPN+nypYCH6OiWXIx4ybCPOBKv+MwYyRZpaE7RVNnmqZOZUpviNp4mOpkmEiReuNlUZKxDPFYmkP1WToaoL0hREdDlFTRqxc7cA35dvIrIEsuZkJEuBxY4HccZpKiqQz1vWmm9GapiyuxfqhKKdGUUJUWoikhmhZCuecJlZxflXQYMuEsqaiSjEIqqiSqoKdO6KkN01sbJh4LU4iODhNzlyotfgdR7iy5mAkRYQbwZr/jMGacDqjyS7+DqAQl1HBpgkSVNrCeNqbkPOl3AJXCkouZjCfBetuYkrFXFVsMrEgsuZgJU6ULeMnvOIwZAwWe8DuISmLJxUzWH4GE30EYM4oXVNnvdxCVxJKLmRRVEsDjfsdhzAj6sLaWorPkYiZNlZfA6rJNYP3BptQvPksuJl8ewa2NYUyQbFdlm99BVCJLLiYvVDkEPO13HMbkSAF/8DuISmXJxeTTc0Cb30EY43lKlR6/g6hUllxM3nhzNT2EjX0x/tsNrPc7iEpmycXklTdy33qPGT/FgfttYkp/WXIxeafKemxqGOMPxU2n3+d3IJXOkosplIdwCzIZU0xrVdntdxDGkospEG9cwT1Axu9YTMXYi/VYDAxLLqZgvPaXR/2Ow1SEBK46zNpZAsKSiykoVTaCLSdrCioL3KtKr9+BmMMsuZhieBBo9TsIU7YesnaW4LHkYgpOlQzwe6DT51BM+XnKm9vOBIwlF1MUqiSBO8G6iJq8eVGVZ/wOwgzNkospGm8qjt+CzVBrJm0XbrJUE1CWXExRqXIQuAubIsZM3EHgblX7HwoySy6m6FTZA9wP1m3UjFsn8FtVUn4HYkYW8TsAU5lU2SqCABdjFzlmbNqBO7zVT03AiapdPBr/iLAEeC0VlWAemgbvuhF6prq/X/cQ/Ow+2FgLl38AOpug8SD87j9gpXWAcNqA33gdQ0wJsORifCfCAuBSKqYk/UQDvNAAN+6A7TE45bPwH9+Gb58HDb3wf7+Dq/8EDtXCfb/wO9oA2A/caUsVl5YKulo0QaXKLlwvsgqpRz+7yyUWgMVJmLUXXmqEdSfDpx5z93/qMXj2FL8iDJBWXInFEkuJseRiAkGVvcBvoNKqPe5pgj0L4Z3boHeqSzwAZ3S5vyvablyJpUIuOsqLJRcTGKrsB26HSlmadmcMrvsgfPBWOGZQI3UYkEqus34B1yss7XcgZmIsuZhAUaUd+CWwz+9YCqs7DBd9EC58Av5xrbuv7pBrjwH3s7bbv/h8kwUeUeURG8dS2iy5mMBRJY4rwZTpnFEZ4ILrYf5e+Pk9h+8/eR18+Vz3+5fPdX9XlCSuGuwFvwMxk2e9xUygiXACcA5ldSH0zeXwkY/DjN2Hq74+9Eu4Zhtc8QHomg6N7XDnd+H4SumK3AH8XtVWLy0XllxM4IkwBzcWptbvWExBbMct9GUN92XEkospCSLU4kbzL/A5FJM/GeAJVdb7HYjJP0supqSIsApXTVYhAy7LVgdu9ch2vwMxhWHJxZQcEaYClwCz/Y7FTMjzwJPWzbi8WXIxJcmb9PJk4AzKqrG/rB0CHvQGzJoyZ8nFlDQRpuPaYmb4HIoZXhbYgFuS2BrtK4QlF1PyvFLMccCZQLXP4Zgj7QIeVaXT70BMcVlyMWVDhCrgNOAErKrMb4eAx1TZ7ncgxh+WXEzZ8Rr8zwGW+BxKJUoBa4HnVcn4HYzxjyUXU7ZEmIdLMtYeU3gZYBPwjCqVMquAGYElF1P2RFgInALM9TmUcpTGzWD8nCUVk8uSi6kYIszGJZnFPodSDpK4HmDrbU17MxRLLqbiiDANN0ZmOdbwP159uEGQL1i3YjMSSy6mYolQBxzr3Rp8DifIssAO4EVgpyp20jCjsuRiDK9Uma0ElgJVPocTFB24RvrN3ho7xoyZJRdjcogQwXVhPhaYD4ivARVfH24K/E3estPGTIglF2OG4Q3KXAAs8n6W63oy+3HVXjtUafM7GFMeLLkYM0YizAAWerfZlG6pJgnsxiWUnVblZQrBkosxE+CVamZ4t5nezyB2CkgDB4EDuBLKAVW6/A3JVAJLLsbkyaCE0wRMAepx1WmF7vLcD3R7tx5cY/x+oEOVbIGPbcxRLLkYU2DerM11uEQzcKsFooNuoZyb4ObpSns/B34f+DvJ4WTSrUqyeK/ImNFZcjHGGJN3NjrZGGNM3llyMcYYk3eWXIwxxuSdJRdjjDF5Z8nFGGNM3llyMcYYk3eWXIwxxuSdJRdjjDF5Z8nFGGNM3llyMcYYk3eWXIwJKBG5SURSItKTc1ua8/gpIvK0iPR5P0/JeUxE5KsictC7fVVExrREgHdcFZG/HnT/X3v335Sv12jKlyUXY/LIO6lP+HslIrMH3XWLqtbn3F72tqsCfgXcDEwDfgj8yrsf4APA1cDJwEnAVcCfjSOUl4DrB913g3e/MaOy5GIqkoh8UkR2i0i3iGwSkdeISFhEPi0iW737nxaRhd7254nIkyLS5f08L2dfD4jIF0XkD7hlgpeKyHEicreItHv7v2aEWKIi8iYR+TWwZYwv4WIgAvyzqiZV9Vu4mZRf7T1+A/B1Vd2lqruBrwPv9o63xCuB3CAiO0SkTUQ+M2j/TwK1IrLae85qoNq735hRWXIxFUdEVgIfBs5U1SnAZUAL8FHgHcDlwFTgPUCfiEwHfgN8C7dOyzeA34hIU85u34UrLUzBLcx1N/ATYBZwLfBtEVk1KI4TReQbuFUhP+EdY+GgcK/yEtQGEflQzv2rgef0yGnNn/PuH3h8Xc5j63IeG3A+sBJ4DfA5ETl+0OM/5nDp5Qbvb2PGxJKLqUQZIAasEpGoqrao6lbgfcBnVXWTOutU9SBwBbBZVX+sqmlV/SnwIq6qacAPVHWDqqaBPwFaVPX73vZrgf8F3gYgIq8WkaeAO4EEcL6qnquq31XVzpx93gocj1vp8v24BPAO77F6OGpFyS5cchvq8S6gflC7yxpVjavqOlzyOXnQ/m4G3iEiUVyCvHmE99SYI1hyMRVHVbcAHwFuAvaLyM9EZB6u1LB1iKfMA7YPum87MD/n7505vy8GzhaRzoEbcB0wx3t8FrAcWI87qe8YJs4XVHWPqmZU9VHgm8BbvYd7cKWrXFNxi4cN9fhUoGdQSac15/c+XELKPf4OXDXdl3DJNfc1GjMiSy6mIqnqT1T1fFwiUOCruASxbIjN93jb5VqEq856ZZc5v+8EHlTVxpxbvap+yDv2z3CJ5sfAe4E9IvI9ETl/tLBx7SoAG4CTBpVETvLuH3g8tyRycs5j4/Ej4GPeT2PGzJKLqTgistKrmorhqqXiQBb4T+DzIrLC6/V1kteucidwrIi8U0QiIvJ2YBVwxzCHuMPb/l1eY31URM7MbdNQ1YSX4C7FnfhbgO+LyCslJxF5o4hM82I5C/grXA8xgAdw1Xt/JSIxEfmwd/993s8fAR8VkfleqexjwA8m8HbdAlyKq6IzZswsuZhKFAO+ArThqoZmAZ/CNdTfCtwFHAL+C6jx2l2uxJ2gD+Ia369U1bahdq6q3bgT8rW4Uk8rrmQUG2b7nar6RVVdgWs4H3AtrlqqG5csvqqqP/Se04/ranw90InrfHC1dz/Ad4Hbgedx1W+/8e4bF69N5h5VjY/3uaayyZFVsMYYY8zkWcnFGGNM3llyMcYYk3eWXIwxxuSdJRdjjDF5Z8nFGGNM3llyMcYYk3eWXIwxxuSdJRdjjDF5Z8nFGGNM3llyMcYYk3f/H8/gwXAQKXSZAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "rb = P.promiscuous_binders(n=3, cutoff_method='rank', cutoff=10)\n",
    "pb = P.promiscuous_binders(n=3, cutoff=.95)\n",
    "pb2 = P.promiscuous_binders(n=3, cutoff_method='score', cutoff=500)\n",
    "fig,ax=plt.subplots(1,1,figsize=(6,6))\n",
    "ep.utilities.venndiagram([pb.peptide, rb.peptide, pb2.peptide],\n",
    "                      ['global percentile <5%','rank<10', 'score<500nM'],ax=ax)\n",
    "txt=ax.set_title('promiscuous binders derived \\nusing the three cutoff methods')\n",
    "plt.savefig('promisc_comparison.png',dpi=150)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute clusters of promiscuous binders\n",
    "\n",
    "Here we have detected clusters of MHC II binders. These clusters vary in length and for peptides synthesis we must create a list of 20-mers covering shorter sequences or splitting longer sequences into 2. For the 20mers hydrophobicity and net charge are calculated."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>name</th>\n",
       "      <th>start</th>\n",
       "      <th>end</th>\n",
       "      <th>binders</th>\n",
       "      <th>length</th>\n",
       "      <th>locus_tag</th>\n",
       "      <th>20mer</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>313</td>\n",
       "      <td>342</td>\n",
       "      <td>7</td>\n",
       "      <td>29</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>QDFLALIPIDQILAIAADYL</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>313</td>\n",
       "      <td>342</td>\n",
       "      <td>7</td>\n",
       "      <td>29</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>DQILAIAADYLANDAEVQAA</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>148</td>\n",
       "      <td>168</td>\n",
       "      <td>7</td>\n",
       "      <td>20</td>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>KLTWADFYFVAILDYLNHMA</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>Bla-g-4</td>\n",
       "      <td>123</td>\n",
       "      <td>142</td>\n",
       "      <td>6</td>\n",
       "      <td>19</td>\n",
       "      <td>Bla-g-4</td>\n",
       "      <td>NGHVIYVQIRFSVRRFHPKL</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>174</td>\n",
       "      <td>201</td>\n",
       "      <td>5</td>\n",
       "      <td>27</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>EFKNFLNFLQTNGLNAIEFL</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>174</td>\n",
       "      <td>201</td>\n",
       "      <td>5</td>\n",
       "      <td>27</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>FLQTNGLNAIEFLNNIHDLL</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>408</td>\n",
       "      <td>432</td>\n",
       "      <td>5</td>\n",
       "      <td>24</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>NGLIDDVIAILPVDELYALF</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>408</td>\n",
       "      <td>432</td>\n",
       "      <td>5</td>\n",
       "      <td>24</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>DDVIAILPVDELYALFQEKL</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>9</td>\n",
       "      <td>42</td>\n",
       "      <td>4</td>\n",
       "      <td>33</td>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>LFAVATITHAAELQRVPLYK</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>9</td>\n",
       "      <td>42</td>\n",
       "      <td>4</td>\n",
       "      <td>33</td>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>QRVPLYKLVHVFINTQYAGI</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>220</td>\n",
       "      <td>243</td>\n",
       "      <td>4</td>\n",
       "      <td>23</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>TGLIDDIIAILPVDDLYALF</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>220</td>\n",
       "      <td>243</td>\n",
       "      <td>4</td>\n",
       "      <td>23</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>IDDIIAILPVDDLYALFQEK</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>126</td>\n",
       "      <td>147</td>\n",
       "      <td>4</td>\n",
       "      <td>21</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>FLALIPTDQVLAIAADYLAN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>Bla-g-6.0201</td>\n",
       "      <td>61</td>\n",
       "      <td>84</td>\n",
       "      <td>3</td>\n",
       "      <td>23</td>\n",
       "      <td>Bla-g-6.0201</td>\n",
       "      <td>ELEFEEFCTLASRFLVEEDA</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>Bla-g-6.0201</td>\n",
       "      <td>61</td>\n",
       "      <td>84</td>\n",
       "      <td>3</td>\n",
       "      <td>23</td>\n",
       "      <td>Bla-g-6.0201</td>\n",
       "      <td>FEEFCTLASRFLVEEDAEAM</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>319</td>\n",
       "      <td>336</td>\n",
       "      <td>3</td>\n",
       "      <td>17</td>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>HSDHFFIGDFFVDHYYSEFN</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>Bla-g-6.0101</td>\n",
       "      <td>62</td>\n",
       "      <td>77</td>\n",
       "      <td>3</td>\n",
       "      <td>15</td>\n",
       "      <td>Bla-g-6.0101</td>\n",
       "      <td>SGELEFEEFVSLASRFLVEE</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>103</td>\n",
       "      <td>117</td>\n",
       "      <td>3</td>\n",
       "      <td>14</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>VDHIIELIHQIFNIVRDTRG</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>17</td>\n",
       "      <td>30</td>\n",
       "      <td>3</td>\n",
       "      <td>13</td>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>ALGEPIRFLLSYGEKDFEDY</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>292</td>\n",
       "      <td>310</td>\n",
       "      <td>2</td>\n",
       "      <td>18</td>\n",
       "      <td>Bla-g-2</td>\n",
       "      <td>RNFNISSQYYIQQNGNLCYS</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>Bla-g-6.0301</td>\n",
       "      <td>32</td>\n",
       "      <td>48</td>\n",
       "      <td>2</td>\n",
       "      <td>16</td>\n",
       "      <td>Bla-g-6.0301</td>\n",
       "      <td>ISTNMVEEILRLMGQPFNRR</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>477</td>\n",
       "      <td>492</td>\n",
       "      <td>2</td>\n",
       "      <td>15</td>\n",
       "      <td>Bla-g-1</td>\n",
       "      <td>VDVDHFIELIKKLFGLSH</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>Bla-g-6.0101</td>\n",
       "      <td>28</td>\n",
       "      <td>43</td>\n",
       "      <td>2</td>\n",
       "      <td>15</td>\n",
       "      <td>Bla-g-6.0101</td>\n",
       "      <td>CISTEMVGTILEMLGHRLDD</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>Bla-g-6.0301</td>\n",
       "      <td>66</td>\n",
       "      <td>81</td>\n",
       "      <td>2</td>\n",
       "      <td>15</td>\n",
       "      <td>Bla-g-6.0301</td>\n",
       "      <td>SGRLEFDEFVTLAAKFIIEE</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>Bla-g-4</td>\n",
       "      <td>34</td>\n",
       "      <td>47</td>\n",
       "      <td>2</td>\n",
       "      <td>13</td>\n",
       "      <td>Bla-g-4</td>\n",
       "      <td>FRGSWIIAAGTSEALTQYKC</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>97</td>\n",
       "      <td>109</td>\n",
       "      <td>2</td>\n",
       "      <td>12</td>\n",
       "      <td>Bla-g-5</td>\n",
       "      <td>DTISDFRAAIANYHYDADEN</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "            name  start  end  binders  length     locus_tag                 20mer\n",
       "0        Bla-g-1    313  342        7      29       Bla-g-1  QDFLALIPIDQILAIAADYL\n",
       "1        Bla-g-1    313  342        7      29       Bla-g-1  DQILAIAADYLANDAEVQAA\n",
       "2        Bla-g-5    148  168        7      20       Bla-g-5  KLTWADFYFVAILDYLNHMA\n",
       "3        Bla-g-4    123  142        6      19       Bla-g-4  NGHVIYVQIRFSVRRFHPKL\n",
       "4        Bla-g-1    174  201        5      27       Bla-g-1  EFKNFLNFLQTNGLNAIEFL\n",
       "5        Bla-g-1    174  201        5      27       Bla-g-1  FLQTNGLNAIEFLNNIHDLL\n",
       "6        Bla-g-1    408  432        5      24       Bla-g-1  NGLIDDVIAILPVDELYALF\n",
       "7        Bla-g-1    408  432        5      24       Bla-g-1  DDVIAILPVDELYALFQEKL\n",
       "8        Bla-g-2      9   42        4      33       Bla-g-2  LFAVATITHAAELQRVPLYK\n",
       "9        Bla-g-2      9   42        4      33       Bla-g-2  QRVPLYKLVHVFINTQYAGI\n",
       "10       Bla-g-1    220  243        4      23       Bla-g-1  TGLIDDIIAILPVDDLYALF\n",
       "11       Bla-g-1    220  243        4      23       Bla-g-1  IDDIIAILPVDDLYALFQEK\n",
       "12       Bla-g-1    126  147        4      21       Bla-g-1  FLALIPTDQVLAIAADYLAN\n",
       "13  Bla-g-6.0201     61   84        3      23  Bla-g-6.0201  ELEFEEFCTLASRFLVEEDA\n",
       "14  Bla-g-6.0201     61   84        3      23  Bla-g-6.0201  FEEFCTLASRFLVEEDAEAM\n",
       "15       Bla-g-2    319  336        3      17       Bla-g-2  HSDHFFIGDFFVDHYYSEFN\n",
       "16  Bla-g-6.0101     62   77        3      15  Bla-g-6.0101  SGELEFEEFVSLASRFLVEE\n",
       "17       Bla-g-1    103  117        3      14       Bla-g-1  VDHIIELIHQIFNIVRDTRG\n",
       "18       Bla-g-5     17   30        3      13       Bla-g-5  ALGEPIRFLLSYGEKDFEDY\n",
       "19       Bla-g-2    292  310        2      18       Bla-g-2  RNFNISSQYYIQQNGNLCYS\n",
       "20  Bla-g-6.0301     32   48        2      16  Bla-g-6.0301  ISTNMVEEILRLMGQPFNRR\n",
       "21       Bla-g-1    477  492        2      15       Bla-g-1    VDVDHFIELIKKLFGLSH\n",
       "22  Bla-g-6.0101     28   43        2      15  Bla-g-6.0101  CISTEMVGTILEMLGHRLDD\n",
       "23  Bla-g-6.0301     66   81        2      15  Bla-g-6.0301  SGRLEFDEFVTLAAKFIIEE\n",
       "24       Bla-g-4     34   47        2      13       Bla-g-4  FRGSWIIAAGTSEALTQYKC\n",
       "25       Bla-g-5     97  109        2      12       Bla-g-5  DTISDFRAAIANYHYDADEN"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pb=P.promiscuous_binders(n=3)#, cutoff_method='score', cutoff=500)\n",
    "#print (pb[:20])\n",
    "reload(analysis)\n",
    "cl = analysis.find_clusters(pb, min_binders=2)\n",
    "final = analysis.create_nmers(cl, prots, length=20, margin=2, key='20mer')\n",
    "final"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "reload(plotting)\n",
    "bc = ep.get_coords(pb)\n",
    "r1 = plotting.binders_to_coords(bc)\n",
    "c = ep.plotting.binders_to_coords(cl)\n",
    "f = ep.plotting.plot_overview(prots, coords={'clusters':c,'MHC-II 11mer':r1},\n",
    "                              cols=2, figsize=(13,4))\n",
    "plt.savefig('cockroach_clusters.png',dpi=150)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "reload(plotting)\n",
    "#print (P.get_names())\n",
    "prot='Bla-g-6.0101'\n",
    "ax = plotting.plot_multiple([P],names=prots.locus_tag,cutoff=.95, cutoff_method='default',n=1, legend=True, figsize=(12,5),regions=cl)\n",
    "#plt.show()\n",
    "plt.savefig('cockroach_%s.png' %prot, dpi=150)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "reload(plotting)\n",
    "from bokeh.io import show, output_notebook\n",
    "output_notebook()\n",
    "p = plotting.bokeh_plot_tracks([P],name='Bla-g-5',cutoff=.95,n=2,width=800)\n",
    "show(p)\n",
    "p = plotting.bokeh_plot_tracks([P],name='Bla-g-5',cutoff=500,cutoff_method='score',n=2,width=800)\n",
    "#p= plotting.bokeh_plot_grid(P, name='Bla-g-5', width=800)\n",
    "show(p)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Test for effect of cutoffs/method on epitope discovery rate\n",
    "\n",
    "195 peptides ranking in the top 40% of predicted affinities for 7 or more of 20 HLA class II alleles were selected in Sette study. These were arranged into 13 pools which were were tested with PBMC cultures for production of IL-5 (Th2) and IFN-γ (Th1). 32 unique peptides were identified that elicited a positive response in at least one donor. These were joined into overlapping contiguous epitopes defining 25 distinct antigenic regions 15–20 in length. Their results indicated that a small number of epitopes that elicit the most dominant and prevalent responses encompass a significant fraction of the total response in this population of donors.\n",
    "Based on previous timothy grass studies, which had shown that predictions of this level of stringency would identify ∼75% of the total response detected with complete sets of overlapping peptides\n",
    "Here we test the use of "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEKCAYAAADpfBXhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvDW2N/gAAIABJREFUeJzt3Xd0VNXax/HvnmRSSAVCD5AAoYcamiC92rChYAP1iqBysWBv2K+A7Vp4RUVA9CoqotIxdKQkIL2GmlBCTwipk9nvH2eAEFImZJIpPJ+1sjJz5syZZ85KftnZs8/eSmuNEEIIz2JydgFCCCEcT8JdCCE8kIS7EEJ4IAl3IYTwQBLuQgjhgSTchRDCA0m4CyGEB5JwF0IIDyThLoQQHsjbWS8cFhamIyIinPXyQgjhltavX39Sa12luP2cFu4RERHEx8c76+WFEMItKaUO2rOfdMsIIYQHknAXQggPJOEuhBAeyGl97gXJyckhKSmJzMxMZ5fiVH5+foSHh2M2m51dihDCTblUuCclJREUFERERARKKWeX4xRaa06dOkVSUhKRkZHOLkcI4aZcqlsmMzOTypUrX7PBDqCUonLlytf8fy9CiNJxqXAHrulgv0DOgRAezJIN5bACnsuFuxBCeKzD6+HL62HzjDJ/KQl3Ow0bNoxffvmlyH0iIiI4efKk3cecMmUKTzzxRGlLE0K4OksWxL4FX/eBEzsh7qsyb7271AeqQgjhcY5shFkj4fh2QEGnJ6DnK1DG3a/Sci/AW2+9RaNGjejSpQtDhgxhwoQJlz0eGxtL69atiY6O5qGHHiIrK+viY+PGjSM6Opr27duTkJAAwJ9//kmHDh1o3bo1vXv3Jjk5uVzfjxDCCSzZsORd+KqnEeyV6sFD86HfO2D2L/OXd9mWe8QLc8rkuAf+c2ORj8fFxfHrr7+yadMmcnJyaNOmDW3btr34eGZmJsOGDSM2NpaGDRvywAMPMHHiRJ588kkAQkJC2LJlC9OmTePJJ59k9uzZdOnShTVr1qCU4uuvv2bcuHF88MEHZfL+hBAu4NgW+G0kJG8x7ncYCb1eA58K5VaCtNzzWbVqFQMHDsTPz4+goCBuvvnmyx7ftWsXkZGRNGzYEIChQ4eyfPnyi48PGTLk4vfVq1cDxvj9fv36ER0dzfjx49m2bVs5vRshRLnKzYFl42BSdyPYK0bAsDkw4D/lGuzgwi334lrYrirvMMYLt0eNGsXTTz/NLbfcwtKlSxk7dqyTqhNClJnkbUbf+tFNxv32w6H3WPAJcEo50nLPp3Pnzvz5559kZmaSlpbG7NmzL3u8UaNGHDhw4GJ/+nfffUe3bt0uPv7TTz9d/N6pUycAUlJSqFWrFgBTp04tj7chhCgvuRZYPgG+7GYEe2gdGPon3DDeacEOLtxyd5Z27dpxyy230KJFC6pVq0Z0dDQhISEXH/fz8+Pbb79l0KBBWCwW2rVrx4gRIy4+fubMGVq0aIGvry//+9//ABg7diyDBg2iYsWK9OzZk/3795f7+xJClIHjO2HWCDjyj3E/5iHo8yb4Bjm3LkDpcrhSqiAxMTE6/2IdO3bsoEmTJk6pJ6+0tDQCAwNJT0+na9euTJo0iTZt2pRrDa5yLoQQBcg+D39/Cis+gNxsCA6HgZ9C/Z5l/tJKqfVa65ji9pOWewGGDx/O9u3byczMZOjQoeUe7EIIF2W1wuYfIfZNOHfU2NbmAej7DvgFO7e2fCTcC/DDDz84uwQhhKvZvxwWvAzHNhv3a7SEfu9CRBfn1lUICXchhCjKid2w6DXYPc+4H1zLGLMefReYXHdMioS7EEIU5PwpWPoexE8GnQs+gdDlSej4eLmPWb8aEu5CCJFXTias+xKWfwBZKaBM0GYo9HgZgqo5uzq7SbgLIQQYszRumwl/jYWzh4xt9XtB37ehWlOnlnY1JNwdbOzYsQQGBjJmzBhnlyKEsFfiOljwEiTFGferNoW+b0GD3s6tqxQk3IugtUZrjcmFPzQRQpTCmQNGS33bb8b9gKrQ82VodR94uXc8unf1ZeDAgQP069ePDh06sH79etq3b8+WLVvIyMjgzjvv5I033gCMhTmGDh3Kn3/+SU5ODj///DONGze+7FhfffUVM2fOZObMmfj7l/0Un0LY7cQuWD4eTu+z+ynZaL43ZRBrysJahqWVq5x0ozumZg2jPz2wOhxbCPMXFvvUQQ0HcVvUbeVQ5NVx3XAfG1L8Pld13JRid9mzZw9Tp06lY8eOnD59mkqVKpGbm0uvXr3YvHkzLVq0ACAsLIwNGzbwxRdfMGHCBL7++uuLx/jss89YtGgRs2bNwtfXt2zeixAllXrUGAHyz3eg7Y/olf5+vF+pIge8zWVYnBP4+ly6nX0aTp+2+6ndwrsVv5MTuW64O1HdunXp2LEjADNmzGDSpElYLBaOHj3K9u3bL4b77bffDkDbtm2ZOXPmxedPmzaN2rVrM2vWLMxmD/tlEO4pMxVWfQKrPwdLBigvYx6UFoPB5FXo05IyTjBu9/9YctKYOyWiQnX+Xf8OqvlWLK/Ky5Z/RQiqflVPrVbBtUfOuG6429HCLisBAcZMbvv372fChAnExcVRsWJFhg0bRmZm5sX9LrTIvby8sFgsF7dHR0ezceNGkpKSiIyMLN/ihcjLkg3rv4Vl70P6KWNb45uMqWjDogp9WoYlg8lbJzN5y2SyrdlU8K7AiJYjuK/JfZi9pMHiDlw33F1AamoqAQEBhISEkJyczLx58+jevXuxz2vdujUjR47klltuYcGCBdSsWbPsixUiL62NDwlj34QztllIa3c0Ziys06GIp2liD8UyPm48R84fAeDGejfydNunqVqhanlULhxEwr0ILVu2pHXr1jRu3JjatWvTuXNnu5/bpUsXJkyYwI033siiRYsICwsrw0qFyGP/CuNy+SMbjPuVo6DPG9DohiIXZd6Xso//rP0Pq48aK4g1rNiQlzq8RNtqbQt9jnBdMuWvi5JzIUoseTv89TrssY30CKwG3V+E1vcXOazvfM55/m/T/zF9+3Qs2kKQTxCjWo9iUMNBeJuk/edqHDrlr1KqP/AJ4AV8rbX+T77H6wBTgVDbPi9oreeWuGohRMmlHIYl78KmH4wRMD6B0Hk0dHq8yJWAtNbM2T+HD+M/5ETGCRSKO6Lu4N9t/k0lv0rl+AZEWSg23JVSXsDnQB8gCYhTSv2htd6eZ7dXgBla64lKqabAXCCiDOoVQlyQcRZWfQxrJoIlE0ze0O5f0PU5CKxS5FN3nt7Je2vfY8Nxo+umRVgLXurwEs3CmpVH5aIc2NNybw8kaK33ASilfgQGAnnDXQMXZqoPAY44skghPM2u07v47J/PSE5Pvopnazh/EtKSwZoLVUPBLwSCa0DuQVjyeDHP1uw+sxurtlLJrxJPtnmSgQ0GYlJyJbYnsSfcawGJee4nAfk/bh8LLFRKjQICgAInZFBKDQeGA9SpU6ektQrh9jIsGUzcNJFp26aRq3NLdzCzF0YvKKAzIMX+q029lBf3NbmPka1GEuzjWisICcdw1KclQ4ApWusPlFKdgO+UUs21vvwSOK31JGASGB+oOui1hXALqw6v4q01b3E47TAKxZDGQxhYfyAUPoDlkqR4WPslnNxt3K8YCR1GQJ2O9j0/nyr+VWRoo4ezJ9wPA7Xz3A+3bcvrYaA/gNZ6tVLKDwgDjjuiSCHc2cmMk4yLG8e8/cZKPg0rNuT1Tq/TokqL4p98bAsseh32xhr3g2pCj5eg1T1FXlkqhD3hHgdEKaUiMUJ9MHBPvn0OAb2AKUqpJoAfcMKRhQrhbqzaym97fuOD9R9wLvscfl5+jGw1kvub3o/ZVMxVnmcPweJ3YPNPgAbfYGMVoA4j3WIVIOF8xYa71tqilHoCWIDRwTdZa71NKfUmEK+1/gN4BvhKKfUUxoerw7SzBtC7EIvFgre3jBO+Fu07u483Vr9xcTRK55qdeaXjK4QHhRf9xPTTsOIDWPcV5GaByQztH4Hrx0BA5XKoXHgKu5LHNmZ9br5tr+W5vR2w//JNF3b+/HnuuusukpKSyM3N5dVXX6VevXqMHj2a8+fP4+vrS2xsLGazmZEjRxIfH4+3tzcffvghPXr0YMqUKcycOZO0tDRyc3NZtmwZ48ePZ8aMGWRlZXHbbbddnDZYeJ6s3Cy+3vI1X2/5GovVQiW/Sjzf7nkGRA5AFXF16MWl3VZ8AJm2eZWa3wk9X4FKMj+RKDmXbVZGT40uk+NuGbqlyMfnz59PzZo1mTNnDgApKSm0bt2an376iXbt2pGamoq/vz+ffPIJSim2bNnCzp076du3L7t3Gx92bdiwgc2bN1OpUiUWLlzInj17WLduHVprbrnlFpYvX07Xrl3L5P0J54k7Fsebq9/kQOoBAO6IuoOn2j5FiG8R01dbc2HzDFj8NqQmGdsiuxpzwNRsXfZFC4/lsuHuLNHR0TzzzDM8//zz3HTTTYSGhlKjRg3atWsHQHCwMWxs5cqVjBo1CoDGjRtTt27di+Hep08fKlUyrvBbuHAhCxcupHVr4xc1LS2NPXv2SLh7kLOZZ/lg/QfMSpgFQGRIJK93er3oOVm0hoRYY7qA5K3GtmrNofcb0KBXkXPACGEPlw334lrYZaVhw4Zs2LCBuXPn8sorr9CzZ88SH+PClMFgXOL94osv8uijjzqyTOEi5u1Zy6urnyJLnwPtTVBmP7LO9OHlA+eB5VfsH2Y9wXU5a+iWvYJmuTsAOK7CmOZ/H4szu2Od4wWsKHEd3RpW4cUbZC4icYnLhruzHDlyhEqVKnHfffcRGhrKF198wdGjR4mLi6Ndu3acO3cOf39/rr/+er7//nt69uzJ7t27OXToEI0aNWLDhg2XHa9fv368+uqr3HvvvQQGBnL48GHMZjNVq8oYY3eVnm1h9uaj/C9uB7vNb2Iyn8NyPpLMY7dzLrsKkGn7AtBEqcP0NcXTzyuOFqb9F4+ToivwuWUgU3P7kZXhA6RfdU0NqgaW5i0JDyThns+WLVt49tlnMZlMmM1mJk6ciNaaUaNGkZGRgb+/P3/99RePPfYYI0eOJDo6Gm9vb6ZMmVLgcnp9+/Zlx44ddOrUCYDAwECmT58u4e5mtNZsTDzLjPhE/tx0lLSsbPzDp+FdIYWKpije6fMFlQNsQxS1Ff/j/xC0fz5BBxfgm3Ip0K3e/qSFd+dcZH/O1e3NrT7B3OqA+oL95VdZXE6m/HVRci5cw+nz2czckMSM+ER2J6dd3B5Zbx0nfWcS5BPMrzf/Qg2/ynBgOeyYDbvmGvO+XOBfyZhLvfGNUL8HmGWxdHH1HDrlrxDXklyrZmXCSX6KO8Si7cnk5BoNoMoBPtzephYtG6Tw6to/QMM7tW+kxvxXjDnUs1IvHSSkjhHmTW4yVkAqYj51IcqC/MQJYZN4Op2f1yfxS3wiR1KMPnOTgu6NqjC4XW16Nq5GRup+7po/Bou2cH/qeXr89f6lA1RtdinQq7eQES/CqVwu3LXWRV/scQ2Qi3vLT2ZOLgu3JzMjLpFVe09y4dTXruTPXW1rc2dMODVyj8HOX9HT/mRMZgJHAioQnZnFU6dOG63yJjcZ3S6V6zv3zQiRh0uFu5+fH6dOnaJy5crXbMBrrTl16hR+fn7OLsWj7Tiayk9xiczaeJiz6TkA+HibGNC8One3DadjwBFMu36GH+ZcHIf+fXAQSypXJAgT41o8jvneIRAoH4wL1+RS4R4eHk5SUhInTlzbc475+fkRHl7MHCSixFIzc/hj4xFmxCeyOSnl4vZmNYMZHFOT2yolEnhgBsyZbUzcdYFPENvqd+aDrF2gc3mz+weE1y1wyQIhXIZLhbvZbCYyUubREI51+GwGXyxJ4NcNSWTmGEsMBPl5c2urWtzdrjbNTy+CeUMh/dSlJwVUhcY3QOObOVerNWPm3YclM5chjYfQW4JduAGXCnchHOnI2Qw+X5LAjPjEiyNeOtWrzN3tatO/eXX8zF6QuA5+GwHWHKhUDxrfZHyFtwOTCa01ry97hqS0JJpUasKYmDFOfldC2EfCXXicoykZfLFkLz/FJZKda0UpGNiqJqN6NqBB1aBLO547Bj/dbwR7+0dhwPtXjHCZsWsGiw4uIsAcwIRuE/Dx8inndyPE1ZFwFx7jWEomE5cm8L91l0L95pY1Gd0rX6gDWLJhxlBIOwZ1O0O/d64I9p2ndzIubhwAYzuNpU6wrPsr3IeEu3B7yamZTFy6lx/WHSLbYoT6jS1qMLpXFA2rBRX8pAUvQeIaY9m6QVPA6/KVkc7nnGfMsjFkW7MZ1HAQ/SP7l/0bEcKBJNyF2zqemsnEZXv5Ye0hsizGB6U3Rtfg372iaFS9kFAH+Od7iPsKvHzg7u+uGM6oteaN1W9wMPUgDSs25Ll2z5Xl2xCiTEi4C7dz/FwmXy7bx/Q1By+G+oDm1RndO4rG1YOLfvLhDTD7KeP2DRMg/MopOmbumcm8/fPw9/ZnQrcJ+HnLNQfC/Ui4C7dx4lwWXy7by/S1By8OaezXrBqjezWkac1iQh3g/EnjA9TcLGj7ILQdesUuu8/s5r117wHwasdXiQyRobnCPUm4C5d3Mi2LScv3MW31gYuh3rdpNUb3jqJZzSKWsMsr1wI/DzOWsgtvZ4yMySc9J50xy8aQlZvFbQ1u4+b6NzvuTQhRziTchUubvuYg78zZQUZOLgC9m1Tjyd5RNKpegR93/cjr8b+TlZtV/IHOn4ScM1A7HELN8OedV+ySnpPO8Yzj1A+pz4sdXnT0WxGiXEm4C5ekteajv/bw39g9APRqXJUnezekea1gFh9azPO/f8ihc4eKOUo+PrYRMWlJhe4S6hvKhG4T8PeWOdeFe5NwFy4n16p5/Y+tTF9zCJOC/9zegrva1WbbyW0Mmz+ODceNpQwjgiMY3WY09ULrFX6wE7vhl4chNxO6PgctBhX52tUrVKeCuYIj344QTiHhLlxKliWXp3/axJwtR/HxNvHZkNZE19W8sOIF5uybA0BF34qMbDWSOxveidlkLvxg6adh9nOQmQat7oUuz8kc6+KaIeEuXEZaloVHv4tnVcIpgny9+e+9Tdic9hMvzfqOrNwszCYz9zW9j0eiHyHIp4hx7ADWXPj1YTh7EGq0ghs/lGAX1xQJd+ESTqVl8eCUODYnpVA50JsH+h5j7D/vcjrzNAD9I/ozus1owoPsnAp58duwdzFUqAx3TwezjFUX1xYJd+F0SWfSeeCbdew7mUaN6gepVHsBX+/YD0DLKi15tt2ztKzS0v4Dbv8dVn4IysuYWiC0dtkULoQLk3AXTrU7+Rz3f7OWE1kHCWswnzTzTtLSoFZgLZ5q+xR96/Yt2apcx3fCrMeM233fgsiuZVO4EC5Owl04zfqDZ3jwu1iyAucSUDOeLKUJMgcxvMVw7mlyT8mn181MgR/vgew0iB4EHR8rm8KFcAMS7uKq7T27l6nbprI+eT2aki3qnZGdy4lzWVArFR9TDl7Km7sb3cWIliOo6Fex5MVYrTBzOJzeC9Wi4eb/ygeo4pom4S5KRGtNfHI8U7ZNYXnS8lIdS9ka5t3De/BMzNNEhERc/cGWj4Pd88Ev1Jjp0UfGqotrm4S7sEuuNZfYQ7FM2TaFLSe3AODr5cutDW7l9qjbCTQH2nWcn9cn8tniBADuaV+Hf/dsSrWAaqUrbtc8WPoeKBPcORkqyWRfQki4iyJlWDL4PeF3pm2fRuK5RMC4RH9I4yEMbjyYSn6V7DqO1poJC3fx+ZKzQBgv39CER7oWcWWpvU4mGN0xAD1fhQa9Sn9MITyAhLso0OnM0/y480d+3PkjZ7LOABAeGM4DzR7g1ga3lmjuFUuulVdmbeXHuES8TIpxd7TgjrZ2jlcvStY54wPUrFRocgt0ear0xxTCQ0i4i8scSj3EtO3T+D3hdzJzMwFoXrk5w5oPo3ed3niZvEp0vMycXEb/+A8LtiXjZzbxxb1t6Nm4lN0wAFobQx5P7oIqjeHWL+QDVCHysCvclVL9gU8AL+BrrfV/CtjnLmAsoIFNWut7HFinKGNbTmzh223fEnsoFqs25kzvGt6VYc2GEVMtpmRjzW3OZebwyLR41uw7TbCfN5OHtSMmwr5unGKt/Ah2/AG+wXD39+BbzHQEQlxjig13pZQX8DnQB0gC4pRSf2itt+fZJwp4EeistT6jlKpa8NGEK7FqKyuSVvDttm9Zn7weAG+TNzfXu5lhzYbRoGKDqz72iXNZDPt2HduOpFI1yJdpD7cvfgk8eyX8BbFvGrdv/wrCrr5OITyVPS339kCC1nofgFLqR2AgsD3PPo8An2utzwBorY87ulDhONm52czZN4cp26awL2UfAIHmQAY1GsS9je8t9eiVQ6fSuX/yWg6eSicyLIBpD7WndiUHDU08vd+YwhcN3V+ERv0dc1whPIw94V4LSMxzPwnokG+fhgBKqVUYXTdjtdbzHVKhsE/GGWMZOZM3RPWDhv2gYt0rdltzdA0vr3iZ4xnG399qFapxf9P7uSPqDgJ97BvOWJSNiWd5ZFo8J85l0bxWMFMebE9YoG+pj8vZQ7BtFsR9DZlnoeEAY352IUSBHPWBqjcQBXQHwoHlSqlorfXZvDsppYYDwwHq1KnjoJcWACyfAPuWGrcT/oJ5z0KVJkbIN+wP4e2YkTCTd9e+S67OJapiFA82e5D+Ef0xexUxJ3oJ/ByfyMuztpJtsdKpXmUmPdCWIL9SHPtsojEJ2Lbf4HD8pe1VmsDtX4LJVPqihfBQ9oT7YSDvtHrhtm15JQFrtdY5wH6l1G6MsI/Lu5PWehIwCSAmJqZk16uLwp3eD2u/BBT0HgtHNkDCYjixA07swLLqYz6oWp3pAcYloQ81vofR7Z/HpBwTjjm5Vt6Zs4Mpfx8A4P6OdXn1pqb4eF/F8VOSLgV6Up4fH3MF449Us1shqi+YZRk8IYpiT7jHAVFKqUiMUB8M5B8JMwsYAnyrlArD6KbZ58hCRRFi3wBrDrQcAl2eNLZZsuHQ35zbOZvnjixkpVnjrTWvnzzNrQfGw/YVl1r1YQ2vehjhqbQsnvjhH1bvO4XZS/HWwOYMbl/C/8pSj1wK9MS1l7Z7+xs1NrvNCHSZUkAIuxUb7lpri1LqCWABRn/6ZK31NqXUm0C81voP22N9lVLbgVzgWa31qbIsXNgkrjNC0dvPuELzAm8fEqvUZ9TmBPaaNaHmID6u1oO2bIL01XBwlfG16DUIrWuEfMN+ENEFvO3rI996OIVHv1vP4bMZVA3yZeJ9bWlb185Jv1KPGkMZt/0Gh1bnqdvPCPJmtxn1+ASU4GQIIS5QWjundyQmJkbHx8cXv6MonNbwTV9IWgfXj4Fel8J9ffJ6nlryFGeyzlA/pD6f9vqU2kG23rWMs8YqRbsXQMIiSM/zd9gcAPV7GMEa1ReCqhf40r9vPMzzv24mM8dKq9qhfHl/W6oFF7Pa0bljsP0P2D4LDv4NF2aS9PKFqD62QO8PvqX/YFcIT6WUWq+1jiluP7lC1Z1tn2UEe0CVS90xwO8JvzN29VgsVguda3VmfNfxl6856h8KzW83vqy5cHi9MaPi7gWQvBV2zja+AGq2vtSqr96SXBTj5u/ky+VGr9tdMeG8dWtzfL0LuXI17bity2WW8Z9CgYHeTy5CEsLBpOXurixZ8Hl7OHMAbvoIYh7Cqq18suETJm+dDMC9Te5lTMwYvE0l+Bt+NhH2LDSCfv8ysGRefMgaUI3ltOaHM01YQwvG3NyG+zvWvfLq1bQTl7pcDq4C2xWvePlAg96XWuh+DrqoSYhriL0tdwl3d/X3Z7DwZWNelRGrSLdm8+KKF1mcuBgv5cVLHV7irkZ3le41stNh/3LYs4CcHfMwnz968SGryQdT5PW2Vn1f8AmEHX8agX5gxaVAN5mNmRqb3QaNBoBfSOlqEuIaJ+HuydJPw39bGcvK3fMzx2q1YNTiUew8vZMgnyA+7P4hHWt0dNjLzdtylGd+3kjdnP0MCd3O4JDt+BzbAJetvqQu3TeZoX5PY9hioxuMbiAhhENIn7snWz7BCPbIbmypWJ1/zxnCyYyT1A2uy6c9PyUyxDGLVVitmg8X7eazJcbiGo1bX8ddtz+Kj9kLzp+EPYuMvvq9iyEnHer1MFrojW8A/6tYKk8I4TAS7u7m1F5YNwlQzIu+kVcXPERWbhbtq7fnw+4fEuLrmG6P1MwcnvxxI4t3Hsek4KUbmvBwl8hL/esBYdBqiPGVmwNWi1xYJIQLkXB3N7FvoK05TGzchYlbPgPgjqg7eLnjy5hNjplGIOF4GsOnxbPv5HlCK5j5bEgbukSFFf4EL7PxJYRwGRLu7uTQGjJ3/MGr1aoyP+sQJmViTMwY7mty31XNt16Qv7Yn8+RPG0nLstC4ehBfPRDjuBkdhRDlRsLdXWjNiQUvMLp6Vbb4+RJgDmBc13F0De/qkMNbrZrPliTw4aLdANzYogbj72xBBR/5ERHCHclvrpvYGfc5T5iSSfbxpWaF6nzW+wuiKkY55NhpWRbGzNjE/G3HUAqe69eYEd3qOey/ASFE+ZNwdwOL9y/ghe1fkuHtTasKNfn4ph+o7F/ZIcc+cPI8w7+LZ3dyGkF+3vx3SGt6NJKFtIRwdxLuLkxrzeStk/lkw8doBTfneDN24G/4OGh2xGW7TzDqhw2kZlpoUDWQrx6IITJMJuoSwhNIuLsoi9XCm6vf5LeE3wAYffosD984GeWgYJ8Rl8gLMzdj1dCnaTU+ursVgb7y4yCEp5DfZheUacnk2eXPsjRxKX6YeC85md7VO0JUb4ccf/XeU7z02xasGkb3imJ0ryhMJulfF8KTSLi7mNTsVEbFjmLD8Q0EmwP5/NB+WmVmQt+3r3pBjbwOnUpn5PfrsVg1j3atx1N9GjqgaiGEq5FFKF3IifQTPDj/QTYc30DVClWZqmrRKjMDWt8L1ZuX+vjnMnN4eGocZ9Nz6Nm4Ks/1b+yAqoUQrkjC3UUkpibywLwH2H1mNxHBEXzX4ika7FpkrB3a45VSHz/Xqhn940b2HE8jqmognwxuhZd0xQjhsSTcXcCOUzu4f96dMcuoAAAZxklEQVT9JKUl0axyM6b2+5aayz80Hrzu3xBco9SvMW7+ThbvPE7FCma+GdqOID+ZLkAITybh7mRxx+J4aMFDnMo8RccaHfmm3zdU2rvUWB0psBpcN6rUr/HL+iS+XL4Pb5Pii3vbUqeyTCcghKeTcHei2EOxjFg0grScNPpF9OPzXp8TgBf89YaxQ89XSr2e6PqDp3lp5hYA3hzYnE71HXPxkxDCtUm4O8nMPTN5eunTZFuzubvR3bx//fv4ePnAui8h5RBUbQqt7i3Vaxw+m8Gj360nO9fK0E51uadDHQdVL4RwdTIUspxprflm6zd8suETAB5r+RgjWo4w5nE5fwqWf2Ds2PctMBWy6LQdzmdZ+NfUeE6mZdOlQRiv3tTUEeULIdyEhHs5smorE+In8N3271AoXuzwIkMaD7m0w7L3ISsF6vcyFpK+2texap6ZsYkdR1OJDAvg83va4O0l/6QJcS2RcC8nOdYcXlv1GrP3zcbb5M17Xd6jf2T/SzucTID4b0CZjFZ7KXwcu4f5244R5OfNVw/EEFJBRsYIca2RcC8HGZYMnln6DCsOr8Df25+Pe3zMdTWvu3ynv143lqpr8wBUa3bVr/XnpiP8N3YPJgWf3dOGBlVL94GsEMI9SbiXsZSsFJ6IfYKNJzYS6hvKxN4TaR6W72rTA6tg52wwB0CPl6/6tTYnnWXMz5sAeOXGpnRrWKU0pQsh3JiEexlKPp/MiL9GkHA2geoB1fmyz5fUC6l3+U5WKyy0BXrn0RBU/epeKzWTR6bFk2WxMrhdbR7sHFG64oUQbk3CvYwcSDnAo4se5cj5I9QLqceXfb6kekABwb31VzjyDwRWh+ueuKrXyszJZfi0eJJTs2gfUYk3BzaXVZSEuMZJuJeBbSe3MfKvkZzJOkOLKi34vOfnhPqFXrljTgbE5rlgyafkC2VorXn+181sSkohvKI/E+9rg4+3jIwR4lonKeBga46u4aEFD3Em6wyda3Xmqz5fFRzsAGv/D1ISoVpzaHXPVb3eF0v38vvGIwT4ePH10BgqB/qWonohhKeQlrsDbTqxicf+eowcaw43RN7A213exmwqZBji+ZOwwjY52FVesLRw2zHGL9iFUvDx4NY0rh5ciuqFEJ5EWu4O9P2O78mx5nBrg1t57/r3Cg92gKX/gaxUaNAH6vcs8WvtOJrKkz9tBOC5fo3p07Ta1ZYthPBAEu4Okp6TztLEpQA82uJRTKqIU3tyD8RPNi5Y6vNmiV/rVFoW/5oaT3p2Lre1rsWIbvWKf5IQ4poi4e4gy5KWkWHJoEWVFoQHhRe986LXQOdC6/uhWsnmfMm2WBkxfT2Hz2bQsnYo790eLSNjhBBXkHB3kLn75wJwQ+QNRe+4fwXsmntVFyxprXll1hbiDpyherAfX93fFj/z1U8uJoTwXBLuDpCSlcLKwysxKRP9IvoVvqPVCgttS+Z1eRKCStZPPnnVAWbEJ+FnNvHVAzFUDfYrRdVCCE9mV7grpforpXYppRKUUi8Usd8dSimtlIpxXImuL/ZQLBarhXbV2xHmH1b4jlt+hqMbIagGdCrZBUtLdx3nnTnbAfhgUCuiw0NKU7IQwsMVG+5KKS/gc2AA0BQYopS6oqNYKRUEjAbWOrpIV2dXl0xOBsTaPjzt+Sr42L/UXcLxNEb98A9WDaN7RXFji9KvqSqE8Gz2tNzbAwla631a62zgR2BgAfu9BbwPZDqwPpd3MuMkccfi8DZ506tOr8J3XPMFpCZB9WhoOdju459Nz+ZfU+M4l2XhhujqjO4V5YCqhRCezp5wrwUk5rmfZNt2kVKqDVBbaz3HgbW5hQUHFmDVVrrU6kKIbyFdJeeSYcVHxu2+b9t9wVJyaiYPT43nwKl0mtUMZsKglphMMjJGCFG8Ul+hqpQyAR8Cw+zYdzgwHKBOHc9Yz9OuLpn5z0P2OWjYH+p1t+u4y3ef4KmfNnLqfDbVgn356oEYKvjIBcVCCPvY03I/DNTOcz/ctu2CIKA5sFQpdQDoCPxR0IeqWutJWusYrXVMlSruP9d40rkkNp/YjL+3P93CuxW80675sO03MFeAG8YXe0xLrpXxC3Yy9Nt1nDpvrH86e9T11Az1d3D1QghPZk9TMA6IUkpFYoT6YODiLFda6xTg4hARpdRSYIzWOt6xpbqe+QfmA9C9dncqmAv4gDTrHMx52rjd8xUILfq/laMpGfz7f/8Qd+AMJgVP92nIYz0a4CVdMUKIEio23LXWFqXUE8ACwAuYrLXeppR6E4jXWv9R1kW6qmK7ZGLfgtTDULM1dBhR5LGW7DzO0zM2ciY9h2rBvnwyuDUd61V2dMlCiGuEXZ24Wuu5wNx8214rZN/upS/L9SWcSWDPmT0E+wTTuWbnK3dIiod1k0B5wS2fFvohak6ulQkLdvHl8n0AdGtYhQ/vailT9wohSkU+obtKF1rtfer2weyVb/bH3Bz449+AhutGGcMfC3D4bAajftjAhkNn8TIpxvRtxKNd68mIGCFEqUm4XwWtNfP2zwNgQOSAK3dY9Qkc3wYVI6Db8wUeY9H2ZMb8vImUjBxqhPjx6ZDWxERUKsOqhRDXEgn3q7D15FaS0pKo4l+FmGr5BgWd2gvLxhm3b/r4iitRsy1W3p+/k29W7gegZ+OqfDCoJRUDfMqjdCHENULC/Spc6JLpF9EPr7x96VrDn6MhNwta3gP1e1z2vMTT6Tzxv3/YlHgWb5Pi+f6NebhLpHTDCCEcTsK9hHKtuSw4sAAooEvmn+lwYAVUqAz93rnsoflbj/LsL5s5l2mhVqg/n97TmjZ1KpZX2UKIa4yEewmtT17PiYwThAeGEx2W54PStOOXpvPt/x+oYPSfZ1lyeXfODqauPghAn6bVmHBnS0IqFLEEnxBClJKEewld6JIZEDng8hWQ5r8AmWehfi+IHgTAgZPneeJ/G9h6OBWzl+LFAU14sHOErJwkhChzEu4lkJObw6KDi4B8XTK7F8LWX40pBm76CJRi9uYjvPDrFtKyLNSu5M9nQ9rQsnaokyoXQlxrJNxL4O8jf5OanUqD0AZEVbRNvZuVdmmKgR4vkxkYzlu/beH7tYcAGNC8Ov+5owUh/tINI4QoPxLuJVDgdAOL34aURKjRin317+PxL/5mx9FUfLxMvHJTE+7vWFe6YYQQ5U7C3U4ZlgyWJC4BoH9kf2Nj0npY+3+gvFjS8BUe/3wN6dm5RFSuwGf3tKF5LVkKTwjhHBLudlqWuIwMSwYtwlpQO6i2McXAn8YUA8sq38WDC7IBuLllTd69rTlBftINI4RwHgl3O+UdJQPA359C8laOmqrzaFJffLxNjL25GUPa15ZuGCGE00m42yE1O5WVh1diUib6RfSzTTHwPgDPZj5IaHAI3z7YjiY1gp1cqRBCGOxZiemaF3swlhxrDu2qtaOKfxjMfhIsmSz07s5KazRP920owS6EcCkS7na4rEtm4w+wfzlZPhV5Pm0wdSpV4LbWtYo5ghBClC/plinGyYyTrDu2Dm+TN73DWsNvPQH4wDSMMwTzYo8GmL3kb6QQwrVIKhVjwYEFWLWVLjW7ELL4Hcg4w/GqnZl0Nobalfy5rY202oUQrkfCvRgXF+XwD4etv6C9/XkmfSigeLy7tNqFEK5JkqkIh9MOs+nEJvy9/Oi+djoA2xs/wYqTgdQK9ef2NuFOrlAIIQom4V6EC6327t6VqJByCF2jJU8fvA6Ax3s0wMdbTp8QwjVJOhXhYpfMoc2gTKxo/Cq7TmRQK9SfO9tKq10I4bok3Aux9+xedp/ZTZCGzunp6A6P8fYGY53Tx3rUl1a7EMKlSUIV4sLY9j5pafiE1mFhtYfYnZxGzRA/BrWt7eTqhBCiaBLuBdBaMy/hDwAGpJ3HeuPHfLTsMAAjpa9dCOEGJKUKsO3kVhLTjxFmyaVd1EAWZDZl57Fz1Ajx464Y6WsXQrg+CfcCzF33MQD9snJRfd/lk9g9AIzsXh9fby9nliaEEHaRcM8nNy2ZBclrABjQ8hEWHsxl57FzVA/2464Y6WsXQrgHCfd8NswZxXEvE7XwJrrjU/w3T6vdzyytdiGEe5BwzyvhL+YeXwfAgAa38tfOE2w/mkq1YF/ubietdiGE+5BwvyD7PDmzn2JRQAUA+jcZfLGvfUQ3abULIdyLhPsFS95lddZxUry8aBBSn8RjoWw7kkqVIF+GtK/j7OqEEKJEJNwBjvwDa75gbmAAAP0jB/Bx7G5AWu1CCPck4Z5rgT/+TQaaxUHGUnmh1vZsPWy02u/tIK12IYT7kXBf8wUc28yysNpk6FyahzXn+1XnAXi0az1ptQsh3NK1He6n98OSdwGYVycagKgKXdmclEJYoC/3dqjrzOqEEOKq2RXuSqn+SqldSqkEpdQLBTz+tFJqu1Jqs1IqVinl+qmoNcx+CiwZpDa/jRVnd6JQrN9hdMM82rUe/j7SahdCuKdiw10p5QV8DgwAmgJDlFJN8+32DxCjtW4B/AKMc3ShDrf5J9i3BPwrEtu4JznWHBoEt2TbIagc4MO9HaWvXQjhvuxpubcHErTW+7TW2cCPwMC8O2itl2it02131wCuPbvW+VMw/0Xjdt93mHtkBQBnjzcDYHjXelTw8XZWdUIIUWr2hHstIDHP/STbtsI8DMwrTVFlbsFLkHEaIrsRW6k6646tw0t5s+9gfSoF+HB/J9fvVRJCiKI4tHmqlLoPiAG6FfL4cGA4QJ06Tur2yM2B7DQyzH6Mr9OQn5c+CUBIThfOWitIq10I4RHsabkfBvJOrBJu23YZpVRv4GXgFq11VkEH0lpP0lrHaK1jqlSpcjX1lp6XmV19X2Nw4zb8fGgBZpOZOyMeIymhv9Fq7yitdiGE+7Mn3OOAKKVUpFLKBxgM/JF3B6VUa+BLjGA/7vgyHUNrzfTt0xkyZwj70pKIDInk+xu+Z8v2loCJf10fSYCvtNqFEO6v2CTTWluUUk8ACwAvYLLWeptS6k0gXmv9BzAeCAR+VkoBHNJa31KGdZfYyYyTvLrqVVYeXgnAoIaDeLbds2w4cJ74g/sIrWDmgU4Rzi1SCCEcxK5mqtZ6LjA337bX8tzu7eC6HGrl4ZW8vPJlTmeeJsQ3hDc6vUGvur3QWvPxX/8A8Mj19QiUVrsQwkN4dJpl52bz0fqPmL5jOgDtqrfj3S7vUj2gOgCr954i7sAZW6td+tqFEJ7DY8N939l9PLf8OXad2YW38ubx1o/zYLMH8TJduur0Y9t87f/qEkmQn9lZpQohhMN5XLhrrfl598+MjxtPZm4mtYNq8/717xNdJfqy/VbvPcW6/acJ8Tcz9LoI5xQrhBBlxOPC/astX/HpP58CcEv9W3ipw0sEmAOu2O8T23ztD0urXQjhgTwq3LXW/LL7FwDGdhrLHQ3vKHC/tCwLVisE+3kzrHNEOVYohBDlw6PCffeZ3Rw9f5TKfpW5Leq2QvcL9PXmp0c7cjQlk2BptQshPJBHzee+PGk5AF3Du2JSRb81pRQ1Q/3LoywhhCh3HhXuy5KWAdCtdoFT2wghxDXDY8L9VMYpNp/YjNlkplONTs4uRwghnMpjwn3l4ZVoNO2rt6eCuYKzyxFCCKfymHCXLhkhhLjEI8I9JzeHv4/8DRgfpgohxLXOI8I9Pjme8znnaRDagFqBRS0SJYQQ1waPCPcLQyC7hUuXjBBCgAeEu9aapYlLAeheu7tTaxFCCFfh9uG+P3U/SWlJhPqGEh0WXfwThBDiGuD24b4s0Rglc32t6y+bzlcIIa5l7h/uMgRSCCGu4NbhnpKVwsbjG/FW3lxX8zpnlyOEEC7DrcN95eGV5Opc2lZrS5BPkLPLEUIIl+HW4S5dMkIIUTC3DXeL1cLKwysBGd8uhBD5uW24bzy+kXPZ54gIjqBOcB1nlyOEEC7FbcP9QpeMXLgkhBBXcvtwl4nChBDiSm4Z7odSD7E/ZT9BPkG0qtrK2eUIIYTLcctwvzBRWJeaXTCbZIFrIYTIzy3DfWnSUkCGQAohRGHcLtzTstNYf2w9JmWiS60uzi5HCCFcktuF+99H/saiLbSq0ooQ3xBnlyOEEC7J7cJdhkAKIUTx3C7cLVYLZpNZrkoVQogiKK21U144JiZGx8fHX9Vz03PS8ff2Rynl4KqEEMK1KaXWa61jitvPuzyKcbQK5grOLkEIIVya23XLCCGEKJ6EuxBCeCC7wl0p1V8ptUsplaCUeqGAx32VUj/ZHl+rlIpwdKFCCCHsV2y4K6W8gM+BAUBTYIhSqmm+3R4GzmitGwAfAe87ulAhhBD2s6fl3h5I0Frv01pnAz8CA/PtMxCYarv9C9BLyVAWIYRwGnvCvRaQmOd+km1bgftorS1AClDZEQUKIYQouXL9QFUpNVwpFa+Uij9x4kR5vrQQQlxT7Bnnfhioned+uG1bQfskKaW8gRDgVP4Daa0nAZMAlFInlFIHr6LmMODkVTzPWdytXnC/mt2tXpCay4O71Qv21VzXngPZE+5xQJRSKhIjxAcD9+Tb5w9gKLAauBNYrIu59FVrXcWeAvNTSsXbc3WWq3C3esH9ana3ekFqLg/uVi84tuZiw11rbVFKPQEsALyAyVrrbUqpN4F4rfUfwDfAd0qpBOA0xh8AIYQQTmLX9ANa67nA3HzbXstzOxMY5NjShBBCXC13vEJ1krMLKCF3qxfcr2Z3qxek5vLgbvWCA2t22qyQQgghyo47ttyFEEIUw23Cvbj5bVyFUuqAUmqLUmqjUiretq2SUmqRUmqP7XtFJ9Y3WSl1XCm1Nc+2AutThv/azvlmpVQbF6p5rFLqsO08b1RK3ZDnsRdtNe9SSvVzQr21lVJLlFLblVLblFKjbdtd9jwXUbMrn2c/pdQ6pdQmW81v2LZH2ua4SrDNeeVj2+7UObCKqHeKUmp/nnPcyra9dD8XWmuX/8IYpbMXqAf4AJuAps6uq5BaDwBh+baNA16w3X4BeN+J9XUF2gBbi6sPuAGYByigI7DWhWoeC4wpYN+mtp8PXyDS9nPjVc711gDa2G4HAbttdbnseS6iZlc+zwoItN02A2tt528GMNi2/f+AkbbbjwH/Z7s9GPjJReqdAtxZwP6l+rlwl5a7PfPbuLK8c+9MBW51ViFa6+UYw1XzKqy+gcA0bVgDhCqlapRPpZcUUnNhBgI/aq2ztNb7gQSMn59yo7U+qrXeYLt9DtiBMUWHy57nImoujCucZ621TrPdNdu+NNATY44ruPI8O20OrCLqLUypfi7cJdztmd/GVWhgoVJqvVJquG1bNa31UdvtY0A155RWqMLqc/Xz/oTt39XJebq6XKpm27/+rTFaaW5xnvPVDC58npVSXkqpjcBxYBHGfxBntTHHVf66nD4HVv56tdYXzvE7tnP8kVLKN3+9NiU6x+4S7u6ki9a6DcYUyY8rpbrmfVAb/2+57BAlV68vj4lAfaAVcBT4wLnlXEkpFQj8CjyptU7N+5irnucCanbp86y1ztVat8KYFqU90NjJJRUpf71KqebAixh1twMqAc874rXcJdztmd/GJWitD9u+Hwd+w/iBS77w75Tt+3HnVVigwupz2fOutU62/aJYga+41CXgEjUrpcwYIfm91nqmbbNLn+eCanb183yB1vossATohNF9ceECzbx1XaxZFTEHVnnIU29/W5eY1lpnAd/ioHPsLuF+cX4b2yffgzHms3EpSqkApVTQhdtAX2Arl+bewfb9d+dUWKjC6vsDeMD2qX1HICVPt4JT5et7vA3jPINR82DbyIhIIApYV861KYwpOXZorT/M85DLnufCanbx81xFKRVqu+0P9MH4rGAJxhxXcOV5vnD+7ZoDqxzq3ZnnD77C+Hwg7zm++p+L8vy0uDRfGJ8c78boU3vZ2fUUUmM9jBEEm4BtF+rE6NeLBfYAfwGVnFjj/zD+vc7B6MN7uLD6MD6l/9x2zrcAMS5U83e2mjbbfglq5Nn/ZVvNu4ABTqi3C0aXy2Zgo+3rBlc+z0XU7MrnuQXwj622rcBrtu31MP7QJAA/A7627X62+wm2x+u5SL2Lbed4KzCdSyNqSvVzIVeoCiGEB3KXbhkhhBAlIOEuhBAeSMJdCCE8kIS7EEJ4IAl3IYTwQBLuwqMopSJUntkj82z/WinVtITHSit+LyFck13L7Anh7rTW/yrL49suQFHauJJTCKeTlrvwRN5Kqe+VUjuUUr8opSoopZYqpWLAaJErpd6xzau9RilVzbY9Uim1Whnz8b+d94BKqWeVUnG2yZ0uzMMdoYy5zKdhXIBS2zY391bbMZ4q7zcuxAUS7sITNQK+0Fo3AVIx5vHOKwBYo7VuCSwHHrFt/wSYqLWOxrgiFgClVF+My+vbY0yg1TbPhHBRttdqBoQBtbTWzW3H+LZM3p0QdpBwF54oUWu9ynZ7Osal9XllA7Ntt9cDEbbbnTGmOgDjsvsL+tq+/gE2YMzgF2V77KA25toG2AfUU0p9qpTqj/GHRQinkD534Ynyz6mR/36OvjTvRi6X/x4UNB+HAt7TWn952UZj3vPzF5+o9RmlVEugHzACuAt4qKTFC+EI0nIXnqiOUqqT7fY9wEo7n7cKY8ZRgHvzbF8APGSb6xylVC2lVNX8T1ZKhQEmrfWvwCsYSwMK4RQS7sIT7cJYKGUHUBFjwQl7jLY9bwt5VrzRWi8EfgBW2x77BWOd0fxqAUttK+1Mx1iEQQinkFkhhRDCA0nLXQghPJCEuxBCeCAJdyGE8EAS7kII4YEk3IUQwgNJuAshhAeScBdCCA8k4S6EEB7o/wGQvEYuOE1AgQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "reload(analysis)\n",
    "resp = pd.read_csv('cockroach_responses.csv')\n",
    "resp['start'] = resp.pos\n",
    "resp['end'] = resp.sequence.str.len()+resp.pos\n",
    "#print (resp)\n",
    "def overlap(x):\n",
    "    return x.seq\n",
    "res=[]\n",
    "methods = {'rank': range(2,40,3), 'score':range(50,2000,100), 'global':np.arange(.9,.99,.01)}\n",
    "f,ax=plt.subplots(1,1)\n",
    "#for n in range(1,10): \n",
    "for m in methods:\n",
    "    rng = methods[m]\n",
    "    for n in rng:\n",
    "        pb=P.promiscuous_binders(n=3, cutoff_method=m, cutoff=n)\n",
    "        #cl = analysis.find_clusters(pb, min_binders=2)\n",
    "        x = analysis.get_overlaps(resp,pb,how='inside')\n",
    "        found=x[x.overlap>0]\n",
    "        res.append((len(pb),len(found)/len(resp),m))\n",
    "res=pd.DataFrame(res,columns=['binders','positives','method'])\n",
    "\n",
    "for i,df in res.groupby('method'):\n",
    "    df.plot(0,1,kind='line',lw=2,ax=ax, label=i)\n",
    "plt.savefig('promiscuity_tpr.png',dpi=150)\n",
    "#res"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using the command line interface\n",
    "\n",
    "You can re-produce most of this analysis, without using Python, from the command line. The configuration settings are given below. Put this in a file file.conf. You can also place the list of alleles in a text file for convenience. Then run:\n",
    "\n",
    "`epitopepredict -c file.conf -r`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "file.conf:\n",
    "```\n",
    "[base]\n",
    "cpus = 4\n",
    "cutoff = .95\n",
    "cutoff_method = default\n",
    "fasta_header_sep =  \n",
    "genome_analysis = no\n",
    "mhc2_alleles = mhc2alleles.txt\n",
    "mhc2_length = 11\n",
    "n = 3\n",
    "names = \n",
    "overwrite = yes\n",
    "path = results_cockroach\n",
    "predictors = netmhciipan\n",
    "sequence_file = cockroach_allergens.fa\n",
    "verbose = yes\n",
    "```\n",
    "\n",
    "mhc2alleles.txt:\n",
    "```\n",
    "HLA-DPA10103-DPB10201\n",
    "HLA-DPA10201-DPB10101\n",
    "HLA-DPA10201-DPB10501\n",
    "HLA-DPA10301-DPB10402\n",
    "HLA-DQA0101-DQB10501\n",
    "HLA-DQA10301-DQB10302\n",
    "HLA-DQA10401-DQB10402\n",
    "HLA-DQA10501-DQB10301\n",
    "HLA-DRB1*0101\n",
    "HLA-DRB1*0301\n",
    "HLA-DRB1*0401\n",
    "HLA-DRB1*0405\n",
    "HLA-DRB1*0701\n",
    "HLA-DRB1*0802\n",
    "HLA-DRB1*1101\n",
    "HLA-DRB1*1302\n",
    "HLA-DRB1*1501\n",
    "HLA-DRB3*0101\n",
    "HLA-DRB4*0101\n",
    "HLA-DRB5*0101\n",
    "\n",
    "\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}