TrackPain-Bayesian-Preliminary.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "7cd5c239",
   "metadata": {},
   "source": [
    "# A discrete Bayesian model for pain level prediction\n",
    "\n",
    "## Introduction\n",
    "\n",
    "_motivations and assumptions for using the model_\n",
    "\n",
    "### _why Bayesian network_\n",
    "Compared to linear regression and artificial neural network, Bayesian network uses statistical inference to make prediction, which can overcome the difficulty of missing parameters. \n",
    "\n",
    "### _relevance to our solution_\n",
    "There are two steps in this algorithm: \n",
    "\n",
    "(1) create the Bayesian network by collecting pain questionnaire answers and sensor data from *communicative* patients;\n",
    "\n",
    "(2) use the network to infer the level of pains given the biological data, for *non-communicative patients*\n",
    "\n",
    "Because non-communicative children might not demonstrate behaviors that can be well observed or are not able to self report their pain levels, some parameters are missing. Bayesian network can come in handy.\n",
    "\n",
    "### _variable choice_\n",
    "- To show a proof of concept, each variable in this model will be discrete and binary. \n",
    "- \"True\" means \"parameter measured falls in normal range\"\n",
    "- \"False\" means \"parameter measured falls in abnormal range\"\n",
    "\n",
    "### _where do we get the probability?_\n",
    "Probability values used in this example are artificially generated. But we will show how to calculate these probabilities from collected data.\n",
    "\n",
    "### _future consideration_\n",
    "- We plan to assign more categories to each variable. By having multiple categories for each variable (for example, social behaviors can be specified as five types: unhappy, withdrawn, seeking closeness, hard to pacify), the model can make more subtle predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "dae6c259",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline\n",
    "from pomegranate import *"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "07d3f75b",
   "metadata": {},
   "source": [
    "## Initialize the root nodes\n",
    "\n",
    "### _background_\n",
    "- Root nodes are variables for which we don't care about their probabilities conditioned on anything else. We want to know, given these data from the beginning, what is the final prediction\n",
    "\n",
    "- In other words, these should be data from the sensor because they can be simply measured and objectively obtained, and we want to infer about the pain level using just these sensor data. \n",
    "\n",
    "- For, root nodes of the causal graph, each has a 50/50 chance of being TRUE (normal) or FALSE (false). \n",
    "\n",
    "### _what are root nodes in our case?_\n",
    "- The solution in our device will measure heart rate variability, respiration rate, actigraphy, and skin conductance\n",
    "- To make each variable binary, we used a unique numerical range to binarize collected data. That is, data points lie inside of the range are considered normal, and those outside considered abnormal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "a594bccd",
   "metadata": {},
   "outputs": [],
   "source": [
    "HeartRate_var = DiscreteDistribution( { 'True' : 0.5, 'False' : 0.5 } )\n",
    "\n",
    "RespRate = DiscreteDistribution( { 'True' : 0.5, 'False' : 0.5 } )\n",
    "\n",
    "SkinCond = DiscreteDistribution( { 'True' : 0.5, 'False' : 0.5 } )\n",
    "\n",
    "Actigraphy = DiscreteDistribution( { 'True' : 0.5, 'False' : 0.5 } )\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c323856c",
   "metadata": {},
   "source": [
    "## Conditional dependence (single dependence)\n",
    "\n",
    "\n",
    "### _background_\n",
    "- in a causal network, non-root nodes have conditional probability, meaning the probability of them happening is with respect to (conditioned on) their parent nodes \n",
    "- some events depend on only one parent node (single dependence)\n",
    "- more strictly speaking, we put them as non-root because we want to know how these factors play into the final prediction.\n",
    "\n",
    "### _in our case_\n",
    "- here we model parameters that can be collected by the phone app as non-root nodes\n",
    "- recall that our phone allows the third-party to put in their observations of the patients, which includes the 4 parameters: sleeping patterns, eating patterns, social behaviors, and verbal activity.\n",
    "- All 4 parameters here each can depends on some or all 4 types of sensor data (root nodes in the last section), ideally there will be 4*3 =12 total edges in the graph at max. \n",
    "- Here we only show some examples of singular depedence (meaning the node only is conditioned on only one event)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "91316232",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "# P (sleeping = normal or abnormal | heart rate variability = normal or abnormal)\n",
    "sleeping_HR = ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 0.2 ], \n",
    "     # given heart rate var is normal, the probability of sleeping is normal\n",
    "     [ 'True', 'False', 0.8 ],\n",
    "     # given .... normal, the prob ... abnormal\n",
    "     [ 'False', 'True', 0.01 ],\n",
    "     [ 'False', 'False', 0.99 ]], [HeartRate_var])\n",
    "\n",
    "# P (social behaviors = normal or abnormal | respiration rate)\n",
    "social_RP = ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 0.2 ], \n",
    "     [ 'True', 'False', 0.8 ],\n",
    "     [ 'False', 'True', 0.44 ],\n",
    "     [ 'False', 'False', 0.56 ]], [RespRate])\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "017ea964",
   "metadata": {},
   "source": [
    "## More conditional dependence \n",
    "\n",
    "### _background_\n",
    "- Sometimes the probability of an event happening can conditionally depend on more than one event. Here we look at parameters that depend on more than one other parameters. \n",
    "- The parameters that are collected in the same phase might also depend on each other\n",
    "\n",
    "### _what are some examples of multiple dependence in our model?_\n",
    "- eating behaviors depend on both heart rate and skin conductance\n",
    "- verbal activities depend on respiration rate and  and heart rate\n",
    "- facial scales depend on verbal activity and actigraphy\n",
    "- emotional scales depend on eating, sleeping, and actigraphy\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "9ab061d5",
   "metadata": {},
   "outputs": [],
   "source": [
    "# P (eating | respiration rate, skin conductance)\n",
    "eating_RP_SC =  ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 'True', 0.96 ],\n",
    "     [ 'True', 'True', 'False', 0.04 ],\n",
    "     [ 'True', 'False', 'True', 0.89 ],\n",
    "     [ 'True', 'False', 'False', 0.11 ],\n",
    "     [ 'False', 'True', 'True', 0.96 ],\n",
    "     [ 'False', 'True', 'False', 0.04 ],\n",
    "     [ 'False', 'False', 'True', 0.89 ],\n",
    "     [ 'False', 'False', 'False', 0.11 ]], [RespRate, SkinCond])\n",
    "\n",
    "verbal_HR_RP =  ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 'True', 0.96 ],\n",
    "     [ 'True', 'True', 'False', 0.04 ],\n",
    "     [ 'True', 'False', 'True', 0.89 ],\n",
    "     [ 'True', 'False', 'False', 0.11 ],\n",
    "     [ 'False', 'True', 'True', 0.96 ],\n",
    "     [ 'False', 'True', 'False', 0.04 ],\n",
    "     [ 'False', 'False', 'True', 0.89 ],\n",
    "     [ 'False', 'False', 'False', 0.11 ]], [HeartRate_var, RespRate])\n",
    "\n",
    "facial_verbal_SC =  ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 'True', 0.46 ],\n",
    "     [ 'True', 'True', 'False', 0.54 ],\n",
    "     [ 'True', 'False', 'True', 0.39 ],\n",
    "     [ 'True', 'False', 'False', 0.61 ],\n",
    "     [ 'False', 'True', 'True', 0.96 ],\n",
    "     [ 'False', 'True', 'False', 0.04 ],\n",
    "     [ 'False', 'False', 'True', 0.03 ],\n",
    "     [ 'False', 'False', 'False', 0.97 ]], [verbal_HR_RP, Actigraphy])\n",
    "\n",
    "emotional_eat_sleep_ag =  ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 'True','True', 0.44 ],\n",
    "     [ 'True', 'True', 'True','False', 0.56 ],\n",
    "     [ 'True', 'True', 'False','True', 0.90 ],\n",
    "     [ 'True', 'True', 'False','False', 0.10 ],\n",
    "     [ 'True', 'False', 'True','True', 0.38],\n",
    "     [ 'True', 'False', 'True','False', 0.62 ],\n",
    "     [ 'True', 'False', 'False','True', 0.29 ],\n",
    "     [ 'True', 'False', 'False','False', 0.71 ],\n",
    "     [ 'False', 'True', 'True','True', 0.96 ],\n",
    "     [ 'False', 'True', 'True','False', 0.04 ],\n",
    "     [ 'False', 'True', 'False','True', 0.52 ],\n",
    "     [ 'False', 'True', 'False','False', 0.48 ],\n",
    "     [ 'False', 'False', 'True','True', 0.77 ],\n",
    "     [ 'False', 'False', 'True','False', 0.23 ],\n",
    "     [ 'False', 'False', 'False','True', 0.19 ],\n",
    "     [ 'False', 'False', 'False','False', 0.81 ]],[eating_RP_SC, sleeping_HR, Actigraphy])\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "08e51e42",
   "metadata": {},
   "source": [
    "## Last layer (leaf nodes)\n",
    "\n",
    "Now we construct the leaf nodes which are the end of this network. Since we would like to know the probability of pain or not pain given all parameters provided above, we will set pain level as the child of facial and emotional scale and verbal activity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "id": "c574b042",
   "metadata": {},
   "outputs": [],
   "source": [
    "pain_facial_emo_verbal = ConditionalProbabilityTable(\n",
    "    [[ 'True', 'True', 'True','True', 0.44 ],\n",
    "     [ 'True', 'True', 'True','False', 0.56 ],\n",
    "     [ 'True', 'True', 'False','True', 0.90 ],\n",
    "     [ 'True', 'True', 'False','False', 0.10 ],\n",
    "     [ 'True', 'False', 'True','True', 0.88],\n",
    "     [ 'True', 'False', 'True','False', 0.12 ],\n",
    "     [ 'True', 'False', 'False','True', 0.29 ],\n",
    "     [ 'True', 'False', 'False','False', 0.71 ],\n",
    "     [ 'False', 'True', 'True','True', 0.76 ],\n",
    "     [ 'False', 'True', 'True','False', 0.24 ],\n",
    "     [ 'False', 'True', 'False','True', 0.52 ],\n",
    "     [ 'False', 'True', 'False','False', 0.48 ],\n",
    "     [ 'False', 'False', 'True','True', 0.77 ],\n",
    "     [ 'False', 'False', 'True','False', 0.23 ],\n",
    "     [ 'False', 'False', 'False','True', 0.19 ],\n",
    "     [ 'False', 'False', 'False','False', 0.81 ]],[emotional_eat_sleep_ag, facial_verbal_SC, verbal_HR_RP])\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "61c3432f",
   "metadata": {},
   "source": [
    "## Create and plot the network\n",
    "After we defined all events and their conditional probabilities, we can construct and plot our network!!\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "9ad7e0f1",
   "metadata": {},
   "outputs": [],
   "source": [
    "s0 = State( HeartRate_var, name=\"Heart Rate Variability\" )\n",
    "s1 = State( RespRate, name=\"Respiration Rate\" )\n",
    "s2 = State( SkinCond, name=\"Skin Conductance\" )\n",
    "s3 = State( Actigraphy, name=\"Actigraphy\" )\n",
    "s4 = State( sleeping_HR, name=\"Sleeping Behaviors\" )\n",
    "s5 = State( social_RP, name=\"Social Behaviors\" )\n",
    "s6 = State( eating_RP_SC, name=\"Eating Behaviors\")\n",
    "s7 = State( verbal_HR_RP, name=\"Verbal Activities\" )\n",
    "s8 = State( facial_verbal_SC, name='Facial Scale' )\n",
    "s9 = State( emotional_eat_sleep_ag,name=\"Emotional Scale\")\n",
    "s10 = State( pain_facial_emo_verbal,name=\"Pain Level\")\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "id": "368b2ee6",
   "metadata": {},
   "outputs": [],
   "source": [
    "network = BayesianNetwork( \"Pain-Track\" )\n",
    "network.add_nodes(s0, s1, s2, s3, s4, s5, s6, s7,s8,s9,s10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "id": "b1b016b5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAV0AAAB7CAYAAADXClMoAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABQO0lEQVR4nO29d3Bc15Xn/2l0N7oBNHLORCCYQTCJmVSWFSxZkpM8ljzyOMxu7c7W2Oup9dTUr9ahalxe1Xi2XDs7O+PRjCxbybIkK1JUoCiRIkgwAwQJIhC5gQYa6IyO7/3+uLwPDRKJJAgC5PtWdTXIBvrde+6533vuueeca1BVFR06dOjQMT9IuNEN0KFDh45bCTrp6tChQ8c8QiddHTp06JhH6KSrQ4cOHfMInXR16NChYx6hk64OHTp0zCNMM3x+a8STXU3YnMEw9+3QMTvMdrz0Mbo6XMl80GU8FaYUzEyku/igquKlKBCLQSAAwaB4BQLg8UAoBOGweEWjE1+KIr4nFhPvBgMkJIy/m0zjL7MZEhPFKzkZ0tLAahWvpCTxLv9WvnRMDjnR5bj5/WK8/H5wucR7KCRe0ej42MmxVtVxWcuxkS85JhkZYLOJsbLZxGfx43MzI35ehEJCnmNj4POB2y3mh5wX0ShEIuPzQco3fi5IGZtMQv8tFvFKS4PUVKH/KSnjc0D+nQ4MMyRHLExLVypQJCIm5sgIDA7C8DB4vUKZ5CSUE85qnUiMkjDjSdRkAqNRPEMqiVQ4VRVkEE/QkYhQ0khEKLF8dig0/g7iWTabeHZODmRlQXa2aIfJdGsqo6II2blc0N8PfX3i52BQyCMpSYxXSoogy5SU8YktydRkmjihJUHEYmJM5PjI8XC5BMlIMo/FxPdkZIhxKSyE3FzxrMVKElJXQyHR38FBIV+XS/RbVUWfk5OFjG02SE8XP0sD4lL5SrKVpK0oE4k5EhHPCwaFUSPnQSAg/t9gEM9LTxfyzcsT8k5KWrxynhlTdmpxkK4kWLcbenqgs1MoUSwmlCYzUwxmVpZYZW22iQoz34Mab7Upimi71ytedjs4naL9kYhQxpwcKC2F4mLRdkn8NxNUVRDgwAC0tAiSjcXEmOXnQ0GB+Nlmu/6WUbzVJ4l/cBB6e2FoSHyWnQ2VlbBkidCphUwMiiIWkZ4eOH9eGB8Gg5Bnbq6QbXq6WPTlvIDr26dL54DfL+av3S5eo6PidzIyoLwcKirEzzcPCS8y0pVt8vmgtVUokscjFKekRBBUdrZYKeUALaaBildIn09M+AsXhDJGIsLiWrFC9NNsXlx9i4cktt5eOHZMkEF+PtTUjC8wC2lrL8clHBbke+4cdHWJMVizRoyJxbIw2quqgsQaG8UcMRgEeVVWCkvSal1Yso1HvJydTmhrE4ZUKCQWubVrhSGyUNs/OywS0pWTtKMDjhwR25Xqali2TAyCtAAX70BMDTkOwaDYDp45I0i4pAQ2bxaLzGLpt7QkW1vhwAFhYW3aJPpiNC6+foyOwsmTghyWLYMtW24c+aqqWKT37xcL9po1sHy5sMZh8cg2HlL3x8agvV3IWlFg505BwouTfBcB6aqqIJu9e4VFu3WrsIoSLo9qU1WVgYEBRkZGKC8vp7W1lerqamw2G4ZJBkf20WAwMDw8THd3NxkZGRQVFWG1Wi/7nXg4HA4cDgfLly+nr6+PtLQ0MjMztc+dTicGg4GsrCxUVaW/v5+UlBQyMjJwOBwYjUYGBwcpKChAURSys7Mnfc6k8ohEBHEdPiy2iXfdNdG6X6gIBuGdd8TEufNOsc29pM2KotDe3k4oFCI9PZ3i4mISJhlriUAggMvlorCwcILs5Lj19vaSl5eHxWKZ9O9DoRDNzc2kpKSQk5NDZmam9j1Tjf0ESPfI0aPQ1AQPPQRFRfM7FpEIfPyxmCd33AFlZZPOD0VR6OjoIBwOk5ubi8fjoaqqilAohMPhoKSk5DIZjoyM4HA4MBgMlJeXk5SUNGkTRkZGCIVCFBYWXvbZrOQ4G6gqOBywb59YpB98cHHo/URM2diFEaerqsKF8O67cP/98OijYos9zSQcHBzkzTffJCEhgeeee45IJILH4yEcDhMMBnG73YTDYdxuN42NjbS2tqKqKh6Ph1deeYX29nZ++9vfEo1GGR0dJRAI8MknnxAOh3G5XEQiEQCi0SjPPvss0WiUs2fPEggEJnx3KBQiFosxNjaG1+uloaGB06dP4/P5iMViRCIR3nzzTcLhMKFQiE8++QSXy4XX6yUQCBCNRifvoMEgDjVWrYJvfUtYib//vfCNLWTEYvDaa2Kr+/jjU1roBoOBI0eO0NTUxMsvv0xTUxOhUAi3200sFsPlchEKhfB4PIyMjKAoCtFolLGxMdxut/be09PD6dOniUajKIqCz+fD5/MRDofxeDx4vV5UVcVgMPD888/j8/n45S9/STAYxOPxEAgEOHToECMjI/h8Pvx+P5MaIgaDsG63bRP6+eabYms8X4jFxDOtVnjySWEBTjE/HA4Hr732Gj6fj9bWVt566y28Xi/BYHCCDEMXD3r7+vp49tlnSU1N5fz58zgcDk3+Ut9DoRChUIjGxkYOHjyIx+PR5oGU+eHDhxkeHtbkGD9WkUgEl8tFLBbD4/EQDAYnlzMIWefnw1e/KnYWL70kFrybBAsjZCwahU8+gW98Y3ybNA0MBgNWqxWHw8GRI0cIBoOcOXOGgYEBfD4fxcXF9PX1kZ2dTVtbG+FwmPXr17N06VIsFgsej4fm5mZqa2s5ffo0PT09GI1GDh8+DAhCN5lMPP744xQWFlJeXs7Ro0exWq2cOXOG7u5uMjMzsdvt5OTkkJOTg6IodHd3k56eTlNTE/v372fZsmUsWbIEs9lMW1sbQ0NDHDlyhJSUFJqbm1mxYgUrV64kdaY+G41QVycm2aefioVpoa763d2CnDZsmLaNBoMBi8VCS0sLRqMRm83Gu+++S2dnJytXriQWi5GZmcnHH39Mamoq5eXluN1u+vv7WblyJSaTiYGBAWw2G9FolGAwyPLly3G5XNjtdqxWq0YGf/VXf4XRaMRgMHDs2DEqKipwu92cOnUKu91OT08PAHa7na6uLv7zf/7P2g5okoaLXcc99wjXySOPzM9YSH//zp3TGiMAGRkZhEIhPv/8c7785S/z8ccfc/LkSQoKCti/fz/BYFDbcT3xxBPU19dTVlZGcXGx9jvJyckcP35cs2qj0SgpKSmYTCY8Hg+/+93vKCwsJDk5WZtTJ0+e5OGHH8bhcNDd3Y3BYCApKQmfz0dpaSk2mw2/34/f72dgYIC//Mu/xDjdoXFCgnCfDA2J3cWGDXMp0RuGhWHpRqNCwMnJV/RnRUVFbNiwgbS0NE6cOEFlZSU7d+7E7/czMjKCqqosW7aMVatWkZWVpf1deno6Pp8Pm81GIBBgbGyM0dFRMjMzOXPmDGvWrOG2227Tfv/BBx/kd7/7HcXFxfj9fs2lsGrVKioqKjRFGhoaIhqNctttt5Gbm6tZtwBWq5VIJEJmZiZr1qzB4XBobZgVDAZxQOL1XpGM5h1+vziFniVycnI0i6m7u5svfelLlJWVcerUKUZHR8nJyWHdunV4PB5CoRBJSUls376dcDjMyMgIBoOBnJwcEhISOHfuHKqqkpeXx8DAAEVFRaiqSuxizLXZbMZisWAwGIhEIgSDQex2O2lpaQwPD5OQkMCXvvSlad0cgBiLjAwREjVfCASE220WBO/3+/ne976HwWDgwIEDdHV10d/fT2ZmJsFgkMTERFasWMHo6CiqqpKTk8PQ0BCqqqIoCsePH9d0PRKJsGLFCrq7uzl//jxlZWUkJiYSDAZRFIXTp08TCoW44447yMzMJCUlhUAgwMDAACaTiZqaGhwOBydPnqSurg6Xy0VWVhb333//1JZuPAwGsVta6Hp/BVgYpGu1irCWzz4T26gZBkNVVcbGxkhPTycSiZCTk8PatWv57LPP8Hg8OJ1OcnNz8V/cii9ZsgSHw4GqqgQCAfLz8/nzP/9zjh49Sn9/PwaDQfPDrly5ko8//phgMAgIi6yiooJHHnmEgoIChoeHyc/PZ2xsjFgshqqqGI1GhoaGKCgo0BS4rq6O/Px8IpEIqampBINBkpOTKS4uZnR0lM2bN1NUVDQ7/5eqCjLbswc2bly4Vi6Iba+MNplmHBVFwWg0UlZWxiOPPEJjYyM2m42GhgZCoRA7d+4kJSWF/v5+AoEAa9euJSkpCZvNhs/nY3h4mLy8PNLT0xkdHSUxMZG1a9cSDoex2Wxs2LABs9ms6Yi07nbt2qXtWJxOJ3l5eZSWlpKenk5XVxednZ3Tk66M196/H1avnnv5TYXSUhFJMTw84/wYGxvj888/p6ioiJqaGrZu3Uo4HObs2bOkpKRgs9kIh8OkpaWhqipbt24lPz+f9957j+PHj3P33XfT0tLCli1bKCsrIxwOU1JSQmlpKc3NzWRlZZGbm8vY2Bi33XYbzc3NnD59mqVLl9LZ2anNhYSEBCKRCAUFBVRVVfH6669TW1vLiRMncDgcMy9uUu8PHxZutpsEC+cgLRKBDz4QW4k77xQhRdOcWk7Wbum7u1bI77n0oGWu0NzczNmzZ3nkkUcwm83TNUTIpalJKN6uXbBy5cImXVUV4W9798K99055+jyVTOPHsKurizfeeIP/8l/+y/Tb0Em+92r0IP7Zk/69qgqL6913xeHg3XfPuNWfM6iqcDH86U8imqW2dspIECmD2c4Hg8Gg7QiknBVFmZkUrxKXzq9JfmFcjz78UOj98uULW+8vxyKIXoDxCIZPPxUn4GvWCEd6SspiDRuZFJFIBEVRSExMvFzxpCU1PAynTwvFq6gQBzhSDgsdMsxqzx7x7+3bxUHgFQa+Swt1qqiUeYOMiT1yRIzHrl1CL+eLcOPbEQiI8w+7HdatE7HDi+9kf3JII6O9Xcg6KUksbJmZi7F/i4R0tadenLRNTSI2MiFBbK8qKoQbYiEHfl8pJMn6fGLBaWkR1n5qqti+VleP93exQVGgrw/10CEiTif+igrSamsx5uYu7PRnGS/u9wsCaGoSZLBundjm3uiEFbkIHD8u2peeLizB8vLxjL6FKtt4SIt2bEwk0Jw5I3S/pETEdefmLo5+TI5FRrra09Vxv05Xl1Cw4WHxWVqaiJOMz+NOTFy4ZCzJNRwW/s7hYWGtOBzCerFaRX8qK0W43ELJfLpCqKpKNBrVTqg7OzsZHBggODpKjqryUH4+Fr9fjF95uZhg2dmivzeKLGQMrscjxqSzU0x+k0ks9CtWjIe+LaQxid8VnT0r5kg0KizDsjKhR5mZEwsv3ci2ypoQbrdIB+/uHp/PxcVi91BcvLAX5NljkZLupZADF18oxeEQ8ZKyyI3FIrbhmZnCWkxNnVjk5mJhG/WiP+xahlYFUBQMshCOLLLi94sJ7POJd1lsJBoVfrjUVLFQFBSIeMS0tEWXbRfvN/T5fAwMDNDb24vdbsfn85GcnExhYSElJSUUFBSQmpqq+QgN0agoUtTVJWowjIwIyzIpSYxbdraIDpCEEV+Y6Gog61/IwiweD4yOorpctHd0kOnzkWWzgc2GobBQ+KHz88WzYXGMiZzHoZCYD93dgthGRwUxywpgWVnCMpYV12SBm4vFnuSEn22Ptd+PL4ITDgt9Hx0Vsh4ZEUQ7NiZkmZ4u5FtWtvBTlq8eNwnpTob49kciYmD9fjHQPp84+PB4hI/4YnUwVVHoHhxkzGRieV7exEpKMLGeg/x+SfhxlZY8wSBnLlxgU3U1xsREDLLcXUqKUHBZWSwjQ4TDWa0Ti9ksIiVTVVWLGhkaGqKzs1OLizYajeTn51NWVkZRURFpaWnagcysMu9gfLEaGREvl0tM2mBQfCZLbhqNEwviGAxEFIVAOEy61XrZGE0oSygry10cEzUzk/ZgkIMNDSSmprJl61bKyspISEi4sT7kuYKUbTQ6XtbU6Rw3BHy+8cVIVVGBYDSKIxikLCsLQ7x1PMlcUBWFjsFBStLSsEjrVM6B5GSh9+np40SfkiI+k7gZZDw1bmLSvUKoqsqZxkY+P3iQrzz+OJk22/jElO/iF8V7vNLJ18VJHwP2ffYZgw4HX3r0UZKTk2+KySoJNhQKMTw8TG9vLz09PYyMjGC1WsnOzqakpISSkhLS09Mxm83Xt99ybGIx8ZI/X5z8/XY7R48e5eEvfvGyMcJoHH9N0UZFUbDb7Rw6dAi328369etZuXLl5AedNyMuyleJRHjt1VepWbaMNcuXC+s1fj5cMhfUhASOnjrFwMAADz78MAk3h1tgrqCTLojJVV9fT3t7O48//jjJV5iMMRlUVeXMmTMcOHCARx99lLy8vEU1USXBRiIRRkZG6Onpoa+vj9HRURRFISsri/LycoqKisjMzCQxMfG6hRJdLfr6+jh8+DCPPfbYNX2PTBM/evQobW1tVFdXs3HjRtLS0hbVmF4NVFXl+PHjdHd388gjj8x6jGOxGH/4wx9Yu3Yty5cvv+nldAXQSTcajbJ3715CoRAPPvjgnFpnqqpqtSB27NjBihUrFqzyqaqq1Y2w2+10dnbidDpRVZWMjAwKCgooLy8nOztby95aqH2RmCvSlZBWfnNzM8ePHyczM5OtW7dqAf83IxwOB6+//jpPPfXUlMVupoLH4+GFF17gm9/85uwzLG9+3NqkGwqF+NOf/kR2dja33377rAPtrwQy2+3111+nqKiIXbt2Ybrag585bFM0GsXj8dDX10dfXx8Oh4NgMEhaWpqWZZSTk6O5RhY6wU6GuSZdCZkW29XVRX19PeFwmM2bN1NdXT19UssiQygU4vnnn+f++++ffZZkHFRVpampidbWVh599NFFqUPXAbcm6cqT9T/84Q/U1tayYcOG66oQMqvno48+YmRkhEceeYSkpKTrroRyDBVFwePxMDg4SHd3N3a7nUgkgtVq1Qr35OXlTUg2uBkmyPUiXQkpX6fTSX19PX19faxevZp169bNy/heT6iqynvvvUdGRgZbt2696r4oisLLL7/Mxo0bqa6uXtQymSPceqSrqipOp5NXX32VO+64g5qamnlTBFVVOX36NPX19Tz22GPk5OTM6bOlHzYQCDA0NMSFCxfo7+/XQrWKioomD9W6SSfC9SbdeEi5nzx5kjNnzlBUVMTWrVvJyspadPJVVZWWlhaOHz/O1772tWveAY6OjvLKK6/w9NNPk5iYOEetXLS4tUhXVVW6urp47733ePjhh69qyzQXbbDb7bz11lvs3r2bZcuWXXU9AOljdDqdWsUol8uFyWQiNzeX0tJSSkpKSE1Nvf6RBAsQ80m6EtJ109raSkNDA2azmS1btlBeXr4oQs5UVcXtdvPiiy/yzW9+c+byorP8zv3792M0GtmxY8eCl8F1xq1DujKa4PPPP+crX/kKGRkZN2zwVVXF7/fz2muvUVpayq5du6a1JuRYhMNhRkdH6evro6urSyv+nJubq4VqZWZmYjabb9qDnSvBjSDdeMiQs/r6elwuF+vWrWPVqlULOuQsGo3ywgsvsG3bNqqqquasnaFQiH//93/XD9VuFdKVIWEdHR089thjC8LfJi2iDz/8EI/Hw8MPP4zVatUqO0UiEbxe74RIgmAwSGZmJsXFxZSWlpKbm7toIgluBG406UrIkLOGhgY6OjqoqqpakCFnqqry6aefEolEuOuuu+bc9XXixAmGhoa49957F1S/5xk3P+lez5Cwa4V0EZw4cYIDBw6wa9cuhoeHGRwcxO/3k5ycTElJCWVlZeTl5ZGSkqIT7Cwh76U7fPgwjz76KHDjfdcyLK+pqYmTJ0+SkZHBli1btPvdbkT74ud5V1cX+/bt45vf/OZ1icKIRCL8+7//O0888cScuC0WKW4+0pUpqSaTCVVVeeONN8jJybluIWFX0i75HggEGB4e1vywXq+XsbEx2tvbufPOO9m9ezdWq1VLOLjRZLEY8e677/Lcc88xNDTEnXfeyd/8zd8smEMcGXLW3d3NoUOHiEQibNq0iaVLl2IymeZ1vBVF4e2336aqqopPP/2Ur33taxNuU5lLqKrKsWPH8Hg83HHHHbeqXk/Z6YVxR9pVQFEUfvGLX2CxWMjPz6euru66h4RNhksPunp6eujt7cXj8WAwGLSDrpUrV2ops16vl9dee41jx46xd+9ebr/9du6///55bffNgrS0NN566y3GxsbYsGHDDY+NjofBYMBoNFJRUcGSJUu0kLMDBw6watUq6urqtPjoeOPneuhwIBDgpz/9KQMDA/z4xz8mLS1tzp8hYTAYWL16Nc899xw7duxYMIvgQsGitHRVVeXUqVPce++9uN1ufv7zn/ODH/xgXixc6YcdHR3V/LCjo6PEYjEyMjK0lNns7GztIGWySRQMBvnBD37Ab37zG+rq6njnnXfmPLTsVoDX6+XOO+/kzJkzvP/+++zcufNGN2layB3QqVOnaGpqori4mM2bN5OVlcU777zD+vXrr0u0zbFjx7j99tvx+/2sX7+eP/7xj5SXl8/pM+Khqirvv/++ZnDcgnp9c1m64XCY//W//hfDw8PYbDZOnz6N3++f89VbHoLFH3QNDw8TDodJTU2lpKSEdevWaQddV+IicDgcHDt2DBAT4p/+6Z/4u7/7uxvqGlmMsNlsfPGLXwSgrq7uxjZmFpD38W3dupVNmzbR1tbG22+/zdjYGM888wzl5eX867/+K5WVlXOapr5v3z6i0Shf//rX+elPf0pZWdmcfPdUMBgMrFu3jo8++oiVK1de12ctNlxf0o2vUBQKiTqbweB4qcVgUPyf/P+LJeaIxUS1IlnGLzFR1ANNTES1WPi0qYkzn37KD59+mq/82Z9Ru3EjltTUyyuDTdksleHhYWKxGAUFBRMyuvx+P3a7ne7ubgYGBgiHwyQlJVFQUMCyZcvYtWsXycnJV+eDjdtVlOTl8daLL3Jg3z7++Pvf88azz/LwunXUFRVhCIfHa79GIuOVtS4tb2gyCbmYzeJlsYiX1SrKSsoawmbz4igpGa8vUgbBoCjPOTYmfpYykfqkKDxkNmOrrMT20UdCLlIWiYniZbWKV2qqKDkoP48Pt5tnmRgMBsxmM8uXL2fZsmX8wz/8AxcuXKCjo4M/+7M/49lnnxU1PKQswuHxWsBjY+JnKYv48pey+pqsqmY2EzGZOLt3L79+6in+7OmnscZiGPr7x0stynrFcyyLnJwcAoEAgUCAlJSUOfnOGSF1SFGEXGRdX59PyE3qUTQqfsdoFH2XpVeTk4WeZGSI+s7yMoE51I+5cS/Ed9TrhcFB8XI4xmt2GgxicJOSRKeSksYLVEtSNZvHy/FJ8pW3LVxUMDUY5OzwMNmBALmKgiESwSC/32IRBbBzc0WB8MxMIURZiu7iqfLevXv56U9/ysMPP8yf//mfa35Yt9uN2WwmNzeX8vJyiouLSUtLu7KMrnh5ykkyOCjqmDqd4/JQVVFM3WwmZrXSB4STkqjOzsYgJ4GUiSyZFy8XVR0vGC0JSsppbEw8R05KqVwWi6htmpUlZJSVJSaenHDzRTzx+uLziVsapL54PKIfsv6txTJOlpPpi9FIWFUZi8VIN5vHC5ZLUo6XiSTvcFjIUNZ9lbeP5OeLuq/SBzkf8rh42Pbsv/4r+995B1dbG16nk9KVK/nF175GkTxwM5vHF46p5k58CUu5UEejhP1+Br1eSsxmDKHQ+EIWCIzLwmgU35mTM15gXxLPVcpCVVX27t1LWVkZK1asmFu5jT9EvAeD4lKDzk7xHgyKPiUnCx6w2URfpOwuFm2XNbY1Qg4EhGxcLvFvEDIvKRHF7XNzx2sCTy+T6xC9ICf/0JC4crunR3Q0KUkMWm7u+K0IspMTnnz1Cq2q6sQeSWsgGBSr2tCQuJHA6RSf5eaiLFlCi6ry97/6FX/84x8JBAJs2LCBH//4x1RUVFBcXExGRobm9L8iK1bKYmQEOjrEfU8ul1hAbLZxgsvJEQMoJ8klz1Cv9LnTtSce8QuXyyXaOTAgrkoJBsX4FBSISv7FxaLNc13JX8poeBhaW4W++P3i2dnZgvjy8sREl1ZqPKZpy6xuAJ5MJqGQmGTDw4L0+/vFIpCYKK5OqqgQk20ub5CIL9putwtZ9PejhkKQkkIsK4toTg7RggIS09IwW60YLg3ruoqCNOLPDBPbIKEoQhbBoJgzQ0NChz0e8XluriCc8nIxn69ANy5cuEBjYyNf/OIX59avq6qivWfPirvVQiHBN/IKKJvt6gumT2Y49faK+xpHR4W+1tYK/Zi6hvAcka4sJt3XBydOCEXNyhL3epWVCSthIRUylm6N4WH8587x2337OLxnD909PdjNZgpqanjzzTevLpZQDnp7OzQ3CzLLyBCDXlw8fu/XXJPXXEISYSAgxvLCBTG2iiLGc9UqochXe7+W1JfeXjh5UjwjIwOqqsb15UZf8hgP2d5AQCxK588LIk5OFrJYtuzqLwmVhkF7OzQ2iomcmyvmTmmpWIynKbR+QyDJeGhIGBMXLoj/q64WN3XP4pZev9/Pyy+/zNNPPz13BoXHAwcPCr1aulSMTXb29b8HTlWFVTw4KPivv19cHrtx4/j9jOOYA9KNxcSqcviwWO3WrxfKspAmzQxQFQXF5yPW1ETwyBF8RiN5X/kKpvz82fdBVYUSfv652A4vWSIGPS9vYS04VwtFEduqCxfELbg+H6xdK1b2yxVrckgyP3MGGhoEuS5CfdFI2OkURNnWJsZ7+3bhlpmtLMbGxLw5f14sNmvXCn1ZaCQ7E+QNyefPw+nTwqjYtUsYGVP0Q1EU/u3f/o2nnnoKi8Vybc+PxeDYMUF4W7eKC0Nv1JyT43r0qODFu+8WC+h4W66BdFVVXCr35pvCStm9e3yLEYdQKITD4SAvL49wOEw4HJ5QeSkQCBCJREhPT9f+RlEUBgcHycvLw2g0arUKZEbWxD6KqmFZWVmaj9Xn82mnwcFgkFAoNOH7pxWYogj/z759YuXesWPmiw8DAfjoI7Ed3bFDbC8mmTjyoC4cDpORkYHb7Z6xALbH48FsNk9ZQFrKSlEUbDabllqqqioul4ukpCSscht8CUKhEGNjY2RkZEzfv0uhqoJ0jx4VE23nTqHoM23jnU546y1hJe/YIay4SWTk9Xpxu91kZGRMutsIBALEYjHtM5/PRygUIjs7W0s8sVgsZGRkaBM6FovhdDrJzc2d1LJSVZWRkRHS09OvLKZX+tDPnBEEum2bsHKmq32hqtDSAp98AnV14gr3SRYuj8eDy+XCbDZjsVjIzMzU2j6b9kYiEQYHBzEajaSlpU2I/R0eHiYrK2vKqJj4OXRFUBRh6e3bJ3jhnnvG3TATRKDy29/+lkcfffTaootiMXjvPfHce++d9Lbs6fhDURScTueswjJlMaDExMSZb5eR/PjGG8L6X79etmvKh8xcLcXrhZdeEqv7F78orJZJGh2NRvnVr37F6OgoXV1dtLW1EQ6HtbjWpqYmPvroI+3/gsEgkUiEF154Ab/fD4Db7eYf/uEfCAaDWrhWOBymr6+PgYEBHA4HoVBII3U5aRsbGwkEAlq8bCgU0hIWYrHY5X2Sh1KVlfCtbwk3wfvvj0cHTIZAAF54Qfj6nnpKEPUUq6zL5eJ3v/sd58+fp6mpieeff55IJEIoFEJRFK1fiqJofXG73drCJNsvP5ML46uvvsqRI0f4zW9+Q3d3N9FolGg0yp49e+jo6CAajaIoCsFgkGg0qslYykZ+JmUbiUQIh8NEIpHJ+2wwCMK8/XZ44gnhIjh48HKfYDxGRuCVV+Cuu+D++yddoKW+PPvss/T19VFfXz+hbbFYjHA4jN/vx+v1Eo1GCYVCXLhwgb179wIizvmXv/wlra2t/L//9/+0cff7/Tz//PMTZBiJRDSZhMNhnE4nkUhEeymKQiQS0WQ2qSEiD4Lr6sT4NzWJxWgqWaiqsAYPH4ZvfAM2b56UKAAGBgb4p3/6J7xeL+fOndP6IqvKOZ1OrR9SH0KhENFoVGvrr3/9a9ra2vjf//t/4/V6CYfDRKNRXnzxRUZGRojFYtr3yu8Jh8N4PB58Pt+EeSPlMuX8AbHYFBeLvuXlwR/+IBalSWAymab+ntlAVcdl/dBDU7p4wuEwv/rVrxgZGZnQj7a2NtxuN0NDQxr3SHmEw2Gt71IPwuEw+/fvp7m5WZOD1J3LdMNgEJwo50d///Tzg5lCxlQVPv1UWCsTTefLkJyczO23386nn35KYWEhkUiEl19+mezsbM6cOUNeXh5dXV20t7ezY8cO3G43KSkpJCYmalldHR0dpKenc/z4cdavX8/hw4fxer1cuHCBqqoqHA4HgUCA9vZ2lixZgtfrRVEUPvzwQ7761a9qsbTnz59n1apVjI2NsXTp0qmDwOWp8L33wu9/L3w1hYWTy2H/ftiwQUy6GVZKi8XC4OAg/f39PPDAAxw4cIAzZ85w4cIFjEYjKSkptLa2kpubi9frxeVyUV1dTVZWFvX19aSkpLBr1y7a29vp6OjgiSeeID8/H6vVisvlwuVyEQgE+Oijj2hvbyclJYVDhw5x9OhRqqurGRgYoLy8nKNHj7Jt2zbcbjdOp5PW1lZMJhM+n4/m5mZWrFhBYmIiJSUlrFmzZuoOGQyCPL/yFfjtb4W1m509uZw++kiQbWnptHIyGo0YjUYOHDjAN77xDT744AOMRiOKopCYmIjX69WqVJ07dw6Hw8HKlStRLi6MSUlJhMNhHA4HiqJw5swZrX+jo6M8//zzrFu3jtHRUUZHR7FYLBQXF+P3++nv7ycYDGK32xkeHiY5OZn+/n4tqWXNmjVkZmZOLYukJHj8cfiP/xDW7mTWUDgMhw6JRX0GP3BKSgp+v5/e3l5ycnL49NNPcTgcGI1GLdXd6XRy9OhRVFXl7rvv5sKFC7S2tvLd735XC2F0Op0EAgEcDgctLS04HA7C4TB79uwhKyuLpKQkHA4HFRUVdHd3U1paitfrxWq1cvToURRFwWw2c+rUKZYtW0ZKSgpVVVVUVlZOLQujEW67TVh7jY3Cmr+kr5dm3F0xFAVOnYInn5xyZyHLuRYXF7Nv3z7uvPNOmpqaGBkZ4fTp09x33300NzcTjUbp6emhs7MTVVXJzc0lLS0Nr9dLbm4uZ86cIRQKkZubi91u58SJE6Snp7N8+XI8Hg/btm273FKWC/I99wi345e/PG13ZrZ0nc4ZJ5B4roHdu3dz4MABgsEgx44dY+PGjRQWFpKamsqmTZuora2lsrKSzs5OgsEgra2t2t+HQiFaWlooLS3l3XffpbOzk97eXu666y4KCwspKSnB6/VSXl5OZ2cn0WiU1NRUAoEAubm5ZGVl0dnZyYkTJ9i2bRtpaWnU19czOjo6YxcxGMSqPTw89e8MDc248EgoisJf//Vf09PTw6effqpVdaqqqiIpKYmOjg527NhBQkICNTU1JCcnMzo6qq22JSUl2O12+vr6KCoqIiMjQ1NaGcr2zjvv4HA42L59uxaIPjAwQCwWY3h4WLOCe3t7KS4upre3lyNHjlBbW0tDQwMGg4G7776b1tZWurq6ZjcpTCYR5TCVTFVV7IwKCmaUUywW46tf/Srl5eX88z//M0ePHqW2tpZXXnmFzs5O7rrrLtLS0hgZGSEajdLV1aXtiCSSkpKoq6ujs7OT999/n6qqKioqKkhLS6Ompoa+vj6i0Sjnz5+ntraWDz/8kMLCQkKhkCaDlJQUzp8/T01NDVVVVezfv59gMDi9HOQkS08X7pfJICN5Znnwlp6ezurVq7HZbEQiEXp6egAoKirSLCzpQhoeHsbpdFJeXq65o4xGIzU1NRiNRt566y3MZjMbNmzAaDSyceNGOjs7icVi9PT0kJeXR0dHB6qqYrPZOH/+PB0dHSxfvpyGhgZSU1PZtm0bR44c0azGGeVRXi7myCSIxWLXlvQTjQqyncJ9BmLONTU1UVhYyP79+6mvrycajXL33Xdrt1c7nU7Nsq+srMRkMrF7924sFgt2ux2Px0NeXh4VFRX09vZSW1vL6OgodXV17Nmzh+zs7KldhAaDCCq4REcnw8ykW109/TYqDmlpadx5553U1NRQWFjI559/rpnvKSkpjI2NkZqaisViYWRkZILZ3tzcTEVFhVYhLBQK0d3dzYcffkhOTg4dHR0oioKqqmzevJmCggJNGa1WK729vSQlJWEymTh9+jThcJht27YRCARmbDfBoAjdmS5LZ9kysbWezgVxES6Xi/3797NkyRJyc3MBqKysZP/+/WRlZeFyuWhvb8doNNLa2kpBQQEZGRkEAgFMJpO2NXS5XIRCIQKBgOaKsNvtmEwmdu3apcUXGwwGjZBOnjxJQkICfr+fXbt2kZ2dTSQSwWKxUF1dTUtLi7ZaezweVq9ePTtLRPquenom3w2AULySEmGVzPB9sipcWloa69evp7KyktbWVr72ta9ht9v54IMPCAQCBINBGhsbSU5O1nza8loio9FIe3s7q1evZs2aNRphGgwGgsEgTqeTlpYWrFYrGRkZVFZWUlhYiKIoVFdXaz7h6upqxsbGiEQibN68+TJyn1QWIyPiFH0qP3lKithuDw7OKAu3200sFuP8+fMcPHiQU6dOkZSUhNlspr+/n3A4jO8iuSuKQiwWo7e3V3OfhEIhDAYDFy5cID8/n+3bt3P06FFGRkZISEjA6/USiUS07w2FQmzatEmbQ1lZWaSmpmK327ntttuIxWL4fD62bt3KmIxVnX4w4fhxwRWTLDChUOja6i8kJgqL2umc8ld6enpITk7mvvvuo7KykmAwyMGDB2loaCAvL4/z588DwiLu7OzU3Hwej4fGxkYsFgsej4fu7m4MBgOVlZWMjY1hMBgoKioiLS2N4uLiqduoqiLCo6Bgxu7MfJAWicDrr4vwkN27Zzx9VhQFg8Gg+UcSExOJxWKYTCYikQhGoxGDwTDBxyMP0QASEhKIRqMYjUbtdy79GdCeIV0T8SZ/LBbDbDZrz53ScS4Pff70JxH2UVs7dd+iUeH3DQSEO2IKX6WUgfSvyv4bjUbC4bD2b5mbbjQaue+++7T+y76cPXuWxsZG8vPzqa2tJScnh2g0qvXdaDQSiUSQtR2kPOLvPktISJiQbWcwGCaMiZSr0Wicvhh6LCYWpX37hOugvHxqOQWDwr9XViYOm6YoHSj9yrFYTHMxycVB9lO2/9LxlfoSjUY1WQAT5Hvp38q/kbon9SxeTlKGM+pMV5c41HnoIbHITPW7AwNi7txxB9TUTLk1lv5FOY7x/VUU5bKx2b9/P16vl/T0dDZt2kRKSsoEmSUkJBAKhbSDN/md8nsv1Rk5HvHzZla6IQ9a33tPuJvuuOOyPgYCAV588UW+/e1vX33ImKqOH3p//euTunPkYizbLcdZ6obs5zvvvIPFYiElJYWNGzditVq13z958iSnTp3iqaee0voci8U4ePAgJpNp6tswZETTH/8I3/ymOAe55pCxaFT4Ks6dEwcCK1YsrtCfSyGttvp6Ybnde68giZn6oyjiNPqzz8SB2qZNItbyKmJx5ak0MOn9WtFoFIfDQXJyMunp6XMbWD67Bgq/ZFsbHDkiFpm77552sdEQiYizgI4O2LJF7BIWu76oqiDRAweEXO6/X2wnZ4rk8HoFKUWj4mykpOSa40lDoRBDQ0NkZmZqkQrzismiWpYvn3RRaWlpob29nfvvv//a2qmq4vDy0CF48EEx/65izo2NjTEyMkJ2drZ2mYCEx+MhEAiQn58/IXrE6XSSmpo6+U0giiLi9A8cgEceiXetzUGcrgxKrq8Xq05pqQiRKChY8PGpWkZOICCI4PRpEfS9YQOsXHll7ZcxqOfPiy1VKCR8vcuWiYyz2aUILjxIPQgEhCXX3Cx2AWVlYnHJzr7yrB63W0yS7m6hL6tXC31ZDDKSROt2C1m0tAif4tat45b+bON0QSSdfP65mENLlwqSyskZJ6qFLgsQcakdHeLAzO8XO8O1a6eM31ZVlZdffpkdO3ZQUlIyN+1wOGDPHrDZULdvh7w8QYTzKb/4kNPPPhO+/XvuuTR2ew5IN/6B4bCYmKdOCd+WzMQqLxerv7RqboQiyf7IzCKHA7Wjg2MnTmAdGaF00yZSt2zBUFCA4VozWOSK39EhAqTdbpF+WFQkyCo7W/xbHiIshIkl5SPjTl0uAr29jJw9S47Hg8VqhaIiDMuXi8PFa11Qpb50dorFbmREKGlZmUg0yMoan7Q3Sj5SJtGoIMW+PtFeh0MQbVWV2N1NES55Rc/x+8XuoalJ6Gd2tpg3paXChSflfaNlIetiOBzjmYogdGLVKrF4zjB/BgcHeffdd/nWt741Z3f5qRcJz9/YSOsLL5BdUEDZtm3C6ElLu35ZadLYGh4WsdodHWLstm0T8eiXP3MOSffShkSjwiLq6hIWjds9HmKUlTVeRETWYJDJBFe7wksLRK420Sj4/QRdLoxOJyanE8PoqFDoxETIyUGtqKDXbOZ8Vxd9/f2EQiFycnKorKykrKyM1NRUzdd81ZCy8HhETn13t5BLMCgmUlqaWJxyc8XPNtt4xSspD/n8a/F9XSqfYFBMdJ9PtGdkRKQsywPGtDTcOTl80t/P0MgIRpuNwqIiqqqqKCws1Lavc5bCGY2KNvT0CJ1xucT/p6YK0snPF3KSRW5kRTVJRFejL/JdVmoLhcS23+MRpOJ0inZEo8JaKSoSC0J+/tWn/c6mXeGw8AV2dYmUVq934tzJyxNzR8oiIWHuZREOC93wekVbnE4RnRIKibji3FyRBFRUJHR2luQpY4R37dpFWVnZNemP5Ci/309raytnzpwhGAxSXVXF+qVLSevpEQtZMCjkVloqDnszM0UfJOfM1Ib4BScWE5a90ykWnL4+oS/p6WIBrqoaL6Y1Oa4T6U7WaFnlSRZWGRwUROzzCaGo6nh8bGLieHUkk2lcoRISJhKHJNdodLxylDyIu1hJqDUUYt+xY6zfsYPa3bsxZ2aOFwqJE4wMLHc4HHR0dNDT08PY2BhZWVlUVFRQXl6uZf9cM9HI1TEUEjJwu4ViezzjpeZkxTApk/hSjPIlFT1eLlI2shJbNCrkHomIn2XbLRahHDabUMjsbKE4NtsEmcurZQKBAP39/bS3t2vRArm5uVRVVVFcXIzNZrv2BSpePtLi9niEvjgcE0vxyWgRk2lcX0ym8bbHk5CiTJSLlImswiblnJgo+p+WJogtK+vyCTrfkO0Oh4WeSFnEz51LZZGYiGo0MhKLkZGUhHEmWcj5IyvPGQziu2w2Qeyywlhmpgh3u8JdjoywcTqdnDx5EpPJdE2XU8oMM0m0oVCIqqoqVq1apWWmTijic7HOCr29wv/uco1XOLRaJ1btk7oj+UXOHVk+FMTfZGUJ676oaOJOZGbME+lO+S1xK6xUhPh6oFIp4hUFxieUJGWpbPFVui6+VMRKeOjQITo6OtiwYQO1tbXTXlAp+x6JRBgaGqKjo4Pu7m78fj/p6enaNSuZmZlzf6fVZFZpfJnGeJlIuUg5xi9OkpjjJuKlsgGu6tBBZu8MDg7S3t5OX18f4XCYzMxMbZcwZwvUpbKJl1F8uUa5qESj4wtWvL7ARJnE19aNJ2lYGO6emTCDLAIeDy+++ipPffnLmOMXZTn28bKQ5UJlCNY17DgVRSEUCjEyMkJ/fz8DAwMMDw8TiUTIyspi2bJlrFix4oqvT5Lp4TKbMxqNUl1dzapVq8jOzp79rmsyq1XW747nGqkT0tiRsdXy/69CNhdxg0l3HiFXx/r6etrb21m3bh11dXWzuh1YyiIajeJ0Ouno6KCrqwuv10taWhpLliyhoqKCrKysCeE4twKkbMLhMMPDwxMWqNTUVMrLy6moqCA7O1u7YfZWkc2NgqqqfPbZZ5hMJrZu3Trn8o4PY/T7/TidTvr7++nv78ftdmMwGEhPTyc/P5/y8vIJEQGzbYt8hs/no6WlRUu9Xbp0KStXrtQiexahLt06pCsh76Kqr6+nra2Nuro66urqJg/7mOY7QJDw6OgonZ2dXLhwQUthLi8vp7KykpycnFuOaC6VTVdXlyYbi8VCaWkplZWV5OXlaYV4bhXZzBeCwSDPPfcc3/rWt6YsdjRbxI+nx+PB4XBoNU8CgQBGo1HL7JK1p+OfeSVjG0+0586d4+zZs0SjUY1oZcGfRa4vtx7pSkjyPXz4MOfPn9cs3ysh30u/T2aLdXV10dHRgdvtxmq1UlZWRlVVFbm5uVdXDH2RQwbYe71eenp66OjoYHh4GKPRSFFRkZYRNqeHc7coVFXl0KFDRKNRdu7cecWkJ/2vo6Oj9Pf3a3UootEoycnJ5OXlUVRUREFBATab7arnS/wzpcUsLdpLiXauIhwWCG5d0pWQ5NvQ0MC5c+c0y9disVyzMsViMS2FsL29nZGRESwWCyUlJVRWVpKfn69dXHkrYarDuVgsRl5eHpWVlZSWlpKSkjJ3h3O3CEKhEP/xH//Bk08+OW35QWkkBAIBhoaGNIL1er2ae6C4uJjCwkJycnKwWq1zOhaKomjV01paWm52oo2HTroSMivlyJEjtLS0UFtby7p1666ZfOO/XypaT08P7e3tDA0NYTabKS4u1qy9pKSkW45kpLUTCoUYGBigo6ODvr4+QqEQmZmZVFRUUFZWRkZGxtwfzt1EUFWVI0eOEAgEuP3224lP85XlGgcGBujt7cXpdBIMBklMTNSs1/z8fDIzM7VzjrmWs6yN0dLSopWqrK6uZvXq1WRkZNwquxyddC+FPJk/cuQIZ8+epba2lvXr188Z+cY/R1VVfD4fvb29mrWXkJBA0cV42KKiIs1auQWUUcOl0SMXLlzQqomlpqayZMkSysvLb0mf+VSQNSp+85vf8NhjjxEIBLDb7fT39+P1egFxLX1eXh4lJSXk5OSQmpp6ZZerXkWbQKTRSh+tqqrU1NSwcuVKLY39Fhs7nXSngiTfhoYGmpubWb16NRs3bpxz8o1/Hojwtr6+Pjo6OrDb7QAUFBRQWVlJcXGxVv3+VlJUKZtYLMbIyAhdXV10dnbicrluycO5+IgRj8eD3W5nYGCA48ePY7fbWbduHVlZWRQUFFBSUkJGRsaEm0eup2ziifbs2bOcO3cOVVW1UDF5S8nNPD4zQCfdmSDJ9+jRozQ3N7Ny5UqtCtH1Vhzpb+7v79e23IqikJ+fT2VlJSUlJaSmpt5yJAzj7hqPxzPhcE7uFCorK7WdwmKVT3wluEAggNPppLe3l8HBQdxuN4qikJycTFFRETk5OezZs4fvfOc7814ISe7apEU7GdEuRvlfJ+ikO1tI8j127BhNTU2sWrVq3shXQl5dMzAwoCUlRCIRLTOstLRU2zLeakoefzjX19fHhQsXsNvtKIpCXl6eljm3UOUjiSsSieB2uxkcHKSvrw+n04nf79fuSCsuLiY/P5+srKwJu67Gxkb6+/v5whe+MC99k+11u92cPXuWs2fPYjAYWL58OStWrCA9Pf1mPgy7Fuike6WQBz7Hjh2jsbGRFStWsGnTphtyACbDe2RmWG9vL8FgkOzs7AmZYbdiBMBCPpyTkS1jY2MMDQ0xMDCA3W7H5XKhKArp6enk5ORQUlJCbm4uNpttyjaeOHGCV199lZSUFP7Tf/pPU18nNEdQFAW3201zczPnz58nISGB5cuXs2zZMu1S1FtN164QOuleLeSEPn78OI2NjSxfvvyGka9sD4jDJ1k/oru7m0AgQEZGxoTU5VuVhEHIJz5zzufzYbPZWLJkCUuWLCE7O/uaY08vfW4wGMTr9WrW69DQEMFgEJPJRG5uLgUFBRQVFZGZmalZr7N9/ksvvcSTTz5JWloaP/zhD/nRj36kHS7OFRRFweVycfbs2QlEu3z58lvWvXUN0En3WiHDcY4fP65d3Ld58+YbHvoVn0kkSUZGAMjUZUky8bdu3CqIP5yTWYWdnZ2Mjo5itVonxFLPdDgX73v1+XxaWqzdbtd8r6mpqRQXF2uXXMoKdtN972z68Mwzz/A3f/M3FBcX8+yzz3LPPfdcVRZYfDvkLsHlcmkWrdFo1FwHqeIGhFtKX+YQOunOFST5njhxgpMnT1JTU8PmzZtvTAX/KdoHgoRHRka4cOECFy5c0G7XjU9dvtXqR8BE4pzN4Vw4HMbtdmO327Hb7Zr1arFYyM7OprCwkKKiItLT069b5ICiKPzwhz/k1Vdf5dlnn+Xuu+++YsIdGxtjz5493HfffSQlJTE6OkpTUxNtbW2YTCbNok1LS5vz9t+i0El3riHJV96rVFVVxZYtWxYM+cYjPnVZ1o9wuVwkJydTVlZGZWXlpKnL8u9m7Qu9tGJaICCqO4VCosLT2Nh4dTBZIc1qFZWdLBbxc3LyxCpP8xA5AuPumsbGRvbv309PTw+FhYUMDQ1RVlbGbbfdxtKlSyksLCQxMZFPPvmELVu2zP7EPl42supVfKW9cHi8CpZ6se4wgMlEDPjRM8/w4B13cMd995EgK+1ZLEJ2l9YcvqR//f39/O3f/i1vvfUWv/jFLzCZTBrRxvtodcwpdNK9XpDke/r0aY4fP051dTWbN2/W4mwXIuQBj9vtpru7m46ODi11WcbC5ufnMzo6yk9+8hO+973vUVdXd/kptSySPjgoapjK+q+xmCgjmJQ0kVBljVajcZx8JBlLYg4EBPEYDKLOa3a2uGWioEAUGJ+h3q28dHIy2V9ac2BgYIC+vj5GR0cJBoMkJSWRlZVFaWkpeXl52Gw2zWXT19fH2NgYmZmZZGRk8P3vf5+ysjJ+8pOfsH379nH/uSRMWTR+cFC8e72in7KcoKzvGl+O02qdWEtZdAg1FsPudlNgMJAQX1NaykzWxpV1k3NzIT8fJSuLI2fP8lf/7b9x9OhRVFXl+9//Pr/85S81H62O6waddK83ZBjQ6dOnOXbsGJWVlWzdunVBk6/EpbGwbW1tOBwOuru7eeaZZ8jOzuYHP/gB3/3ud0lPTcUwNCTuDevuFhM+M1NcOVNQIG59kLVrJ7G8ZmiIeJdF6z0eQVi9vaJyfzAoikmvXCluB4i77FJRFDo6Onj22Wf5i7/4CyorK4lEIvj9fhwOB/39/QwODuL3+1FVlaysLPLy8iguLiY7O5vk5OQpDx7jkxSGhoZ44403+MEPfkAkEiEzM5Pvf+97/OCJJ8gZGsIwMDB+60JGhrh9Ii9PFAmXNy9cjWymklW89ez3C3J3OsHhYLCri//zySf0NjcTGRsjmp3Nst27+buf/QzTHB/C6bgMOunOFy4l34qKCrZu3YrNZlvw5CshSfjnP/85P/nJTzCbzSQnJnL3li388qGHWJKbi2HFCqiuHr+M73r2TepoKAT9/eKOqv5+WLIEdfNmnNEozz33HL/+9a/p7+/nv//3/051dbWWVJCTk0NhYSH5+fmkp6dfU80BVVV5/fXX+Zd//Efy3G5KzGaKV63i9gcfZOXKlRhKSgTB3ugLJ1UVFVBjMUHEfX3Q3g79/RjS0jCsWSMuyLxe1xHp0El3viH9oY2NjTQ0NLBkyRK2bdu2aMhXVVX279+PY2CAHLud7M5Osu65h5ydO7GmpQmNuoGEQjgM587R/uGH/O0bb7CvsZFRv59YLMb/+B//gx//+McTdhlzdb2QOjpKdN8+DA4HxuXLxQ3HF29KXvDjKq1ip1NcEtreLq6E37792i/d1HEpdNK9UZCW75kzZzhy5Ajl5eVs3759cZBvOAzvvCN+vvfeyy7ikyUDpdV4qStFnppbrVYSEhK0eqpJSUkYjUbNpWGxWC7L+JO+clVVpy3QHQmF8NfX43j3XXo2baJ9ZIS8vDweeeSRuZVvNAoHDkBrK+zcKaz8S66iURRFc19IecyUrRWLxQiHw9NmPMpYcYPBgMVi0QooJSQkkJSUdNkzZELGtK4t6XtuaYGDB2HNGti8edyXrONaoZPujUY8+R49epSSkhK2b9++cA80VBXefltYQDt2THoLrM/n42c/+xmPP/44TqeT3bt3a3WDFUXBYDBw+vTpCXdlPfvsszzwwAMUFxfz0UcfMTY2htfr5fHHH8doNGrkDFBfX4/X6+Xee+/ViGxSWamq8C/v3QtPPim2zHOJaBTeeEP4aHfvFr7kSX8tyjPPPMOGDRsAqKurIzs7W/tctl1RFBISEujv78dms9Hb28vKlSu1/5f9lyUbAfbt20diYiI7d+6kvb2dgwcPYrVa2b17N7m5udrvK4qCw+HgzTff5Lvf/S6qqk6fDi13DR99JPr50EOzvvFXx7SYclJf2a1xOq4aBoOBxMRE6urqWL16Nc3Nzbz00kuUlpaydevWeS9eMiPkde0PPDDlJLRarYRCIUwmE6mpqZw4cYLR0VFWr16tXfMi75fr6urCbDZjNpuJXbzJORwO09DQwJNPPklvby8DAwMAWlH4/Px8AoEAR48epa2tjQceeID09PTLG2IwiAiH5cvhxAnYunVuZXH2rDgYu/POaQlJJkFIkuvt7WX//v34/X4yMjIoKirCbDbjdruJxWLs3buXJ554ArvdTigU0gqLt7e3Y7VaWbVqFS6XC5/PR1paGsFg8GJ3DTQ0NPD0009rtXU9Hg8ZGRkMDw+zfPlywuEwZ8+epaWlhTVr1rB06dLJG20wiL594Qvwxz8Kl8NUv6tjTqAvafMMg8GA2WymtraWp59+muLiYl555RXeffdd3G43M+w85g8jI1BYOON2MyEhgZSUFGw2m2bJv/baaxQUFLB27Vr6+voIBoOEQiEaGhq0v1NVla1bt7Jr1y7+8R//kTfffJOamhqqq6uJRCKcPHlSq1nQ0NBAUVERkUhk6oYYDLBkiQjRmmt0dcGqVbOyAA0GA0lJSaSmpmoLktfrpaqqisbGRt577z3WrFnDxx9/jM1mIycnh5aWFvbs2cPq1avZu3cvfX19lJWV0dLSQiQS4fjx44TDYUDILT8/n7/8y7/kX/7lX3jjjTdITExkw4YNRKNROjs76e/vJxKJ8P7771NUVIQib0ueDgkJIiqkp+dapaVjBuike4MgyXfNmjV8+9vfprS0lFdeeYV33nkHl8t148k3MxMGBsavN58EoVCISCRCeno6NTU1fP755xoJ/+EPf6Czs5NgMMjZs2fp6OggGAzi8XgIhUIAfPzxx6Snp1NaWkpOTg4vvfQSnZ2d1NfXYzAYGBwcxGQy4XQ68Xq909cakC6GvLy5loSwos+eHY+imALSx+1yuRgZGaG5uZlwOEw4HMbv9xMOhyksLKS1tZXa2lqSkpLo6ekhGo2Sl5dHe3s7a9asAWBsbAyHw0FDQwOqquJwODRLt7GxEYfDwapVqyguLubtt9/m3LlzHDx4EFVVsdvtRCIRzGYzHR0dEzLlpoSiwLlz4mBNx3WF7tNdIJAJC2fPnqW+vp7CwkJ27Nhx49wOigJvvSUC7bdsmdTKk24AWR3L6/WSkJCAyWTSbn+QV8WEQiHNX2m1WklMTMTv9xMKhbBYLFgsFtxuN+np6fj9fu1wSPqHo9Ho1NlfqipCot59V/h0Z0MyV4JoVGy98/LEIdoUCRryYFA0SdVC02KxGGazmWg0itVq1YrvyMOxWCyG1WrF7/eTkpKiFclRFEVL9gDhvkhKSiIcDhMIBEhISCA1NRW3201SUhLRi1lssqZCYmIiXq+X7OzsqQ/05IHaxx+LxJSHH9YP0+YG+kHaYoEk33PnznHo0CEKCgrYsWPHjSkQHQ7Dm28Kn99ddwkyW0h+ZxCEcfIkHD8Ojz8OWVnXp42RCOzfL1wNu3cLV8YM2XELGqoqEira2kRUxtKlInTMpB/zzBF00l1skOTb0tJCfX09OTk57Ny5k8zMzPkl31hMkFpDg4hJrasbT4i4UVBVQYKtrVBfL3zPd9xx/QP9VRWGh+Gzz4TPe/ly4QdNT5+XOhHXDBmn63JBY6OQX36+iE7JzFz47V9c0El3sUKSb2trKwcPHiQ3N5cdO3aQlZU1f+SrqiIF9/Rp8UpLEwS8ZMllsbvXtQ3hsPAznzkjUoOLimDbtvknDFUVKbfnzol06FBI+EIrK0WbkpMXBglLazYQEO6Xjg4hP4tFLBgrVszf+N160El3sUOm5p4/f56DBw+SnZ3Nrl275p98FUVM3KYmcdKdkCBqLpSUCIvTZhOTOr5Ns63CFY9IRBDb0NB47YWxMcjJEZEES5aIGg832uIG0a7eXkFqdjtqNEokMZGRjAzyystJyMsTcpEFbS7FlfZhsjkrq5T5fKLw0NCQeAUC4plFRWJRKC4ej2PWyfZ6QifdmwWSfKXlm5WVNf/kKxoi3v1+QcK9veLd7xfEbDYLH3Bqqni/WElLNRqJGQwYVRWDLGwTTxh+/3h5Q6tVHOQVF4vXAt7Gy3k0OjrKsfp62k6cYFl5ObeXlpIwPCz6djHsC5NJyCc5WbwsFvFvWYEtvo/xBW0ikfEykH6/IPtIZLwMpMUyXmUsL0+8p6aO+2kXoNxuYuikezMiFovR1tbGgQMHyMrKYufOnWRnZ9/YJAvpb5VlB2VJw4t1Y9VYjHfef5977rgDiyz1KOvDpqYKf7Esd7gIMqNkpmFHRwcNDQ1EIhHWrVvH8uXLJ6b2ynq6kiQjEWGFBgLjNXVlrWFFGQ/Vk1XJjEZBzGazkE1KiljM5P+ZTNe/8JCOK4FOujczYrEY7e3tHDhwgPT0dHbt2kVOTs7CynC7CEVReO655/j6178+u/jRBQp5n9jx48dpa2ujuLiYjRs3kp+fr9+OqwP0NOCbG0ajkZqaGqqqqujo6ODtt98mNTWVXbt2kZubuyDJdzFCWrXt7e0cPXpUs2p37Ngx4Zp0HTqmg066NxGMRiPV1dVUVlZy4cIF3nnnHZ185wCqqgpf7bFjtLe3U1xczD333ENeXp5+Q66OK4ZOujcZDAYDRqORqqoqKioq6Ozs5N133yU5OZndu3eTdzFNVieK6SHrIbe1tU2wanfu3KlbtTquCTrp3qSQ5FtZWUlFRQUXLlxgz549JCUlsWvXLvLz87Xf0yEwIQJBt2p1XCfopHuTQxJFVVUVS5Ysoauri/fff5+kpCR27txJQUHBLU8k0/lqpysurkPH1UAn3VsI0vItLy+np6eHDz74gMTERHbv3n1LnrrLCATdqtUxn9BJ9xaE0WhkyZIllJWV0d3dzYcffojJZGL37t0UFhbe1GQzlVWr+2p1zBd00r2FkZCQQHl5OWVlZfT09PDhhx9iNpvZtWsXRUVFwM3h873UV9vW1kZJSQn33nuvFtVxM/RTx+KAnhyhAxivwSqvmDEajezevXtOyTcQCPDRRx/x/vvvs3PnTrZv307JdSyafWkEQjgcZv369axYsUK3anVcb+gZaTpmD0VR6O3t5bPPPgPQyPdaLUK/388999xDfX09GRkZfPLJJ9TW1s5VszVMFle7adMm8vLybjm/tY4bBj0jTcfskZCQQFlZGU888QR9fX3s378fgF27dlFcXExCQgLhcFi7JWK2SE5O5oEHHqC+vp5169ZRU1MzZ23WfbU6Fgt0S1fHjFAUBbvdzv79+1EUhZ07d/L6669jtVr5zne+M2viVVWVkydPctddd/E//+f/5L/+1/96VWQYDAYZHh6muLgYgJGREY4fP057ezslJSVs2rSJ3Nxc3arVcSOhW7o6rh4JCQkUFRXx9a9/HbvdzhtvvMHf//3fa/d0ffvb3x4n3vhFXJZtvFiS0BAMstxi4c7aWu6rqcHQ3S3KN14s+6hV0JK4hJBVVcXn8/Hzn/+cc+fO8bOf/YxTp07pVq2ORQWddHXMCtKfW1RUhKqquFwuwuEwP/rRjzCbTDz12GMYHQ5RbHxwUNR7jcXGSzeazWC1YjWZ+P++8x0qo1FxC0UwOIGYMRhEmcfcXFEUvaQE0tNRjUaGhob44Q9/yEsvvURiYiJ33303X/3qV/UIBB2LCjrp6rhiPPTQQ6yoqaHv2DG69uyh6YUX6PB4qK6pwVBaKq7ySUkR5HnJFt8ArJnuy0MhUX93aEgURW9uRnW76bNY+OuXX+bThgays7Mxm82Ew2EtkUGHjsUC3aerY/ZQVUGKp06Jiw3T0lDXrEEtKwOrFUNCwtwToKqihsM4u7sZPH6cxJMnsSxZQuLtt2MrLiYlJUUnXR0LEXrImI5rhKKISxgPHhSXGm7YMP+3AsubF1paxC3AlZWwc6dwXejEq2NhQSddHdeAWAzef1/4ab/wBXEP10WSk7cVK4qC0Wgk4aK1K//faDROsERlaFdCQsJln10KeR8ciNTluA/EId3hw3D+PHz1q6JNOnQsHEyp2HpMjY7poarw6aeCZB97TNxjdglRHjt2jP/7f/8vL774Ih0dHVra7enTpxkbG9Oy3WSG2D//8z9z5syZi1+vTvkCOHDgAA0NDRPbZDAI63b7dnEF+x//OH45ow4dCxz6QZqO6eHzQVsbPP20uBxxEiQnJ5OQkEBKSgoffPABdXV12Gw2xsbG6OnpobW1FUVR+MIXvoDZbCYUCpGcnAyA2+2mubkZi8VCb28vW7ZswW63YzQacTgcJCYm4vV6J2+bwQDLlkF7O5w7Jw7wdOhY4NBJV8f0GBqCoqIpCVe6B/r6+li7di2rV6+mvb2dhoYGLBYLK1asoKenh0gkgsfjuey24vfee4+srCxcLhfBYJCGhgasViu9vb0kJyeTkZEx/UGZwQArV8KZMzrp6lgU0N0LOqaHzQYu15Qfq6qK3+8nKyuLHTt2MDAwQE9PDz6fj7GxMc29EI1GCYfDxGIxvF4vfX19HDlyhK6uLhwOB2VlZezcuZPz589TUVFBf38/paWlRCIRgsEgU549qCo4HJCZeX36r0PHHEM/SNMxPWIxeOkl2LQJli6dNEtsdHQUn89HSUkJwWCQoaEhjEYjqqpisViIRqMoikJ2djaJiYn09PRc/OoY2dnZ2O12SkpKSElJwe/3k5KSQldXF4mJiVgsFlRVnfpKeY8Hfv97+OY3hb9Zh46FAT16Qcc1wOOBF18UB1crVwrivdEhWqoqXB+vvw533QVVVTe+TTp0jEMnXR3XAFWFQADeeUf8fOedkJMjPptvolNVkbF26BB0dMADDwifs064OhYWdNLVMQdQFEF0Bw+KmgobN8KSJaK2wvUkPVUVbo7BQTh+XNR3qKuD9esnFsjRoWPhQCddHXMIRREEeOIE9PYKX2pFBZSVQVaWiKFNSLg6IlbVcZL1eATBtrcLV0JGBqxdKzLRTCbdutWxkKGTro7rAEmOIyPQ2Qnd3eORDikpgowzMiA9XRS/SUwUZJyQIP5WUURabygkst1GRsDrFbHB0aj4jqIiQbL5+dffotahY+6gk66OeYC0UiMRQZ5eryBhj0cQazgsSFpRxslX1tC12UTYV2oqpKWNVyjTSVbH4sRVk27P3LdFhw4dOm56lE71wUykq0OHDh065hB6RpoOHTp0zCN00tWhQ4eOeYROujp06NAxj9BJV4cOHTrmETrp6tChQ8c8QiddHTp06JhH/P/xqHAzWG03JQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "network.add_edge( s0, s4 )\n",
    "network.add_edge( s1, s5 )\n",
    "network.add_edge( s1, s6 )\n",
    "network.add_edge( s2, s6 )\n",
    "network.add_edge( s0, s7 )\n",
    "network.add_edge( s1, s7 )\n",
    "network.add_edge( s2, s8 )\n",
    "network.add_edge( s7, s8 )\n",
    "network.add_edge( s4, s9 )\n",
    "network.add_edge( s6, s9 )\n",
    "network.add_edge( s3, s9 )\n",
    "network.add_edge( s8, s10 )\n",
    "network.add_edge( s9, s10 )\n",
    "network.add_edge( s7, s10)\n",
    "\n",
    "try:\n",
    "    import tempfile\n",
    "    import pygraphviz\n",
    "    import matplotlib.pyplot as plt\n",
    "    import matplotlib.image\n",
    "except ImportError:\n",
    "    pygraphviz = None\n",
    "    \n",
    "network.bake()\n",
    "network.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b16a50f1",
   "metadata": {},
   "source": [
    "## Make Inference using our constructed model!\n",
    "Now given the network constructed above, let's test if a person who has normal skin conductance, normal facial scale, and abnormal sleeping behaviors is likely to experience pain\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "id": "6583adce",
   "metadata": {},
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "None is not in list",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-59-fff41d39ff45>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0mobservations\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m{\u001b[0m \u001b[0;34m'verbal_HR_RP'\u001b[0m \u001b[0;34m:\u001b[0m \u001b[0;34m'False'\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'facial_verbal_SC'\u001b[0m \u001b[0;34m:\u001b[0m \u001b[0;34m'True'\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'HeartRate_var'\u001b[0m \u001b[0;34m:\u001b[0m \u001b[0;34m'True'\u001b[0m \u001b[0;34m}\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0mnetwork\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpredict_proba\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mobservations\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m \u001b[0;31m#beliefs = map( str, network.predict_proba( observations ) )\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0;31m#print(\"\\n\".join( \"{}\\t\\t{}\".format( state.name, belief ) for state, belief in zip( network.states, beliefs ) ))\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/opt/anaconda3/lib/python3.8/site-packages/pomegranate/BayesianNetwork.pyx\u001b[0m in \u001b[0;36mpomegranate.BayesianNetwork.BayesianNetwork.predict_proba\u001b[0;34m()\u001b[0m\n",
      "\u001b[0;32m~/opt/anaconda3/lib/python3.8/site-packages/pomegranate/FactorGraph.pyx\u001b[0m in \u001b[0;36mpomegranate.FactorGraph.FactorGraph.predict_proba\u001b[0;34m()\u001b[0m\n",
      "\u001b[0;32m~/opt/anaconda3/lib/python3.8/site-packages/pomegranate/distributions/JointProbabilityTable.pyx\u001b[0m in \u001b[0;36mpomegranate.distributions.JointProbabilityTable.JointProbabilityTable.marginal\u001b[0;34m()\u001b[0m\n",
      "\u001b[0;31mValueError\u001b[0m: None is not in list"
     ]
    }
   ],
   "source": [
    "observations = { 'verbal_HR_RP' : 'False', 'facial_verbal_SC' : 'True', 'HeartRate_var' : 'True' }\n",
    "network.predict_proba(observations)\n",
    "#beliefs = map( str, network.predict_proba( observations ) )\n",
    "#print(\"\\n\".join( \"{}\\t\\t{}\".format( state.name, belief ) for state, belief in zip( network.states, beliefs ) ))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f99741d3",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}