{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# ipyrad-analysis toolkit: BPP\n",
"\n",
"The program *bpp* by Rannala & Yang (2010; 2015) is a powerful tool for inferring species tree parameters and testing species delimitation hypotheses. It is *relatively* easy to use, and best of all, it's *quite fast*, although not highly parallelizable. This notebook describes a streamlined approach to easily setup input files for testing different hypthotheses in *bpp*, and to do so in a clear programmatic way that makes it easy to perform many tests over many different parameter settings. We also show how to distribute many separate jobs to run in parallel on a cluster. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Required software\n",
"All software required for this notebook can be installed using conda. "
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"## conda install -c ipyrad ipyrad\n",
"## conda install -c ipyrad bpp\n",
"## conda install -c eaton-lab toytree"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"import ipyrad.analysis as ipa ## ipyrad analysis tools\n",
"import pandas as pd ## DataFrames\n",
"import numpy as np ## data generation\n",
"import toytree ## tree plotting\n",
"import toyplot ## data plotting"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Analysis setup\n",
"\n",
"You must define a tree with the \"species\" names in your analysis which will act either as a fixed-tree or as a guide-tree for your analysis. You must also define an `IMAP` dictionary which maps sample names to \"species\" names. And you can optionally also define a `MINMAP` dictionary which is used to filter RAD loci to include only loci that have at least N samples with data from each species in a locus. "
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"## set the location of our input .loci file\n",
"locifile = \"./analysis-ipyrad/pedicularis_outfiles/pedicularis.alleles.loci\"\n",
"\n",
"## a tree hypothesis (guidetree) (here based on tetrad results)\n",
"newick = \"((((((rex, lip),rck),tha),cup),(cys,(cya, sup))),prz);\"\n",
"\n",
"## a dictionary mapping sample names to 'species' names\n",
"imap = {\n",
" \"prz\": [\"32082_przewalskii\", \"33588_przewalskii\"],\n",
" \"cys\": [\"41478_cyathophylloides\", \"41954_cyathophylloides\"],\n",
" \"cya\": [\"30686_cyathophylla\"],\n",
" \"sup\": [\"29154_superba\"],\n",
" \"cup\": [\"33413_thamno\"],\n",
" \"tha\": [\"30556_thamno\"],\n",
" \"rck\": [\"35236_rex\"],\n",
" \"rex\": [\"35855_rex\", \"40578_rex\"],\n",
" \"lip\": [\"39618_rex\", \"38362_rex\"], \n",
" }\n",
"\n",
"## optional: loci will be filtered if they do not have data for at\n",
"## least N samples/individuals in each species.\n",
"minmap = {\n",
" \"prz\": 2,\n",
" \"cys\": 2,\n",
" \"cya\": 1,\n",
" \"sup\": 1,\n",
" \"cup\": 1,\n",
" \"tha\": 1, \n",
" \"rck\": 1,\n",
" \"rex\": 2,\n",
" \"lip\": 2,\n",
" }"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"## check your (starting) tree hypothesis\n",
"toytree.tree(newick).draw(use_edge_lengths=False, width=250);\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## The *bpp* Class object\n",
"\n",
"To simplify the creation of input files for *bpp* analyses we've created a bpp job generator object that can be accessed from `ipa.bpp()`. Running *bpp* requires three input files (.ctl, .imap, and .seq) of which the .ctl file is the most important since it contains the parameters for a run and points to the location of the other two files. The `ipa.bpp()` object can be used to easily modify parameter settings for a run, to generate the input files, and if desired, to submit the bpp jobs to run on a cluster (your ipyclient cluster). "
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"## create a bpp object to run algorithm 00\n",
"b = ipa.bpp(\n",
" name=\"test\",\n",
" data=locifile,\n",
" guidetree=newick, \n",
" imap=imap,\n",
" minmap=minmap, \n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"binary bpp \n",
"burnin 2000 \n",
"cleandata 0 \n",
"delimit_alg (0, 5) \n",
"finetune (0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01, 0.01)\n",
"infer_delimit 0 \n",
"infer_sptree 0 \n",
"nsample 20000 \n",
"sampfreq 2 \n",
"seed 33333 \n",
"tauprior (2, 200, 1) \n",
"thetaprior (2, 2000) \n",
"usedata 1 "
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## set some optional params, leaving others at their defaults\n",
"b.params.burnin = 2000\n",
"b.params.nsample = 20000\n",
"b.params.sampfreq = 2\n",
"b.params.seed = 33333\n",
"b.params.tauprior = (2, 200, 1)\n",
"\n",
"## print params\n",
"b.params"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"maxloci 100 \n",
"minmap {'cys': 4, 'rex': 4, 'cup': 2, 'rck': 2, 'cya': 2, 'lip': 4, 'sup': 2, 'tha': 2, 'prz': 4}\n",
"minsnps 0 "
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## set some optional filters leaving others at their defaults\n",
"b.filters.maxloci=100\n",
"b.filters.minsnps=0\n",
"\n",
"## print filters\n",
"b.filters"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Generating files &/or submitting jobs\n",
"When you create a *bpp* object you save it with a variable name (in this example named `test`) which is simply the name you will use to reference the object. You must also provide a *name* (prefix for output files) that will be used by either of the two main functions of the bpp object: **write_bpp_files()** or **run()**. Both functions make it easy to sample different distributions of loci to include in different replicate bpp analyses. Each replicate will start from a different random seed after the initial `seed`. If you used a `maxloci` argument to limit the number of loci that will used in the analysis then you can also use the `randomize_order` argument to select a different random number of N loci in each rep, otherwise just the first N loci that pass filtering will be sampled. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### `write_bpp_files()`\n",
"This writes the .ctl, .seq, and .imap files for the specified run. "
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"input files created for job test (100 loci)\n"
]
}
],
"source": [
"## write files (filtered based on test.filters and test.params)\n",
"b.write_bpp_files()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### `run()`\n",
"This writes the files for each job and submits the bpp jobs to run on the cluster designated by the *ipyclient* object. You can efficiently submit many replicate jobs in this way. "
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"submitted 2 bpp jobs [test] (100 loci)\n"
]
}
],
"source": [
"## or, write files and submit job to run on ipyclient \n",
"b.run(\n",
" nreps=2, \n",
" ipyclient=ipyclient, \n",
" randomize_order=True,\n",
" )"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Accessing job results\n",
"When you submit jobs the results files will be stored in the bpp objects `.files` attribute. In addition, the asychronous result objects (a representation of the running job) from each submitted job will be accessible from the `.asyncs` attribute of the bpp object. You can view these objects to see if your job has finished or use them to trace errors if an error arises, as we do below. "
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"job finished\n",
"job finished\n"
]
}
],
"source": [
"## check async objects from the bpp object\n",
"for job in b.asyncs:\n",
" if job.ready():\n",
" print 'job finished'\n",
" else:\n",
" print 'job running'"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"## uncomment and run this to block until all running jobs are finished\n",
"#ipyclient.wait()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Results"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"## set options to print prettier tables\n",
"## (this suppresses scientific notation)\n",
"pd.set_option('display.float_format', lambda x: '%.4f' % x)"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
count
\n",
"
mean
\n",
"
std
\n",
"
min
\n",
"
25%
\n",
"
50%
\n",
"
75%
\n",
"
max
\n",
"
\n",
" \n",
" \n",
"
\n",
"
theta_1cup
\n",
"
40000.0000
\n",
"
0.0016
\n",
"
0.0005
\n",
"
0.0004
\n",
"
0.0012
\n",
"
0.0016
\n",
"
0.0019
\n",
"
0.0046
\n",
"
\n",
"
\n",
"
theta_2cya
\n",
"
40000.0000
\n",
"
0.0005
\n",
"
0.0003
\n",
"
0.0000
\n",
"
0.0002
\n",
"
0.0004
\n",
"
0.0007
\n",
"
0.0027
\n",
"
\n",
"
\n",
"
theta_3cys
\n",
"
40000.0000
\n",
"
0.0004
\n",
"
0.0003
\n",
"
0.0000
\n",
"
0.0001
\n",
"
0.0004
\n",
"
0.0007
\n",
"
0.0016
\n",
"
\n",
"
\n",
"
theta_4lip
\n",
"
40000.0000
\n",
"
0.0010
\n",
"
0.0003
\n",
"
0.0004
\n",
"
0.0008
\n",
"
0.0010
\n",
"
0.0012
\n",
"
0.0023
\n",
"
\n",
"
\n",
"
theta_5prz
\n",
"
40000.0000
\n",
"
0.0039
\n",
"
0.0007
\n",
"
0.0019
\n",
"
0.0034
\n",
"
0.0038
\n",
"
0.0043
\n",
"
0.0070
\n",
"
\n",
"
\n",
"
theta_6rck
\n",
"
40000.0000
\n",
"
0.0014
\n",
"
0.0005
\n",
"
0.0004
\n",
"
0.0011
\n",
"
0.0014
\n",
"
0.0017
\n",
"
0.0042
\n",
"
\n",
"
\n",
"
theta_7rex
\n",
"
40000.0000
\n",
"
0.0028
\n",
"
0.0008
\n",
"
0.0009
\n",
"
0.0022
\n",
"
0.0028
\n",
"
0.0033
\n",
"
0.0063
\n",
"
\n",
"
\n",
"
theta_8sup
\n",
"
40000.0000
\n",
"
0.0011
\n",
"
0.0004
\n",
"
0.0002
\n",
"
0.0008
\n",
"
0.0010
\n",
"
0.0013
\n",
"
0.0039
\n",
"
\n",
"
\n",
"
theta_9tha
\n",
"
40000.0000
\n",
"
0.0015
\n",
"
0.0005
\n",
"
0.0005
\n",
"
0.0012
\n",
"
0.0014
\n",
"
0.0017
\n",
"
0.0052
\n",
"
\n",
"
\n",
"
theta_10rexliprckthacupcyscyasupprz
\n",
"
40000.0000
\n",
"
0.0039
\n",
"
0.0032
\n",
"
0.0002
\n",
"
0.0011
\n",
"
0.0024
\n",
"
0.0068
\n",
"
0.0144
\n",
"
\n",
"
\n",
"
theta_11rexliprckthacupcyscyasup
\n",
"
40000.0000
\n",
"
0.0029
\n",
"
0.0011
\n",
"
0.0004
\n",
"
0.0021
\n",
"
0.0030
\n",
"
0.0038
\n",
"
0.0084
\n",
"
\n",
"
\n",
"
theta_12rexliprckthacup
\n",
"
40000.0000
\n",
"
0.0044
\n",
"
0.0013
\n",
"
0.0009
\n",
"
0.0035
\n",
"
0.0044
\n",
"
0.0052
\n",
"
0.0116
\n",
"
\n",
"
\n",
"
theta_13rexliprcktha
\n",
"
40000.0000
\n",
"
0.0026
\n",
"
0.0017
\n",
"
0.0001
\n",
"
0.0013
\n",
"
0.0022
\n",
"
0.0037
\n",
"
0.0104
\n",
"
\n",
"
\n",
"
theta_14rexliprck
\n",
"
40000.0000
\n",
"
0.0017
\n",
"
0.0010
\n",
"
0.0000
\n",
"
0.0010
\n",
"
0.0015
\n",
"
0.0022
\n",
"
0.0088
\n",
"
\n",
"
\n",
"
theta_15rexlip
\n",
"
40000.0000
\n",
"
0.0018
\n",
"
0.0010
\n",
"
0.0001
\n",
"
0.0010
\n",
"
0.0016
\n",
"
0.0023
\n",
"
0.0079
\n",
"
\n",
"
\n",
"
theta_16cyscyasup
\n",
"
40000.0000
\n",
"
0.0023
\n",
"
0.0012
\n",
"
0.0001
\n",
"
0.0014
\n",
"
0.0021
\n",
"
0.0030
\n",
"
0.0094
\n",
"
\n",
"
\n",
"
theta_17cyasup
\n",
"
40000.0000
\n",
"
0.0026
\n",
"
0.0011
\n",
"
0.0003
\n",
"
0.0018
\n",
"
0.0025
\n",
"
0.0033
\n",
"
0.0085
\n",
"
\n",
"
\n",
"
tau_10rexliprckthacupcyscyasupprz
\n",
"
40000.0000
\n",
"
0.0111
\n",
"
0.0034
\n",
"
0.0043
\n",
"
0.0075
\n",
"
0.0122
\n",
"
0.0141
\n",
"
0.0191
\n",
"
\n",
"
\n",
"
tau_11rexliprckthacupcyscyasup
\n",
"
40000.0000
\n",
"
0.0041
\n",
"
0.0005
\n",
"
0.0020
\n",
"
0.0037
\n",
"
0.0040
\n",
"
0.0044
\n",
"
0.0064
\n",
"
\n",
"
\n",
"
tau_12rexliprckthacup
\n",
"
40000.0000
\n",
"
0.0017
\n",
"
0.0004
\n",
"
0.0007
\n",
"
0.0014
\n",
"
0.0018
\n",
"
0.0021
\n",
"
0.0033
\n",
"
\n",
"
\n",
"
tau_13rexliprcktha
\n",
"
40000.0000
\n",
"
0.0014
\n",
"
0.0005
\n",
"
0.0005
\n",
"
0.0009
\n",
"
0.0013
\n",
"
0.0018
\n",
"
0.0027
\n",
"
\n",
"
\n",
"
tau_14rexliprck
\n",
"
40000.0000
\n",
"
0.0013
\n",
"
0.0005
\n",
"
0.0005
\n",
"
0.0009
\n",
"
0.0012
\n",
"
0.0017
\n",
"
0.0025
\n",
"
\n",
"
\n",
"
tau_15rexlip
\n",
"
40000.0000
\n",
"
0.0012
\n",
"
0.0005
\n",
"
0.0004
\n",
"
0.0008
\n",
"
0.0011
\n",
"
0.0016
\n",
"
0.0025
\n",
"
\n",
"
\n",
"
tau_16cyscyasup
\n",
"
40000.0000
\n",
"
0.0032
\n",
"
0.0006
\n",
"
0.0010
\n",
"
0.0028
\n",
"
0.0032
\n",
"
0.0036
\n",
"
0.0053
\n",
"
\n",
"
\n",
"
tau_17cyasup
\n",
"
40000.0000
\n",
"
0.0015
\n",
"
0.0008
\n",
"
0.0001
\n",
"
0.0008
\n",
"
0.0014
\n",
"
0.0021
\n",
"
0.0041
\n",
"
\n",
"
\n",
"
lnL
\n",
"
40000.0000
\n",
"
-12151.9644
\n",
"
271.9991
\n",
"
-12524.5230
\n",
"
-12422.8427
\n",
"
-12161.9050
\n",
"
-11879.7943
\n",
"
-11807.7240
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" count mean std \\\n",
"theta_1cup 40000.0000 0.0016 0.0005 \n",
"theta_2cya 40000.0000 0.0005 0.0003 \n",
"theta_3cys 40000.0000 0.0004 0.0003 \n",
"theta_4lip 40000.0000 0.0010 0.0003 \n",
"theta_5prz 40000.0000 0.0039 0.0007 \n",
"theta_6rck 40000.0000 0.0014 0.0005 \n",
"theta_7rex 40000.0000 0.0028 0.0008 \n",
"theta_8sup 40000.0000 0.0011 0.0004 \n",
"theta_9tha 40000.0000 0.0015 0.0005 \n",
"theta_10rexliprckthacupcyscyasupprz 40000.0000 0.0039 0.0032 \n",
"theta_11rexliprckthacupcyscyasup 40000.0000 0.0029 0.0011 \n",
"theta_12rexliprckthacup 40000.0000 0.0044 0.0013 \n",
"theta_13rexliprcktha 40000.0000 0.0026 0.0017 \n",
"theta_14rexliprck 40000.0000 0.0017 0.0010 \n",
"theta_15rexlip 40000.0000 0.0018 0.0010 \n",
"theta_16cyscyasup 40000.0000 0.0023 0.0012 \n",
"theta_17cyasup 40000.0000 0.0026 0.0011 \n",
"tau_10rexliprckthacupcyscyasupprz 40000.0000 0.0111 0.0034 \n",
"tau_11rexliprckthacupcyscyasup 40000.0000 0.0041 0.0005 \n",
"tau_12rexliprckthacup 40000.0000 0.0017 0.0004 \n",
"tau_13rexliprcktha 40000.0000 0.0014 0.0005 \n",
"tau_14rexliprck 40000.0000 0.0013 0.0005 \n",
"tau_15rexlip 40000.0000 0.0012 0.0005 \n",
"tau_16cyscyasup 40000.0000 0.0032 0.0006 \n",
"tau_17cyasup 40000.0000 0.0015 0.0008 \n",
"lnL 40000.0000 -12151.9644 271.9991 \n",
"\n",
" min 25% 50% \\\n",
"theta_1cup 0.0004 0.0012 0.0016 \n",
"theta_2cya 0.0000 0.0002 0.0004 \n",
"theta_3cys 0.0000 0.0001 0.0004 \n",
"theta_4lip 0.0004 0.0008 0.0010 \n",
"theta_5prz 0.0019 0.0034 0.0038 \n",
"theta_6rck 0.0004 0.0011 0.0014 \n",
"theta_7rex 0.0009 0.0022 0.0028 \n",
"theta_8sup 0.0002 0.0008 0.0010 \n",
"theta_9tha 0.0005 0.0012 0.0014 \n",
"theta_10rexliprckthacupcyscyasupprz 0.0002 0.0011 0.0024 \n",
"theta_11rexliprckthacupcyscyasup 0.0004 0.0021 0.0030 \n",
"theta_12rexliprckthacup 0.0009 0.0035 0.0044 \n",
"theta_13rexliprcktha 0.0001 0.0013 0.0022 \n",
"theta_14rexliprck 0.0000 0.0010 0.0015 \n",
"theta_15rexlip 0.0001 0.0010 0.0016 \n",
"theta_16cyscyasup 0.0001 0.0014 0.0021 \n",
"theta_17cyasup 0.0003 0.0018 0.0025 \n",
"tau_10rexliprckthacupcyscyasupprz 0.0043 0.0075 0.0122 \n",
"tau_11rexliprckthacupcyscyasup 0.0020 0.0037 0.0040 \n",
"tau_12rexliprckthacup 0.0007 0.0014 0.0018 \n",
"tau_13rexliprcktha 0.0005 0.0009 0.0013 \n",
"tau_14rexliprck 0.0005 0.0009 0.0012 \n",
"tau_15rexlip 0.0004 0.0008 0.0011 \n",
"tau_16cyscyasup 0.0010 0.0028 0.0032 \n",
"tau_17cyasup 0.0001 0.0008 0.0014 \n",
"lnL -12524.5230 -12422.8427 -12161.9050 \n",
"\n",
" 75% max \n",
"theta_1cup 0.0019 0.0046 \n",
"theta_2cya 0.0007 0.0027 \n",
"theta_3cys 0.0007 0.0016 \n",
"theta_4lip 0.0012 0.0023 \n",
"theta_5prz 0.0043 0.0070 \n",
"theta_6rck 0.0017 0.0042 \n",
"theta_7rex 0.0033 0.0063 \n",
"theta_8sup 0.0013 0.0039 \n",
"theta_9tha 0.0017 0.0052 \n",
"theta_10rexliprckthacupcyscyasupprz 0.0068 0.0144 \n",
"theta_11rexliprckthacupcyscyasup 0.0038 0.0084 \n",
"theta_12rexliprckthacup 0.0052 0.0116 \n",
"theta_13rexliprcktha 0.0037 0.0104 \n",
"theta_14rexliprck 0.0022 0.0088 \n",
"theta_15rexlip 0.0023 0.0079 \n",
"theta_16cyscyasup 0.0030 0.0094 \n",
"theta_17cyasup 0.0033 0.0085 \n",
"tau_10rexliprckthacupcyscyasupprz 0.0141 0.0191 \n",
"tau_11rexliprckthacupcyscyasup 0.0044 0.0064 \n",
"tau_12rexliprckthacup 0.0021 0.0033 \n",
"tau_13rexliprcktha 0.0018 0.0027 \n",
"tau_14rexliprck 0.0017 0.0025 \n",
"tau_15rexlip 0.0016 0.0025 \n",
"tau_16cyscyasup 0.0036 0.0053 \n",
"tau_17cyasup 0.0021 0.0041 \n",
"lnL -11879.7943 -11807.7240 "
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## uncomment to print a subset of the table above\n",
"b.summarize_results()"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {},
"outputs": [],
"source": [
"## uncomment to print a subset of the table above\n",
"#b.summarize_results()[[\"mean\", \"std\", \"min\", \"max\"]]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"--------------\n",
"## Examples \n",
"--------------"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Algorithm 00 - fixed tree parameter inference\n",
"\n",
"The 00 algorithm means `'infer_sptree=0'` and `'infer_delimit=0'`, thus the tree that you enter will be treated as the fixed species tree and the analysis will infer parameters for the tree under the multispecies coalescent model. This will yield values of $\\Theta$ for each branch of the tree, and divergence times ($\\tau$) for each split in the tree. "
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"## create two new copies of the bpp object above, but with new names\n",
"A00 = b.copy(name=\"A00\")\n",
"A00n = b.copy(name=\"A00-nodata\")\n",
"\n",
"## set params on the 'nodata' object to not use seq data\n",
"A00n.params.usedata = 0"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"submitted 3 bpp jobs [A00] (100 loci)\n",
"submitted 1 bpp jobs [A00-nodata] (100 loci)\n"
]
}
],
"source": [
"## submit a few replicate jobs from different random subsets of loci \n",
"A00.run(nreps=3, randomize_order=True, ipyclient=ipyclient, force=True)\n",
"A00n.run(nreps=1, randomize_order=True, ipyclient=ipyclient, force=True)"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"True"
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## wait for jobs to finish\n",
"ipyclient.wait()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Summarize results tables for algorithm 00\n",
"Different bpp algorithms produce different types of results files. For algorithm 00 the mcmc results file is simply a table of $\\Theta$ and $\\tau$ values so we can just parse it as a CSV file to summarize results. The same results will be available in the `.out.txt` file, but I find that parsing the results this way is a bit easier to view. Take note that these results tables are showing percentiles which is similar but not identical to posterior densities. You can get those from the `.out` file produced by bpp, or using the `individual_results=True` argument to `.summarize_results()` as shown below. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### The 'no-data' results\n",
"We expect that when we use no data the priors will be returned. In this case we see no significant differences in theta values among samples, which is good. The tau values will be distributed according to the prior on the root divergence, with the other nodes distributed according to a dirichlet process. It will be of interest to compare these divergence time estimates to the ones produced by the analysis *with sequence data*, to assess whether the sequence data causing the divergence times estimates to diverge significantly from their priors. "
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
\n",
"\n",
"
\n",
" \n",
"
\n",
"
\n",
"
count
\n",
"
mean
\n",
"
std
\n",
"
min
\n",
"
25%
\n",
"
50%
\n",
"
75%
\n",
"
max
\n",
"
\n",
" \n",
" \n",
"
\n",
"
theta_1cup
\n",
"
20000.0000
\n",
"
0.0008
\n",
"
0.0007
\n",
"
0.0000
\n",
"
0.0003
\n",
"
0.0007
\n",
"
0.0011
\n",
"
0.0061
\n",
"
\n",
"
\n",
"
theta_2cya
\n",
"
20000.0000
\n",
"
0.0011
\n",
"
0.0007
\n",
"
0.0001
\n",
"
0.0006
\n",
"
0.0009
\n",
"
0.0014
\n",
"
0.0054
\n",
"
\n",
"
\n",
"
theta_3cys
\n",
"
20000.0000
\n",
"
0.0011
\n",
"
0.0007
\n",
"
0.0002
\n",
"
0.0006
\n",
"
0.0009
\n",
"
0.0014
\n",
"
0.0052
\n",
"
\n",
"
\n",
"
theta_4lip
\n",
"
20000.0000
\n",
"
0.0009
\n",
"
0.0007
\n",
"
0.0001
\n",
"
0.0003
\n",
"
0.0007
\n",
"
0.0013
\n",
"
0.0049
\n",
"
\n",
"
\n",
"
theta_5prz
\n",
"
20000.0000
\n",
"
0.0008
\n",
"
0.0006
\n",
"
0.0001
\n",
"
0.0003
\n",
"
0.0006
\n",
"
0.0011
\n",
"
0.0045
\n",
"
\n",
"
\n",
"
theta_6rck
\n",
"
20000.0000
\n",
"
0.0010
\n",
"
0.0007
\n",
"
0.0001
\n",
"
0.0005
\n",
"
0.0009
\n",
"
0.0013
\n",
"
0.0071
\n",
"
\n",
"
\n",
"
theta_7rex
\n",
"
20000.0000
\n",
"
0.0011
\n",
"
0.0007
\n",
"
0.0001
\n",
"
0.0006
\n",
"
0.0009
\n",
"
0.0014
\n",
"
0.0072
\n",
"
\n",
"
\n",
"
theta_8sup
\n",
"
20000.0000
\n",
"
0.0011
\n",
"
0.0007
\n",
"
0.0001
\n",
"
0.0006
\n",
"
0.0010
\n",
"
0.0015
\n",
"
0.0059
\n",
"
\n",
" \n",
"
\n",
"
"
],
"text/plain": [
" count mean std min 25% 50% 75% max\n",
"theta_1cup 20000.0000 0.0008 0.0007 0.0000 0.0003 0.0007 0.0011 0.0061\n",
"theta_2cya 20000.0000 0.0011 0.0007 0.0001 0.0006 0.0009 0.0014 0.0054\n",
"theta_3cys 20000.0000 0.0011 0.0007 0.0002 0.0006 0.0009 0.0014 0.0052\n",
"theta_4lip 20000.0000 0.0009 0.0007 0.0001 0.0003 0.0007 0.0013 0.0049\n",
"theta_5prz 20000.0000 0.0008 0.0006 0.0001 0.0003 0.0006 0.0011 0.0045\n",
"theta_6rck 20000.0000 0.0010 0.0007 0.0001 0.0005 0.0009 0.0013 0.0071\n",
"theta_7rex 20000.0000 0.0011 0.0007 0.0001 0.0006 0.0009 0.0014 0.0072\n",
"theta_8sup 20000.0000 0.0011 0.0007 0.0001 0.0006 0.0010 0.0015 0.0059"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## summary for the first 8 params of the no-data run\n",
"A00n.summarize_results().head(8)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Results with data\n",
"You can see here that the theta and tau values are quite different from the priors, which is good and tells us that our sequence data is informative. If you check the `.out` file produced by bpp you can also see \"effective-sample-size (ESS)\" which is indicative of whether or not we've run the mcmc long enough. Ideally we would like each ESS value to be >100. "
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"
"
],
"text/plain": [
" mean S.D. ESS*\n",
"theta_1cup 0.001350 0.000420 40.3\n",
"theta_2cya 0.000692 0.000324 43.3\n",
"theta_3cys 0.000717 0.000177 60.4\n",
"theta_4lip 0.000920 0.000204 105.9\n",
"theta_5prz 0.003490 0.000489 377.8\n",
"theta_6rck 0.001130 0.000321 121.2\n",
"theta_7rex 0.002355 0.000554 98.2\n",
"theta_8sup 0.000954 0.000396 51.8\n",
"theta_9tha 0.001655 0.000485 95.3\n",
"theta_10rexliprckthacupcyscyasupprz 0.006313 0.002876 7.6\n",
"theta_11rexliprckthacupcyscyasup 0.002176 0.000815 22.9\n",
"theta_12rexliprckthacup 0.004195 0.001267 36.9\n",
"theta_13rexliprcktha 0.003357 0.001846 14.0\n",
"theta_14rexliprck 0.001622 0.000920 116.8\n",
"theta_15rexlip 0.001853 0.001100 95.7\n",
"theta_16cyscyasup 0.002449 0.001300 172.2\n",
"theta_17cyasup 0.003107 0.001055 110.5\n",
"tau_10rexliprckthacupcyscyasupprz 0.008021 0.002025 7.6\n",
"tau_11rexliprckthacupcyscyasup 0.004015 0.000519 107.8\n",
"tau_12rexliprckthacup 0.001575 0.000461 15.9\n",
"tau_13rexliprcktha 0.000965 0.000194 22.1\n",
"tau_14rexliprck 0.000911 0.000181 21.9\n",
"tau_15rexlip 0.000822 0.000176 27.5\n",
"tau_16cyscyasup 0.002851 0.000577 124.8\n",
"tau_17cyasup 0.000876 0.000438 30.4\n",
"lnL -11880.819322 21.881173 22.1"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## get individual results for each rep\n",
"reps = A00.summarize_results(individual_results=True)\n",
"\n",
"## check convergence (ESS*) for each param in rep 1\n",
"reps[1][[\"mean\", \"S.D.\", \"ESS*\"]]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Draw some histograms\n",
"Below we draw the posterior parameter distributions for the analyses with and without sequence data. In a real analysis you will probably want to run the analyses for longer, but you can see here that the results depart significantly from the prior expectation. "
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {
"scrolled": false
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAyAAAAMgCAYAAADbcAZoAAAJNmlDQ1BkZWZhdWx0X3JnYi5pY2MA\nAHiclZFnUJSHFobP933bCwvssnRYepMqZQHpvUmvogJL7yxLEbEhYgQiiog0RZCggAGjUiRWRLEQ\nFBSxoFkkCCgxGEVUUPLDOxPn3vHHfX49884755yZA0ARBQBARQFSUgV8Pxd7TkhoGAe+IZKXmW7n\n4+MJ3+X9KCAAAPdWfb/zXSjRMZk8AFgGgHxeOl8AgOQCgGaOIF0AgBwFAFZUUroAADkLACx+SGgY\nAHIDAFhxX30cAFhRX30eAFj8AD8HABQHQKLFfeNR3/h/9gIAKNvxBQmxMbkc/7RYQU4kP4aT6edi\nz3FzcOD48NNiE5Jjvjn4/yp/B0FMrgAAwCEtfRM/IS5ewPmfoUYGhobw7y/e+gICAAh78L//AwDf\n9NIaAbgLANi+f7OoaoDuXQBSj//NVI8CMAoBuu7wsvjZXzMcAAAeKMAAFkiDAqiAJuiCEZiBJdiC\nE7iDNwRAKGwAHsRDCvAhB/JhBxRBCeyDg1AD9dAELdAOp6EbzsMVuA634S6MwhMQwhS8gnl4D0sI\nghAROsJEpBFFRA3RQYwQLmKNOCGeiB8SikQgcUgqkoXkIzuREqQcqUEakBbkF+QccgW5iQwjj5AJ\nZBb5G/mEYigNZaHyqDqqj3JRO9QDDUDXo3FoBpqHFqJ70Sq0ET2JdqFX0NvoKCpEX6ELGGBUjI0p\nYboYF3PAvLEwLBbjY1uxYqwSa8TasV5sALuHCbE57COOgGPiODhdnCXOFReI4+EycFtxpbga3Alc\nF64fdw83gZvHfcHT8XJ4HbwF3g0fgo/D5+CL8JX4Znwn/hp+FD+Ff08gENgEDYIZwZUQSkgkbCaU\nEg4TOgiXCcOEScICkUiUJuoQrYjexEiigFhErCaeJF4ijhCniB9IVJIiyYjkTAojpZIKSJWkVtJF\n0ghpmrREFiWrkS3I3uRo8iZyGbmJ3Eu+Q54iL1HEKBoUK0oAJZGyg1JFaadco4xT3lKpVGWqOdWX\nmkDdTq2inqLeoE5QP9LEado0B1o4LYu2l3acdpn2iPaWTqer023pYXQBfS+9hX6V/oz+QYQpoifi\nJhItsk2kVqRLZETkNYPMUGPYMTYw8hiVjDOMO4w5UbKouqiDaKToVtFa0XOiY6ILYkwxQzFvsRSx\nUrFWsZtiM+JEcXVxJ/Fo8ULxY+JXxSeZGFOF6cDkMXcym5jXmFMsAkuD5cZKZJWwfmYNseYlxCWM\nJYIkciVqJS5ICNkYW53txk5ml7FPsx+wP0nKS9pJxkjukWyXHJFclJKVspWKkSqW6pAalfokzZF2\nkk6S3i/dLf1UBiejLeMrkyNzROaazJwsS9ZSlidbLHta9rEcKqct5ye3We6Y3KDcgryCvIt8uny1\n/FX5OQW2gq1CokKFwkWFWUWmorVigmKF4iXFlxwJjh0nmVPF6efMK8kpuSplKTUoDSktKWsoByoX\nKHcoP1WhqHBVYlUqVPpU5lUVVb1U81XbVB+rkdW4avFqh9QG1BbVNdSD1Xerd6vPaEhpuGnkabRp\njGvSNW00MzQbNe9rEbS4Wklah7XuaqPaJtrx2rXad3RQHVOdBJ3DOsOr8KvMV6Wualw1pkvTtdPN\n1m3TndBj63nqFeh1673WV9UP09+vP6D/xcDEINmgyeCJobihu2GBYa/h30baRjyjWqP7q+mrnVdv\nW92z+o2xjnGM8RHjhyZMEy+T3SZ9Jp9NzUz5pu2ms2aqZhFmdWZjXBbXh1vKvWGON7c332Z+3vyj\nhamFwOK0xV+WupZJlq2WM2s01sSsaVozaaVsFWnVYCW05lhHWB+1Ftoo2UTaNNo8t1WxjbZttp22\n07JLtDtp99rewJ5v32m/6GDhsMXhsiPm6OJY7DjkJO4U6FTj9MxZ2TnOuc153sXEZbPLZVe8q4fr\nftcxN3k3nluL27y7mfsW934Pmoe/R43Hc09tT75nrxfq5e51wGt8rdra1LXd3uDt5n3A+6mPhk+G\nz6++BF8f31rfF36Gfvl+A/5M/43+rf7vA+wDygKeBGoGZgX2BTGCwoNaghaDHYPLg4Uh+iFbQm6H\nyoQmhPaEEcOCwprDFtY5rTu4bircJLwo/MF6jfW5629ukNmQvOHCRsbGyI1nIvARwRGtEcuR3pGN\nkQtRblF1UfM8B94h3qto2+iK6NkYq5jymOlYq9jy2Jk4q7gDcbPxNvGV8XMJDgk1CW8SXRPrExeT\nvJOOJ60kByd3pJBSIlLOpYqnJqX2pymk5aYNp+ukF6ULMywyDmbM8z34zZlI5vrMHgFLkC4YzNLM\n2pU1kW2dXZv9ISco50yuWG5q7uAm7U17Nk3nOef9tBm3mbe5L18pf0f+xBa7LQ1bka1RW/u2qWwr\n3Da13WX7iR2UHUk7fiswKCgveLczeGdvoXzh9sLJXS672opEivhFY7std9f/gPsh4YehPav3VO/5\nUhxdfKvEoKSyZLmUV3rrR8Mfq35c2Ru7d6jMtOzIPsK+1H0P9tvsP1EuVp5XPnnA60BXBaeiuOLd\nwY0Hb1YaV9YfohzKOiSs8qzqqVat3le9XBNfM1prX9tRJ1e3p27xcPThkSO2R9rr5etL6j8dTTj6\nsMGloatRvbHyGOFY9rEXTUFNAz9xf2pplmkuaf58PPW48ITfif4Ws5aWVrnWsja0Latt9mT4ybs/\nO/7c067b3tDB7ig5BaeyTr38JeKXB6c9Tved4Z5pP6t2tq6T2VnchXRt6prvju8W9oT2DJ9zP9fX\na9nb+aver8fPK52vvSBxoewi5WLhxZVLeZcWLqdfnrsSd2Wyb2Pfk6shV+/3+/YPXfO4duO68/Wr\nA3YDl25Y3Th/0+LmuVvcW923TW93DZoMdv5m8lvnkOlQ1x2zOz13ze/2Dq8ZvjhiM3LlnuO96/fd\n7t8eXTs6/CDwwcOx8DHhw+iHM4+SH715nP146cn2cfx48VPRp5XP5J41/q71e4fQVHhhwnFi8Ln/\n8yeTvMlXf2T+sTxV+IL+onJacbplxmjm/Kzz7N2X615OvUp/tTRX9KfYn3WvNV+f/cv2r8H5kPmp\nN/w3K3+XvpV+e/yd8bu+BZ+FZ+9T3i8tFn+Q/nDiI/fjwKfgT9NLOcvE5arPWp97v3h8GV9JWVn5\nBy6ikLxSF1/9AAAABmJLR0QA/wD/AP+gvaeTAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAHXRFWHRT\nb2Z0d2FyZQBHUEwgR2hvc3RzY3JpcHQgOS4xOJQFEHMAACAASURBVHic7N17eBvXeeD/FzMYDDC4\nEgBJURRFStbFkuyETlZuY8cR7dTpxrGbxIlcx2lqp1m5a7vZ2E2d/dVOKrmb+GnsZOP04nbt3VTe\ntNnUTjbKk9q7rZsErhM1tTcNc7EdXSxBpiiJNxAgiQEGAwx+f5CQKZo3USQAkt/P8/CRcM7BwUvw\nYGbemTMHIgAAAAAAAAAAAAAAAAAAAAAAAKh3aq0DWAzhcOjaUND/700z96KIyLatm/r9hu+6zMjo\nExfSz4VQFcW/ft3apy2r8PNiqXTmQvtDfWuMR+/UNC1iWdYxw/Bt37Sxo9/lEjHN3HML7WehsaiK\n4m+IhD8cCBjv9er6lmKpdMxxnPxC+0P9q8Pxd1sgYNzg1fUthULhULlcthfaH5aXehqLFdGGyJ6m\nxtgXvV69I5s1n5sa10KPGVA/6mncGYZvZyQc2uM3fF2VH9suvuw4jrnQPlcid60DWAyRUHCPYfh2\nDwymHq11P5qmdUTCwevCodAeTXN3KqoSuZCYsDzEY9G/MM3cU5nMyLO17mf9+tbve3V9p4gMiEhj\nPB7dezzZ8yu2bfdeSGyoX/Uy/lRF0TZ0rP9XTXPvkNfH372vvpq8rOQ42QuJDctDvYzFCk3TGpub\n4l8SEb+Mj8k3GBxKPVywi9PWYXmop3EXCYV2h8PBeyeXZc3cU7ZtM8YmWfZXQMLh0O5wOLRHcblC\nfsPXmhkZ/U5jPHqvbRcHxOUaCAX9t3h1vTWXz/+88pyJ7PQ2v+H7VdsunnAcZ3S6fjRNi0TCoY8E\nAsYNfsN3eaXtbPEEg4H3BQP+W1RVjamqEsuMjD5h28XkpHivCwX9t0zuT9O0jlg08gkRMW27eErT\ntMZYNHKviLhsu5hsjEf3aprW5PXqO0JB/25N0xosy/rlkr2pOC/rWlse0j2eKxVVjbjdqlMo2K9G\nwqG78pZ1UPd4LgoEjBtURbELBbtHZPwgrSESvjkQMG6cfIZ4aj/ZrPkjw/Btnxir13p1vbVQKPyy\nXC47M8ViGL6d8WjDvtRw+sHkayff7TjOkWDAf5vjOGdMM3dwYmzd5jd813l1fUOhUHjJ7XaHYtHI\nf9Y0rdWyrJ+LjH+uQkH/bitvHXS73aHz/Rygeupp/Pl83jdFGyL/cXAo9akTPb03qaoifsP3gTHT\nfNq2iz2LMP522HbxFOOvPtXTWKxoa215PG8V/k33eN5q28WXMyOjT2ma1hgJh+4yc7mEaeae8xu+\nK0uOY1qW1V0ZeyKSiYRDv+PV9Yu5ilff6m3cRRsid+etwrPHjr92+eBQ6oHBodQDtl0cqBzrqYri\nCvj91wUCxrUikrft4ikRkcZ49F5N07aoqqIZhvEur1dfO3G8ePZKikwcFy71e1oNy/4KiG3bA07J\nSauK0mrmckcq5Ybhu9IwfDsnHvo1zb25r3/wgeam+L3RhshDtl3sVlQlEo9F957sPX391H7Gz+S1\n/URVlA4ZP2sSiceie48eO7Fhtiw2kxnZn8mM7G+MR++Mx6J/Mblu7Zrmx8Lh4B4RSYuIPx6L7j3R\n03u5iDTFY9F9IqkB08y9qGnuxonH+0wzlxj/v1Q2fgUR8Xt1z4N9/YP3L9LbiAtgWVYyGPCLU3LS\neatw9ipDtCFyr4hkZfzM276TvaevHh3LJtavb/1Hr65faZq5HxqGb+fEGeLLp/YTDPivWNfakpDx\nv7ktIpFgwP++Ez29N80YS946cqKn92rbLh4REbFtOzlRZWua1rqho+1fVUVpkvEx2BgOB/ccT/Zc\nGQ6FPqxp7sax0bGnSo5jNzfFH3JKzkAqlX5wIZ8DVE89jT/TzHW/cuhok4iIYfi2G+MHdkkrb/10\ngePvCxs62n6hKkrTRLyd8Vj0oaPHTmzlil79qaexKCISDPi7dK9+7fFkz45gwL9npnbxWPRe08y9\nlMmM7I+EgrsNw7c7HpP7J2KORBsinziefO2ykuOQhNSheht3huHbbtvFtdu2biqISDY1nH64r3/w\nQU1zd0w6niuIiCcek30ne0+/b3Qs+0w8Fr235DiiKkrENHM/LDnOwUnj1hj/PcaPCxf5LawJpdYB\nXCjTzCVs235ZRGRgMPVwpbzkOP2HjxxrPnzk2DoRyXp1vUtVFC3aELnPtosvnjrT9+HTp/tuEhFp\nbmr83NR+FFX1Dw+nHz7Ze/ryVw4dbcpkRh+W8URmx0Li9Op6azgc3GOauQOvHDracDzZs8M0cy96\ndf3K+Ty/5Djpo8dOtB4+cqzBtosvRhsi92qaxvSuOlCZsmfb9suZzMhTlXLTzD39yqGjDSd7T79H\nRMTr1XeFw6HrvLreNTqWfXhgKHXX4FDq91RF2RyNRm6b2k/JcezBodTdR4+d2HD4yLF1JcdJGuNn\nQGZUcpz0xFju9ep6R3NT42O2XXwplUrvb4xFP6EqSuupM/3veeXQ0aa+/sHbnZJT8Or69tGxscdF\nJBIIBq41DF+nqigdo2NjTy325wCLr57GX4VX1yPtba0/9Or6tVbeSohIYSHjT/fqa1VFaTXNXGJg\nKPWp48meywaHUp/SNLex2O8jLly9jcXmpsYvDQ+n9y3kZElf/+BNrxw62pAaTu/TNPeOcDi0+3z7\nQHXU07ibOC7rKDklc3AodXfesl6ONkQ+Fw6Hrp0UV+LwkWMNx5M9G0Sk0NzUuK9SpyqKcbL39LtO\nnem/6WTv6ftfOXS06WTv6feJiOQtqzuVSn958d652lr2V0BmYuWtI5PmHJsiIrpX7xCRiKa5d7a3\ntb5Uaatp7s1Tn++USrbm1jY3tETuXqcohm0XL+jmIUVVNouImLnxG6Ly1viZahGR+WxIrbz1YmUj\nmhkZeToei+7UNPda27bTFxIXlo6Zy3WLiJQc5+wZGY/m7hARCQb89wUD/vsq5Ybv7NW6s5ySM6Dr\netemjdH7Js68ReT1K2GzmrgC99DoWPaJ06f7PlVynKw28dqVua2p4fTjqeH04yIikhEz2hB5KBjw\n77Zt+8h4feZri/05QPXUcvzlLSv9yqGjDdGGyG3NTfG/DodDLy5k/Nm23Ts4lLonHArd0t7W+oKI\nZEfHss+kM6NPzfDSqEO1GIvRhsgtE/v2xsZ4dK+IiKZp28Ph0HWTrgzPKDWcPiAiMjqWfSbaENmn\nTuzDsXzUYtzZtp1+5dBRV+Vx1sy92N7W+oJHc++0bfvgRFzPlhzHLllWr2nmXjQM39kTeqaZe3F0\nLJuoPPbqemdLS/PTJcc5dfp0//Ulx1kxx3wrNgGZjpW3kiJi5y3r2ePJnvGMWNc7FFXRprYdvyck\neHdf/+BHU8Pp/WvXND809aai8+GUnKSIiOHz7RKRR7y63hoM+v9DPm8drFzWVRSldSKmy6c+X/fq\nOzVNa3RKpXQ4FHqPiNiVeYNYPgoTczcHh1IfHRhM7VcVRdO9+mZ7mhsg4/Ho3mDAv/t4smeDbdsD\nGzrW/0RV5l7UYGKq322nzvTflMmMHKiU23axV0Qk2hC5LjWcfiYY8F/r9epXpDOj/z1vWcm8ZSUM\nw/cep6Qn85b1rG3bvdGGyG2L+TlAbS31+AsG/FdEGyKfGx3LPp4aTn+tUq6qin8h409VFH/WzP1j\nKpV+XFFVIxIO3haPRR+yLKt7YDD14CK/PaiipR6Lqqp4RMSMx6J3Vco0zb0lGPBfmRpOJ+eKb9I4\nvU5EpFRyjsz1HNS/pR53huHr9Bu+92bN3LdNM9ft1fVOEZFSyTl7vGb4fNeqivKIpmlNhuHbadvF\nl6fry6vrnevXt35fRAZee6336rxlrahppysiAalc6Whva/3WiZ7e98/Szk4Np/dFGyKf27Sx4yci\nIprm7sxkRu83zdyDk/tJj4w+KiLS3BR/bGIFjQua7pS3rKRp5r5mGL5btm3dNCwimojIiZ7et0/c\nUJSNNkTuizZE9ohI49Tnq4rSuGlje69M3AOSGr+svGIy4RUgbRi+6xrj0Tuzs8zPHBsdezbfED4Y\nj0UfCwT8ezRNW6sqSuP4fRv2wOR+nIkN1oaOtl9MPH3GVVwqvLq+feI+I1m7pulba9c0iYjI4FDq\nrtRw+suBoP/DzU3xp5ub4gMi0pi3rIOpVPpBEZFMZvSJ5qb4X6uK0pjqT39U5PUzR4v1OcCSqYvx\nZ5q5n7a0NLcahu9vm5vij4hIY8lxetOZ0a+piqKd7/jTvXp7e1vrCyXHMa289ZLu1XeIiJ01cxe0\n0g2WVF2MxYHB1P6BwdT+yuNtWzeVTTN34GTv6fsNw7d9rl9iYpymRSSSt6wXJ0/tQV2qi3Fn28WB\nhobInngser9M3OtWGT+6V98pImIYvmu3bN44LCIeEZG+/oFp7+dtbop/qZLwbOho+4mIiGnmDpzo\n6b39PN6XurUiEpC+/sFPWZZ1RCYujU1dUm/y477+wQfzVuGgR3NfIeM7so+a5vhlusn9TEwTeNdE\nO7NgF494NPebVeWNV0umkzVzL4qk9k1ereBET++Hw+HQUx7N/WYRMdOZ0a9VbqQ80dN7ud/wvU9E\npGAXD3o0967spPWr85Z1IJMZfUZVlbUFu/hiJjPyzAW9aVhUJ3tPv8fr1a8t2MVTtl0cGBxK7av8\n/SY/LjmO/dprvV3hcOjDqqq0i8hAOjN6dnm+yf2kUunH85Z1RFWV9omzb35VVdZ6db11pjMhiqrI\n4FBq39TyrJl7MW9ZyePJnh2RcHC3iDSWSs6RyjxXEZHRseyB5qb4Y5X/i5ydLrPgzwGqo17GX8lx\nsq++mrwsHA7doqrK2lLJOTU6ln3Ktu20LSLnO/5MM/fy0WMntkbCwfeKSKOZyz2TNXNPm2Zu2jOG\nqL16GYtTDQ6l9hUmzjRPjWu6ZXhPnel/l0dzX1EqOacymZGvcQN6fauXcWfbdu/xZM9lM23nREQy\nmdF9dtE2RUTLmrlnKsegU8dhemT0CTN3bjJVmOFqCbAktm3dVG5va32y1nFg5du2ddPJ9rbWb9U6\nDqxOjD/UWntb65Pbtm4q1zoOrDyG4evatnVTuTEevbPWsdSDFXEFpJoMw3dFe1vrgRmq+185dPSS\nqgaEVae9rfUxY+Jq2VSjY9nHT/aePu/lmQ3Dt3Ptmua/FpHW9Mjoiri8i6XB+EO9WNfa8rmZltdd\nSVNVUF8Yd1g1GuPROycv4QYsNk3TOhrj0b3hcOi6WseC1Yfxh3oRDoeu5Qw1loKmaa2N8eid87kH\nCQAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADmw1XrABZbR0eH1tHR4ZmuLp1O\nF7u7u61qx4TaikQiSmdnp09EpLu720yn0+Vax4TVIxKJKB0dHZ7u7u58rWMBFsM0+9lyIpEwaxYQ\n6lpXV5d/uvJ0Om13d3cXqh0P6oN7oU9sjEf3iogMDKYeWLxwLtxfP/THO98hp2+Yru65Afvn19zV\n/bVqx4TFdz7jr7OzM/7dO999j4jIOx+VLyUSif6ljg8r2/mOvz/a/Wu/xrYHi6XW+99br7vm3/1R\n19bfeD2g9Vk1kfhsLWJB9Z3v+Lvjf3zh93vG0sbU8sIrr/28++bfYbu4Si04AalblqmV+w+1T18X\nfq3K0aBOlHtmGBNANWSGwrUOAVg0lqlP3qa6RIZqGQ7qW89Yuulwpi82tTxojXFMtooptQ4AAAAA\nwOpBAgIAAACgakhAAAAAAFQNCQgAAACAqiEBAQAAAFA1JCAAAAAAqoYEBAAAAEDVkIAAAAAAqBoS\nEAAAAABVQwICAAAAoGpIQAAAAABUjbsaL9LZ2RmIRCLGfNomk8mxZDJpLnVMAAAAAKqvKgnIHfv+\n8K0/jhZ/dz5tW/73Dx5LPvKXiSUOCQAAAEANMAULAAAAQNWQgAAAAACoGhIQAAAAAFVDAgIAAACg\nakhAAAAAAFQNCQgAAACAqqnKMrxAvfiv/+GDt6TfdZn1Pw8PP7l///6hWscDAACw2pCAYFV5syd7\nWVnPaP/c0fHNWscCAACwGjEFCwAAAEDVkIAAAAAAqBoSEAAAAABVQwICAAAAoGpIQAAAAABUDQkI\nAAAAgKohAQGAJdbR0eHv6uoK1ToOAADqAd8DAgBLrN2dv3TXRWvekkhIotaxAOfr1vdff0lHNBCv\nPN71lku2Sv5YLUMCsMytqgSkvaWpae/Hbu6atjIUH3zgS3/+i+pGBABAfbt15+ar3qFnrqo8djVq\ndrmnlhEBWO5WVQLS4Vc2fma7sXW6OuVt1z5PAgIAAAAsLe4BAQAAAFA1JCAAAAAAqoYEBAAAAEDV\nkIAAAAAAqJpVdRM6AABYdOHSk5//s8kF73z0//yXRCLRX6uAANQ3EhAAAHAhyuWeQ7FaBwFg+WAK\nFgAAAICqIQEBAAAAUDUkIAAAAACqhntAAAAAgHno6Ojw3vrOX/3VaSt1o/Tcy8cOzfDUciKRGFi6\nyJYXEhAAAABgHjo6OkKf2W7cPl2dq22rXb5ItBnqhtRE4uNLG93ywRQsrEq7du2K79u3r6mrqytU\n61gAAABWExIQrEq79JH3fzr42iO7LlrzllrHAgAAsJqQgAAAAACoGhIQAAAAAFVDAgIAAACgalgF\nCwAAAFW17uLNmz/5xJ/eMV3dd7/0lb/p7u4erXZMqB4SEAAAAFRVyavFRjfE1kxXF4lEvikiJCAr\n2LJMQL70yf947Ztjvk3T1b15wzqR/heqHRIAAACAeai7BOSK979719rLtm2brc2vhLV3XNR3eFOo\nUOjXCrnM5DpFKb3kLG2IAAAAABao7hIQK2JsHN0Q2zpbm76xzKawojT6FCUz7ddNAgAAERH52B98\n4sptH7nh5vm0HfvWDz63b9++M0sdE4DVre4SEAAAsHjMYkE7nOmLzaftWhEmEQBYcizDCwAAAKBq\nSEAAAAAAVA0JCAAAAICqIQEBAAAAUDUkIAAAAACqhgQEAAAAQNWQgAAAAACoGhIQAAAAAFVDAgIA\nAACgakhAAAAAAFQNCQgAAACAqiEBAQAAAFA17loHAAAAAFS884Pv/dWWzm0jU8vHzgz0fefr33il\nFjFhcZGAAAAAoG74rrr0A8E3NWlTy9ceH3peSEBWBKZgAQAAAKgaEhAAAAAAVcMULAAAAGBpeb//\n1b/aNW2NbpSuvum3f1DleGqKBAQAAABYWp6r+l/43ekqXG1bh0RkVSUgTMECAAAAUDUkIAAAAACq\nhilYAABARETa29ujXV1d55RForGgFHLeswWqpoqma2Jb+Zn66ezsNETEX3mcTCatZDJZXPyIASxH\nJCCvM/Z+7OauaWt0w3rg0a/8S3XDAQCguuLXvPUDWy7y7JhcdmIss1nvHW6sPHYUlxMwAnZbxvrx\nTP188b2/cpdc0XI24fjjxKFvP/DoV360NFEDWG5IQF7n/cx24/bpKlxtW4ceeFRIQAAAmI+B15rK\nPYde/yI5y/TO0hrAKsM9IABQBbvecsmWrq6uUK3jAACg1rgCAgBVsKtRu7Krq+tAIpEYqXUsAICZ\ndXZ2ev/rng/eKJmh8NS6yPqLdOn7t1qEtaKQgAAAAAATIpGIuqtRu7Scz7RPrXN5sna5FkGtMEzB\nAgAAAFA1yzsBcblUUVT3OT8ul+cNZePlrlqHi/rzR7uv/c3Sk5//s+/9xWdvq3UsAAAAq8GynoI1\n4jXWDrndTZPLUvnR5qg/GJ/atsku9HpFTlcvOiwLI4NGueeQJlbYV+tQAKAeNOjGOStWecysWqtY\nAKxMyzoBsUU00+XSJ5eli4Wi1+XyT21bEtfyvtoDYNn7o4/d/OuJROLpRCKRqnUswEw+2HNkhzHc\nd3bf6sTXObWMB8DKs6wTEGCxfewPPnGVWSzMerYv5vWn/vxPvvCzasWEFaTn0HphuwsAWOXYEQKT\nbPvIDTcdzvTFZmvz6+2XPk8CAmBVU1StJxx7q4iIo/scJRw7Z5bB5LJQqZQSMWsRJYA6xbQkAKii\nzs7OaK1jABaDWXZ0s+zoplPSz/5/mrKilDnZCeAcJCAAUEX/9aarb691DAAA1BIJCAAAAICq4bIo\nVge3RxeXS0TVVNF0tyiqXzRdrzzuaL0ouG/fvibhMwEAQF2KRCJGV1dX03R13d3dqXQ6Xax2TFgY\nDrawKgx5fesHFaXRUVyO4jOUM6MDm9b4jHDl8eGQsv7UNZs7t4jYtY4VK5731vdff+kT3/r7n9c6\nEABYTnwda9+05Y/3PDJt5R89fncikeivckhYIKZgAUB1eW7dufnttQ4CAIBaIQEBAAAAUDUkIAAA\nAACqhntAAAAAgNoJfO/Be+6YtiYcy11z16f3VzecpUcCApy/wIfuvqNrrkav/r+fHn3hBwdPViEe\nAKhru95yycV75eazj5OpsYEnvvX3L9UwJKCeKO/QM1dNV+FqXDMkIvurG87SIwEBzp8neOPb5/wy\nuYtEHiMBAbDa5T3eaPua0Cdu09vOTvvu6dhx6tJ9d/9iatuxb/3gs/v27euraoAAqo4EZH6Umdad\nFhFh2TcA56O9pan51g/dtP6J//Xka7WOBVhqTrnsMp2SrpSdswnIUMEMHM70xaa2XStSrm50AGph\n1SQgOa+v2fR4A65AWJ2ufqxUCAQC4RYRkVDB6tcK+cyk6sA/3bBx2nWnXW1bh9RE4uNLEDKAFarD\nr2y47bprNpCAAABWo1WTgJguxT9WdvyKoky78teZ/JhnjaKERUR8iprRqhseAAAAsCqsmgQEAGpK\n9XiOlwq/Yrg9MZfHqxbWNL/3k0/86fapzT7w6+9Rh6xsabau/uEv/+fzf/4nX3jD/HkAWK3e+cH3\nvq2lc1tmavnYmYH+73z9Gy/XIibMjAQEEJGwS/FeWnKa1pSK9mGR4VrHgxXI5VKS+bH2uC/gV8qO\nknEVftW1tXX9SCFfmNzsYN9x+3Cmb9aLsKP57CtLGywALC++qy69MfimpjdsO9ceH3peSEDqTt0m\nINcWSxs7UmcapqvzeP11GzeWp8Z8zn/NySObJL4u92ysaVESkHfd8ds3v/WjH/jAXO3+/hOf/Wwi\nkWDVl1Umljf9mt7gqnUcQD1pb2+PdXV1zaepi+0mLkQkEnF1dnY2TlfX2dlpCF/WvaTq9kDenR9T\njeE+fbo6J77OqXY8wPkaymd9hzN9oXk0ZdUXnJeOjo7QbCvziYik0+lCd3d3uloxAYshfs1bb9xy\nkWfHXO22hJuHEm9mARgsXGdnp/HdO999T7nnUPvUOldbC6ubLrG6TUCAWnArinppyWkSETns8QxZ\npeKsc/GBWohd/dYPbHnLuptnaxM8PvR8963/6S+rFROWpwbd8H6w58jZA/6W4aENw6q25OuwqC7F\n5XN7zh6DlMvlcr5ks70F3kjp6uryz1DnJBKJXDWDWSwkINNR3R7RPN7XH2uqaJ7pN8gu1TPbmcju\n7u50Op0uzFSP+uJWVOWao92bREROvuntI9VIQGa6BDxFOZFIDCx1LKiuN5WcNS+49VPZomXXOhas\nXpNnG2gu1ZBYy5K/ZlvRjnyw5/DZxKcQjBb/LhDim9Gx6NZdvHnLJ5/40zumlsejMV1yJ4O1iOk8\nBb9757s/PV2Fa/uVLymXvP1vqh3QYiABmcaQR1/XoyrrK48dxeUoPv+0cwEvVt3/OtN3hIiI/JrI\n3XxRIWaz7qPXf8rM7Jz1jCPTDVamprFMwBOKqNmiLHoCsu7izZum2+lOFtMD2ftu/p3/udivDczF\nXSwq0ZHU2bO6pqpZEpjPjNWzPB+6+46u+TQcOzPQ952vf4OFG1apkleLjm6IXTW1vOj2uOWkLIdp\nquXppomJiLjWXbxsv0uKBORCKapHNM+Mc/gnn93masjyskVcMbPk2K/q3pRZLBRrHM68dralwczR\nJ//mayerEA/qXMmrxUc3zH4quyXcnOnq6vq/c/XV3d09kk6n84sXHXDBvMEb3377fBpOrIJEAoI3\nyHu8a0pevzG13OPxFbRCbjkkJ8tWVRIQw63ZW8KxE/Np2xaIuEWkGPf4db8jb1jPWUTECTe6FU0v\n6orb7XKK5xwYNjSssf1W7g1TVXTF7S4HG0TR9GkPJCc/b2q/ldeb7nnpSHNw1D3tl6uLiMjt79jy\nSDl7VUFEZJseOCplZ6xS91P3mh+l0+lpV1z69re/Xeru7p729xeRwsRPhS0i850qlBORed3E393d\nbabT6RVxg7TWctFrfkXJVP6Wlb935XHlX0+0xae6lJyIyNUF261kR+Vf2i62ckU7myvZWq5ol2K6\nX9sSbh6T8fdxxhWzKmN5rtjaAhFV5vj7tQUibvnoB359rr4uV6Ph7Zu2zDl/u7u725NOp+eaN+qS\n8fEyl5KIzHlwmkwmrWQyWetErjYa1/eJqgUbDP+pgMcXdJWl2OAL5EREbglEgs9a2UOVpvMZN/Ns\nM59xpV7/5U/fM1f4t/5bz4ETJ07MeiU3nU57u7u7s3N0Nd8xNZ/tlEvmWMAhnU6Xuru7V33iZLg9\n+S3h5mn3wUGP7ve3XHT20oPmUmx/KOabvM+bug90wo1ur+6Xyn5yun3k5LLp9r9Tn+MNx90fCIRd\nFxWsNk8w7hIRGQg0FAql4tnxYpYKhXyxmBWZ3/iu2Np16fquEz/7g/m0Hfzej5+fa6xP0BKJxHzH\nliUiY3O2EpHu7u58Op1eUffCtAUi/TLN7z/bdmymuoU8Z5cR3vxmx3nDyV+XxysjLs1rud+4ywzH\n2/xaIT/kcusjb6iMrQ26FPcb//ZlxyXxdUXX+N97qmFpXO9yzTRmG9e7XTNt02erC8Xc3//qX71n\n+rp49oEv/fm0yw+n02m7u7u7pl85sOAlIDVN84qI2LY95wews7PTE4lE5ntTW2WnosvMS6BV2ky3\nA5ppp1T5XWfaYU1+3tQ+ZtvRKTL7jnIp+rXl3ME45454IW3rOQE5n/EXiUSUzs7OiIy/p1PHztR/\nJ7/vc7Upy+wH3vN9r+fTbjH7qslrrqQETn76uQAAIABJREFUZAHjzxARVUQ88sb3qyznHpQv1t+m\n2n/junu9lZqAnM/4ExHp6OjQOjo6PDNUKzK+v62Ybl853b5rcpu59sXz7XPqfn1qm4K8fgC3JPu9\n82i7VPvduk9Aznf8zXID9Wzvy2zHcuf7HGORX2em47PZjjPzs9TN9XqLXrdSt40AAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADAarPgLyI8H+Fw6N97NPevVOO1sHgGBlMP1DqGxaBp\nWnMkHPyPtY5jtcuaub8zzdwvax1HtTH+aqdgF3+WyYx8q9Zx1Jph+C7zG77fqHUcmFs6M/pXtm33\n1TqOxdQYj+6tdQw4V8EunshkRvbXMgZ3NV7Eq3sMr643VOO1gKlURXEbPh/jr8YKdtFz7hd+rw6M\nv9pRFMvI1DoIYJVj+1d/FMUaYtsIAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAMBKpNY6gMXUGI/e\nqWlaxLKsY4bh275pY0e/yyVimrnnFtrPhcRjGL6dkXDoNr/h+1XbLp5wHGf0fPvYtnVT2W/4dmRG\nRp+6kFiw9Opl/BmGb3skHLrTb/i6Jv9omtZkWdbLC+kT9SscDl0bCvr/vWnmXhQR2bZ1U7/f8F2X\nGRl94kL6WahJ272rRGTItosDF9Iflod62f6JiKiKojVEwr8dCBjv9ep6U6FQOFIul53FjBNLp962\naSIimqY1trW2PG1ZheeKpVJ6ymt8xKvrGwqFwkuVcbYYpo7P9rbW769taX54cCj18GK9Ri25ax3A\nYorHon9hmrmnMpmRZ2vdT2M8uiceiz4mImkR0eKx6N6jx05cYtt28kJiQ/2ql/HnN3w74rHovqnl\no2PZ/ZnMCInsChMJBfcYhm/3wGDq0Vr30xiP3haPRf9aRLIiIvGY7DvZe/p9o2PZZy4kNtS/etn+\niYisX9/6j15d7xKRARFpDIeDzxxP9rxnMePE0qmnbZpX1zcHg/5rw6HQnZrm3qGoilGpa26K7402\nRPbJ+HFeJBwO7jme7LnyQmKeTXpk9Akzl2tcqv6rbcVcAVnX2vKQ7vFcqahqxO1WnULBfjUSDt2V\nt6yDusdzUSBg3KAqil0o2D0iZ8+Q3BwIGDd6dX1LoVA4VC6X7an9ZLPmjybOKN/mN3zXenW9tVAo\n/HK2LFdVFG3t2jVP53P5Z48eO3HZ8HD6y2Om+ZTjOAOqqq6JRSOfEBEzGAh0uRSX17aLpya9Rpem\nad7K2Z/GeHSfbRdfzoyMPjVxZnGPqiieQsG+oKszWFz1NP5MM/fy4FDqgcqPR9P8mkfbcvp0/2+5\nFCUyw/g7e7VORNK2XRzw6vrmhobwxzVNa7Us6+fBgL8rHAreqiqKq1Agka4H4XBodzgc2qO4XCG/\n4WvNjIx+pzEevde2iwPicg2Egv5bvLremsvnf155znRXZqfrR9O0SCQc+kggYNzgN3yXz+cqbmMs\n+lDJcVLHjr+2ZSg1/PlQKLhb93haMyOjT03pb4dtF085jjOqaVpHZUzadvGUpmmNsWjkXhFx2XYx\n2RiP7tU0rcnr1XeEgv7dmqY1WJb1y6V+bzF/9bT9C4dD1zZEwp8ZHErdc6Kn9wMul2RCweDH7WLx\nxcZ4dM/5xLmQzwAuTL1t02Kxht/2eb03qqrapKpKJDMy+qhtFwc0TYu0rl3znbxlfe/Iq8lLVFWx\ng4HAbXax+KLjOPZi7GfL5XI2Eg7dZeZyCdPMPef16p2qqvhNM3cwHA7tDgX9u0UkEwmHfser6xdX\nPkeTt6mTX39p/3Lnb8VcAbEsKxkM+MUpOem8VeitlEcbIvfK+Nk4v4jsO9l7+urRsWxi4gzJlaaZ\n+6Fh+HbG49F7X301efnUfoIB/xXrWlsSIlIQEVtEIsGA/30nenpvmikW3avvUBWlUUTMTRs7fqKo\nin94OP3IwGDqUcPwdcZj0X0NDc4eVVFaB4dSd2ma1rh2TdOBya/h1T339/UPPljp06vr29e1tvyj\nbdvdqVT6C0vyJmLB6mn8TRYM+K8Nh4P3nuw9/Z68ZSUNw9c1dfz5x8fkY7ZdfElEPPFY9HOnzvS/\nL5MZeaalpanLq+v7nFIp2dLS/Dci4jme7Lmgs1JYPLZtDzglJ60qSquZyx2plBuG70rD8O2ceOjX\nNPfmvv7BB5qb4vdGGyIP2XaxW1GVSDwW3Xuy9/T1U/tRFUXb0NH2E1VROmT8LHJk4iruBtu2Z5xS\nZeZyB+Kx6CPNTY0PTbzulszIyBdURfFv6Gj7haooTRNjvjMeiz509NiJrZrm7hi/YpcaMM3ci5rm\nbpx4vM80c4mJq3n2xEsURMTv1T0P9vUP3r8EbykWoJ62fx7N3SEiks9bL0382z1RvuN84jTN3A8X\n8hnAham3bVpf/+AjIvJIYzy6d/LMAk1zbxYRbWwsmxARGR3LJqINEfFo7h225jYXYz+rae5zrnZE\nQsFbDcO3Y2Aw9XAkFNxtGL7d8ZjcL+NjNxJtiHziePK1yyrb1MmvvxjT0BabUusAFkvlEptt2y9P\nnmZimrmnXzl0tOFk7+n3iIh4vfqucDh0nVfXu0bHsg8PDKXuGhxK/Z6qKJuj0chtU/spOY49OJS6\n++ixExsOHzm2ruQ4ScPwdc0RTkRERPfq142OjT3jlJxsPBb9i3A4dG2lgW3bJ0709F6WSqUfbxwf\n1NnDR46te+XQ0YZMZvRxr653Vtpqmtaxfn3rPzol58hrr/VeX3Kc7GK9b1gcdTb+RGT8LGNLS/Nj\nppk7MHUKzOTx19AQ2WvbxSOnzvTd1Nc/cJOIpBtj0YdERE72nvmtkuMMrGttSaiK0nr6dN9N7Hzr\nh2nmErZtvywiMjD4+rzgkuP0Hz5yrPnwkWPrRCTr1fUuVVG0aEPkPtsuvnjqTN+HT5/uu0lEpLmp\n8XNT+1FU1T88nH74ZO/py185dLQpkxl9WMZ3+jtmi6dgF5Mlx8nqXk+n7vVsLjlOtmAX+3WvvlZV\nlFbTzCUGhlKfOp7suWxwKPUpTXMbs/U36fdJHz12ovXwkWMNtl18MdoQuVfTtMhC3zcsrnra/o2O\nZp8VEYnHo59rjEfvjMfHt2UiYpxPnAv9DODC1Ns2bRb+iX/TU/49u02rxn62r3/wplcOHW1IDaf3\naZp7Rzgc2j3d6y/wd1xSK+YKyEzMXK5bRKTkOGfPdlTOkAQD/vuCAf99lXLDdza7PsspOQO6rndt\n2hi9r+Q4townF/bUdlNkRUSGh9MPDAymHhkdyz7V3tb6E4/m3mnb9kERkbGx7FOmOR6bprk7TDPX\nXXKctIjIqTN9t0/uTNPclbg8iqoaJCDLR43Gn4iIhMOhD6uK0pEaTn94al1l/Hl1PaIqSquqKNLe\n1vpSpV7T3BEREdu2e8dGs4+Hw8H7bLuYGJ0424P6ZuWtI5O2E6aIiO7VO0QkomnunVP+1punPt8p\nlWzNrW1uaIncvU5RDNsumvN53eam+H7btl8+nuy5WkRkQ0fbD5ub4vsPHznWNDiUuiccCt3S3tb6\ngohkR8eyz6Qzo09pmrt1Hr/Pi5UdcmZk5Ol4LLpT09xrbdtOz/Vc1E4ttn95y0r29Q9+NB6PPhSP\nRR8xzdwBEdkpIjNOQZkuzoV+BrA0arVNm0UllkoiUjkhcrbfauxnU8PpAyIio2PZZ6INkX2qqmyW\nibE++TizHq34BGQ6BbuYFBEZHEp9dGAwtV9VFE336punW60lHo/uDQb8u48nezbYtj2woWP9T1RF\nmfXMm20Xj4iIrev6ZhGRytWMUsmZvAEsTGrfaxi+Tq+ut+YtqzfaELlTVRV/Jfu37WL3wFDqgbVr\nmr7V3BR/6GTv6Y9e8JuAmlnq8VcRDgd3i8jA6Fj24HRhiIjkLSstIum8ZXVXDhq9ut6hqIpW+X8g\n6L9TRLKa5u4Kh0Pvy2RGDizoF0dNWXkrKSJ23rKerdyQO/lvPdn4/Ong3X39gx9NDaf3r13T/FA4\nHLx3rtdQFUVzFDWsKop//LEaVhVFUxXFnzVz/5hKpR9XVNWIhIO3xWPRhyzL6s6auYMiIoqitE7E\ndPnUfnWvvlPTtEanVEqHQ6H3iIhdj3OaMbel3v6pihLJW1byeLLnMtu2exvj0fsMw7d7dDT7w/OJ\nc6GfAVRPNbZpM6kc5wUC/msHBlMPBwP+LhGRwvgUq4ol389GGyLXpYbTzwQD/utEREol58ik6sIM\nT6sLKy0BSRuG77rGePTOrJlLzNRobHTs2XxD+GA8Fn0sEPDv0TRtraoojSd6eq+eOMt2th9nImnY\n0NH2i4mn+2V8/uCMbNtOp4bTn4s2RPZt27rpFhGJlBznSCYz8pTu1d9wlmciuXhyQ0fbcZmYy5ca\nTu+b1N+RTGbkQCQUPBAM+G8LBvxPcCa6LtXF+Kvw6voV5sTB3WwGh1L3x2PRv9i0seNwySmlvbq+\nM5MZfdA0c/e3tDT9raoo/pO9p7taWpr/trkp/php5n7INKz6UTkr2N7W+q0TPb3vn6WdnRpO74s2\nRD63aWPHT0RENM3dmcmM3m+auQcn95MeGX1URKS5Kf5Yc1P8S/L62b1ZpYbTD0YbIg9t2byxb6LI\nnxpOf0r36u3tba0vlBzHtPLWS7pX3yEidtbMPTuxI89GGyL3RRsie0TkDau8qIrSuGlje69M3AOS\nGk7v4+pH3amL7Z+iqrKuteVvVEVpkvFpMY2mmfta/vUlyOcVZ+VqyPl+BnDh6mmbNpOJ47yHow2R\n+7Zt3TQsIpG8ZR3MZEaemW6a4EL2s/OJo7kp/nRzUzw98fovThxnLtlKXItpxayCJSJiWdbzBds+\nU7CLRy2r8EvHcbJZM/ecbReTIi6pPC4U7GOjI2NfLZZKJ4rF4lA+n3/21Jn+2y3LOjK1n+Hh9OPF\nUulYPp9/aXQs+0Q2a/5DPp//pVNyjhZLpRlXT8hmzefsYvGgZVm9o2PZr50+3fd74x8GlziOk56I\n69TE6/3SzOW+WSwWe8xc7kfpkdHPD6WGvyIi4nKJZHP55yzLetmyCv9aKpXSLpfiyuXzdXtZbbWq\np/GnKopWlrKWzeWfPfe7P944/kwz96KZyz3jOM5ooWD/cmAo9Yep4fRXvbq+3eUSb3pk9L+NjI79\nQy6fP1iaWP+8ULCPTP/KqLZc3jroOE7aLhZPmGbu4MQ246Bljd98O/lxNms+bxeLzxWL9qhlWT8d\nGEp9KjWcfmpqP6nh9Fcntl9JM5f7Xnpk9FHLsn7plJxTs63Al82aB81c7plisXjKzOWeHxhK/WE6\nPfKUbRcHMiOjT0i5fMYuFgfy+fxzA0OpT5pm7qeO4+TNXO7AxHN+lB4Z/axlWcnKZ6cxHt2Xt6wD\nQ6n0l/P5/M/TI6MPD6WG/3u13l/MT71s/xzHyY+Mjn3VcZweM5d7KT0y+uXJC7rMN86xsez3FvIZ\nwIWrp23aJK5isZg0zfxzjuOYIiLZrPm9Scd5T/T3D35yfIW2xdvPnjM+XS6XZVndppk7GAmHdmua\ntuPUmf53VY4z+/sH/6DkOPnpXh8AAJyHbVs3ldvbWp+sdRwAUC/a21qf3LZ1U7nWcVyIlTYFq2ra\n21ofMwzf+6arGx3LPn6y9zRLRGLJMP5QC4bhu6K9rXWmucn9rxw6eklVA8KqtK615XPBgH/PdHWm\nmTtwoqf39unqgKnYpgEAgDdojEfvnLyEOQCsduFw6NrGePTOWscBAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAoLZctQ5gITo6OrSOjg7P5LJ0Ol3q7u7O1yomYDF0dXX5RESZrq67\nuzufTqdLVQ4JAABgUbkX+sTGeHSviMjAYOqBxQtnfm697pp/90ddW39jctlzA/Yvrrmr+2+rHQtq\n43zH36f/5LNvyrWErpyrXUwPjDz5J3/6zVols9/d+7uflIHX/NPVvfNReSSRSPRVOyYAAIDFtOAE\npKYsUy/3HGo/tyz8Wo2iwTLQdvHmhh9Hi1fN1a4l3HwiEokcqEZM0xp4rbnccyg2Q225qrEAqEuG\n4btYRMQ0c7+s1mve+v7rt3dEA01nC3Sj8MCjXzlYrdcHsLIszwQEAIBVym/4flNExDRzVZuBcOvO\nzbveoWfOnsRxtW0deuBRIQEBsCDTzjUHAAAAgKXAFRAAAADMS2M8eo+IyMBg6ktztf3YH3zibds+\ncsMt8+m3uWfkyY9cf+PzFxoflgcSEAA4T5qmNYiI2LY9vNSvFYlE3J2dndEpxeVEIjGw1K8NANMI\nzbehWSzohzN9M93XeI5g0c2snFWEBAQAzlMkHPxPItVZBbCzszP6TzdsfGRymatt65CaSHx8qV8b\nAIClQLYJAAAAoGpIQAAAAABUDQkIAAAAgKohAQEAAABQNSQgAAAAAKqGBAQAAABA1ZCAAAAAAKga\nEhAAAAAAVUMCAgAAAKBq+CZ04Fzud37wvb/S0rnNmavhq//vp0de+MHB3moEBQAAsFKQgADn0nxX\nXfobwTc1Nc3V8CKRx0hAAAAAzg9TsAAAAABUDQkIAAAAgKohAQEAAABQNSQgAAAAAKqGBAQAAABA\n1ZCAAAAAAKgaluEFAADny/e9B++5Y3LBT4dyh+754l99r1YBAVg+VkwCsmvHxreWnvz8n51T2Lg+\no179oU/XKCRgUd16662NXV1d09bt379/NJlM5qobEYBVTHuHnrnqnJJYWESEBATAnFZMAiKWqZV7\nDhmTi1y1igWYRen7/+sBGXgtOkN1aKbn/ba//1Pl4GvadHWJjo67SUAAAMBysHISEGC5GHgtWu45\nFJuuytW21a52OAAAANXETegAAAAAqoYEBAAAAEDVMAULAAAAdeHWD93Udtt112ycT9svf/P//vTA\ngQPppY4Ji48EBACWH6Wrq6tpSpmTSCQGaxEMACyWDkO56Kr+F26fT9vnOjvvPnDgwFKHhCVAAgIA\ny0/gn27Y+MjkAlfb1iE1kfh4rQICAGC+uAcEAAAAQNVwBQQAAAA11bRu7doP3X1H1+YdF3WahRNT\np5iew1ssjChFO1+t2LD4SEAAAABQU+7G8MbgjW+/YbjkNPUcPb1ptrZtIkcNEpBljSlYAAAAAKqG\nBAQAAABA1TAFCwDqRFdXV9N373z3Z84pbFxfLv/bszWKCACAxUcCAgB1pNxzKDb5sUvErlUsAAAs\nhaolIDfc/MHtgTWNs65qICIydmbgzHe+/o1fViMmAAAAANVVtQRky7vfsWt0Q+yqudqtPT70vJCA\nAAAAACsSN6EDAAAAqBruAQEAAGd1ve3yll3bN26dXNbe0tQkqUytQgKwwpCAAACAs3Zt37j1M9uN\n2yeXufyKXU7VKiIAKw1TsAAAAABUDQkIAAAAgKohAQEAAABQNSQgAAAAAKqGm9ABAFihPv0nn700\n1xJ6+3zaHv4///z8d77+jV8sdUwAQAICAMAK1Xbx5uiPo8U5vwRYRCSwpvGVpY4HAERIQACgrpg+\nf3NOdYcqj8seb8kVCKuT21TKYrlsUkpFu/pRAgCwcCQgAFBHcqo7NKgojZXHTtlxFEU55369SllM\nUXqkJCQgAIBlpe4SkC2XbG/5bwf+btdsbXb41Tdbh59r1PPZgWrFBQAAAODC1V0CMuqW9sPR4u9W\nHl+keiJbMgONk9sUXD69X/fqbSQgAAAAq1J7e3tDV1fXvNomEon+pY0G56PuEpCp/Lbl2dL76jkJ\niBNf54gw7QAAAKDeaYqqbBFXVC3aM379Q9PIsPkmp9zU6lJCM7WZ6tat0Q9+5IaNO+Zq52rbOqQm\nEh+fb79YenWfgAAAAGD5ciuq8paew63RkZR/pjbxoz9rDJRL4YmTzFjh+CJCAAAAAFVDAgIAAACg\nakhAAAAAAFQNCQgAAACAqiEBAQAAAFA1rIIFAADkive/+x1rL9u2bXtDcMvQ6Z9tnlxX9nhLrkBY\nrTwOe3y9bhGz+lECWAlIQAAAgFgR46LRDbGL+0tO02Cfcu73b5UdR1GUs7MmDEUZJAEBsFAkIACw\nxLq6upqv//KnPz1Xu6BH9+V/9ExUioVqhAUAQE2QgADA0isfzvTF5mrUoBtep1x2VSMgAABqhQQE\nWAE6OzujM9V1d3en0+k0p9QBLKlIJKJ1dXWd803X3d3d+XQ6XapVTADqEwkIsAJ88YqWPyy3jWjT\n1f2ayN2JRKK/2jEBWF3e3Bq95Lt3vvucqYbvfFQeSSQSfbWKCUB9IgEBgJUh4Bz89h3nlIRieeWS\nt/91jeLBamOZernnUPuU0nJNYgFQ10hAAGBlUJx/+fZVkwtcbVuHRIQEBMCKktN9TTm3OzxWKgQC\ngXDLbG1juWyySmHhPJCAAAAAYNkwFSVglpXQmfyYZ42ihGdrG1OUnmrFhfnjm9ABAAAAVM1KvwIS\nKj35+T87p6RxfVq9+kOfqVE8AAAAwKq20hMQKfccOmftfRbYBwAAAGqHKVgAAAAAqoYEBAAAAEDV\nkIAAAAAAqBoSEAAAAABVs+JvQgeWyg03f/Dirq6uOb/l13Br5Y9cf+M/VyMmrC4jRrDdLhUdEZGy\nx1tyBcLq5HrbrbV88ok/vWNLuPnY777vN/+hNlECAHAuEhBggU75ym87HC2+Y652W8KxIREhAcGi\nyyhqg1l2FBERp+w4iqKcc1XbdIqh0Q3NV61t31abAAGgPuh7P3Zz13waPvfysVcT//ICX164xEhA\ngCUQ8wa8IuWyiEjI4/Xt27evaVK1p0ZhAQCwGumf2W7cPr+mGx8jAVl6JCDAErjxtUM7jOE+XUTE\ncClWW3Dokddrt9q1igsAAKDWuAkdAAAAQNWQgAAAAACoGhIQAAAAAFVDAgIAAACgargJHQAAnB+X\n4hHN4z2nTNVU0Tza5KLOzs7GZDKZTSaT2arGB6CukYAAAIDz0i9Oe97nP+cYwlFcjuLznzOzQv/Q\n2//L28T60+Qjf5moaoAA6tryTUBcimp6/a9/t4LbUxKv/5xvARa3p6R6DZdu5YakXJ7zG6sBAAAA\nLK1lnIC4lB6PZ1PloeMSR/F4zjnz4rjECWieXFvBSku5VKx+kAAAAAAm4yZ0AAAAAFVDAgIAAACg\napbvFCwAWMZ01a0GNP2cbXDY49NrFQ+wFLaUXbFLd1zUueVjN48X6Ib1wKNf+ZfaRoVVRfXoomr2\n1BXapmUX8lWICEICAqx4nZ2dsZnquru70+l0ulDNeDDOcHu03/rp82+ZXObE1zm1igdYCpefPt7W\nXjhxvbHduFxExNW2deiBR4UEBFXT4zO2jSkumbpC21SGS7HaMkM/3vWWS7bslZvn1fcT3/3RC8lk\n0lyUQFcZEhBghfviFS3/X7ltZNozP78mcncikeivdkyoDreiqJeWnKamzNDmvR+7uUtERHSj+MCj\nX/lBbSMDgPq0q1G78h3bja75tH3u1Y6XSUAWhgQEAFYot0tVrjnavSl++N9Cn9lurBE5ewaaBAQA\nUDMkIACw0rlc6tlvrVY1X1dXV9NMTfnWagDAUluNCYj/ew/ec8c5JeHYyO8//o1vdnd3c/MRgJVH\ncbsP+fxvERExy0Vryx/veWSmpi3/+weP8a3Vq4vPrblFXKKr7rLP7XF7xFbnfhYALNxqTEDUd+iZ\nqyYXuBrXnIhEIgdqFRBWPNfZs88iIqqmzrwah6tcpZgAQEREfiOd2uoZHXKHtNecTjuv6MHYajw2\nAFBFbGSApaaoZ88+i4g4isuZaTWO1rJzOCAyWL3gAKx2IXPEa4ykdEPVHHfJVhyPwWpsAJYUCQgA\nAABwnr503++/Lf2xm9Nztes+dnL4nr2f/Vk1YlouSEAAAMDSUVXtPBZBGEkmk9yPiWXhTelXbiz3\nH5rzCw7LWvh5ESEBmaSuEpD1ihraaJmaXnLO3gAXU92+WsYEAAAWbsjjbetR1XYRFkEAMK6uEpAd\nI6nmzjPdsa0l++z8+DMdl2QurFeXSzSPLk5p/Hed7gZgVfN3dnY2JpPJEmdeAAAAVhlFcQ8FwpvL\nHm/JFQjPuhKcz3HGDHP0dLVCW4nqKgFZEoqi9HiNrWbZ8YpMfwNwqpjbWr7pyjW3RSK/t2/fPhIQ\nAKvWLffc0XXLPXdsm0/bL976n/4ukUikljomrBxeVdP29Pe+dXLZD5rXn3ilXGLxDdSYyzWoKI1O\n2XEURZl2oZiKuEjZqFZYK9TKT0CAJXLNQO9lu9J9+nR1XlWbc05oPfjen3/2wzIyNG3Sfc19X3oq\nkUhwULDCTHcAmA802F81jJ+JiBxK9284nOnbMs/uvrnoAWLFM4bP3W66Yy2KqLMe7wHL2q4dG99a\nevLzfzaftn+cOPT1Bx79yg+XOqZaIwEREdWlKNtEib3r1664YtdFfzVcKb/ni3/14+7u7rFaxob6\npY+mNBmePgFx4usWtIzliBFsHw7H2qftU/c5SjimiIjoLpfdlB684Bvayi//sLPcM+MNdBxcrlBT\nDwDHCzmfBwBLwjK1cs+h+W1kLXNZnMC8UCQgIqK6FNdlvUdb159OfmS9bVmV8kgkcreIkICgamwp\nq/myM+3n0nFKjlJ2Jk4TcrYQtdfZ2WmIiH8eTZ1EIpFb6niw/F3x/ndfvfaybfOaAjj0/X97av/+\n/VylRfUpqtv0+pvE7SmJ1z/r/SKaU8p7RIbm23VHR0dktpXiJkeRSCTOzLffekMCMknGH9owUnbO\n7kw/9JcfffD6gvWGneYrX/3O3/2PL3z5B9WNDphk4ma5ysPZbpoLieu0RiKNWWiK6n73WGaziEjc\nLpYusrKq7Q0U/8mtHp/teW/7+G/dtW7s+uJc/V/RvOGlK9Zc9DeLFS9WlktKpea23FhIRCSulzY2\nrgs7B72+nnzJLs32vI7jHVylRW19GcZoAAAgAElEQVS4FLXH49nkuMRRPJ7Z7xdxnIH4eSQgv721\n4f0fMTbeNGcIbVuH1ETi4/Ptt95ULQGJ6YFMS7j5xGxt4h6/xyvurMspnt2hNTSssf1WbmByOyfc\n6Pbqfpnczgk3uhVNL76hneYtiIi4pKzM1E4JRBwREZcRLrnKztlLX9tEiVuOY05u+6rudV1/8wev\nufztV1w8qTgnIpacyx49MzD291//xpHZfucJqojMuqGdZERE7Hm0c4lIea5G6XS62N3dPTX2Fcdw\na/aWcGzW8Sci0haIVA7U5zxba7g8a1yaPu1Uq8njTFfc7rnGaoUn2uJVXcq092RM7TM3tU+3Nm2f\nrujacNkfmjY5KUea3C5/uCgioomo4Wymr1L3vT//nftlZGjaPn//ye9/s7u7u3+6OhFJyeufB7eI\nTDtNbapkMplOJpOF+bRdhlxbpmz/gh7d72+5KDS57Hy2bYqmFyePrdnaaS7F9pedgdnaXZYZGB9b\nluVuc6RoqR5lrGn9yJTfQ0tZ5uRtVUREHBn/m88m/N8O/N175mgjIiKv/r+fDrzwg4PzPatny/y2\nnY7M4zMtIpJMJgvJZHI+29hlYT7bvmBZiXu9AaUynmYbS5XH022rZtwPTxrT07Vp1/1uJW9OjL/8\n/8/e/Ye5cdeHvv9opNFIo9VKq5W0tjfrlR07xnZ62UANJSlYSU9oSaCkFKchlCaUEwpOzyW0DadN\nUuz0NjyHhJZwTknb5AHiUnLa5PYkPIVc2hSyLo2hpC1L28Q4dmJt1mt7f2ml1Wqk0Uij+8euzHqt\n3ZW9u6P98X49T55YM19956PVR9/RZ+Y7I0+rTy9Y8e68VamcGw9KdsU3aZmumc97z13/5QNv/vD7\nGpquMvqdf/3H/v7+kYVbivT19ZUymUyjOdDQvnZaQaZysZEYjEwm02i/K5bu8RavCHX0a26PNxQz\nWwKBNu9cbWvj1Hz7yJpartb7jlivbTXYJo322cj2LyZWt8ttG6GoX3yBBccqrWyV3W0bSq4Lv1Ne\nKLZZff6rf/ahBduJuFJ5u//QXz052EDbal9f30Amk1n2fbFr4Sb1qarqExGxLKuhu0b19PT4wuHw\nvKepRMQrU19YZn7o6n24a3E30q6RZXOtW8y2azvGRgaQixnAGm27pguQS8g/bzgcbmRHVe/9nct8\nv1Ez8+/fSF7VKDL3Dmol9TnfOnPG9hrO7VQqVUylUo0W4k11sfkXDofdPT09vlmLFbmwOLvYsW2+\n9+9S2s31nLliqT1e6Mv9coxxF9P2YvJwxRcgF5N/DY59tbHsYt77euNKI/vMRvNvtrn2swuegVug\n38W2vRgN97uSC5CLyb9EIuFJJBKaTL322WPgbPPl32LbyjL0udTbv9i28+3b6/XbSJ/Vvr6+8Uwm\nc0nXsQIAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABw1iX/EOHFCoVaf8Wret6w\ncEs4qWSV/zmbnfhWs+NYbrruTwZ0/95mx4G5ZbK5Q5ZlpZodx3JQVTURDgVva3YcWFjeKHzLMAr/\n3Ow4lkoo1PoLXtXz1mbHgYszMpq+v9kxLAVd978hoPt/pdlxoHGZbO7PLMsaWu7teJZ7AzU+zdvi\n07Q2p7aHxiiKqWebHQSwxrkVxaP7/Yx/q0DJKvsW/mH31cOneXX2vWgWVVW9jH2rSy6X91jNDgIA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAKw1oVDr9bFoZH/t8c4d24a7uzqfX2w/i6Gqaqy7q/N5\nn6YlZq8LtgSu7+7qfHoptjOTrvt37dyxrRqLRg6IiHR3dT6/c8e24aXeDs63WvLPrSiBSFv49lg0\nciDSFr5DVdXwUmyrhvxrjlWWf3fU8s+tKIGl2FYN+dccqyX/aiJt4Tu6uzqf74hHDyzFtmrIv+bS\ndf+uWDSyX1XVThGRWDRyYOeObVVd9+9aTD+LtWlDxxdDodbba4/diqLGopEDM/8LhVr3LcW2amZ+\nBi/177DSeJodwFzCrcE7dN2/b2Q0/Uiz+/Fp2vZgMHB9qLV1v6p6dituRa+tC7YErtd1/9WhUOtd\nbkWxFhNrIzITuUNGoRBb7u2sd6sl/zZv7nzep2l7RGRERGLRaOTAydTAWy3LGlxM3HMh/5yxGvLP\nrSjqlsTmf1ZVz275Sf7d/eqrqasqtp1fTNxzIf+csRryr0ZV1VhHPPp5EQnIVB4uG/LPWQHdn4y2\nR76YNwovL2afthT9uBUl0BJsuSHcGtyn6/59Vtk6WFun+bTt0fbIwZntDaPwVDY78dSlxjyfvFE4\nLJI+aFnlZc335bYiC5BQqHWf5tOuFhHp7up8tH9g8KMz1t3gVT17KhX7eHo880Rtua779wR0//Ui\nYmWyuScsyxqs14+qquFgS2Cf261sEhGj1nb+eII3+jTtvTI1wJ0n0hbeLyIXHHVWVTURDgVvyxuF\nZ32atr1omscNo/DizDjzRuGbhlF4eXqAvbVklY9nsxNPBFsCSZ9P21ssmocrtj3n0ZZQqHWfV/Xs\nyhuFrwd0/w2Vij2SzU48UbHt/FzbX/CPj1WTf7ru3+PTtD3p8cxnhoZH7420hW/tiEe/Fg4Fb81k\nc0+Rf6vTask/zaftVlVPfHQsfefIaPqRjnj0QKQtfFDzaW+0rPJp8m91Wi35V7NpQ/zh3GT+iWBL\n4I7aMva/q59P0zpDra13iIh0xKMHzpwZPvdeqKraGYv6bxQRmZlDqqrGwqHgB0UkVLLKR7LZiefq\n9VM0zZdDodabvKrnjTKVC88t9P6oqhoLtwb3K+4Lz/L6NG23iOSPHjvRMntdLBq5u2SVRyzLeklV\n1d3Z7MTj9eKstRURPZ3OfE5EJBIJ/06lYmfT45mH54krEQ4FbysWzcOqqu5yu5VwI69nJViRBYhl\nWSN2xc64FaXTKBSO15bruv8aXffvmX4YUFXP9qHh0fs74tG7I23hBy2r3Ke4lXC0PXLg1OCZd8/u\nZ+qIXdcP3YqSkKkjJeFoe+TAidf6t1iWNWclOTQ8+rCIPByLRg7MrnL7BwZ/SWTq1Kyu+3fXlquq\nJxFtjxxsa7PvcCtK5+hY+s6A7u+JtkcetazySyLijbZHHjh9dvimbHbi2Y0b40mfph20K5XUxo0d\nfyki3pOpgUdU1XPe0ZZwa/A2XffvHhlNP1SrxKPtcq+I5EUkHGkLf+Jk6vWr6m1/NSTkSrBa8s8s\nmsf7Bwavtazy8em4U7WXQP6tXqsl/wyj0Hf02In4dGy7dN1/TcW2U2bR/JHm0/aQf6vTask/EZFg\nSyCp+bTrT6YGdp9fgDD+rXYV2y5ZlnVaVT09ZrGUqti2UVu3aUP872X6rGtbW/juV19NbVFVNbZ5\nc+cPRKb2hT5NOxjw+x8ZGUsfnN3P1BSq4H4RyYiIGm2XB06fHb4xm514dq54iqaZ6h8YvFbX/bu6\nuzpfmrlOVT27RMS6YvvWk25F6TSMwnNDw6MfLZrmYLQ9cnfFtsWtKGHDKLxgFs3eenGePjt0Z6Vi\nj3TEo19RPVNTqUOh4F1Dw6MfnB1LQPfvjbZHDuaNwlMiEp/+XFgiUhIRb7RdDp4aPHNTbjI/5+tZ\nCZRmB1CPYRR6Lct6WURkZDT9UG15xbaHXzn+Wscrx1+7TETyPk1LuhVFjbSF77Gs8ounzw598MyZ\noZtFRDrisQdm96O43YHx8cxDpwbPvOXosRPxbDb3kEwNpLvrxbEULMvq7x8YvCqdzjzW1hY+YFnl\n46fPDt08NDxys4hkYu2RB0VETg2e/dWKbY9c1rmx160onWfODN0836A809Dw6M1Hj51oS49nDqqq\nZ/fMuYczt79ML3HNWS35V7HtzPQ2Bn2aluiIxx61rPJL6XTm8Vob8m/1WS35V+PTtHB3V+cLPk27\n3iyavTK1ExQR8m81Wk351xGPfX58PHNwrlwh/1Yvy7JGjELhmyIimYmJQzMOsMnQ8OiHjx47Ec9N\n5h9zK0pM82nd0WjkHreixIaGR28eGh693TAKT4VCwf1uRQnM7scqWy8MDY9+8OixE22nzw7/koiI\nV/XsqRtIA6bPgKjZ7MRj2WzuYV3339ARj/5Jbb1bUfRTg2feefrs8M1zxenTtER6PPP49OM7Q6Hg\nnYZReGrmmcb5GEah95Xjr7WdTA1sEZFSRzx28FJfj1NW5BmQuZhF8/iMucWGiIjm0xIiElZVz56Z\nVamqerbPfr5dqViqR93etjF812WKoltW2ZjdZqlNTuafMoxCn0/Twm5F6XQrisyKMywiYlnW4GQu\n/1goFLzHssq9ucl8b6PbSI9nnhERyU3mn420hQ+63cp2ETk9c/tL+6rWp5Waf7FoZH+0PfJgbjJ/\n6MyZoU/NnH9P/q0dKzX/iqaZOXrsRFukLXx7Rzz6lVCo9cWiab4sQv6tJSst/yJt4VuntxOrXSSu\nququUKj1BsuyDBHyb60qmuYPRERM0xwMtkzNiHIrynYRkU0b4ufdLEHX/T2zn1+yyplIW+jejnj0\n4YptZxYbT//A4M0zH2s+7zW67r+m9tgwCi/Wcmq+OIummRoaHr13S6Jrn4jI0PDopxqNwSgUnqvY\ntlUxzcHp6YbLdmB9qayqAqQes2imRMQqmuZzJ1MDN4qI+DQtobgVdXbbUKh13/QprQ+nxzOPb9rQ\n8WAoFLx7mUMsiUztpEUkUzTNvpOpgWtnx+nTtERLMLBfRPKq6kmGQq03ZbMTzzSygUhb+Ib0eObZ\nYEvgBhGRSsU+PmN1aY6nYQk0O/82beh4NBQK3n767PDNc+QL+beGNTP/gi2BqyNt4Qdyk/nHZh6l\nc58/R5r8W8OamX9ut+IVESPaHrmztkxVPVcEWwLXpMczz00vIv/WiYptHxeRa145/lq8YtsZVVXD\nqurZZBbN/kgkfNvMtps2xL9mWeXjR4+diIdCrTds2hD/5mK2HYtG7qlU7Ex6PPOIqqoxVVW7Latc\n9/qh+eIUEemIRx+UqelU0hGPPji7uJmL7vdf71aUh1VVjeu6f49llV9ezGtywoqcgiUiUjvSstCt\nbSu2baXHMwd9mnbDtq2JH27bmvjhlkTXyXDr1GnQmf1UbHtQRKQjHn10545t4w4UH+cZHUvf69O0\n5LatiVe2JLp+MB3n7SIiGzfGv+ZWlMCpwTPvrNh2qiMefVRV1YbuttERj35z545t45G28MGiab64\nXHdeWE9WQ/75NG1XKBS8Q0TUTRviT+/csa06fcvIure9JP9Wj9WQf4ZR+JHm0zo74tGv7dyxbbgj\nHv1KxbYHM9lc3SkD5N/qsRryb2Q0/fjRYyfitf9ERAyj8MypwTP31mtP/q0+lcrUdR+bNnQ8utAt\ndEdH0w9XbDt/+eWJV7q7Op/ftrX7VEc8+ni9fiq2PaKqnj07d2wbXmzxISKi+/27O+LRL+7csW18\n29buk25FiY+MpeuevZgvzlg0sl/X/TelxzMPpMczD+m6f1+jt7HWdf/1V2zfOr4l0XVSRLxDwyN1\nPwcribvZAcylUDSP2LadscrlfsMoHHG5RPKF4hHTNPtERGY+zueN71rl8uFy2cqZpvmjkbH0p9Lj\nmadm95Mez3zVKpePmKaZMgqF72Qmco+Ypvlju2KfLpWs1xoIy1Uul1OGUTxsz7ggamqNy2WaZp9h\nFI5MLxDbtjN5o3DYssqnRaZOwxmFwrO2bedKJevHI2Pp30uPZ77q07RdLpf4MhO5P5/ITf5doVg8\nUqlUMiIipZJ13Lbt/HQ/qZnbCYda96mquvv02eF3mqY5mJvMPzE8PPo7Fdsu1ts+Grca8s+reWPl\ncnnEKBR6Z/5Xe8/Jv9VrNeRftVq1Mpnsl8uVSqpYLP44N5l/Ymh49Den5s4z/q1mqyH/Llg5FdNh\n0zRfJv/WhlKp9FK5UhkoWaVU0Sz9W8W2c7Ny4FxOlEql/onc5F9KtZqzyuXBzETuT0dGxn63YtvF\n2f1kMhOHbNseMwqF72cmcn9smubxSsU+XSgW/23hqFxyXk6ISHYi93+scvlF0zQHjELhH4aGRz+e\nzxsvipz/WRERKVcqI3PFGdD9SaNQeGFsbPxzhlH4blWqJRFx1fkMnnvdbrcSD4dab89mcwcn8/l/\nMAqF746MpT+VzxtHlut9AaS7q/PJnTu2VZsdB9Yn8g/NRP6hmcg/rAS67k/ON/NhJVv114AsBV33\nX93d1TnXfM/ho8dOXOloQFhXyD80E/mHZiL/sFLs3LHtP0UkXm9d/8DgTT+Z4QKsc6FQ6/WrserF\n2kD+oZnIPzQT+YeVYOpHGSP7dd2/q9mxAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAHCeq9kBAM2QTCYDMx9nMhmrr6+v1Kx4AGA5JBIJNZFIeOutS6VSpVQqZTkdEwB4mh0AcCli\n0cgBEZGR0fT9F/vcZDIZ+PiXPvexgclMtLasY2Dibz/07vcdWcoYAWA5XMz493Pvf+9P7/zQe36x\n3rqjX/3br3/pc1/4/lLHBwALoQDBujQwmYm+kh3qrj0Olj1qM+MBgOVglEvazLFu1jqf0/EAgAgF\nCAAAABoUi0Y+KSIyMpr+/GL6SSaT+se/9LnfGpjMdNRbf/Srf/s1ztCtXRQgAAAAaFTrEvXjGpjM\nBF7JDrXXW8kZurVNaXYAAAAAANYPChAAAAAAjrnkKVixaORjIiIjo+k/u9Q+kslk/N1fuO/3Zy/v\naglnf+Xyn77vUvsFAAAAsDIt5hqQuhcNXay55v4BAIDlk0wmQzvCHfF66x5//PFcKpUqOB0TgPWB\ni9ABAFiHclva33c6sv1X6q1L9CbuogDBYr3nlvfvvuJd73hHvXWXd3RahojmdExYGShAAAAAsORa\nNsRiuS3tb6+3zgiFhp2OBysHBQgAXKRQqPUXRESy2YlvXWof77nl/W+44l3vuHb28natZfKeW379\nq4uJDwCAlYwCBAAuklf1vHX6n5dcgLRsiG2od2RwY6hjTEQoQAAAaxYFCNa899zy/iuveNc7zn3R\nu7yjs8S8UwAAgOagAMGa17IhFp15pJl5pwAAAM3DDxECAAAAcAwFCAAAAADHUIAAAAAAcAwFCAAA\nAADHUIAAAAAAcAwFCAAAAADHUIAAAAAAcAwFCAAAAADHUIAAAAAAcAwFCAAAAADHUIAAAAAAcAwF\nCAAAAADHUIAAAAAAcAwFCAAAAADHUIAAAAAAcAwFCAAAAADHUIAAAAAAcIyn2QEAAAAA7YrHv8kq\nBUVE2nZf3nPFR26Zs23KsF899L+fHHAsOCwpChAAAAA03SarFLzuRN82EZGuUund+i79LXO1/W78\nLX9OAbJ6MQULAAAAgGMoQAAAAAA4xpEpWIlEwpNIJLTZy3t6evyGEwEAAIB5dSue0K6JsbiIyMdu\n+S+3Z955VW6utodePH740NPfeNm56ACsJY4UID/3/vf+9M4Pvee9s5d3tYTd3x485kQIAABgHq1W\nSbti8NWYiMiOUv5nRCsV52p7ONJyVEQoQLBsyqrWOtdBarVqlxwNBkvOkQLEKJd8r2SHuuusspzY\nPgAAAFaPrMcTNVzVeL11UdsecToeLK11fxesRCLhSyQSrQs0q/T29o45EhCANaWnp6clHA7rs5cn\nEolWBhUAwHq07guQ22+/vfW+4OsPz9dGedt7v6tc3funTsUEYO34+MHfe/O/Rsq/MXt5e6jDGssO\nNSMkAACairtgAQAAAHAMBQgAAAAAx1CAAAAAAHAMBQgAAAAAx6z7i9ABEZGuN2zb+duH/ueu2uN2\nrSX75P/4n3/T19dnNjMuAACAtYYCBBCRil9tz21p3117vDHU0R8Oh78uIhQgAAAAS4gpWAAAAAAc\nQwECAAAAwDErdQqW+oG7Pp6st+J7z/x/30+lUkWH4wGAZdeuePxbi0bngY/cklyobSo9OXLo6W+8\n5EBYAAAsqZVagPiD7/vZj9Zbkeg7+jIFCIC1aJNVCl7143/p+sVdet3xb6Z/NDd+99DTQgECAFh1\nVmoBAgAAANQVCoeDyWQy3kjbVCo1mUqljOWOCY2jAAGAlcTlconq9dVdV7ZMqVarDkcEACtOT7v/\njf/wnq23NtL2/3l5w6P3fynVu8wh4SJQgADASuJSPMf8gTfVW7WjkP83sUpMQQUArGrcBQsAAACA\nYzgDAgAAzudyu0Vx1/+OYFfKDkcDYI2hAAEAAOcZ1lu2mxWr7vVGXUXjmAjX8wK4dBQgAADgPKbL\n5TNcLqZpA1gWDC4AAAAAHEMBAgAAAMAxTMECAADAsgh5/do16aGu2cvbS9YGEZFtpUJLbVnY6/c7\nGRuahwKkMd6Ffm0zk8lYfX19404FhLUhmUzq3z7wG78lI693zNfuD3qP/dX9j3z5BafiAgBgKSgu\nl2v74IkLvkPpLqVVRCRStc/98Kodvcx2MjY0DwVIY1r+4T1bH56vwT+aoe9e19f3p04FhJXvtg/c\nvPn2G67bMl+bcMcmj4wcbasOHGuftzPTUJc0OAAAgCahAAGWSUJXtr59+Acfna+NS9sx7FQ8ANYf\nn0d1/8rkxO5WI+ebva4lPWrvLZvnrgVVAmGuCwXgiDVdgITDYV9PT0/rfG26u7vbJf26UyEBAOAY\nRVzizaU9+kRam73O51Ztu2KdKzpst9rw9JcDH/nAtQc+8oGdCzZsbTeVK3/2yw0HDGBdWNMFSE9P\nT+tCU6eUHZGX7O85FREAAKtf9dSPE9WBY9sXaufq2jHmRDwAVhdOtwIAAABwDAUIAAAAAMdQgAAA\nAABwDAUIAAAAAMdQgAAAAABwDAUIAAAAAMdQgAAAAABwzJr+HRAAWKvC4bAnmUwGGmmbyWSsvr6+\n0nLHBABAIyhAgPq0m3/3v33gzR9+n1pb4D8z8cIf/u59/97MoICaN3ZGfurb+991XyNtD49Y/3Hd\nnX1PLHdMAOAEy6O2nKyU3trSEto4e52nWi2FjMkBqdp2M2JDYyhAgPoUY2PoypxejNcWXBHqeLmZ\nAQHnMQ2tOnCsu7G2odeXORoAcIzlUrRUcbJ7g6KEZq/Tq9V8yOUalKpQgKxgXAMCAAAAwDEUIAAA\nAAAcQwECAAAAwDGr7hqQm3/3v9365g+/z5y9/Nuf//Jf9fX1jTcjJgCYS5um++4YPvrm2ctb0qP2\n3rJ53kEgd0uYg0IAsMT2vunKNxyQWxpqe+jb3/+XVCo1ucwhrXurrgAxNoauyulFdfbycDj8N82I\nBwAWoo8PabOX+dyqbVes8woO263Oe9HkcDBypVmxqiIitua3lVC7IiLSZUy+LJZZWMqYAWCt2BtT\n3/aOXfo7Gml7+NXEyxQgy2/VFSAAsF6ZVdtjVG1FRMS2K7Yy/W8R7vYCAFg9KEAAAEDDhoORKwua\nv1I7A1ejuRQrnhnht5IALIgCBAAANMys2h7DrnhmnIEDgIvC4AEAAADAMZwBAVaB237ube/Ye/nG\nnfM2CrUXrrvzvsediQgAAODSUIAAq0C3p3j5Zi37hvnauGIbxkTkcWciAgAAuDRMwQIAAADgGAoQ\nAAAAAI6hAAEAAADgGK4BAS7BwYMHN8gCP/629/INIRn+gUMRAcCKFKw8+dn/1UjDwyPWi9fded9f\nLHdAAJqPAgS4BL+/y39vdeBY+3xtXJrfqjoVEACsTK6FxspzzFBgmWMBsEIwBQsAAACAYzgDskS6\nN8bjBz5yS3K+Nqn05Mihp7/xkkMhYYklEonQwYMH4yIi4lb9omo+KZdMqVY50YEVLRwOB5LJZLyR\ntplMptjX1zex3DEBwLJwuRTTF2iv2JWf7Js9Wll8AY+IiK9iTSpWyWhafBARCpAlkwgoW39/l75j\nvjb/aG787qGnhQJklSq1+X/q9HXbbxER+bfhwTdHff5yV778klQr5WbHBsznjVHfT/3De7Y+3Ejb\n78bf8ufXfuhjh5c7JgBYFi7FNawqnUbV7astsl1VW/F6FRGRrpKc0ClAmo4pWAAAAAAcs2bOgNx2\n222xZDJ53rLu7u6oZI54xOYINYC1Kxto3TJRtS+4gDduFvu1Yn60GTEBADCXNVOAjL2p61OnL/eq\nM5cVNN1n/IuvVS/mh5sVFwAst3K16jaqtjZ7eYWz3ACAFWjNFCAAAKCJXIrb8AXOv9mBx1sRX8Ct\nVu2SahYyTYoMwApDAQIAABbP5VIGvN5tMxfZLrEVr1eJ2vZIOwUIVoE//q/vvzXzzqvMBRuG2ovX\n3XnfVxwIaU1atQWI4nK5NLfqrj3W3J6q3+M97/X4ZqwHAAAA5vNGb/6qqpZVF2rnim0YExEKkEu0\naguQsFfXbky9dO62t63q63aPVTxvvrMS7nCJyGnHgwOaI1h58rP/a6FGf9B77Gv3P/Ll7zsREAAA\nwGyrtgCpSlUiE+lzd33R3artqVjnFSC2V7edjwxoGld14Fj7gq1Mw7dgG6xb3TN/cHMBvb29hd7e\n3txyx4T14Y07r9j4/Ff/bG8jbb/wN9/60TPPPMOULmCVWrUFCABg6SXc5k/dF3z9lkbaVi/f8Ghv\nr/Quc0hYJ8LWRPfbh4/9RiNtD/f03PXMM88sd0hYgyb0lu4xzddVb11rpZwOi7zidEzrEQUIcKlc\ninu4NbLLrNrn5orGS+brWmFypJlhATUFzR8veDyhqtdXcbWEzl0T569Ucnph8mwzYwOAZrCmblte\n9xoPne/Fjlnzf+iy6m016iz3la2cUi4VHA8Ia4fLJWbV9sz8/YVKtepqZkjATIaitBhVpdWu2rai\nKOemqEZFRBehAIGDXC5R3PW/c1SrtlRtpkwD68iaL0CyHjVquOSC+cxdIid0ChAAAJZdUfWGBwLB\n3fXWxS1rUCvmR52OCUDzrMb3pV4AACAASURBVPkCZCXZu3vrmxe6S9HhEetfr7vzvscdCglLzPZ4\ndcM3XfB6vBW3T3dpZjG9ko7u3fZzb0vuvXzjzvnauLp2vHbthz72d07FhNXptp97296Fcqnm0IvH\nDx96+hsvL3dMWJlsEcVwuQL11lVElHrLAce5FLe4VZ+o3vq34bVKxRmP9O985pMfb6hb9qkXoABx\nkmmo1YFj+vxtQn6HokEdPrfq9ntUNejVdBGRNk33zVgXqD1WXC6XiFRnP990uwOjineTyNQPcLWo\n3kJXycxIVVZMAdLtKW7ZrGWvmK+NcvlGp8JBE1TcHt8Fv1gtIm4R+2IGoG5PcetmLbtj4ZYihyMt\nR0WEAgRL4tMfueXnP/3Oq1oWbNjaXlSu/Fl+qwENMTxq+GXFtUfxBy4oivVqNd9VqbwkdqU8vcjz\nDi379kb6ZZ96oZVagHj/L7saq7ei0yyUfXbVExBlwR+JmU/t4szJSqmlpSV0LjP8lcqEXpgcWkzf\nWL22l8z2617+5226SymKiBhV+1wBEv2P712+qVoJiYjY0ctsESkt3KPLJapXE7vyk8+abZelUi7P\n86SVQE0mkwveirWvr288k8lYTgSEpWMq7sD4rF+sFhHRXYrZvUzb3PumK99wQBq6uZYcfvXsv/X2\n9k4sUyhYYUyfHh3TfF225reVUPt5X/ziZrG/7vSsgWOb7e99ve6UrplcXTvSyWTyyUbiyGQyVl9f\nXwPjOtanWftzt+o+70xJpVKeUZzMph/4yC3JRraSSk+ePfT0N3682GhXOkcKEN3jLV4R6uifvbyr\nJewRkQverK6WcPCtA6/ULUA0w/CZdrlsu7WqsvHyE+eWKx6Pyz7/S50dinl8WkBmLxcRsRWPx7TL\n5ay/RXV3bs/Wlrtdiu2qlM+7kNhqjXTpl79RUUtFs05IU0fCY5vFJTLzi1hVRM6c1zK22e0SqdR7\nXTU97Vtan//qjhvna5NKT2YPPf2NhW4Tp4jMf9Q9k8kYfX19kwv0s+rpHq85M/+6WsK113zuGqB2\nLTB5Raijv9O2zcDGy0UTly0i4pLquZ2h6lKsQNUeEZnOLdVXaqTNiIjfnNEmbtsZrVye+nIV7ayI\niNslLkOkOvP9KojIT+6fENvscdX5rJyngfxquK/Wdte397/rkwt19Rf5257u7++f965fmUzG29fX\nt+D1Vn19fWfWYjHjVdx2YMZYVXOxY9bs59ihmEdRtbLI+XlXr/3MtnNtW0REc3vUVCD4psAb3nLe\nb8WoItVQyRyWavX8L2fRTsXlUtxStevl3YSInBtf9sY2X743pjZU3/xj7C1GMpls6G5yfX19WiaT\nqXevkXosaeQzMtXvSCaTueAs52qjuT1WKNbdHwi0eS9YNysPZubJXDlSW1cNtsnMnJrrObU+5+vP\n7VIsV9VWXaGYx+X2nNem7FIse9Y+WWRqvxzq3i1il+caMzIytW9WvvOZT868tW9VRIr1ntAf6Hyt\ngX2riIj09fVZmUymkWKl7tnyelKpVD6VSuUbabta6B5vcWtr9FRg4+Wts9fV289eTA62tW0oBMzC\nhePeIvutl9u1fmfuz+1QzKNovnPt4nY1o5Wt+gdO2uLBT//8VTfVXTfl3D6yX+8c2PLGn673ffMC\nqVTKnUqlGi2aV9T4d8l37FFV1SciYllW3Q/yTIlEwpNIJLQ6q+b6YLpF5ILBctZzZj+3Xl+111dv\nG8vdj8iMhJqn77n6W/Y2q/loz0Xmn5pIJGbmU733s/b3aiT3Zv97MW3mamfJ+UXCUuXFiu2rr6/P\nWC1f+C4m/3p6enzhcLjeVJGLHWtmP+dicq2RnGtk+xfT/lJyeMW0Xcn5eJHjnzeRSFzw5W/afHkx\n399qrvd9rpyeax/ZSJvF5urs589ZgCwQ46W2vZgCpJRKpVb8QZhL+P7nF5F63wHrWWwOOtlvo+Oq\nyMIHhWd+X2T8AwAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAANCQS/4hwouh6/6r\nArr/F53YFhpXssr/nM1OfKvZcSw3XfcnA7p/b7PjwPwy2dxnG/lhq9UmFGr9Ba/qeWuz48DFGRlN\n39/sGJaCrvvfGtD9v9DsOLCwvFE4bBiF3mbHsZT4/rc6OTH+eZZ7AyIiqqp6db+/zYltoXGKYurZ\nZgcBrHE+zav7NI3xD02hqqqP/e/qULLKvvN/EHv14/sfAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAACAk9zNDkBEJBRqvb41GPgFwyi8KCKyc8e24YDuvyE7kTu0mH4WQ1XVWFfnxm+aZulwuVLJ1Jbr\nun9PONR6e0D3/4xllftt284tdltYOWLRyH5VVcOmab6m6/5d27Ymhl0uEcMoHL7UfhYbU7AlcP3G\njvgfZSdyfy0iEgq17msNBvYFdH9y5n8iMmJZ5ZHFbg/NtZS5M5+5xrjpcfRDPk3bUiqVXqpWq/ZS\nbXP2Z+pSx/qFLFe/60Gz80/kwjFvqZB/K18z88+tKIG2cOiDLS36e32adkW5UnnNtu3iUm2T/Duf\np9kBiIiEW4N36Lp/38ho+pFm9+PTtO3BYOD6UGvrflX17Fbcil5bF4tG7oi2Rx4VkYyIqNH2yIET\nr/VfaVlWajFxY+WItke+aBiFp7LZieea3U+wJXC9rvuvDoVa73IrilVbHm4N7tN1/77Z7SuV0R8Z\nRuHlS90eVoalysG5zDfGdcSjByJt4YMyNcaFQ6HgHSdTA9csRxwiIqNj6YdKFM0rSjPzb64xb7mQ\nfytPM/Nv8+bO532atkdERkQkFo1GDpxMDbzVsqzB5Yhlvedf0wuQUKh1n+bTrhYR6e7qfLR/YPCj\nM9bd4FU9eyoV+3h6PPNEbbmu+/cEdP/1ImJlsrknLMsarNePqqrhYEtgn9utbBIRo9Z2/niCN/o0\n7b0iEpi53K0oaltb+AHDKDzTPzD4S25FCWg+rVumdtSiqmosHAruE5FYpWL3Z7MTX1Pc7kA4FPxE\nySofz2Ynnqi9Xq/q2ZVOZz6juN2BGfGNZLK5r9fiu5T+RMQbCrXeWu/1+jQtEQwG9omIXqnYx7PZ\niacqtr3sA/xqclnnxgdFRDSfdnVHPHpXbjL/9yIiiqJ4I23hO9xuZVOxaD6Xm8wfEZnKiVCodZ/b\nrWyvVOzT2ezEExXbzs/uZ2h49GFd9+8K6P4b5SL+/pG28H4RCc9e3j8weHPt325FCWze3Pm8iIyk\nxzPPiEwdZZnellqyyi9msxPP6br/6oDuvz5vFL5uGIU+t6KokUj4nlouzRdfKNR6vVf17BERKRbN\n3trrj0Ujd5es8kg2O/H4dLvbvaonNjKafqiWl3mj8PWA7r+hUrFHan+fxbxHa93F5E7tbzwymr5f\nRESfOhO2N5PNPWJZ1pw7tbnGOFVVw5G28L1F03zuZGrgnR3x6D2RtvADoVDrDYZReDkcCt6WNwrP\n+jRte9E0jxtG4cWZY3HeKHzTMAovT+/gb63lVrAlkPT5tL3Fonm4YtvDc8c1d86oqpqYY/sX5Prs\nfn2a1hMMBt5bLJpHcpP5ZflSs1Y0M/9E5h7z5nn/yb81pJn5p+v+PT5N25Mez3xmaHj03khb+NaO\nePRr4VDw1kw29xT5t/SaXoBYljViV+yMW1E6jULheG25rvuv0XX/numHAVX1bB8aHr2/Ix69O9IW\nftCyyn2KWwlH2yMHTg2eeffsftyKom5JdP3QrSgJmapmw9NnLLbMl5xDw6MPi8jDsWjkQLQ9crC2\nXPNpu92KEhMRY9vWxA8VtxIYH888PDKafkRV1c4tia5/ditKXKYKkljtyGGotfWDquqJTeYmn6rY\nttURjz5oV+yRdDrzuS2Jrv90K0rcMAov6Lq/J9oeefDEa/07REQutr+R0fT927Ymfqiqnt2zX6+q\nehLdXZ0vVGw7Y1nWcZ+m3RsKBe9cziObq5FpmqlgS0Dsip0pmqVzhWqkLXy3iORlasA6eGrwzLW5\nyXzv5s2df+/TtGum37890Wjk7ldfTb1ldj/BlsDVl3Vu7BWRkohYIhIOtgRumllI1NM/MPhLIiLd\nXZ3P67p/d702Gzd2PKiqauJkauBGkamifdOG+DMzt+XTvPcaRuG5aHvkoOrJdRpG4aO67k9G2yMH\nh4ZH75wvvpmfNxGRaLvngfR45lNDw6MPRdsjdxtG4aVaARJuDd6m6/7dI6Pph2pnaaLtcu/03y4c\naQt/4mTq9asofOd2MblT+xvXdsAB3b832h45mDcKT13KGKeqnu0iok5O5ntFRHKT+d5IW1i8qme3\npXqMaHvkYFubfYdbUTpHx9J3BqbGrEctq/ySiHij7ZEHTp8dvimbnXh248Z40qdpB+1KJbVxY8df\nioj3ZGrgEVX1xGbGMjOH5ssZVfUkZm9fVdVYvVwfGh79TK1/n6Zt37y58xuWZR1PpzOfW4K3aE1r\nZv6JzD3m1Xv/yb+1p5n5ZxbN4/0Dg9daVvm4iMiMmS0W+bc8lGYHYBiFXsuyXhYRGRlNP1RbXrHt\n4VeOv9bxyvHXLhORvE/Tkm5FUSNt4Xssq/zi6bNDHzxzZuhmEZGOeOyB2f0obndgfDzz0KnBM285\neuxEPJvNPSRThUzdL3INCIuIaD7thtzk5LN2xc5H2yNfDIVar4+1Rz7hVpTO02eHbzx67ER8aHj0\no3bFLvk0bVducvIxEQm3BFuu13V/j1tRErnJyac0n7bJrSidhlHoHRlLf+pkauCq0bH0p1TVo19K\nf6qqBrITE4+cGjxz7ezXG9D9PSKiZrMTj585M3z7qcEz75yczH/NrSjqJb9xa1Bt6p5lWS9nsxNP\n1ZYbRuGbR4+daDs1eOZGERGfT9sbCrXe4NO0ZG4y/9DIWPrO0bH0b7oVZXskEr59dj8V27ZGx9J3\nnXitf8srx1+7rGLbKX3qmo1F0XX/nmBLYP/oaPpTtQE3NjWg5l85/tplR4+daMtmc4/5NK0nN5l/\nsWLbqZZg4AYRkWBLy/UiYuUm8/PG59O060VERsbSn+kfGHzn6Fj6o0WzdLx+RBcaGh69+eixE23p\n8cxBVfXsDoVaL5g6hp9wKnfmUDsimJn1/3NTFCzL6u8fGLwqnc481tYWPmBZ5eOnzw7dPDQ8crOI\nZGLtkQdFRE4Nnv3Vim2PXNa5sdetKJ1nzgzdPN+Xgpnmy5mZ258r12ttVVWNbd7c+bxdsTOvvz74\nbs6+LazJ+bcg8m9ta2b+VWw7M/09ctCnaYmOeOxRyyq/lE5nHq+1If+WVtPPgMzFLJrHZ/zBDBER\nzaclRCSsqp493V2dL9XaTh+5O49dqViqR93etjF812WKoltW2VhkSHkRkfHxzP0jo+mHc5P5p7q7\nOn/oVT17VNWTEBGpnf5Kj2ceS49nHhMRkawYkbbwg8GWwD7Lso5Prc8+YVnW4OhY+pOh1tZbu7s6\nfyAi+dxk/tlMNvfUJfaXVxQl0BGP/cllnRvjM5Mtnc48oXrUXaFQ675IW/juim0PTubyTMFqkFEo\n9ImIVGz73FkR7/R7FGwJ3BNsCdxTW677z521O8eu2COapiW3bY3cM/03D8vUEYtFibVH7qnY9kg2\nO/G12jJV9SQMo9BXse2MiMjps0PnpjRmsxOPRdrCDwRbAntagoF9RdN8zrKsEbeiBOaKb2h49K5o\nNHJw04b4V0QkYFnll0bG0p+RBtWmheUm889G2sIH3W7lgs8q5rZcuTOH2phRK0RqU2HOjZ2Tk/mn\nDKPQ59O0sFtROt2KIrPG4rCIiGVZg5O5/GOhUPAeyyr35qbPqjRijpw5PXP709uaM9en1+8WEXEr\nSlhxu8MrcQe80jmcfwsi/9aXZuRfLBrZH22PPJibzB86c2boUzPfN/JvaTX9DMjFMItmSkSsomk+\ne/TYCdfRYydcJ1MDW/oHBi+YShQKte4LhYJ3jY6mP3P02InLDKPwzGK2PX1aztI0bbvI1Lw6EZFK\nxT5tWeVBEZFIW3j66HLg+lg0ckBV1c6iaaaKptmr6/4bgy0tN01/6Rt0K0ogbxT+/mTq9WtPvNYf\nHx1L3x9sCewLh4J3XEp/uu7vibSFHzSMwrNHj52IT+by547gq6oazkxMPPbK8de2nEwNXGEYhcOh\nUPAuXfdfvZi/yXpWssopEZHRsfSHjx474Xrl+Gve/oHB3afPDt81u200GjkQbAnsO5ka+NlXX01d\naVfsRV90pqpqWNf9NxpG4ZszC0nLKg/qur/Hp2mdIlNzqmPRyN0iItls7onpeO5xK0oim809tVB8\niluxRkfTnzx67ETLqcEz14iId9OG+IMz4oi5FUWdviYqMTvOGTl8g4hIpWI3fPYEjeVO7b3W/fWn\n6TWqNsa1tASuFxEJtgSSIiKlqSkGNSURkaJpZkQkUzTN3llj8e7pmBItwcB+EcmrqicZCrXe1Ggc\nC+RMaUa8c+Z67fWcGjxzrYgENm2If/4i/hSY5mT+NYj8W0eczr9NGzoejbZHHj59dvhXTw2eubPO\nl3bybwmtiDMgtTe5u6vz6doc0DnaWenxzMFIW/iBbVsTPxQRUVVPTzabu9cwCp+Z2U9mIveIiEhH\nPPpoRzz6ealzYdvFsCwrkx7PPBBpCx/cuWPbrSISrthTF0QZqtrbEgx8sCMe/WZHPDoiIrGiaR6Z\nvjBcstncoY549CtuRYmlhzMfFhHRfFp3d1fnDyq2bZhF8yXNp+2WqQuZnrMr+ZGL7c+u2BkRsUKh\n4F2hUPCOma83GAzsi7ZHPm9Z5ZcsyxrRdf81FdseNIvmjxbzN1mjMrruvyEWjezPG4XeuRpN5iaf\nK7aFjkTbI4+2tATuUFV1k1tRYlNzSK2Rmf3YFfu0iMiWRNd/Tj89IFPX6Vwyn+btERHVNM0jM5eP\njKXv37Qh/uSWRNdJmZ5Hmh7PHBQRKZpmyrLKL/o07SaZmn71jIjIfPF1xKNf9GlasmiaL7oVt66q\nnu2GUXhcRGT6Irgbrti+dVhE1HqvazqHMyISLprmizOntmFODeWOUSi8qOv+fVsSXcdkeg7wYjY6\nPcY9FGkL37Nzx7ZxmXrPjmSzE8/Wm/YwOpa+N9oe+eK2rYlXKnYl49O0Pdls7jOGUbh348b419yK\nEjg1eCa5cWPH1zri0UcNo/BCI3HUyxnNp11wkGm+XJ9+PYO5yXxvbjL/eLAlcHso1HpTNjuxqANR\n60RT8u9ikX9rVlPyz6dpu6a/O8mmDfGnN22Ii4jI6Fj6znydu0uSf4u3In4HpFA0j9i2nbHK5X7D\nKBxxuUTyheIR0zT7RERmPs7nje9a5fLhctnKmab5o5Gx9KfS45mnZveTHs981SqXj5immTIKhe9k\nJnKPmKb5Y7tiny6VrEbuL+0ql8spwygetm3bEBHJ543D030O5ibzT5w5M/SbFdvOlyuVzERu8nHb\ntlNGofBSbjL/2PDw6KdqR6bLlUoq0hb+bRGxzw6NfMS27aJllUeyE7lDUq2etcrlkWKxeHhkLP3b\nhlH40aX0V65UMkah8K1yuXx29usdz0w8YZrmc5VKJW2Vy4O5yfxfDw2PfnL2vdchYprmd0uWdbZk\nlU+YZunHtm3n80bhsGWVUyIuqT0ulazXchOTXy1XKv3lcnmsWCw+d/rs8EdN0zw+u5/x8cxj5Url\ntWKx+FJuMn8onzf+rlgs/tiu2CfKlcrCvyPjcrlM0+wzjMK5YsPj8YTL5fLZvFH41szf/jBN88dG\nofA35XJ5wCgUvp+ZyH12LD3+5Rm9WS0B/abcZP6vM5nsEyIihlH47lzxjY9nHjMt6+VKpTJhmuaP\nMxO5LwwNj35WRGRyMv/1cqWSmvW8fzOMwpFwqHWfqqq7T58dfmft8zI8PPo7lSW8p/pa1WjuTOQm\nv2WVyz8yTfNVo1D428xE7pBpmi/NHLMWUG+M+86MMe7Q8PDob0/9DohLbNvOTH8WTotMFaBGofCs\nbdu5Usn68chY+vfS45mv+jRtl8slvsxE7s8ncpN/VygWj1Smx5pSyTo+8zM1c2yfP2cu3P58uT6z\n3+n9Ql5E3IVicdG/EbXWNTP/frJm9phH/q0Xzco/r+aNlcvlEaNQ6J35X+09J/+wau3cse1Ud1fn\n0yu1P6wPwZZAcueObdVQqPWG5dxOd1fnkzt3bKsu5zawtpAzaCbyD820HvNvRUzBcpKu+6/u7uqc\n6zTU8NFjJ65c4u3t2bSh4ysi0pmZyH10wSc43B+ap7ur81Fd99edG5qbzD92avDMvUu5vU0bOh5o\nCQbuqNj2yGRuckXdDxxLx+kxDpiJ/EMzkX/ANFVVE7Fo5MBSHXFe6v6wfoRCrfti0cjdPu3Ci8WX\nYVvXx6KR/cu9Hawd5AyaifxDM5F/AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAC6Fq9kB1CSTycB86/v6+oxMJrOufqYezujp6fGFw2H3XOv7+vqKmUym4mRMAAAAa5XnUp8Yi0YO\niIiMjKbvX2wQyWQy8PEvfe5jA5OZ6JyNPvGHn+/t7R1e7LawNixl/t38u//3L3h3bn7rnA0+8YcP\n9/b2Di12OwAAAFhEAbLUBiYz0VeyQ93NjgPrz5g52ZqbP/c48wYAALBElGYHAAAAgNVBVdWEqqqJ\nxfSRSCR8d911V2yJQsIqRAECAACAhoRDwdvCoeBti+kjkUi0/vHN175/qWLC6kMBAgDAKhKLRj4W\ni0Y+tlz9V57/33+QTCZbl6t/AFgx14AAAICGdCxr7yOvt4kId/4DsGw4AwIAAADAMRQgAABAREQO\n7P/1a0TE2+w4AKxtFCAAAEBERD6d3HGLiPiaHQeAtY0CBAAAAIBjHLkIPZlMtr/7C/cdnGt9V0s4\nOzCZ4YgLgFUhFo38dxGRkdH0Z5sdy1I5sP/X3/bp5I5bRUT+oPfYX9//yJf/qdkxAQDWJqfuguV+\nJTvUPs967rYBYDVZewdMTEOrDhxrn/43d0gEACwbpmABAAAAcAwFCAAAAADHUIAAAAAAcAwFCAAA\nAADHUIAAAAAAcAwFCAAAAADHUIAAAADAEeFw2N3T0xNtdhxoLgoQAAAAOKKnp6f9j67e+N+bHQea\nix+bwprX09MTCofD2lzrw+GwP+dkQAAAAOsYBQjWvI8f/L2ef42Uf2Ou9f5Qx7Bkh5wMCQAAYN1i\nChYAAAAAx1CAAAAAAHAMBQgAAAAAx1CAAAAAAHAMBQgAAAAAx1CAAAAAAHAMBQgAAAAAx1CAAAAA\nAHAMBQgAAAAAx1CAAAAAAHCMp9kBACtVyOvXrkkPdf3Z//id28X+ZG6+tvc/9dwL93/+T/7dqdgA\nAABWKwoQYA6Ky+XaPngibp/40c+IVSrO23jCOOpQWAAAAKsaU7AAAAAAOIYzIACwTvX09IT++OZr\nbxURSSQScTnzo2aHBABYByhAAGCdCofD2ju07NtFRFyeolVtdkAAgHWBKVgAAAAAHEMBAgAAAKd5\ne3p62podBJqDAgQAAABOa/njm6+9pdlBoDkoQAAAAAA4hgIEAAAAgGO4CxawBG77ube9Y+/lG3cu\n2DDUnv2tx/7f/9PX1zf/DxsCAACsUSuyAPEobiWoat6Zy3p6emL12mYyGbOvry/rTGRAfd2e4uWb\ntewbFmrnim3oD4fDX3ciJgAAgJVoRRYgQVXz/tq//9ObZi7bsd3+rGx9Q2l223/3Xvb93/qr8ON9\nfX2ZTCZzwXoAAAAAK8eKLEDqGfDru4yKesE1K5nLott/+cBH3nnbC0f/JHXytcx8fRx++bUTvd/7\nwanlixIAVr9EItGaTCbjIiKpVGoilUoxZXCd6enp8adSqWIqlao0OxYAa8+qKUDmEisWAtedOr5t\nR2v+12SXvsBOcuujFCC4WMPBtt1mpTzn+q587j9FxHIuImB5/dqOtl/+kL71FhGRP8wl7zp48CAF\nyDrzR+9960dDpYnv3P/Il7/f7FgArD2rvgABlptZrapG1Z77jnEucTkYDgAsv5HX42IavmaHAWBt\n4ja8AAAAABxDAQIAAADAMRQgAAAAABxDAQIAAADAMRQgAAAAABxDAQIAAADAMRQgAAAAABxDAQIA\nwDrX09Pjff7Jv7haRPjtDwDLjgIEAIB1LhwOq++QM+8REW+zYwGw9lGAAAAAAHAMBQgAAAAAx1CA\nAAAAAHAMBQgAAAAAx1CAAAAAAHAMBQgAAAAAx1CAAAAAAHAMBQgAAAAAx1CAAAAAAHCMp9kBOCmZ\nTIZcXW+IN9L28ccfz6VSqcJyxwQAALBe9PT06M2OAc23rgqQd2jZ9709+PqvNNK2N5G4iwIEAABg\n6fzxf33/7dWXX2h2GGgypmABAADAGRNjrc0OAc1HAQIAAADAMRQgAAAAABxDAQIAAADAMRQgAAAA\nABxDAQIAAADAMRQgAAAAABxDAQIAAADAMevqhwiBFcD/7QO/cY/sf9eC90E/PGL963V33ve4AzEB\nAAA4hgIEcNrI64HqwLH2BduZIb8D0QAAADiKKVgAAAAAHMMZEKw7CcUTClolrfa4s5j3ahU7PLud\nLi7V2cgAZySTybiISE9PT0xkotnhAADWGQoQrDs7J8biVwy+Gqs91l1KcUfV9s1uZ0cvs52NDHDG\nt/e/6/erA8faXV0breoABQgAwFlrpgAxNb294latudb7KtakWyTrZEwAAAAAzrdmCpBh1dtlKK45\nr2npKsmJAAUIloGp6e0VjzcgvoB79jq3SFUr5keaERcAAMBK5GgB4vd4625Pc3vUmet9bvWCL3LA\nSjWsersmXSKK13tBAay7FLOLAgQAAOAcxwoQv8fr+fmxs9uCxoQ2e12r+6QhItJTsXQRESXc4XIq\nLgAAADhv7+6tb3r+GLz6IAAAIABJREFUyb/42Wtv/rV/anYscJajZ0CCxoQWmUgHZi/XXYpbRMQz\nfSGw7dWbfvHvbTdc95a9l2+YXLChphcPPfudf02lUnNefwIAAIBZTMMrJrNe1qM1cw3IUvu1hLav\n6tYXvA2rq2tH/+GXX/sPChAASymZTMbf/YX7fr+Rtke/+rdPfOlzX/jecsTx6Y/c8vOffudVLdLa\nnrvuN+/7m97e3sJybAcAsH5QgADACvVKdqi9kXZGuXTB1NYlM3Bss/29r+92de3oF5Gnl207AIB1\ng19CBwAAF9j7pit3JN/2lk3NjgPA2kMBAgAALrA3pl69d9fWK5odB4C1hwIEAAAAgGMoQAAAALDs\nvvPFP7xdRC64GyrWHy5CB4BVLpFItCaTyXgDTV29vb1Dyx4QANSxN6a+WUS47S4oQABgtWu/9s2/\nfMWbLrtloXZXhDrGet/Y+9+ciAkAgLkwBQsAAACAYyhAAABYx8LhsNrT0xMXpsYAcAhTsAAAWMd6\nenra/ujqjZ9odhwA1g8KEGCF2vumK99gH/n6xxtpe909n3+yt7d3bLljwqqnfeCujydzbrVT8QVa\nxeOtiC8wddR71r9Vr6+slorZZgYLYO34/G9/7OdFpEVEKs2OBc23fgoQt8crbtUnqlb/NdsVWypl\ny+GosNa5XIqomnbusVvVRUTOLatWRcols+5zJ0bD9ksvvL3BLf3NouKEY5LJpF8amP7a09PjN5Z+\n81rwfT/70ePDg7t1r1ezXWIrXq8iIjL733GvNtZOAbJ+udyevNu9YfvP/HTPB4Jt8zatjGZPPPmX\nT5xyKDKsUm9s92+VqbGPAgTrpwAZ82qX9SuuyxS/XnfH31UqndAr5WGn48Ia51Lcx/z6m2sPJ22r\nKCLS4td9IiJR2x5unyydaFZ4cN5n/upLv3xk6OTuhdp1tYTd3x485kRIwIXcHk+/q7p9fFuHO7hj\n41vma/rmtOfPKUCwVI6cffXDR4ZOvmmhdpNP/9MDBw8ePOtETFh666YAWUb+bx/4jXtk/7taF2wZ\n25x1X/uB+xyICcAKNWbmg69kh7obaNrcM7KK2yOqponiDoiqaeJWA7fddlvs/2/v7qPjqs8Djz9z\n79y586q58yZZloTG+A3bJIgQp7wlHkqcNkBOgBQSkm5Nsoe8kG4LyyG7C7QWu4R2Ay2QNqSHHLbm\nZEPSkLakadLTksQiJKSNk41Igh2/YI8tZCNpNJ6RNDMa3Zmr/UMaR8iSNbJHc2XN93OOj6WZn36/\nR3ce3Xuf+/K7iUTitIeI7dq1aySZTI7bESaAlWW4mHMfyA5EFmq3WsSqRzxYGhQgtTB0zDfZt3/B\nPxZHPWIBgBrIq85gyuO97I3RoXWrPN5gulS4qHdjx6rC2vWl2W3jPfG7KEAaV3P76vbb7vpUopq2\nP37+X36STCaX4OpCAOcTChAAaAABze0KuHQ1pHs1B8dDUCWf5tI8isNZKJmnFZ4VzlgwHrj56uur\n6S/eu28vBUhjMsIRr6iaWxTVKarmFEX1i6q5g6FwKJFINM9o6rItSNQNBQgANICrTw50du3bE8mX\nTcWKtnPpQoN64M8femvHRevfdFd5NBwOjQ79sl1ERHEofhERUdWiWCXZOnS8I2dEzV+qCvdI4py0\nb926fu+xvVutUsFSFIcyODp0UbPiCOxvKrdf+dk7B1NlszDd9NVq+rvixvdee8+asLFQuwP/8oPd\n3/raN359TsGj5ihAAAALcojDEXJ5dLfqdIrI5Mlifu7Z27CsFVqbrvpZuPSm2fVC+pi73TI3i4h4\nna6po89lsyiO2p8ou+ehP73ittSnMgu16/v1weGH/vsDv6p5AFh2movj/lans5BSlcLCrX9jIuS5\ncHRNZMEJPfyrYvtEhAJkmaEAqS9Xd3d388LNRJLJZHnXrl0812GlcyiqaC73nO+pmiqqU6lmeugd\nO3Y0JxKJBYdLJpPWrl27UouOEw1PU1TlvUde3ejPDLrzoZbil5rbflZ571NPP/rADWOZBfdWJ/Yd\n+9l9H/rYriUNFMva/oB18wGrpC3ULtDa9JKIUIDgnL3l8q0b5a6F273201f2/+SHL59Y+oggQgFy\nyojX3zmsuzsq31u6x1KCkVNT9jZPFPv0wti5noJ2PxA49ng1DX+wIfTSLpEvnuN4WO4U1bnf45tz\nukFLcVjNHt9wZCx7cKFu/sA3eO9k4NiCG3XyCkuhbywTPJAdWHinsjjmqUc8OP+1X7R+/T3PfL6q\nB7F+77H/87Xe3t6TSx0Tzk+uTRdcGVitb1uo3VqRpyhA6ocCZJo5OanmJ61TG1DLKlvKpHWqAClP\nTtoTGAAADabs1iKjayKrqmlrGAYPYsU5u/Km916z+tJNmxZqd+BffvCDb33tG1Xdp4L5UYBUyXJq\n3rzbd9rlU5aquUREdJfbo9XwqcGXXNi+7vsP313V0Z9n9hx88Zl//Oe9tRoby8iZLtESEbGscjWX\naFUYhuGZNdvIvJLJ5Biz1aCW2i9av67ao9rDu//f17kMdem5VU3dII4zTiO/pVxqzrn8E4fLEwve\nu2GHW//7f/nQZR+9ecH14IZgy+FP3PjBf61HTFiczaVSc1bzjB+zyiN2xVA0vPHRNZH1C7W74sbr\nLP+qWGyhdq/99JUDP/nhy8drE93KU7cCZF1TrC8YK/p9vtBp06vp4rBERBwyqYiIWMGYU9H0N035\npytOp8MqnTYNYKXtfO/P/PnJQEhm9ztf/7NjsByKmZlxhqSi4DeaRERaXLrpKo7PPw/+5KRTom3i\nmP/hYpMicmrlbsRWa9uU0sbZw831g451l65fc8nbq7p5q7e315nJZKq6ebS3t7cvk8mc96d+vE7N\n3BCMHK18H3X5dJ8lp4pFXRxWJfdmqkVuzf5ZxW9YIiIeT0AREdEciumbtIbm6tcKxpyW7pPXrdK8\nD7lstiYzeskckWib4nDqp+eAVfKJVT4149ElsQsu+d6d7ktmtDBlnoc5/UBa//3FvYd7ZCo3z8TV\n09NT7eUPjir6ExGR3t7eoZWQf7NFdN/ohmDL0YXadfgNp4jMm3ez2qoiUp79+hrVueqi7HCziIgW\n8o26xZlzWKVSJbdnrudmf+3WfeKwSqVQaJXpKxaGZuauJxjT/lhRgiIiA76gNeQ30jKVS6PT/+a0\n2m84y13aBdX8Trf91m9fecWN711wti6vU5t8+tHPV3uvQEFEqiqqe3t78ysz//wjrTPyL+DSfdek\n3mgpxjqOiYjo0wfVZHLSEodDERHxiUjWFyw4HcqIiJz2tx7RfWPV5LRI9Xk9X07PJdJmdORbgwv2\n+ZaWNcru3bsPVNPniy++6O7p6clV0dQh82yb52CKyEQ1DVdi/mmBcMYXbT822RTTHJNihjz+gs/t\n94Q8/oJRGPNYRnMxKJOD1eZTte0WkXPVtWsz1rk2XRBfqN11794eeu2nv+hdqF0ymdSqfI5STkSq\naVfVdjaZTE4kk0nbHnh71lNcaJrmFhExTXPBhWEYhqOrq8snIvMfyT09rtkLb74FWnl9oQVe+V3n\nazP75xf6fqGfr/X4kzJ/4lW9U7eYtst5BbiY/Ovq6nIZhjGzeNRF5LSCYw61yK1zyZuFcqaaGM8l\nrysb1oVygPxbhEQi4ZHF5V815muryW8ONM3Mp7nyZvbXs9vON97MtqaceQO+FLlC/i1CV1eX2zAM\ndcZLikxtm2d/nnMpydw70OfL53q+9Lki8y+RSBgytU4603IYn36vlsuednOwuwABAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAcO7O+kGEi6FpWsgIBv6oHmOhdoZS6QftjqEWgsGm\n33Vpzt+yOw4szkrJv9li0fAnRaTF7jhwZis1/4LBpg+6NOdFdseBM1uJ+RcMNt3u0pyddseBKbl8\n4Z/y+cLP7RrfuXCTc6cqiuL1eEL1GAuYza27vG5dJ/+wLHg9niYRIR9hC7fu8rM+hB3IveVlwiy5\nRAp2hwEAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAGpKtTuAucSi4Ts1TTOKxeLhpRxH07RYR1vr\nt4vFiRdL5XJGRERVFF/ICN7u93vf59b1DRMTE/snJyfNWo3p9Xo2r7swPuhwiOTzhRc3bVw36PN6\nrsuOjD5TqzFERJaq30ZgZ/5VhEPGHc2xyF+43Xo8l8u/WKsxyb+Vwc4c9Xo9W41g0x0+rydR+Wea\npb2WZeVrMebsHO3saNu9urXlkdRw+pFa9I8pwWDT9qaA73fz+cKepR5r9aqWLyiq2lwsFnsrr2ma\n1hYJG5/0eT0JTdPcS53LOD/ZmaexaPhen9ezfea6zuf1JPL5Qs22yY3MaXcAc4lGwl/I5wvPZbMj\nLyxF/25dXx8I+LYHm5ru1DTnFkVVvCIiqqJoa+IX/IemObeIyJCIxKLR8L2vvZa8tGxZuaWIJTWc\nfmTCLA0tRd84O3blX4WmabGW5uhjIuKTqTxcMuTf+cnOHDWamm4JBgP3zmyfyxeeM01zSfIoMzL6\nTL5QiC1F343MaArc4fV6bhlKpZ9civ5VRfH5A/7rjKbALV6v5xazZHZX3vN6PV2dHW0/nP7WFBHD\npTnvHkqlH1+KWHD+sjNPo5HwvSIye91jDqXSDy5FLI1m2RUg7W2tnxMR0d36lS3N0bsGBlOPe72e\nzT6v53oR8ZbL1sFsduS5smWZwWDTLS7NubmSDN6p6nRbJjv65Jk2hsFg4Hq3rr9fpnbwTtHd+hZN\nczanhtOfHkqln2xpju4Mh4xu3a1fYpql40YwsCOXL3zHrevrx4vFg/l8YY/X69nq83q2i4iZyxe+\nnc8X9k5vvD88YZYOZrMjzwb8voTbrW8bHy++WLaswfnjmvp9cvnCN31ez3XlsjWUzY48W7asnKZp\n8XnGrywbbcIs7Zlrh8St612BgO/94+PFl0fHckuyw7JS2Jl/FatXNT8+OpZ7NuD33VF57QyfP/nX\nYOzOUU1zbhkdyz35ev+JT7/5dXtz9CwXZ0MKBptu0d36lSIinR1tTx3t6/+4pmlGwO+7RVWV1SKS\nz2RHnzVNs3/m55rPF/ZomhYzgoE7c/nCi/l8oWe+MTRNixlNgTsVVTkth2KR8M6yZeWPJPu2mKY5\n5PV6NltlKzsjvutcmnPrzDjCIeMOVVVWp9OZh8uWZXq9ni6f1/P+XL7wzXy+0BsMNm2f/hkZHy/2\njI7lXj7H/io/Y85ct6mKogWDTbeoqrJ+McspFg3fO32wZ9ClObdOmKVXstmR52vyga5Qdufpvv2H\nmitfe72erZ0dbbsHBlN/OCM+8uocLLsCpFgsJgN+n1hlKzNenOgP+H1Xtre19ojIhEwfKQn4fTce\n7eu/tVKxVjauPq9nWzQS7l7oaNzAYOpxEXk8Fg3vjEbC3ZXX8/lCbyXhvF7PZq/Xc1XZspLF8eIr\nulvfGo2Eu0Mh6w5VUdpSw+lP+7yermgk/JRpll4VEVc0Ev7s8TcGb8xmR77T2tqccOt6t1UuJ1tb\nW/6viLiOJPue1DTnm6rpaCR8bz5feDWbHdlV+X2iEblfRHIiYoRDxh8fSR67VNOc8dnja5oWW72q\n+fmZy8atu+4fGEw9XOnfrevrL7ig7Z9N0zyYTmcercmHtILZmX8iIgG/L6G79e1Hkn1b3lyAnP75\nk3+Nye4c9Xo9m02ztHrTxnUTIpJLn8w8MjCYenipctRoCuzwej1bhlLpR86Uo2XLqtmlsiudaZpD\nVtnKqIrSli8UDk6d/e/4uaoocZk662pEI+Gdhw4fXVP5XEXSQ1M7QM7Y9PfdZ9qxGy8Wk0f7+q/x\nej2bOzvaXp35nnfqsr3k6lXNT3m9nqtGx3LfPnFi4G4RkdWrWp4KBgN3iEhGRHzRSHjn0b7+d7h1\nfX0wGLh3eqftO0ZT045gMHDX6GjfMy3N0XvDIeNzplnqFRGJRpyfTZ/MfGZgMPXIWfb3uXDIuHf6\nZzQR8bk056eHUuknL7ig7Xm3rm8fLxb3aJq2utKfiDSfaTlNH003phfBhIj4An7fl17vP/HxGn2s\nK47deVqhKoqvva316/l84Tvpk5ldImedp+TVDIrdAcxWOc1mmubeylG81HD6rkOHj645cPBwe9my\nkl6vJ7GUMbh13ejsaPuRW9e3F8eLPTL1ocp0XEeP9vVfmk5nvhQKGTtNs3Tw+BsDtw4MDt0qIplY\nJPw5EZHX+9/4/bJlDbW3tfaoitJ24sTArdVeojAwmLp13/5DofTJTLemObcEg023zDV+bGrHIHfg\n4OH2ffsPhbLZ0S+5db2r0lbTtNgFF7TttspW5tix/hs4Srgwu/OvpTn22MmTme75coX8g505qmma\nISLxslXOp4bTd40Xi3vDIeOzwWDT9kobu3MUC8vnCz2mae4VERlKpR9RVNV38mTmkdf7T7xj3/5D\nzdns6CMi4pu+HHkpGJrm7DLN0lA+X3gh4Pfd3tra8jm3rrcFg4E78vnC8/v2HwodSfZtyecLe9y6\nflVlxy/g990iIuIP+G4cLxZfHi8Wk25d3y4iMjScfvhoX/97UsPpj48XJw6ebX+mWdozMJj6yL79\nh0LTO4Hi9Uz9Tbl1fbtplpKpVLr72LH+a1LD6c8s4veeOJLsW3Pg4OHQ9O99h1vX22q1UFeaZZCn\nIiLS0hzbqSpK2/E3Bu8WESGvamPZnQGZzSpbQ7quJ9ZdGL5v+giXIVNH+ZbMeLGY2bf/UCgcMm5v\naY7+bTDYtGe8WNwrIjI2lnsuny/0unXdUBWlTVUUmVk1a5rTEBExTbN/bDT3pWAwcJ9plnpGx3I9\n1Y6fPpl5XkRkdCz3nXDI6J4+JXd85vjTY8Xz+UJv2bIyIiLH3xh4U8Vb+aNUFcVQVNVgB3Dx6pl/\n4ZDxYU1zrheRWCwa3ikiomna5mCw6TrTNPMi5B9OV88cNU0zs2//IUfl+1y+sKezo+0nLs251TTN\nl0VszVGcJatcNjWntj7UatzVrihe0yzVZEKBM8jl84UXKuuMdRfGNwf8vuvSJzPPiojkC1M3+Y4X\niweP9vVfU/kh0yzt8Xo913u9ni5VUeKp7OiDIiIDg6m7otFw9+pVzX8rIj7TLL06NJx+WJnOi8X2\nZ5pmJhoNd7c0Rx83zVL/zMBf7z9xYzQavr+9rfXb033uGR3NvayoCx/PzecLe8aLxf7pmF7wej3b\np2PsX+BHIbbkqaiK4gsGA3eOjuW+Yppmv4gIeVUby+4MyGzRaHhnwO+75Uiy7+rXXktebJWt046Q\nVSo9r8dzTlVwwO+7srOjbXc4ZHx45uvqm68NnBCZKlJEJDNeLPbs23/IsW//IceRZN+ao339W6Zj\nivsDvjtFJKdpzkQw2HRjtXGEQ8Z10/FcJyJSLlsHZ48vImKapX6v19NV+f3DIePOWDR874z3D77e\nf+IaEfGtXtX8WLXj4zfqmX+qqrhEJB+NhD8djYQ/LSKiac4NAb/vqhnNyD+8ST1z1Ov1dMWi4Z1e\nr6drut8uEZFy2To+o5ndOYpFCgabbgkGA3elUumH9+0/1J7PF067hlxRlDYREbeuv+NcxzPN0l5N\n09ariuLTNC2mqErMNEtHrbKVFBHxejzbpsdqi03l93YRkezIyLOqosRamqP3iYg5OpZ7XkREURUz\nlUrfvW//If/r/SeuEhHX6lXNnzvb/trbWr+uKqq2b/+h5uNvDHykEreqKL6yZSVf73/j+gMHD4eO\nvzF4o1vXr2xpjp5a751pOXm9nq1uXW9TFUXzejzbRUQscrdq9c5TERF/wH+diPhGx3LPVV4jr2pj\nuZ4ByXi9nuti0fCd1vSGbU2841fT752aGShfKOzxej23rIl37Jfpa5/PZdB8vvBKa2tLm9fr+UpL\nc/RxEYmVLas/kx19dvrI9JukhtP3RyPhL6y7MH6gbJUzbl3fms2OPpzPF+5vbW3+iqoovtf7TyRa\nW1u+0tIcfSqfL/yomjhamqPfbmmOZkTEGC8W92SzI8/pbv2q2e2GhtMPrl7V/PU18Y4jMn09dPpk\nprvyvmma/aNjuZ7RsdyugN93ezDYdOP5cnOSzWzJv6FUetdQKr2r8v2mjesm8/nC86/3n7h/rktq\nyL+GZkuOmmZpKBQy7ohGwvfL1HXMsRk5snV2+3rm6Ln8Xo2ockays6PtHzMjo0+KiLQ0R5+anoHv\nVJ5M31eRC4eM+8Ih4w45fVagRRsYHLq/va312xvWXzgw/ZLrxImBj44Xi8l8vvCs1+v58KaN607K\n1HXycrSv/2oRkUx29LloJPw5t67fMl4sfsc0zcx03F9w63pivFjcoyqqV9Oc6/P5wq6z7a9sWUOa\n5ty6aeO6wZnLQlFVo7Oj7SciYubzhV7drcdFpo46V7mcfNPrywkR8WWzo49XjlxjbnbmqYiIz+O5\nUkSkOF48NQ0weVUby7IAeb3/xPVut759wiwdT6czXxovFg+qqtI5fZTLp6rKareutw2l0o9MmKWD\nLs15iYgMTZilIZfm3GxWOa1oLl94USTdXWlftqzca68lLw0Gmz6sqsrqctk6PjqWe246ebTUcLo7\nN2Mu6qFU+slcvrBneoYXGRhMfXr6OsDNY2O5f0ufzD4+OpZ7udx/4laf13OdW3dtHS9O7Jnu50WR\nuadBPf7G4HtcmvPKctk6Pj3Di2mapeTs8bPZkedN0+yaaxaimf0ODKY+UywWk6qirD6Xz6VR2JV/\ns6WG090TZmmviMhcnz/517jsylHTNPuPJPsuNYKBW0QkNnPGraXK0bmm4Z0rR2u0aBvG9N/lQREx\np/9u3+PSnFeKSL6SM6qiaHnTzBzt63+Hz+u5UURkwiy97NKc23JVPgvBNEtDMz9PEZHRsdwLR/v6\nT607cvnC8/l8Ya+IyNG+/o8Eg03PTefsqdmApvoy+8eLxRfcun5d+mT2K5X+jh3rv94f8N/o0pyb\np2N8sHKw42z6O5Lsu8oIBm6XqfVar0tzdk01N/sPHT660QgG3i8isXyhYObyhRcqz6hYaDmNF4sv\nZ7OjX1JVpfN8mq3ITnbmqYhIrlB42SyZmfFi8U3rTPIKK0pnR9vXN21cN2l3HGhM5B+WO3IUq1e1\nfGHTxnVjqnL6tKnLob/5bNq4brCzo233Uo6B5YO8WtiyPANyrrxez5WdHW3zVYCD+/YfuriuAaGh\nkH9Y7shR1MKmjet+JSLNc713tK//xny+8PJc752tdRfG/03TnInRsdyuWkxqUev+sDyd73kKLLlg\nsGl7LBq+0+440JjIPyx35Ghji0XD94ZDxh21Oqpc6/4WGOv2xUy0gPMXeQUAAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAQJ0YhuFIJBIeu+MAALslEglvPB532h0HAADLTiwa3hmL\nhnfWoq9EItFc3v3VB2rRFwAstVqu/2Yr7/7qn+y882OXL0XfAFCh2B0AsCwMHWuxOwQAsN3QsWYp\n5t12hwFgZaMAAQAAAFA3Z12ABINNtweDTbfXIojd3/xG1867//AttegLAAAAwPJ11jeauTRnZ82i\nGEkFZSTlqFl/AAAAAJYlZroAAABAVTRNaxERMU1zoFZ9xuNxbzwe92cyGbO3t/dkrfrF8sU9IAAA\nAKiKEQx80ggGPlnLPndce/k7vvu+Cx//y1uv+VAt+8XyRQECAECDSyQSzeWv/++/EpGA3bEAWPko\nQABgkWo5CUdF+pf/ft3Ou//wrbXsE1iMyb79ERHhfkwAS457QABgkWo6Ccc0Y3SgU0ZS+Vr3CwDA\ncsMZEAAAAAB1QwECAAAAoG4oQAAAAADUDQUIAAAAgLqhAAEAmyUSiWYRcdkdBwAA9UABAgA2+96d\n7/0TEfHbHQcAAPVAAQIAAACgbihAAAAAANQNBQgAAACAuqEAAQAAAFA3FCBoaIlEwvv9v37oDrvj\nAAAAaBQUIGh0DhkZdtsdBAAAQKOgAAEAAIDtLtm0YdXuL//NtscefOCtdseCpeW0OwAAAADAMEfi\n7xzc/4lJLfiSiPzC7niwdGw/A7Ljphs2X7KmfZPdcQAAAABYerYXIPGwv9lQyhG74wAAAACw9LgE\nCwCWiR3XXvEuwzC0u//ib16wOxYAWGo77/zY5X+a2PgRaYpqk6/+yO5wUEe2nwGpiMfjTfF4nNmI\nYBd9x003bLE7CDS2Tuf42ksinnV2xwEAdVHMuyf79kdkJOWxOxTU17IpQP5gY+gDO669/HK740DD\n0ndsXf8uu4MAAABY6ZZNAQIAAABg5aMAAQAAAFA3FCAAAAAA6oYCBAAAAEDdUIAAAAAAqBsKEAAA\nAAB1QwECAAAAoG4oQAAAAADUDQUIAAAAgLqhAAEAAKfsuPaKbY/d88lr7Y4DwMpFAQIAAE7pdI5f\neEnEs8HuOACsXBQgAAAAAOqGAgQAAABA3VCAAAAAAKgbChAAAAAAdUMBAgAAAKBuKEAAAAAA1A0F\nCAAAAIC6oQABAAAAUDcUIAAAAADqhgIEAAAAQN1QgAAA0MB23HZrx2P3/dfL7Y4DQOOgAAEAoIHF\nvcrat2b2fcDuOAA0DgoQAAAAAHVDAQIAAACgbihAAAAAANQNBQgAAACAuqEAAQAAAFA3FCAAAAAA\n6oYCBAAAAEDdUIAAAAAAqBsKEAAAAAB147Q7AAAAADSOHTfdsGXH1vXvisfjq+TEK3aHAxtQgAAA\nAKBu4mF/7F169p0O57g5aXcwsAWXYAEAAACoGwoQAAAALBuGYXgSiURzV1dXyO5YsDS4BAsAAADL\nxiVR9yXffd+Fb/9BMfjSb/f2ftHueFB7nAEBAAAAUDcUIAAAAADqhkuwAMAmXV1dob+89ZoPiYjP\n7lgAAKgXzoDDiyplAAAL+klEQVQAgE0Mw9DepWffKSKq3bEAAFAvFCAAAAAA6oYCBJjW2drcsuOm\nGy62Ow4AAICVjAIEmBb3KWt2bF3/TrvjAAAAWMlsvQl9x003bN5x7RXvEpGyiMi2t128Mf69+J5k\nMpmzMy4AAACcm091/4/3eNe0rp39+uZQYMPwiV+sn3S5yw5/UA26PMedIuz7NRBbC5B42N/c6Rxf\nKyIHRES2xbQr4/H4P1KAAGhU27Zc+Lby7q8+qF5z2067YwGAc+Fd07p2dE3ktCsLBstWc2pAiVmT\nlqUoiuJVlBQFSGPhEiwAWE6KeZcMHQvbHQYAAEuFAgQN7ft//dAfiYjb7jgAAAAaBQUIGtvIsFNE\nHHaHAQAA0CgoQAAAAADUDQUIAAAAgLqhAAEAAABQNxQgAAAAAOqGAgQAAABA3VCAAAAAAKgbW5+E\nDgAr3Zf/+R/eOdDRdOtc7wVcuqfvR/90maV7rMxEbsNFqpZ1l82BescIAEA9UYAAwBLKl0zlQHYg\nMtd7Id3rzk9aumWVrYxZ9FoiI/WODwCAeuMSLACwQZPL7bpIlDkLEwAAVjIKEACwgaaojrcfP9xh\ndxwAANQbl2ABAIA3icfj/u7u7ubu7u5Bu2PBypdRlAtTwYhYusdSghGl8r8zvqX10Ve+u7nSbmLf\nsZ/d96GP7bIxVNQIZ0AAAGhQO2664aJtb7t4w+zXO53jF//pey79gB0xofFMTFpqftLS81ZZn/n/\n8ETefyA7EKn8Gy6OeeyOFbXBGRAAABrUjq3rr9kW0y63Ow4AjYUzIAAAAADqhgIEDev7X3jodhHx\n2h0HAABAI7GtALn99tsjO669ImHX+MC2mHaZiGh2xwEAANBIbLsHJB6Pq53O8TV2jQ8AQCNIJBKe\nTz396D19Y5mW2e85h16/OJkbWSUiogQjikv3jq8q5nvrHyUawQWK2rRlJN0SdHncdscCe3ETOgDU\n2WaHGrtoYoKHEKJelL6xjO9AduC0nFs7kfeHrLIuIqJMWoo1aVn1Dw+NImia7g39r8WsaHtVeaY6\nFIfH6XKWJ63JiXKpvNTxoX5sKUASiUTkT//zh35n8sfftGN4NKhEItFc+fr7X/yzD0/u/XGLKGpK\nRERUzS0iYoQjgUQi0dLT0zNgU5hoAGrJdGilkmp3HMBpHA6lsj4URQ3MXG/OlslkJnp7ezN1iw0N\np71kGr/Xd2BLX2R1psflOmp3PKgdu86AqNK3/4K53vjezk/cd63IQyJi3fPM5z9YTWfHf77v8Cdu\n/OC/1jRCrDg3PPHA3QeyA50iIv90pHfdWsURHStkiyIifsWhi4gUO8KX/MnjD1zVumvTY9X0OfbG\n0NC3vvaNV5cuajQURXWJiFtEzO7u7nl3/EREenp6Cj09PaP1CQwNQ1HUvYpjq4jIa+m+NRv+5x1d\n8zUNHBl+qXfHH32xfsGh0WilkhIeSftSgVBeXC67w0ENLb9LsIaOhWTq5njlW0d/+c5qfiSQHV7a\nmNAw/BNFV1MuGwvcfPXHq2m/+sjwS0IBghoZ1rSOpGl1jpdN8/hvr3/8TG1bMwNPSU9PT51CAwDb\ntV+0fv09z3z+Uwu1yx85cfCL3X/23XrEhLOz/AoQAAAAYJayW4uMromsWqhdYOo/CpBljOeAAAAA\nAKgbChAAAAAAdcMlWACwzGiKqm5yqNFDDkmbVplpUbEsrX37JW979JXv/lU1bfd9+Vtfe/rRJ360\n1DEBOD/UvQB5+Y3XPvrLk8ev+vlA3xZD9xgiIpmJ3IbK1yIit33xfz1suDyF7/Xvr6rPiy/fuvFb\nR3+x4E1JIiI/+9t/+Pvu7u7BswoeAOpAdSjK24+/1pbs2JChAMFylS+ZrgPZAV91bSe0pY4HwPmj\n7gXIcDHnPjySCq8zi17X9MOPMjO+FhE5PJIKx9yBqqeXHB7PhX48cLiqGbNWi/z94qMGAADA+aC5\nva3ltrs+lVioXcTtS//1nz/6izqEhFnqXoA0aW7XtomJznqPCyyFRV6C8HdPP/rED5c6JgBYbhKJ\nRHCj0XLGZ9tU7Nq1azSZTBaWOiasXGo0uCZw89UbFmr3O51veYkCxB51L0BURVE6UsdDg37Dlgdo\nveeT/+n3/TddvWahdh1+I/vBtW9/oB4xob7eNulY5XO6XFIcP+e+FnkJAvdcoWbi8bieSCSqyr1k\nMllMJpOlpY4J5xe3qqmqQ3HM936T5tY/ODay5ZvByK/Hy2b5XMYaXRO5+Xh4fVUPF473xO+iAFlZ\nOhxq08Wj6Zagy6Mv3BqNoOF2iIaLOc+B7ECkmrbVbtwzmUy5t7f33PdmURctuUxA0TzzbnQVh8Pp\ncWrOQsms6Q5bPB4PJBKJqo4AJpPJkWQySU5hXhve+67rI9dclqim7b4vf+urTz/6BEf58Cbvzw5f\nFC6X5t3OTYqIazTtdBjRedeXS2HHjh2xRCJRTVNl165dQ8lk8pyKIyw9o2S6N/S/FrOi7Wd1T5uu\namrE7dO9Tk0RkVreF+ev5lItr9M18fSjT7xcw3Eb3rIsQCK6z606HGM2h9F0wxMPVHUGxHsiu+/r\nf/5XVd1bkkwmx1lZLm+u4ePhdf7gyC9VpaaTFUSuuez3Nryt/bZq2ia+f/Cu7u5uCpAV6N2l8pq1\nxfFoxqnlz6Wf4fFc8EB2IFpN2+b21TGKX8zWlB91O1Vt2U3HP/y2js8cX+ta8Kb1DcGWo/GenoeS\nyWSuHnHBPh35MeP3j/36spNru97odSqHa9i1K3Dz1R9fqNGGYMtotz90qIr+lO7u7jdqENeKV/cC\nJKh5cr5o+7GQx1/wuf0eEZGZX4uIvKdsOZNGS7o8OVnVKdgOv6GKSFU79RHdN7Yh2HK0ij6dfWOZ\narqUjk3r19/zzOf/oJq2xX3Hkr/695/+uoqmam9vby6TyVTze5VFpNodBodMHdhaUG9vbz6TyVTV\n9nzQ4TeGRESiLp/uc+pen9vvUfyGJSLicbpObYTdqtPZ5vG9XlSUE1X06RSRqs6ULKbt++569wcu\n++jNVZ2BS33/Zz88evToQBVNPSJS1d9UJpOR3t7eajfqizkiVZQql8FKyT+vUzM3BCOn1jmthTFX\nSPMUHdPrPUXTS1Yw5gx5/AXdLPqskmmJiKhls/hRPeD7Sch3eGSieNoBmcWs9zret/HyyDWXvbua\ntl2RdrNQnhippq3/8PBLL774YlWFem9vr57JZKotukpSZZ7I1Lqv2oM6jbz+S1W+1lXVvcGcWBUS\n7bjicJw6u6FoeskVbvU4yiVTZGrb7NZ05aN6wPcfId9roxPF3Kw+q8/BJVhXdviNsYe/9vQdw8Xc\ngvsyEd3n+Le/+fI3qhk/k8k4e3t7i9W0lanf36yybcOt/yK6P9sabDnaZllFX+tasYIxZ2WdN/N/\nV7jVozqUwnzvV/6X0KqJLpdrvDxpmcVyOT3fuNXm5qL2CW+6+u4q2pV3b9v29ELtRMT1xBNPlDKZ\nTDW5Y0p1f2dVrQuXy4Hwsz6tqmmaW0TENM1FHSlLJBJeEalmx2pcqt+pqXqjsoi2S9HnYseXRbSt\neazLeQV4Nvk3nXsOEdFFRD1DU4dUv6GwO08W09bu8RfVdqXkX1dXl8swjJlHc13y5oM/k/Kbv/eZ\ny6fyf1HmXhcuZrkvhu2f/XJou1LyT0TEMAxHV1eXd8ZLDhFxy+n7AJMydUBhctZrInPnod2fk+3b\nyaVqu1Lyr6ury20YhipT21yX/GYZzP6/ciBrvvdntyuLyMQZhq71dnEp2okdYy+XAgQAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADAuTvrBxGerWCw6XaX5uys97hYvKFU+kG7Y6il\nWDT832Tq4Vs4D6yU/AsGm37XpTl/y+44sDgrJf8qYtHwTrtjQPVWSv4Fg00fdGnOi+yOA7+Ryxf+\nLp8v/NruOJwLN6ktt+7yu3U9VO9xAcAObt3lZZ0Hu3k9HnIQdcc+3/IzYZZcIgW7wwAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA\nnJf+PyNYjsyKh1ZuAAAAAElFTkSuQmCC\n",
"text/html": [
"
"
],
"text/plain": [
""
]
},
"execution_count": 21,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## get full MCMC table from the .files attribute\n",
"table_1 = pd.read_csv(A00n.files.mcmcfiles[0], sep='\\t', index_col=0)\n",
"table_2 = pd.read_csv(A00.files.mcmcfiles[0], sep='\\t', index_col=0)\n",
"\n",
"## draw barplots\n",
"canvas = toyplot.Canvas(width=800, height=800)\n",
"\n",
"for aidx in range(table_1.shape[1]):\n",
" ## get data\n",
" param = table_1.columns[aidx]\n",
" hist1 = np.histogram(table_1[param], density=True)\n",
" hist2 = np.histogram(table_2[param], density=True)\n",
"\n",
" ## build plot\n",
" nx = int(np.ceil(np.sqrt(table_1.shape[1])))\n",
" axes = canvas.cartesian(grid=(nx, nx, aidx), margin=25)\n",
" axes.bars(hist1, opacity=0.65, style={\"stroke-width\": 0})\n",
" axes.bars(hist2, opacity=0.65, style={\"stroke-width\": 0})\n",
"\n",
" ## style axes (hide tick labels for clarity)\n",
" axes.label.text = param[:15]\n",
" axes.label.style[\"font-size\"] = \"12px\"\n",
" axes.x.ticks.labels.show = False\n",
" axes.y.ticks.labels.show = False\n",
" \n",
"## save plot to pdf\n",
"import toyplot.pdf\n",
"toyplot.pdf.render(canvas, \"analysis-bpp/bpp-A00.pdf\")\n",
"canvas"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### The FigTree file\n",
"Only algorithm 00 produces a FigTree file, so be aware that you will not find one when you run the other algorithms. There are other types of results files that are produced instead. "
]
},
{
"cell_type": "code",
"execution_count": 40,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"#NEXUS\r\n",
"BEGIN TREES;\r\n",
"\r\n",
"\tUTREE 1 = ((((((rex: 0.001598, lip: 0.001598) [&height_95%_HPD={0.0009684, 0.002092}, theta=0.00171779]: 0.000096, rck: 0.001695) [&height_95%_HPD={0.0010488, 0.0021839}, theta=0.0017499]: 0.000091, tha: 0.001785) [&height_95%_HPD={0.0010807, 0.0023218}, theta=0.00193679]: 0.000123, cup: 0.001908) [&height_95%_HPD={0.001117, 0.0025796}, theta=0.00467129]: 0.002187, (cys: 0.003522, (cya: 0.002113, sup: 0.002113) [&height_95%_HPD={0.0010666, 0.0031479}, theta=0.00212346]: 0.001409) [&height_95%_HPD={0.0026035, 0.0044186}, theta=0.00208407]: 0.000574) [&height_95%_HPD={0.0033103, 0.0050114}, theta=0.00372337]: 0.010006, prz: 0.014101) [&height_95%_HPD={0.011928, 0.016307}, theta=0.00146896];\r\n",
"\r\n",
"END;\r\n",
"\r\n",
"\r\n",
"[Species tree with tau as branch lengths and theta as labels, for FigTree.\r\n",
"In FigTree, choose 95%HPD for Node Bars and label for Node Labels]\r\n"
]
}
],
"source": [
"## the tree file -- can be input to the program FigTree\n",
"tre = A00n.files.treefiles[0]\n",
"! cat $tre"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"----------------------------------------\n",
"\n",
"## Algorithm 10 - species tree inference\n",
"\n",
"The algorithm 10 aims to infer the correct species tree from the data by implemented a tree search method, thus the input tree is treated only as a starting tree. Based on the results above I'm also going to increase the priors on the divergence times (tau) by about 10X. "
]
},
{
"cell_type": "code",
"execution_count": 55,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"## create new bpp objects\n",
"A10 = A00.copy(\"A10\")\n",
"A10n = A00.copy(\"A10-nodata\")\n",
"\n",
"## set new params\n",
"A10.params.infer_sptree = 1\n",
"A10.params.infer_delimit = 0\n",
"A10.params.tauprior = (2, 200, 1)\n",
"A10n.params.infer_sptree = 1\n",
"A10n.params.infer_delimit = 0\n",
"A10n.params.tauprior = (2, 200, 1)\n",
"\n",
"## also set no-data on the 'n' run\n",
"A10n.params.usedata = 0"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"submitted 2 bpp jobs [A10] (100 loci)\n",
"submitted 1 bpp jobs [A10-nodata] (100 loci)\n"
]
}
],
"source": [
"## submit job reps to the cluster\n",
"A10.run(nreps=2, randomize_order=True, ipyclient=ipyclient, force=True)\n",
"A10n.run(nreps=1, randomize_order=True, ipyclient=ipyclient, force=True)"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"True"
]
},
"execution_count": 27,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"## block until finished\n",
"ipyclient.wait()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### Plot the distribution of species trees from algorithm 10\n",
"Here we use the `multitree()` object of toytree to plot the posterior distribution of trees. The code for drawing 'cloudtrees' here is still in development and a littl clunky, but will improve in the future. From the example here (which has not been run long enough, btw), you can see clearly that the two replicates yield quite different trees. In this case the replicates are starting from different random seeds, but probably more importantly they are sampling a different random subset of 100 RAD loci that met our filtering requirements (e.g., minsamples, minSNPs). "
]
},
{
"cell_type": "code",
"execution_count": 52,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"