{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "INFO:stancache.seed:Setting seed to 1245502385\n" ] } ], "source": [ "%load_ext autoreload\n", "%autoreload 2\n", "%matplotlib inline\n", "import random\n", "random.seed(1100038344)\n", "import survivalstan\n", "import numpy as np\n", "import pandas as pd\n", "from stancache import stancache\n", "from matplotlib import pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## The model\n", "\n", "This style of modeling is often called the \"piecewise exponential model\", or PEM. It is the simplest case where we estimate the *hazard* of an event occurring in a time period as the outcome, rather than estimating the *survival* (ie, time to event) as the outcome.\n", "\n", "Recall that, in the context of survival modeling, we have two models:\n", "\n", "1. A model for **Survival ($S$)**, ie the probability of surviving to time $t$:\n", "\n", " $$S(t)=Pr(Y > t)$$\n", "\n", "2. A model for the **instantaneous *hazard* $\\lambda$**, ie the probability of a failure event occuring in the interval [$t$, $t+\\delta t$], given survival to time $t$:\n", "\n", " $$\\lambda(t) = \\lim_{\\delta t \\rightarrow 0 } \\; \\frac{Pr( t \\le Y \\le t + \\delta t | Y > t)}{\\delta t}$$\n", "\n", "\n", "By definition, these two are related to one another by the following equation:\n", "\n", " $$\\lambda(t) = \\frac{-S'(t)}{S(t)}$$\n", " \n", "Solving this, yields the following:\n", "\n", " $$S(t) = \\exp\\left( -\\int_0^t \\lambda(z) dz \\right)$$\n", "\n", "This model is called the **piecewise exponential model** because of this relationship between the Survival and hazard functions. It's piecewise because we are not estimating the *instantaneous* hazard; we are instead breaking time periods up into pieces and estimating the hazard for each piece.\n", "\n", "There are several variations on the PEM model implemented in survivalstan. In this notebook, we are exploring just one of them." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### A note about data formatting \n", "\n", "When we model *Survival*, we typically operate on data in time-to-event form. In this form, we have one record per Subject (ie, per patient). Each record contains [event_status, time_to_event] as the outcome. This data format is sometimes called *per-subject*.\n", "\n", "When we model the *hazard* by comparison, we typically operate on data that are transformed to include one record per Subject per time_period. This is called *per-timepoint* or *long* form.\n", "\n", "All other things being equal, a model for *Survival* will typically estimate more efficiently (faster & smaller memory footprint) than one for *hazard* simply because the data are larger in the per-timepoint form than the per-subject form. The benefit of the *hazard* models is increased flexibility in terms of specifying the baseline hazard, time-varying effects, and introducing time-varying covariates.\n", "\n", "In this example, we are demonstrating use of the standard **PEM survival model**, which uses data in long form. The stan code expects to recieve data in this structure." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Stan code for the model\n", "\n", "This model is provided in survivalstan.models.pem_survival_model. Let's take a look at the stan code. " ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "/* Variable naming:\n", " // dimensions\n", " N = total number of observations (length of data)\n", " S = number of sample ids\n", " T = max timepoint (number of timepoint ids)\n", " M = number of covariates\n", " \n", " // main data matrix (per observed timepoint*record)\n", " s = sample id for each obs\n", " t = timepoint id for each obs\n", " event = integer indicating if there was an event at time t for sample s\n", " x = matrix of real-valued covariates at time t for sample n [N, X]\n", " \n", " // timepoint-specific data (per timepoint, ordered by timepoint id)\n", " t_obs = observed time since origin for each timepoint id (end of period)\n", " t_dur = duration of each timepoint period (first diff of t_obs)\n", " \n", "*/\n", "// Jacqueline Buros Novik \n", "\n", "data {\n", " // dimensions\n", " int N;\n", " int S;\n", " int T;\n", " int M;\n", " \n", " // data matrix\n", " int s[N]; // sample id\n", " int t[N]; // timepoint id\n", " int event[N]; // 1: event, 0:censor\n", " matrix[N, M] x; // explanatory vars\n", " \n", " // timepoint data\n", " vector[T] t_obs;\n", " vector[T] t_dur;\n", "}\n", "transformed data {\n", " vector[T] log_t_dur; // log-duration for each timepoint\n", " int n_trans[S, T]; \n", " \n", " log_t_dur = log(t_dur);\n", "\n", " // n_trans used to map each sample*timepoint to n (used in gen quantities)\n", " // map each patient/timepoint combination to n values\n", " for (n in 1:N) {\n", " n_trans[s[n], t[n]] = n;\n", " }\n", "\n", " // fill in missing values with n for max t for that patient\n", " // ie assume \"last observed\" state applies forward (may be problematic for TVC)\n", " // this allows us to predict failure times >= observed survival times\n", " for (samp in 1:S) {\n", " int last_value;\n", " last_value = 0;\n", " for (tp in 1:T) {\n", " // manual says ints are initialized to neg values\n", " // so <=0 is a shorthand for \"unassigned\"\n", " if (n_trans[samp, tp] <= 0 && last_value != 0) {\n", " n_trans[samp, tp] = last_value;\n", " } else {\n", " last_value = n_trans[samp, tp];\n", " }\n", " }\n", " } \n", "}\n", "parameters {\n", " vector[T] log_baseline_raw; // unstructured baseline hazard for each timepoint t\n", " vector[M] beta; // beta for each covariate\n", " real baseline_sigma;\n", " real log_baseline_mu;\n", "}\n", "transformed parameters {\n", " vector[N] log_hazard;\n", " vector[T] log_baseline; // unstructured baseline hazard for each timepoint t\n", " \n", " log_baseline = log_baseline_mu + log_baseline_raw + log_t_dur;\n", " \n", " for (n in 1:N) {\n", " log_hazard[n] = log_baseline[t[n]] + x[n,]*beta;\n", " }\n", "}\n", "model {\n", " beta ~ cauchy(0, 2);\n", " event ~ poisson_log(log_hazard);\n", " log_baseline_mu ~ normal(0, 1);\n", " baseline_sigma ~ normal(0, 1);\n", " log_baseline_raw ~ normal(0, baseline_sigma);\n", "}\n", "generated quantities {\n", " real log_lik[N];\n", " vector[T] baseline;\n", " real y_hat_time[S]; // predicted failure time for each sample\n", " int y_hat_event[S]; // predicted event (0:censor, 1:event)\n", " \n", " // compute raw baseline hazard, for summary/plotting\n", " baseline = exp(log_baseline_mu + log_baseline_raw);\n", " \n", " // prepare log_lik for loo-psis\n", " for (n in 1:N) {\n", " log_lik[n] = poisson_log_log(event[n], log_hazard[n]);\n", " }\n", "\n", " // posterior predicted values\n", " for (samp in 1:S) {\n", " int sample_alive;\n", " sample_alive = 1;\n", " for (tp in 1:T) {\n", " if (sample_alive == 1) {\n", " int n;\n", " int pred_y;\n", " real log_haz;\n", " \n", " // determine predicted value of this sample's hazard\n", " n = n_trans[samp, tp];\n", " log_haz = log_baseline[tp] + x[n,] * beta;\n", " \n", " // now, make posterior prediction of an event at this tp\n", " if (log_haz < log(pow(2, 30))) \n", " pred_y = poisson_log_rng(log_haz);\n", " else\n", " pred_y = 9; \n", " \n", " // summarize survival time (observed) for this pt\n", " if (pred_y >= 1) {\n", " // mark this patient as ineligible for future tps\n", " // note: deliberately treat 9s as events \n", " sample_alive = 0;\n", " y_hat_time[samp] = t_obs[tp];\n", " y_hat_event[samp] = 1;\n", " }\n", " \n", " }\n", " } // end per-timepoint loop\n", " \n", " // if patient still alive at max\n", " if (sample_alive == 1) {\n", " y_hat_time[samp] = t_obs[T];\n", " y_hat_event[samp] = 0;\n", " }\n", " } // end per-sample loop \n", "}\n", "\n" ] } ], "source": [ "print(survivalstan.models.pem_survival_model)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulate survival data " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In order to demonstrate the use of this model, we will first simulate some survival data using survivalstan.sim.sim_data_exp_correlated. As the name implies, this function simulates data assuming a constant hazard throughout the follow-up time period, which is consistent with the Exponential survival function.\n", "\n", "This function includes two simulated covariates by default (age and sex). We also simulate a situation where hazard is a function of the simulated value for sex. \n", "\n", "We also center the age variable since this will make it easier to interpret estimates of the baseline hazard.\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "INFO:stancache.stancache:sim_data_exp_correlated: cache_filename set to sim_data_exp_correlated.cached.N_100.censor_time_20.rate_coefs_54462717316.rate_form_1 + sex.pkl\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "INFO:stancache.stancache:sim_data_exp_correlated: Loading result from cache\n" ] } ], "source": [ "d = stancache.cached(\n", " survivalstan.sim.sim_data_exp_correlated,\n", " N=100,\n", " censor_time=20,\n", " rate_form='1 + sex',\n", " rate_coefs=[-3, 0.5],\n", ")\n", "d['age_centered'] = d['age'] - d['age'].mean()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*Aside: In order to make this a more reproducible example, this code is using a file-caching function stancache.cached to wrap a function call to survivalstan.sim.sim_data_exp_correlated. *" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Explore simulated data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here is what these data look like - this is per-subject or time-to-event form:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
sexageratetrue_tteventindexage_centered
0male540.0820851.0138551.013855True0-1.12
1male390.0820854.8905974.890597True1-16.12
2female450.0497874.0934044.093404True2-10.12
3female430.0497877.0362267.036226True3-12.12
4female570.0497875.7122995.712299True41.88
\n", "
" ], "text/plain": [ " sex age rate true_t t event index age_centered\n", "0 male 54 0.082085 1.013855 1.013855 True 0 -1.12\n", "1 male 39 0.082085 4.890597 4.890597 True 1 -16.12\n", "2 female 45 0.049787 4.093404 4.093404 True 2 -10.12\n", "3 female 43 0.049787 7.036226 7.036226 True 3 -12.12\n", "4 female 57 0.049787 5.712299 5.712299 True 4 1.88" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "d.head()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*It's not that obvious from the field names, but in this example \"subjects\" are indexed by the field index.*" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can plot these data using lifelines, or the rudimentary plotting functions provided by survivalstan." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAIABJREFUeJzt3Xd4VFX+x/H3mfTeE0JCSEIPHUK30AREir2uws+C2Nay7i7rqquurrqWtawNFAVFBcXCIqKAIL33DgkBAqRDSALp5/fHHSCEhCRkZm4y+b6eZ55pd+795s7kk5Mz956jtNYIIYRwLhazCxBCCGF7Eu5CCOGEJNyFEMIJSbgLIYQTknAXQggnJOEuhBBOSMJdCCGckIS7EEI4IQl3IYRwQq5mbTg0NFTHxsaatXkhhGiUNmzYkKW1DqtpOdPCPTY2lvXr15u1eSGEaJSUUgdrs5x0ywghhBOScBdCCCck4S6EEE7ItD53IYQAKCkpITU1lcLCQrNLaVA8PT2Jjo7Gzc3tkl4v4S6EMFVqaip+fn7ExsailDK7nAZBa012djapqanExcVd0jpq7JZRSk1VSmUopbZX87xSSr2jlNqvlNqqlOpxSZUIIZqkwsJCQkJCJNgrUEoREhJSr/9matPn/hkw4iLPXw20sV4mAB9ccjVCiCZJgv1C9d0nNXbLaK2XKqViL7LIWGC6NubrW62UClRKRWqtj9WrsmrsWvMLGTuXs7XZ9ZS4eNtjEwD4ebrxfwNicXWR75yFEI2PLfrco4DDFe6nWh+7INyVUhMwWvfExMRc0sZy967gyoPvEHHge+4veYJDRFzSei7mzLSyPVoG0rNlsM3XL4RoeN555x0++OADevTowYwZM2y+/ueeew5fX1+efPJJm6+7Kg79QlVrPRmYDJCYmHhJM3P3vfMFSBpI+2/v5nf9PNz4CbQeatM6dx07ydVvLyPjZJFN1yuEaLjef/99Fi5cSHR0tNml2IQt+hyOAC0q3I+2PmY/rQbDfYvBPwpm3ATL3zrX3LaBMD8PADLyJNyFaAomTpxIcnIyV199NS+99BJ33303vXv3pnv37vz4448AfPbZZ1x77bVcddVVxMbG8t///pc333yT7t2707dvX3JycgCYMmUKvXr1omvXrtxwww2cOnXqgu0lJSUxYsQIevbsyeWXX87u3btt/jPZouU+B3hYKfU10AfItVd/+3mC4+DeBfDjQ7DwH3BsM4x9D9x96r3qIG93XCyKTAl3IRzq+f/tYOfRkzZdZ0Jzf/4xuuNFl/nwww+ZP38+ixcv5s0332Tw4MFMnTqVEydO0Lt3b4YONXoHtm/fzqZNmygsLKR169a8+uqrbNq0iccff5zp06fz2GOPcf3113PfffcB8PTTT/PJJ5/wyCOPnLe9CRMm8OGHH9KmTRvWrFnDgw8+yG+//WbTn7vGcFdKfQUMBEKVUqnAPwA3AK31h8A8YCSwHzgF/J9NK7wYdx+48VOI7AoLn4eM3XDDx9CsU71W62JRhPi4S7gL0QT9+uuvzJkzh9dffx0wDtU8dOgQAIMGDcLPzw8/Pz8CAgIYPXo0AJ07d2br1q2A8Qfg6aef5sSJE+Tn5zN8+PDz1p+fn8/KlSu56aabzj5WVGT7rKnN0TK31fC8Bh6yWUV1pRRc9jg06wLfT4Qpg+Gq56HPROO5SxTq60FmvoS7EI5UUwvbEbTWzJ49m3bt2p33+Jo1a/Dw8Dh732KxnL1vsVgoLS0FYPz48fzwww907dqVzz77jCVLlpy3nvLycgIDA9m8ebNdfw7nOc6v9RB4YCXED4T5k2DGjZCfccmrC/F1J7ug2GblCSEah+HDh/Puu++ird/jbdq0qU6vz8vLIzIykpKSkiqPuvH39ycuLo5vvvkGMP6YbNmypf6FV+I84Q7gGwa3z4SRr0PKcni/H+z95ZJWFerrQba03IVocp555hlKSkro0qULHTt25JlnnqnT6//5z3/Sp08fBgwYQPv27atcZsaMGXzyySd07dqVjh07nv3S1paUtuFRJnWRmJio7TpZR8YumH0vpG+Hfg/D8Jfq9PJ/zt3Jl2sOseufFzs5VwhRX7t27aJDhw5ml9EgVbVvlFIbtNaJNb3WuVruFYV3gHsXQeebYNV/IS+9Ti8P8XXndEmZtN6FEI2S84Y7gJsntLW2vAtP1Oml8aHGIZX9XvmNB2dsYP72NApLymxdoRBC2IXzD/nrGWBcF+bW6WUjOkXy3YP9mbP5KHO3HmXetjT8PF0Z2SmSsd2a0yc+BBeLDHYkhGiYmlC41/3EiB4xQfSICeLpazqwIimbHzcfYe7Wo8xcf5gIfw9Gd2nO2G5RdIryl1HthBANShMI90DjOnsftLm0MWhcXSxc2TaMK9uGcfraMhbtTufHzUeZtiqFj5cfID7Mh7FdoxjTrTlxofU/Q1YIIerL+cM9pDXE9IPF/4L2oyCwRc2vuQgvdxdGdWnOqC7NOXGqmJ+3p/Hj5iO8tWgv/1m4l67RAYztFsUtvVrg4+H8u1cI0TA59xeqABYLXPch6HL44QEoL7fZqgO93bmtdwxfT+jHykmDeWpke0rLNS/M3cm0VSk2244QouFasmQJo0aNMruMCzh/uAMExcKIVyBlGax+3y6biAzwYsIVrfjpj5cT4e9BcmaBXbYjhBC10TTCHaD7H6DdNbDoeeMEJztqGezDoewLh/kUQjRMKSkptG/fnvHjx9O2bVvuuOMOFi5cyIABA2jTpg1r165l7dq19OvXj+7du9O/f3/27NlzwXoKCgqqHC7YDE2nU1gpGP02vNEOtn0LQ+p2SnFdtAj2Zvn+TLutXwin9fMkSNtm23U26wxXv1LjYvv37+ebb75h6tSp9OrViy+//JLly5czZ84c/vWvfzF9+nSWLVuGq6srCxcu5KmnnmL27NnnreOll16qcrhgHx/HH2jRdMIdjLFnAqLhxEG7bqZliDezNxZxqrgUb/emtYuFaKzi4uLo3LkzAB07dmTIkCEopejcuTMpKSnk5uYybtw49u3bh1KKkpKSC9ZR3XDBZgyv0PSSJ6glHE+x6ybahPsC0P+V3xjcPpxhCRFc0TZMgl6ImtSihW0vNQ3n+8wzzzBo0CC+//57UlJSGDhw4AXrqG64YDM0nT73M4Ji4bh9W+7DOzbjozt7Mrh9OIt2ZTDxi410f2EB905bx6x1h8mS8WqEaHRyc3OJiooCjCn3qlLf4YJtqek1JYPjoSADNn8JXW+r14Qe1bFYFMM7NmN4x2aUlpWzNiWHBTvT+XVHOgt3ZaAUJLYMYlhCM65KiCBWTnwSosH7y1/+wrhx43jxxRe55pprqlzmmWee4bHHHqNLly6Ul5cTFxfH3LlzHVypwXmH/K1OQTbMuhMOroCO18Go/4BXkEM2rbVm57GT/LojnQU709l5zBgSoW2EL1clRDAsoRmdowKwyJg1ogmRIX+rV58hf5teuAOUl8GKt2HxS+AbAdd9BHGXO7yMwzmnWLDTCPq1KTmUlWua+XsyNCGcYQnNuKx1qAS9cHoS7tWTcL9URzYaE3rkJMNlj8HAp8DV3ZRSjhcU89vuDBbsTOf3vZmcLiljdNfmvHVLNxl9Ujg1Cffq1Sfcm16fe0VRPWDiMpj/N1j+H0haDHd8axwy6WBBPu7c0DOaG3pGU1hSxpSlybyxYC/uLhZeu7GLtOCFU9Nay8iqldS34d30jpapzN0HxrwDN0+HY5th2yyzK8LTzYVHhrTh8aFtmb0xlWd+3F7vN1qIhsrT05Ps7Gz5jFegtSY7OxtPT89LXkfTbrlXlDAWvEPtPjRBXfxxSGsKS8v4YEkSHq4uPDOqg7RuhNOJjo4mNTWVzEw5q7siT09PoqOjL/n1Eu4VhXeAzN1mV3GWUoq/DG9HYUkZU1ccwNPNwp+Ht5OAF07Fzc2NuLg4s8twOtItU1FYe8jcAw3o30OlFM+OSuD2PjG8vySJtQdyzC5JCNEISLhXFJEARSdhziOQc8Dsas5SSjG+fywA2QXF5hYjhGgUJNwr6nob9LoXts6Cd3vCd/dD5l6zqxJCiDqTcK/IzQuueQMe3QJ9H4Bdc+C93vDNeEjbbnZ1QIPqMRJCNGAS7lXxj4ThL8Fj2+DyJ2DfQvhwAHx1G6RuMLs6IYSokYT7xfiEwpBn4fFtxtmrB1fCx4Ph8+sgfadDS3G1nsS0JfWEHA8shKiRhHtteAXBwL/C49th6PNwdDNMvhKWv2WMU+MAcaE+XNc9islLk3nl590S8EKIi5JwrwsPP2MMmofWQpthsPAf8OlIY2waO1NK8cZNXbmrX0s+WprMpNnbKCuXgBdCVE3C/VL4hsEtX8B1k40zWj8YAOs+tvu3nRaL4vkxHfnj4NbMXH+YR77aSFGpY/5zEEI0LhLul0op6HoLPLgKWvSBn/4EX1wPuUfsvFnFE8Pa8fQ1HZi3LY17p62noKjUrtsUQjQ+Eu71FRAFd35vHEJ5aDVMGQzFBXbf7L2Xx/PvG7uwYn8Wj35t3lReQoiGScLdFpQyTn664xvIT4NNMxyy2ZsTW/DokLYs3JXB/ow8h2xTCNE4SLjbUuxlEN0bVv3XYUfR3NE3BncXC5+vsu+k30KIxkXC3dYG/BFOHDTObnWAUF8PrukSyeyNR8iXvnchhFWtwl0pNUIptUcptV8pNamK52OUUouVUpuUUluVUiNtX2oj0W4kBMfDinccNlbAnf1akl9UyvcbUx2yPSFEw1djuCulXID3gKuBBOA2pVRCpcWeBmZprbsDtwLv27rQRsPiAv0ehqMb4eAKh2yye4tAOkcF8NbCffxr3i7WHsihtKzcIdsWQjRMtWm59wb2a62TtdbFwNfA2ErLaMDfejsAOGq7EhuhbreDT5hxBqsDKKX413WdSWjuz6crDnDzR6vo9dJCnpi1mXnbjkl3jRBNUG1mYooCDle4nwr0qbTMc8CvSqlHAB9gqE2qa6zcvKDP/fDbi8Zoks062X2TnaMD+PyePuQVlrB0bxaLdqXz2+4Mvtt4BHcXC31bhTC0QzhDOkQQFehl93qEEOZSNY1RopS6ERihtb7Xev9OoI/W+uEKyzxhXdcbSql+wCdAJ611eaV1TQAmAMTExPQ8eNCJj/A4fRze7AgdRsH1k00pobSsnA0Hj7NodwYLdqZzIMs4/r5DpD9XdQhnaEIEnZoHYLHItH1CNBZKqQ1a68Qal6tFuPcDntNaD7fe/xuA1vrlCsvswPgDcNh6Pxnoq7XOqG69iYmJev369bX5WRqv+U/Bmg+NAcf8m5tdDUmZ+Szalc7CnRmsP5hDuYZwPw+GdIjgqoRw+rcKxdPNxewyhRAXYctwdwX2AkOAI8A64Hat9Y4Ky/wMzNRaf6aU6gAsAqL0RVbeJMI9Owne7WEMG3z5n8yu5jw5BcUs2ZPBwl3p/L4nk4LiMrzcXHh+bEduTmxhdnlCiGrUNtxr/EJVa10KPAz8AuzCOCpmh1LqBaXUGOtifwLuU0ptAb4Cxl8s2JuMkFbQ8jLY9EWDm0Ip2Med63tE8/4dPdn47FVMv7s3Ph4uLN2baXZpQggbqM0Xqmit5wHzKj32bIXbO4EBti3NSXT/A/ww0ZjoI7Zh7iIPVxeuaBtGdJA3uadLzC5HCGEDcoaqvSWMAXc/o/XewAV5u3HilIS7EM5Awt3e3H2MgN/zE5Q17OPNA73dOX6q2OwyhBA2IOHuCG2GQWEuHGnYXyAHeruRcbKIdxftY11KDsWlcparEI1VrfrcRT3FDwTlAvsWQExfs6up1oiOzViVlM0bC/bCAvB0s5DYMpi+8cH0axVC56hA3F2lPSBEY1DjoZD20iQOhaxo6ggoOQ33/252JTU6XlDMmgM5rE7OZnVyNrvTjLHivdxcSIwNom98CH3jQ+gSHYCbi4S9EI5U20MhpeXuKK0Gw+KX4PQJ8Ao0u5qLCvJxZ0SnZozo1AwwjolfeyCbVUnZrE7O4bVf9gDg7e5CYqy1ZR8fQueoAFwl7IVoECTcHcUv0rguymvw4V5ZsI87IzpFMqKT8TNk5Rex1tqyX5WUzb/nG2HvYw37fq2Mln2n5v4S9kKYRMLdUVw9jevSInPrsIFQXw9Gdo5kZGcj7DPzjLBflZzF6uQcXvl5NwC+Hq5c2S6Mf9/QBR8P+agJ4UjyG+corh7GdWmhuXXYQZifMRvUNV3Ohf3q5GxWJmXx1drDNA/w5O/XVJ4CQAhhTxLujnIm3Msaf8u9JmF+Hozu2pzRXZsDiqkrUhjbLYpOUQFmlyZEkyEdoo7iHWpc5x4xtw4HmzSiPUHe7jz1/TbKyhvW+DpCODMJd0eJ6AgWVzi6yexKHCrA241nRyewNTWX6atSzC5HiCZDumUcxc0TwhOMuVWbmNFdIpm9IZXXf9lDWbmmT1wIHSL95EgaIexIwt2RmneHnT8Yw/+qpjP7kVKKF6/txD3T1vHiT7sA47DJnrHB9IkLpndcMF2iA/BwlYlChLAVCXdHat4dNk6DE4cgqKXZ1ThUi2Bvfn38StJyC1mbksPaA9msO3D87AlR7q4WurUIPBv2PWKC5PBJIepBfnscyTfCuD6d0+TC/YxmAZ6M6dqcMV2NaQePFxSzLiWHtQdyWJuSw/tLknj3t/24WBSdmvvTOy6Y3nEh9IoNItDb3eTqhWg8JNwdycPPuC7KM7eOBiTIx51hHZsxrKMx1EF+USkbDx43wv5ADtNWHWTKsgMAtG/mR6/YYGvgBxPh72lm6UI0aBLujiThXiNfD1euaBvGFW3DACgsKWNrai5rD2Sz5kAO321M5fPVBwGIDfGmf+tQJl3dHn9PNzPLFqLBkXB3pDPhnrXP3DoaEU83l7Mt9YeB0rJydh47ydoDOaw5kMNXaw8R5uvB41e1NbtUIRoUORbNkQKiIaw9LPwHzLgZMveYXVGj4+pioUt0IPdeHs+UuxK5sm0YX609REmZTCwiREUS7o7k6gETfoerXoBDq+H9fjD3ccjPMLuyRuvOvi3JyCtiwc50s0sRokGRcHc0N08Y8Cj8cRP0ugc2Tod3esDS143JPESdDGwXTlSgF5+vOmh2KUI0KBLuZvEJgZGvwYOrIe4K+O2f8G5P2PI1lEsXQ225WBR39I1hVXI2v+1Ox6yZxYRoaCTczRbaBm77Esb/BD5h8P398NEVsPcX40xWUaNbElsQ5ufB3Z+tZ/Abv/PBkiQy8pxvaGUh6kLmUG1Iysth+2xY/CIcT4GYfjDkH9Cyn9mVNXinikuZty2NWesOszYlBxeLYnD7cG5JbMHAdmEyjo1wGrWdQ1XCvSEqLYZN0+H3f0N+OrQZBkOehWadza6sUUjKzGfW+sPM3nCErPwiwv08uLFnNDcntiA21Mfs8oSoFwl3Z1B8CtZ+BMv/A4W50OlGGPQUhLQyu7JGoaSsnN92ZzBr3WEW78mgXEPf+GBu7RXDiE7N8HSTgcpE4yPh7kxOH4cV78CaD6GsGHqMg8FPg3ew2ZU1Gmm5hXy74TCz1qdyKOcU/p6ujO0WxSODWxMuwxiIRkTC3RnlpcPSf8P6T8HT3+iq6TEOLNICra3ycs3qA9nMXHeYn7enEe7nwRf39JHuGtFoSLg7s/SdMO/PcHA5RHaDa96A6Brfa1HJ1tQTjJu6FheLhel39yahub/ZJQlRo9qGuxxC0BhFJMD4uXDDJ5CXBh8PgR8fhoIssytrVLpEB/LNxP64uShumbyKdSk5ZpckhM1IuDdWSkHnG+GR9dD/EdjylXES1NopUF5mdnWNRutwX759oD9hvh7c+ckaFu+WoSCEc5Bwb+w8/GDYizBxhXGo5Lwn4Ze/m11VoxIV6MU3E/vROtyX+6av5+u1hygrlxPIROMm4e4swtvDuP9B4j2w5gNIWWF2RY1KiK8HX93Xl8TYICZ9t40r/r2YD5YkkVNQbHZpQlwS+ULV2RQXwAf9AQUPrAR3b7MralRKy8pZuCudaSsPsio5G3dXC2O6Nueufi3pEh1odnlCyNEyTdqBZTBtFPR9EEa8bHY1jdbe9Dymr0rhu41HOFVcRrcWgYzr35KRnSPxcJXDT4U5JNybup+ehHUfw//9LGPT1NPJwhK+25DK9FUHSc4qIMTHndt6x3B7nxiaB3qZXZ5oYiTcm7qifGMyEDcvmLgcXN3NrqjRKy/XrEjKYvqqgyzalY5Siqs6RDBxYCu6tZAuG+EYcpx7U+fha5zclLUHVr5tdjVOwWJRXN4mjCl3JfL7nwdx3+XxrDmQzS0frSIzr8js8oQ4T63CXSk1Qim1Rym1Xyk1qZplblZK7VRK7VBKfWnbMsUlaTsMEsYaszxlJ5ldjVNpEezNpKvb892DAyguK+fTFQfMLkmI89QY7kopF+A94GogAbhNKZVQaZk2wN+AAVrrjsBjdqhVXIoRr4LFDX76k0z+YQdxoT6M7BTJ56sOcrKwxOxyhDirNi333sB+rXWy1roY+BoYW2mZ+4D3tNbHAbTWcppfQ+EfaQwwlrzYmMpPAt7mHhjYiryiUmasPmR2KUKcVZtwjwIOV7ifan2sorZAW6XUCqXUaqXUiKpWpJSaoJRar5Ran5mZeWkVi7rrda8xeuSyN2D+3yTgbaxTVACXtwnlk+UHyD0lrXfRMLjacD1tgIFANLBUKdVZa32i4kJa68nAZDCOlrHRtkVNLBYY/Ta4+8Dq96GkAEa9JUMF29BDg1pz6+TVdP/nr3SOCqBvqxD6xYfQKzYYHw9b/ZoJUXu1+dQdAVpUuB9tfayiVGCN1roEOKCU2osR9utsUqWoP6Vg+L+MgF/6mjHL03Ufgoub2ZU5hb7xIXz3YH+W7M5gVXI2U5cf4KPfk3G1KLpEB9CvVQj94kPp2TIIL3f5oyrsrzbhvg5oo5SKwwj1W4HbKy3zA3Ab8KlSKhSjmybZloUKG1DKmMHJzRsWPQ+lhXDTNHCRlqUt9IgJokdMEGBM2L3h4HFWJWWzKjmbD39P5r3FSbi5KLq1CKRffAh9W4XQIyZIpvsTdlGrk5iUUiOBtwAXYKrW+iWl1AvAeq31HKWUAt4ARgBlwEta668vtk45iclkqz+A+ZNkiAIHyS8qZV1KDqutYb/9SC7lGtxdLfSICaRffCj9WoXQrUUg7q5y+omonpyhKmr28yRjBMnRb0PP8WZX06Tkni5h3YEcViVnsyopm11pJ9EaPN0sJLYMpl+rEPrGB5MQGSDdOOI8Eu6iZmWl8NUtkLwE7vwB4i43u6Im68SpYlYn57DaGvZ70vMAoyetZbA3bSP8aN/Mj7bNjOvYEB9cXaSF3xRJuIvaKcyFj6+Cggy4dxGEtDK7IgFk5xexLuU4u9NOsjc9j91peaRkFXBmDhF3Fwutwn2NwLcGf7tmfkQGeGL0kgpnJeEuai8nGaYMNlryve+Dfg+BT6jZVYlKCkvK2J+Rz560vLOBvyctj7SThWeX8fN0pV3EuRZ+uwgj9AO9ZeA4ZyHhLuomcy8seRl2fA+unpB4tzE3q3+k2ZWJGuSeKmFPep5xSTvJnjQj+PMKS88uE+HvUaGF70+7CD/aRPjKkTqNkIS7uDRZ+2DZm7B1pnGSU/c7YcCjENTS7MpEHWitSTtZyB5r635PmhH++zLyKS4tB8CioGWIz9nWfbtmfiRE+hMb6mNy9eJiJNxF/RxPgeVvwaYvAA1dboXLn5A++UautKyclOxTFbp1TrI3PZ+U7IKzo1L8dUR7Hhgo73NDJeEubCP3CKx8BzZ8BmXF0HsCXP2q2VUJGztdbPTnf/D7fuZtS2PynT0Z1rGZ2WWJKshkHcI2AqKMMH9sG7QbCWs+hKI8s6sSNubl7kLn6ADevLkbXaMDeGzmZnYePWl2WaIeJNxF7fiGQ4cxxu28dHNrEXbj6ebClLsS8fd0477p62WGqUZMwl3Unp/13/S8o+bWIewq3N+TKXclkl1QxMQvNlBUWmZ2SeISSLiL2vNvblznpZlbh7C7ztEBvHFTNzYcPM7Vby9jytJksvOlFd+YSLiL2vOLBBSkbTO7EuEA13SJ5IM7ehDk7c5L83bR9+VFPDRjI0v3ZlJeLtMxNHRytIyom5l/gANL4fEd4OFndjXCQfam5zFz3WG+25jK8VMlRAV6cUuvFtyUGE1kgJfZ5TUpciiksI/UDfDxYBj2onEGq2hSikrL+HVHOjPXHWb5/iwsCga2C+eWXi0Y3D4cNxnMzO4k3IX9TBttnMn66BZw9TC7GmGSQ9mnmLX+MN9sOEz6ySLC/Dy4sWc0tyS2kLNc7UiOcxf2c9njkHcM1n9qdiXCRDEh3jw5vB0r/jqYT8Yl0q1FIJOXJjPw9SX8uLnyTJzC0STcRd3FDzIuv/wNtsw0uxphMlcXC0M6RDDlrkRWThpMbIg3325INbusJk/CXdSdUnDrDIi9DL6/3zr+jBAQ4e/J0A4RrEnO4VRxac0vEHYj4S4ujbsP3D4LWg2CHx+SLhpx1sB24RSXlbMqKdvsUpo0CXdx6dy84NavoM0wmPsYrJ1idkWiAegVF4SXmwu/7800u5QmTcJd1I+bJ9zyBbS7BuY9CV/fYQwXLJosD1cX+rcK4dcd6eSeKjG7nCZLwl3Un6sH3DwNhjwLSb/Be31g8b+g+JTZlQmT3H9lK7ILirhv+noKS2RsGjNIuAvbcHGDy/8ED6+H9tfA768aIb9zDph0LoUwT++4YN68uRtrU3J47OvNlMlwBQ4n4S5sKyAKbpwK438yhieYdSd8fi1k7jG7MuFgo7s259lRCczfkcZzc3Zg1gmTTZWr2QUIJxV7Gdy/FNZPhcUvwgf9IWEsNO8OzTpDsy7gHWx2lcLO7r4sjvSThXy0NJkIfw8eHtzG7JKaDAl3YT8urtBnAnS63uiD3/MzbJ997nn/KGvQV7gExoJF/qF0Jn8d0Z6MvCJe/3Uvob4e3No7xuySmgQJd2F/PqEw6k3jUpBlDBlc8bJvAWjrl27uftCs0/mBH9bBOCpHNEoWi+LVG7qQU1DM377fho+HK6O7Nje7LKcnA4cJ85Wchoxd5wd++nYozjeeVy4Q1u78wI/oDD4h5tYt6uR0cRnjpq5l46HwHrjiAAASzklEQVTjTL6rJ4PbR5hdUqMko0KKxq28HI4fuLCVX3GKv8rdOq0GyxjzDdzJwhJun7Kafen5TLu7N33j5Q90XUm4C+dUVbdO1h7Q5dBqCNz5ndkVihrkFBRz80erOHbiNF/e15euLQLNLqlRkXAXTUfJaVj8Eqx8Fx7bDoEtzK5I1CAtt5CbPlpJ7qkSBrYLJy7Uh/gwH+JDfYkN9cbP083sEhus2oa7fKEqGj83L0i8xwj3bbOMk6lEg9YswJMZ9/Tlhbk72XT4OP/bevS8c93C/DyItwZ+XKgPcaG+xIf50CLIG3dXOZqqNqTlLpzHJ8Pg9Al4aI0xLLFoNApLyjiUc4rkzAKSs/I5kFnAgawCkrMKyCkoPruci0XRIsiL+DBfa+j7WP8I+BLh74FqAu+7tNxF09PlFvjpCTi6CaJ6mF2NqANPNxfaRvjRNuLCL8RPnCo2gt4a+AeyCkjKzGdlUhaFJeVnl/N2dyE25Ez3jg9x1m6euDAf/JtgN4+03IXzOH0c/tMZWg+Gm6ebXY2ws/JyTdrJQmvw55NsDf7kzAJSj5+i4nA2ob7u1la+EfZxoT60CvOhRbA3Hq4u5v0Ql0Ba7qLp8QqCPvfDstchbbtxMpRwWhaLonmgF80DvRjQOvS854pKyzicc4qkM619a3fPot3pZK0/181jURAd5E2bcF/6xAfTv1UoCZH+WCyNv3tHWu7CuZzKgbe7QvyVxjjzQlSSe7rE2r1j9O0nZRWw69hJkjMLAAj0dqNvXAgDWofQv3Uo8aE+DaovX1ruomnyDoa+DxhDDqcsh5j+MlaNOE+AlxvdWgTSrdLx9Wm5haxKzmLF/mxW7s9i/o40ACL8PRjQKpR+rYywjwr0MqPsOqtVy10pNQJ4G3ABPtZav1LNcjcA3wK9tNYXbZZLy13YzekTRuu98AS4+UB4B4joCBGdICIBwhNkREpxUVprDmafYmVSNiuSsliVlH32qJ3YEG/6tQplQOsQ+sWHEOLr4dDabHYSk1LKBdgLXAWkAuuA27TWOyst5wf8BLgDD0u4C1OdOATJSyB9h/Wy3fjC9Qz/KCPkIzqeu4S0AVd300oWDVd5uWZPeh4rk4xW/ZoDOeQXlQLQvpkfA1qH0r9VCL3jgu1+ApYtw70f8JzWerj1/t8AtNYvV1ruLWAB8GfgSQl30aBoDXlpRtBn7DgX+pl7oNw6z6fFDULbnh/4ER3BL1KOmxfnKS0rZ+uRXFYlZbMyKYt1KccpLi3HxaLoEh3A40PbckXbMLts25Z97lHA4Qr3U4E+lTbWA2ihtf5JKfXnOlUqhCMoBf6RxqXN0HOPl5VA1r7zQ//gSuNM1zO8giD8TNgnGN074R3A3cfxP4doEFxdLPSICaJHTBAPDWpNYUkZGw8dZ+X+bOZsOcojX21i4RNXEubn2C6b82qs7wqUUhbgTWB8LZadAEwAiImRAftFA+DiZg3sBOCmc4+fPg7pOyFjp9Glk74DNs84NwwxCoJioe0IGPw0ePiaULxoKDzdXOjfKpT+rUIZ260517yznBfm7uTd27qbVlNtwv0IUHEkpmjrY2f4AZ2AJdbDhZoBc5RSYyp3zWitJwOTweiWqUfdQtiXVxDEDjAuZ5SXw4mD1lb+Tji2BdZ8CHt/hus+gpi+5tUrGow2EX48NKg1/1m4l2u7NWdIB3PGra9Nn7srxheqQzBCfR1wu9Z6RzXLL0H63EVTcXAlfD/R+AJ3wKMw6ClwNe9fcdEwFJeWM+rdZeQVlrLgiSvx9bDdUee17XOv8QBgrXUp8DDwC7ALmKW13qGUekEpNab+pQrRiLXsDw+sgB53wYq3YPIgY4x50aS5u1p4+foupJ0s5OV5u0ypQc5QFcJW9v4CPz5s9NcPegr6/9GYJFw0WS/9tJMpyw7w3OgExg+Is8k6bdZyF0LUUtvh8OBqaD8SFj0PHw+Go5vNrkqYaNLVHbgqIYLn5+7kp63HHLptCXchbMknBG6aBjd+ahxXP2UQzH8KivJrfq1wOi4Wxbu3dadHTBCPz9zMqqRsh21bwl0IW1MKOl0PD62FHuNg9Xvwfl/YM9/syoQJPN1c+GRcIjEh3kyYvp5dx046ZLsS7kLYi1cgjH4L7v7FOOHpq1tg1l1Gi140KYHe7ky7uzc+Hq6M/3QtqcdP2X2bEu5C2FtMX7h/mXGy05758N9esPlLs6sSDhYV6MW0u3tzuriMZfuy7L49OVpGCEfKToI5j8DBFdDlVrjmdfC4cGo54byy84vqNZKkHC0jREMU0grG/Q8GPmWMX/PRlcaZrqLJcNQQwRLuQjiaxQUG/tUI+ZJT8PFQWPORMXKlEDYi4S6EWWIvg4krIH4Q/PwX+PoOyM8wuyrhJCTchTCTTwjcPhOGvwz7foU3Oxghv+dnKCs1uzrRiMm50UKYTSno9yC0GQYbP4MtM2H3XPAJh663Qvc/QFg7s6sUjYwcLSNEQ1NWAvsXwqYvYO98KC+F6F7Q7Q7j5CjPALMrFCay2TR79iLhLkQt5GfC1plG0GfuAlcvSBhjtOZbXgYW6VltaiTchXAmWsPRjbBpBmz7FopyISAGontCcDwEt7Jex4NvuMz56sRsOYeqEMJsSkFUT+My/CXY/RNsn20cI79zDuiyc8u6+RghHxJ/LvDPBn8zae03ERLuQjQ2bl7Q+UbjAkYffe5hyE6GnAqX9J2wex6Ul5x7rauXNejjLgx+/ygJfici4S5EY+fidi6gKysvM4L/bOgfMK6z98O+BVBWVGE9HpVCP+5cl09AtHHylWg0JNyFcGYWFwiKNS6tBp//XHk5nDxyfmv/zCVpMZSerrAeN2MdFVv64e3lS90GTMJdiKbKYoHAFsYl/srznysvh/y0C0M/OxlSlkNJgbFcZFcY+tyFfziE6STchRAXsljAv7lxib3s/Oe0NoZJSFoEi1+Gz6+DuCuNkI/qYUa1ogry/5QQom6UAr8I6HY7PLIeRrwC6duNKQVnjYOs/WZXKJBwF0LUh6sH9H0AHt0CV04yvqR9rzf87zGZccpkchKTEMJ28jNg6Wuw/lOwuELsAONIG/9o4zogGgKijMMuXR0zrrmzkZOYhBCO5xsOI1+Dvg/C8jfh2FbjRKuCzCqWjbAGfxQEtDgX/AHRxn3vUDkSpx4k3IUQthccB2PePXe/5DScPAq5qRUuh41DMTN3GwOllVSaNNrFw/hC90zYnw3+Cv8JePg69udqRCTchRD25+ZlTDEY0qrq57WG08fPD/+TFW4fWAp5R0GXn/86z8Dqgz8gGvwiwaVpxlzT/KmFEA2LUuAdbFwiu1S9TFkp5B2zBv8Ro+Wfmwq5R4zrQ6uh8ESl9VqMgD8b/NYuoOA4YwYsJw5+5/3JhBDOxcX13ElX1SnKrzr4cw/D0U2w639QVmwsG9beGISt9VDH1O9gEu5CCOfh4WvMWlXdzFXl5XAqCw6uhIXPwRc3QOurjJB3stmu5KtoIUTTYbEYR/R0vBYeWgPDXoTDa+H9fjDvL3Aqx+wKbUbCXQjRNLl6QP9H4I8boed4WDcF3ukOqz8whlFu5CTchRBNm08ojHoTJq6A5t1h/iTjLNvFLxvH6Zt0omd9yRmqQghxhtaw71dY/hYcWgVoCIyB9qOMS0xf08e1lzlUhRCiPvIzYe/PsGsuJC8xJjbxDoF2VxtBHz/QOH7fwSTchRDCVoryjLNod/8Ee3+BopPGXLWthxhB33YYeAU5pBQZW0YIIWzFww86XmdcSoshZZkR9Lt/gl1zrIOkXWYEfddbjeVNJi13IYS4VOXlcHSjcXLU7rnG3LThHeH2mRc/2aoeattyl6NlhBDiUlksEJ0IVz0Pj2yAP8w2zob9eAgc2WhuaaZuXQghnEnroXDPr8aIlp+ONFr0JqlVuCulRiil9iil9iulJlXx/BNKqZ1Kqa1KqUVKqZa2L1UIIRqB8A5w3yKI6Agz74SV75pyrHyN4a6UcgHeA64GEoDblFIJlRbbBCRqrbsA3wL/tnWhQgjRaPiGw/i5kDAGfn0alr3u8BJq03LvDezXWidrrYuBr4GxFRfQWi/WWp8ZaX81EG3bMoUQopFx84IbP4O4K2DL1w7ffG3CPQo4XOF+qvWx6twD/FyfooQQwilYLMbJTtn7jclIHLlpW65MKfUHIBF4rZrnJyil1iul1mdmVjGnohBCOJso61GLRzY4dLO1CfcjQMUDNqOtj51HKTUU+DswRmtdVNWKtNaTtdaJWuvEsLCwS6lXCCEal+bdAQWpjj2vpzbhvg5oo5SKU0q5A7cCcyouoJTqDnyEEewZti9TCCEaKU9/aNYZtn/n0KGEawx3rXUp8DDwC7ALmKW13qGUekEpNca62GuAL/CNUmqzUmpONasTQoimZ9BTkLUHVr/vsE3K8ANCCOEIX91mjC758Dpjsu5LJMMPCCFEQzLiFeNkpvl/c8jmJNyFEMIRglrCFU8ao0juW2j3zcmQv0II4Sj9HzEm5Hb1sPumJNyFEMJRXD3gjlkO2ZR0ywghhBOScBdCCCck4S6EEE5Iwl0IIZyQhLsQQjghCXchhHBCEu5CCOGEJNyFEMIJmTZwmFIqEzh4iS8PBbJsWI6tSF11I3XVXUOtTeqqm/rU1VJrXeOEGKaFe30opdbXZlQ0R5O66kbqqruGWpvUVTeOqEu6ZYQQwglJuAshhBNqrOE+2ewCqiF11Y3UVXcNtTapq27sXlej7HMXQghxcY215S6EEOIiGnS4K6VGKKX2KKX2K6UmVfG8h1JqpvX5NUqpWAfU1EIptVgptVMptUMp9WgVywxUSuVaJwvfrJR61t51WbebopTaZt3mBRPUKsM71v21VSnVwwE1tauwHzYrpU4qpR6rtIzD9pdSaqpSKkMptb3CY8FKqQVKqX3W66BqXjvOusw+pdQ4O9f0mlJqt/V9+l4pFVjNay/6ntuptueUUkcqvF8jq3ntRX9/7VDXzAo1pSilNlfzWrvss+qywbTPl9a6QV4AFyAJiAfcgS1AQqVlHgQ+tN6+FZjpgLoigR7W237A3irqGgjMNWGfpQChF3l+JPAzoIC+wBoT3tM0jON0TdlfwBVAD2B7hcf+DUyy3p4EvFrF64KBZOt1kPV2kB1rGga4Wm+/WlVNtXnP7VTbc8CTtXivL/r7a+u6Kj3/BvCsI/dZddlg1uerIbfcewP7tdbJWuti4GtgbKVlxgLTrLe/BYYopZQ9i9JaH9Nab7TezgN2AVH23KYNjQWma8NqIFApFenA7Q8BkrTWl3ryWr1prZcCOZUervg5mgZcW8VLhwMLtNY5WuvjwAJghL1q0lr/qrUutd5dDUTbYlt1Vc3+qo3a/P7apS5rBtwMfGWr7dWypuqywZTPV0MO9yjgcIX7qVwYomeXsf4i5AIhDqkOsHYDdQfWVPF0P6XUFqXUz0qpjg4qSQO/KqU2KKUmVPF8bfapPd1K9b9wZuyvMyK01sest9OAiCqWMXPf3Y3xH1dVanrP7eVha5fR1Gq6GczcX5cD6VrrfdU8b/d9VikbTPl8NeRwb9CUUr7AbOAxrfXJSk9vxOh66Aq8C/zgoLIu01r3AK4GHlJKXeGg7dZIKeUOjAG+qeJps/bXBbTxP3KDOYRMKfV3oBSYUc0iZrznHwCtgG7AMYwukIbkNi7earfrPrtYNjjy89WQw/0I0KLC/WjrY1Uuo5RyBQKAbHsXppRyw3jzZmitv6v8vNb6pNY633p7HuCmlAq1d11a6yPW6wzge4x/jSuqzT61l6uBjVrr9MpPmLW/Kkg/0z1lvc6oYhmH7zul1HhgFHCHNRQuUIv33Oa01ula6zKtdTkwpZptmvJZs+bA9cDM6pax5z6rJhtM+Xw15HBfB7RRSsVZW323AnMqLTMHOPOt8o3Ab9X9EtiKtT/vE2CX1vrNapZpdqbvXynVG2M/2/WPjlLKRynld+Y2xhdy2ystNge4Sxn6ArkV/l20t2pbU2bsr0oqfo7GAT9WscwvwDClVJC1G2KY9TG7UEqNAP4CjNFan6pmmdq85/aoreL3NNdVs83a/P7aw1Bgt9Y6taon7bnPLpIN5ny+bP2NsS0vGEd37MX41v3v1sdewPjAA3hi/Ju/H1gLxDugpssw/q3aCmy2XkYCE4GJ1mUeBnZgHCGwGujvgLrirdvbYt32mf1VsS4FvGfdn9uARAe9jz4YYR1Q4TFT9hfGH5hjQAlGv+Y9GN/TLAL2AQuBYOuyicDHFV57t/Wzth/4PzvXtB+jD/bMZ+zMUWHNgXkXe88dsL8+t35+tmIEV2Tl2qz3L/j9tWdd1sc/O/O5qrCsQ/bZRbLBlM+XnKEqhBBOqCF3ywghhLhEEu5CCOGEJNyFEMIJSbgLIYQTknAXQggnJOEuhBBOSMJdCCGckIS7EEI4of8HOxEhzr3d4vsAAAAASUVORK5CYII=\n", "text/plain": [ "