notebooks/example_lockdown_social_distance.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# OpenABM-Covid19 - Example Lockdown then Social Distance\n",
    "\n",
    "This examples show how to run a simulation on 50k people and turn on a lockdown after 20 days.\n",
    "\n",
    "First we get the model with baseline parameters but change the population to 50k."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from COVID19.model import Model as abm\n",
    "\n",
    "model = abm( params = { \"n_total\" : 50000 })"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we run for 30 days."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "for t in range(30):\n",
    "    model.one_time_step()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we turn on a lockdown where the number of interactions on the occupation networks and random network are cut to 10% of their pre-lockdown level and the strength of household transmissions is 1.5 times the pre-lockdown level. And run for another 10 days."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.update_running_params(\"lockdown_occupation_multiplier_primary_network\", 0.1)\n",
    "model.update_running_params(\"lockdown_occupation_multiplier_secondary_network\", 0.1)\n",
    "model.update_running_params(\"lockdown_occupation_multiplier_working_network\", 0.1)\n",
    "model.update_running_params(\"lockdown_occupation_multiplier_retired_network\", 0.1)\n",
    "model.update_running_params(\"lockdown_occupation_multiplier_elderly_network\", 0.1)\n",
    "model.update_running_params(\"lockdown_random_network_multiplier\", 0.1)\n",
    "model.update_running_params(\"lockdown_house_interaction_multiplier\", 1.5)\n",
    "model.update_running_params(\"lockdown_on\", 1)\n",
    "for t in range(10):\n",
    "    model.one_time_step()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lockdown is ended, and replaced by social distances measures (e.g. such as wearing masks). These are modelled by a 70% reduction in the rate of transmission in occupation and random interactions (but no change to household transmission rates)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.update_running_params(\"lockdown_on\", 0)\n",
    "model.update_running_params(\"relative_transmission_occupation\", 0.3 )\n",
    "model.update_running_params(\"relative_transmission_random\", 0.3)\n",
    "for t in range(10):\n",
    "    model.one_time_step()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The social distancing measures are then removed and transmission rates return their pre-lockdown levels."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.update_running_params(\"lockdown_on\", 0)\n",
    "model.update_running_params(\"relative_transmission_occupation\", 1 )\n",
    "model.update_running_params(\"relative_transmission_random\", 1)\n",
    "for t in range(30):\n",
    "    model.one_time_step()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot the number of daily new infections over the course of the simulation. Note the increase in infections post-lockdown is much less than that without social-distancing measures."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.axes._subplots.AxesSubplot at 0x7fd4ef93f190>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZQAAAD4CAYAAADLhBA1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOydeXiU1dn/P/fMZCcJCQkQCPu+hbAGkN2KIipWi4gLiFbU4mvr21q1r6226q+ttVapCsUdRUVFcSluWJBF9h0SkLCHBMgC2ZdZzu+PmaRDyDIJkzxPkvO5rrkyc57zPOfOJPN853u2W5RSaDQajUZzqViMDkCj0Wg0zQMtKBqNRqPxC1pQNBqNRuMXtKBoNBqNxi9oQdFoNBqNX7AZ1bDFYlEhISFGNa/RaDRNkqKiIqWUMqUZMExQQkJCKCwsNKp5jUajaZKISLHRMVSHKVVOo9FoNE0PLSgajUaj8QtaUDQajUbjF7SgaDQaTTNCRF4XkbMiss+rbJmI7PI8jonILk95VxEp9jq2yOucYSKyV0RSRWSBiEhtbRs2KK/RaDSaBuFN4EVgSXmBUmpm+XMR+TuQ61X/sFIqsYrrLATmAZuAlcBVwJc1NVyrQ6lK7SodF496pYrIHhEZWts1NRqNRtMwKKXWAjlVHfO4jJuA92q6hojEARFKqY3KvYPwEuD62tr2pcvrTdzKVB1TgV6exzzcqqbRaDQa8zEOOKOUOuRV1k1EdorI9yIyzlPWEUjzqpPmKauRWru8lFJrRaRrDVWmA0s8KrZJRFqLSJxSKqO2a2s0G05s4OvDXxsdBgABlgDuGX4PbcPagt0OL7wAeXlGh+Vm0CCYMcPoKExDcZmTFbtOcU1CHOHBAUaH09jYRGSb1+vFSqnFPp47iwvdSQbQWSmVLSLDgBUiMgCoaryk1lwn/hhD6Qic9HpdrmQXCYqIzMPtYggMDPRD05qmzmOrH2PNsTVIlf+/jYfyfFaiQ6KZP3I+bN8ODz3kPlj7WGTDohSEh2tB8VBc5uTnS7ayITWb71LO8srsYfgwXtyccCilhtf1JBGxATcAw8rLlFKlQKnn+XYROQz0xn0fj/c6PR5Ir60NfwiKz0rmUdHFAGFhYTqzl4YSRwlTekzh69uMdSl5pXlE/iWSUmepJ7AS98/Vq2HiRMPiAuDRR+G554yNwSSU2J3cvWQbPxzOZurA9ny57zSL1x7hngk9jA6tKfAT4IBSqqIrS0RigRyllFNEuuMeujiilMoRkXwRGQVsBmYD/6ytAX9MG04DOnm99knJNBoAu9NOgMX4LovyGOxOu7vA7vkZYHxsBAT8N54WTLmYbDicxd9+NpiXbx3KtEFxPPP1QTYfyTY6PNMgIu8BG4E+IpImInd5Dt3MxYPx44E9IrIb+Ai4VylVPqB/H/AqkAocppYZXuAfh/IZcL+IvA8kAbl6/ETjK3aXnQCr8Tft8hjsLpMKilLgdILVanQ0hqCU4r53trM+NYtnbkzgZ8PcvTF/uXEQKRl53P/eTv79wFjahgcbHKnxKKVmVVN+RxVly4Hl1dTfBgysS9u1CopH7SYCMSKSBjwOBHgaXIR7fvLVuFWsCJhblwA0LRuzOBSruG/UpnUo4I6phQrK8ewiVh/M5NdX9GbG8P92iIQHB7DwtmFc/9IG7l6ynct6tKHU4aLU4aRdeDD3TOhBoE2v324sfJnlVaXaeR1XwHy/RaRpUZjFoYgIAZYA8zoUcMcU3DK/gSdnuGfbTezT9qJjfdqH89efJfDbj3aTnJ5LkM1KoM1CTmEZm4/msPC2oS1xJpgh6JXyGkOxO+3YLOb4N7RZbBc7FJsJYiuPoQWPoySn52G1CL3atary+HWDO3BtQtwFs70+2p7GI8v3cNO/NvHm3BG0i2iZYtyYaC+oMRS7yxxdXuAeRzG9Q2mhJGfk0TO2FcEB1Xf5VZ46/LNh8bx+xwhOZBdyw8s/kHo2v6HDbPFoQdEYilnGUMA908v0YygtlOT0PPp3iKjzeeN7x7LsntGUOlzcvWQ77h56TUOhBUVjKGYZQwHtUMxKdkEpp/NK6B9Xd0EBGNgxkgev6MXRrEIOZxb4OTqNN1pQNIZiOoeiBcV0pGS4u6oG1MOhlFM+mP+fA2f9EpOmarSgaAzFdA5Fd3mZjuQM907r/erpUAA6tg6hb/twVh/I9FdYmirQgqIxFO1QfKClC0p6Hh0ig4kKu7T9/yb2acvWYznklbTM97Ex0IKiMQyny4lCaYdSGy1dUDLqNyBfmcl92+JwKTYcyvJDVJqq0IKiMYxyN2Bqh2KGdSgtWFBK7E4OZxbWe0Dem6GdWxMRbNPjKA2IFhSNYZS7AdM6FJvN+K3roUULyo9n8nG6lF8cis1qYXzvWNb8mInLpacPNwRaUDSGYXqHYobuLmjRgpKc7t5ypX9cpF+uN6lPWzLzS9mfbpLEac0MLSgawzC9Q9GCYjjJGXmEB9mIjwrxy/Um9IlFRE8fbii0oGgMQzsUH2nJgpKeR7+4CCwW/3Q9xrQKIiG+NasPakFpCLSgaAxDOxQfaaGC4nIpUvw0w8ubyX3asjvtPNkFpX69rkYLisZAtEPxkRYqKMdziigsc/plhpc3k/rGohR8/6Ne5OhvtKBoDEM7FB9poYJSMSDvZ4cysEMkMa2C+Db5jF+vq9GCojEQ7VB8pIUISk5hGSdziipeJ2fkYrMIPdtWnQOlvlgswvTEDnybfIazeSV+vXZLxwSrtjQtFe1QfKSFCMr8pTvYeCSbHrFhXN6vHZuPZNOzbc05UOrL7aO68Nr6o7y75QS/+klvv1+/paIdisYwtEPxkRYgKC6XYnfaeYZ1iaJD6xDe2HCU3Wm5DOzon/UnlekaE8bEPrEs3XyCMoerQdowChF5XUTOisg+r7InROSUiOzyPK72OvaoiKSKyEERudKrfJiI7PUcWyCVM5hVgXYoGsMwnUOpnGBLC0qjcTyniKIyJzOHd+KmEZ0oKHWw5Wg2Azs0jKAAzBndlblvbuXr/ae5dnCHBmvHAN4EXgSWVCr/h1LqWe8CEekP3AwMADoAq0Skt1LKCSwE5gGbgJXAVcCXNTWsHYrGMEznUCon2NKC0mikZLgH4Mu3qG8VZGNy33a0bcA88BN6x9I5OpQlG481WBtGoJRaC+T4WH068L5SqlQpdRRIBUaKSBwQoZTaqNxpLpcA19d2MS0oGsPQDsVHWoCgJKfnYbUIvdr5dwC+JiwWYfboLmw9do796bmN1q4fsInINq/HPB/Pu19E9ni6xKI8ZR2Bk1510jxlHT3PK5fXiBYUjWFoh+IjLUBQUjLy6BEb1iAD8DUxY1gnQgKsvL3xeKO2e4k4lFLDvR6LfThnIdADSAQygL97yqsaF1E1lNeIFhSNYWiH4iMiYLU2e0G5lIyM9SUyNIDrh3Rkxa5TnC8qa/T2Gwul1BmllFMp5QJeAUZ6DqUBnbyqxgPpnvL4KsprRAuKxjDM6FAUCqfLaS5BAXcszVRQzheVkZ5bYoigAMwe3YUSu4sPtp2svXITxTMmUs5PgfIZYJ8BN4tIkIh0A3oBW5RSGUC+iIzyzO6aDXxaWztaUDSGYUaHAh6h04LSaCRXGpBvbPrFRTCkc2s+3nHKkPb9jYi8B2wE+ohImojcBTzjmQK8B5gEPAiglNoPfAAkA18B8z0zvADuA17FPVB/mFpmeIGeNqwxEDM6FHALXbAWlEYjJSMfwO97dtWF6xM78vhn+zlwOo++7Y2Lwx8opWZVUfxaDfWfBp6uonwbMLAubWuHojEM7VDqQLMWlDxiWgURGx5kWAzTEuKwWoQVO2sdJtDUgBYUjWGY2aFoQWk83DlPwg2NIaZVEON7xfDZrlM6PfAloAVFYxjaodSBZioodqeL1LMFft9RuD5cP6Qj6bklbD3m65pATWW0oGgMQzuUOtBMBeVwZgFlTpeh4yflXNG/HaGBVlbs0t1e9cUnQRGRqzwbh6WKyCNVHI8Ukc9FZLeI7BeRuf4PVdPc0A6lDjRTQam85YqRhAbamNK/HSv3ZjS7DSMbi1oFRUSswEvAVKA/MMuzoZg384FkpdRgYCLwdxEJ9HOsmmaGdih1oJkKSnJ6HoE2C91jwowOBYDpQzqSW2xnjc45Xy98cSgjgVSl1BGlVBnwPu4NxbxRQLhnAUwr3BuTOfwaqabZYXfasYoVH3bFbhS0Q2l8UjLy6dMuHJvVHL3v43rG0CYskE91t1e98OWvWN3mYd68CPTDvTR/L/BLzxL/CxCReeUbmjkcWm9aOnaX3TTdXeDlUMpKQCktKA2MUsqz5YqxM7y8sVktXJMQx6qUM+SXNK/3uzHwRVB82STsSmAX7v30E4EXReSiTlGl1OLyDc1sNr2msqVjd9pN090F/3UojtJiT4F5YmuOgnI2v5TswjJTjJ94M31IR0odLp1zvh74IijVbR7mzVzgY+UmFTgK9PVPiJrmilkdirPUk2dcC0qDUr7lihlmeHmTGN+ayJAAthzV04frii+CshXoJSLdPAPtN+PeUMybE8DlACLSDugDHPFnoJrmh2kdSpl2KI1B+QyvviYTFItFSOzUmp0nzhsdSpOjVkFRSjmA+4GvgRTgA6XUfhG5V0Tu9VR7EhgjInuB74CHlVJZDRW0pnmgHUodaGaCopTiu5SzdIsJIzLERO+zh6Gdo/jxbD55ehylTvg0kKGUWok7p7B32SKv5+nAFP+Gpmnu2F3mdCjOMi0oDc2mIzlsP36OP00fYHQoVTK0S2uUgt0nzzOuV6zR4TQZzDFXT9MisTu1Q/GZZiYoL61OJaZVEDcN71R7ZQNI7NQaEdhxXHd71QUtKBrDMKtDcZWVegrME1tzEpSdJ86xPjWLeeO7NXrKX18JDw6gd9twdpw4Z3QoTQotKBrD0A6lDjQjQXlpdSqtQwO4NamL0aHUyNAurdl54pzefbgOaEHRGIZpHYpdO5SGYn96LqtSznLnZd0ICzL3WrQhnaPIK3FwJKvA6FCaDFpQNIZhVofiKtWC0lC8vPowrYJszBnd1ehQamVo5yhAj6PUBS0oGsPQDqUONANBST1bwMp9Gcwe3YXIUBO9t9XQPSaMiGCbHkepA1pQNIZhWoeiB+X9Tondya8/2EVogJU7x3YzOhyfsFiEIZ2jmtwCRxF5XUTOisg+r7K/icgBEdkjIp+ISGtPeVcRKRaRXZ7HIq9zhonIXk/akgXiwy6uWlA0hmFWh6LKyjwF5omtKQuKUoqHl+9hz6lc/jEzkZhWxuWOrytNdIHjm8BVlcq+BQYqpRKAH4FHvY4dVkoleh73epUvBOYBvTyPyte8CC0oGsMwq0MxraC4XO5HE2Ph94f5dFc6v5nShykD2hsdTp3wXuDYVFBKrcWdQsS77BvPricAm3DvyVgtIhIHRCilNiqlFLAEuL62trWgaAzDtA7FblJBgSbnUlYln+FvXx/k2sEd+MXEHkaHU2ea6QLHO4EvvV53E5GdIvK9iIzzlHXEvTFwOVWlLbkIc8/b0zRrzOZQRASrWM3rUMAtKEFNo8voZE4Rv3x/JwM7RPLMjQmmSaRWF0y6wNEmItu8Xi9WSi325UQR+T/cyQ+XeooygM5KqWwRGQasEJEB+Ja25OLAfAlCo2kIzOZQwN3tpUpMLihNhP8cOEthmZMFs4YQEmjOFfG+MLRLa/69JwOXS2GxmEIUHUqp4XU9SUTmANcAl3u6sVBKlQKlnufbReQw0Bu3I/HuFqsqbclF6C4vjWGYbft68HR7ld+0taBcEsnpeUSFBtC1TajRoVwS5QscUzOb7gJHEbkKeBi4TilV5FUeKyJWz/PuuAffjyilMoB8ERnlmd01G/i0tna0oGgMw2zb14PHoegxFL+QcjqP/h0immRXlzfjesUQaLPw0upUo0PxCRF5D9gI9BGRNBG5C3ea9nDg20rTg8cDe0RkN/ARcK9SqnxA/z7gVSAVOMyF4y5Voru8NIahHUodaGKC4nC6OHA6n9mjzL1fly/ERYZw74QeLPjuELeM7ExS9zZGh1QjSqlZVRS/Vk3d5cDyao5tAwbWpW3tUDSGYVaHogXl0jmaVUiZw0X/DubKxlhf7pvQg46tQ3j8s/04nE1v6nZjoQVFYxhmdSji8EzX14JSbyryxTcTQQkJtPL7a/px4HQ+SzefMDoc06IFRWMY2qHUgaYmKOl5BFot9IhtZXQofuPKAe0Z1yuGv39zkKyCUqPDMSVaUDSG4FIuXMplTodi1w7lUknOyKNXu1YEWJvPLUZEePzaARSVOfnbVweNDseUNJ+/tqZJYXe6b4xmdChaUC6dlIw8+sU1j+4ub3q2bcUdY7qybNtJ0s8XGx2O6dCCojEEu8sjKGZ0KA4HWK1gpumuTUhQzuaXkFVQRv9mKCgA1w7uAMCuJrS/V2OhBUVjCKZ2KA6HudwJNClBSU53D8g3R4cC0DcunACrsDtNC0pltKBoDMHcDsWpBeUSSMnIB2i2DiXIZqVv+wj2nMw1OhTToQVFYwhmdigWu3Yol0JyRh4dW4c0iayM9SUhPpJ9p3JxuWrdL7FFoQVFYwhmdigW7VAuieY6IO/N4PjW5Jc6OJpdaHQopkILisYQTO1QtKDUmxK7kyOZBc1mQWN1JHSKBGCPHke5AC0oGkMwt0NxaUGpJwdP5+NS0D8u3OhQGpSesa0ICbCyW4+jXIAWFI0hmNmhWJ3aodSXii1X4iINjqRhsVktDOgQwd5TWlC80YKiMQTtUOpIExGUlIw8WgXZiI8KMTqUBichvjX703P1ZpFe6O3rNYZgWodiCcDqcEFgw8Rlt9tJS0ujpKSkbic6nfDllxAdDSkpDRKbP9h2+DRdIm0cPHjA6FAanBhLISV2F19v2kP3aP+nZQ4ODiY+Pp4As325qQGfBMWT7esFwAq8qpT6SxV1JgLPAwFAllJqgh/j1DQzTOtQrAFYnQ3nUNLS0ggPD6dr1651SzzlcEBJCXTqBO3aNUhsl4rLpTj+3nFuHBZPv379jA6nwQmOLeSZdWfJD4imX7/Ofr22Uors7GzS0tLo1q2bX6/dkNTa5eVJD/kSMBXoD8wSkf6V6rQGXsadXnIAMKMBYtU0I0ztUBpQUEpKSmjTpk3dsxiW11fmXfdwPKeIwjJns58yXE6X6FDCg23saYBxFBGhTZs2dXeyBuPLGMpIIFUpdUQpVQa8D0yvVOcW4GOl1AkApdRZ/4apaW6Y2aHYGngMpV4pcZuAoGw6kg3AiK7RBkfSOFgsQkJ8ZINNHW6KqZN9EZSOwEmv12meMm96A1EiskZEtovI7KouJCLzRGSbiGxzlCcx0rRIzO1QlPkG5f0kKK1a1S8/ycSJE9m2bVuNdd75divBrhJ6xIbVq42mSEJ8aw5k5FNidxodiinwRVCqksnK/9U2YBgwDbgS+L2I9L7oJKUWK6WGK6WG22x6PkBLxtQOpRkLSkPhcimOFAUS48xqkt+s68vg+EgcLkWKZ7p0S8cXQUkDOnm9jgfSq6jzlVKqUCmVBawFBvsnRE1zxMwOJcAJKsCEX3hE/CYoSikeeughBg4cyKBBg1i2bFnFsWeeeYZBgwYxePBgHnnkkQvOc7lczJkzh8ceewyAN954g969e3PZtJsodlqIcbq7vY4fP87ll19OQkICl19+OSdOnMDpdNK9e3eUUpw/fx6LxcLatWsBGDduHKmpqTzxxBPceeedTJw4ke7du7NgwQK//L4NRUJ8awBTrUcRkddF5KyI7PMqixaRb0XkkOdnlNexR0UkVUQOisiVXuXDRGSv59gC8eGbgi+fmq1ALxHpBpwCbsY9ZuLNp8CLImIDAoEk4B8+XFvTQjGzQwlwgctmw9rAbf3qq1+x6/Qu308oKHA7p6Dqp6gmtk/k+auer/VSH3/8Mbt27WL37t1kZWUxYsQIxo8fz65du1ixYgWbN28mNDSUnJycinMcDge33norAwcO5P/+7//IyMjg8ccfZ/v27Xy49xzPfHOoQlDuv/9+Zs+ezZw5c3j99dd54IEHWLFiBb179yY5OZmjR48ybNgw1q1bR1JSEmlpafTs2ROAAwcOsHr1avLz8+nTpw/33XefaafOxkUGE9Mq0L1ifrTR0VTwJvAisMSr7BHgO6XUX0TkEc/rhz0TrG4GBgAdgFUi0lsp5QQWAvOATcBK4Crgy5oartWhKKUcwP3A10AK8IFSar+I3Csi93rqpABfAXuALbinFu+r7poajdkdisvW0HJiLOvXr2fWrFlYrVbatWvHhAkT2Lp1K6tWrWLu3LmEhoYCEB393wH2e+65p0JMADZv3szEiROJjY1l49FztLGVEazcudY3btzILbe4v3fefvvtrF+/HnA7kbVr17J27VoeffRR1q9fz9atWxkxYkRFO9OmTSMoKIiYmBjatm3LmTNnGuU9qQ8iQkJ8a7YdzyGnsMzocABQSq0FcioVTwfe8jx/C7jeq/x9pVSpUuookAqMFJE4IEIptVEppXCL0/XUgk++Xim1ErdCeZctqvT6b8DffLmeRmN+h9LwguKLk7iAXbsgKgq6dLnktlU1XWdKqWrHQMaMGcPq1av59a9/TXBwMOC+oZbYnWw5mkP/oKJq2yu/5rhx41i0aBHp6en86U9/4m9/+xtr1qxh/PjxFXWDvByY1WrF7BN4rhrQnt8u30PS/1vFFf3bMWN4J8b3isVqabCxJJuIeM+QWKyUWlzLOe2UUhkASqkMEWnrKe+I24GUUz7pyu55Xrm8RvTWKxpDML9DMeFHw49jKOPHj2fZsmU4nU4yMzNZu3YtI0eOZMqUKbz++usUFbnFwbvL66677uLqq69mxowZOBwOkpKSWLNmDav3HqfU4eLk1m8r6o4ZM4b3338fgKVLlzJ27FgAkpKS+OGHH7BYLAQHB5OYmMi//vUvxo0b55ffywhuGtGJr381ntmju7LpSA5z39jKve9sb8gmHeWTmzyP2sSkJqqbdOXLZKyLMOGnRtMSML1DsZqwy8uPgvLTn/6UhIQEBg8ezOTJk3nmmWdo3749V111Fddddx3Dhw8nMTGRZ5999oLz/vd//5ehQ4dy++23065dO5544gnmP7kAlIuR3f7bPbZgwQLeeOMNEhISePvtt3nhhRcAt/vo1KkTo0aNAtyOJT8/n0GDBvnl9zKKPu3D+f01/dn06OXcktSZ/xw4S36JqfZdO+PpxsLzs3ytYHWTrtI8zyuX14xSypBHaGio0rRcnln/jOIJVH5pvtGhXMDrO15X2cGo3LtnN8j1k5OT63/ynj1KHT7sv2D8xHX/XKdufHmD0WGYhvWHMlWXh79Q/zlw5pKvVdX/C1Coarm/Al2BfV6v/wY84nn+CPCM5/kAYDcQBHQDjgBWz7GtwCjcbuVL4Ora2tUORWMIZncoTqsJPxp+dCj+IrfIzp5TuYztFWN0KKZhaOcoAqxSsXNAYyMi7wEbgT4ikiYidwF/Aa4QkUPAFZ7XKKX2Ax8AybgnVs1X7hleAPcBr+IeqD9MLTO8QO82rDEIs4+hlDbzMRR/sfFIFkrB2J5aUMoJCbQyOL41m49UnmjVOCilZlVz6PJq6j8NPF1F+TZgYF3aNuGnRtMSsLvsWMSCRcz1L6gdSt1Yn5pFWKCVwZ1aGx2KqUjqHs3eU7kUlJp7hpq/MeGnRtMSsDvtpuvuAgjAilWBQzsUn9iQms2o7m0IMKMAG8io7m1wuhTbj58zOpRGRf8XaAzB7rKbrrsLINDlni3ptJpwPyqTCcrJnCKOZhVyme7uuohhXaKwWYwbRzEKLSgaQzCrQwlyuT8SDi0otbI+NQuA8b21oFQmNNBGQnykFhSNpjEwrUNRbiHRglI927Zt44EHHmDdoUziIoPpEXvxlvhr1qzhmmuuqfVa3tviX3311Zw/X31ukeeff75iwWVTYFT3NuxNy6WwBY2jaEHRGIJZHUqAZ8KkFpTqGT58OP94/gU2pGYztmeM37arX7lyJa1bVz+439QEJal7GxwtbBxFC4rGEEzrUDxjKHYzfjL8ICiFhYVMmzaNwYMHM3DgwIpt67/77juGDBnCoEGDuPPOOyktdW/yuHXrVsaMGcPgwYMZOXIk+fn5rFmzhitm3kVusZ0OllzGjBnDkCFDGDNmDAcPHqyx/eLiYm6++WYSEhKYOXMmxcXFFce6du1KVlZWlTEuWLCA9PR0Jk2axKRJkwC47777GD58OAMGDODxxx+/4DqPP/44Q4cOZdCgQRw4cACAgoIC5s6dy6BBg0hISGD58uUAfPPNN4wePZqhQ4cyY8YMCgoKLuk9Lmd4lyisFmHz0ZbT7aXXoWgMwe4yqUNxuX82ikP51a/cGz76SnExuFwQVkNGxMREeL76TSe/+uorOnTowL///W8AcnNzKSkp4Y477uC7776jd+/ezJ49m4ULF/KLX/yCmTNnsmzZMkaMGEFeXh4hISEAFLSKRwRuuGwA/zN9LTabjVWrVvG73/2u4kZdFQsXLiQ0NJQ9e/awZ88ehg4d6lOMkZGRPPfcc6xevZqYGPeYzdNPP010dDROp5PLL7+cPXv2kJCQAEBMTAw7duzg5Zdf5tlnn+XVV1/lySefJDIykr179wJw7tw5srKyeOqpp1i1ahVhYWH89a9/5bnnnuMPf/hDDX8I3wgLKh9HMWY9ihGY8XuYpgVgd2qHUmf80LU0aNAgVq1axcMPP8y6deuIjIzk4MGDdOvWjd693UlW58yZw9q1azl48CBxcXEVW8tHRERQnmm1oFUnBnaIxGIvYsaMGQwcOJAHH3yQ/fv319j+2rVrue222wBISEioEIDaYqyKDz74gKFDhzJkyBD2799PcnJyxbEbbrgBgGHDhnHs2DEAVq1axfz58yvqREVFsWnTJpKTk7nssstITEzkrbfe4vjx4768lT6R1K0Ne9LOU1TWMsZRtEPRGILD5TClQ7E1pkOpwUlUybFjkJcHVdyEfaV3795s376dlStX8uijjzJlyhSuu+66KuuqarayL3YoikLbMa5XDL///e+ZNGkSn3zyCceOHWPixIm1xlDbmEtVMVZ2DEePHuXZZ59l69atRJIMQFYAACAASURBVEVFcccdd1BSUlJxvHwLfO/t76v6fZRSXHHFFbz33nu1xl0fRnWPZtH3h9lx/HyL2J7GjN/DNC0As46hlA/Kl1mNH/y+CD+MoaSnpxMaGsptt93Gb37zG3bs2EHfvn05duwYqampALz99ttMmDCBvn37kp6eztatWwHIz8/H4XBwIMcJYmVsrxhyc3Pp2NGdJuPNN9+stf3x48ezdOlSAPbt28eePXt8ihEgPDyc/Px8APLy8ggLCyMyMpIzZ87w5Ze1bjPFlClTePHFFytenzt3jlGjRrFhw4aK372oqIgff/yx1mv5yvCu0Vhb0HoU7VA0hmD2WV6m7fK6REHZu3cvDz30EBaLhYCAABYuXEhwcDBvvPFGRZ6TESNGcO+99xIYGMiyZcv4n//5H4qLiwkJCWHVqlXsy3IiLjvDukTx29/+ljlz5vDcc88xefLkWtu/7777mDt3LgkJCSQmJjJy5EifYgSYN28eU6dOJS4ujtWrVzNkyBAGDBhA9+7dueyyy2pt+7HHHmP+/PkMHDgQq9XK448/zg033MCbb77JrFmzKiYiPPXUUxXdf5dKqyAbAztGtpiBeVEGTUMMCwtThYWFhrStMZ5Jb03CpVx8f8f3RodyATmrVxI9eRqfv/ALrn3gJb9fPyUlhX79+tXv5JMnISsLhgzxb1B1ZPKza+jSJpQ35l4sBpqL+fPKFN7YcIzdj08hJLBueXaq+n8RkSKlVA0zM4zDjN/DNC0A0zsUEy5DMcM6lJM5RRzJKmRcr1hD42hKJHWPpszpYueJ5r8eRQuKxhDMOoZic7pv2M11DOVS0dut1J3hXaOxCGw62vynD2tB0RiCWR2KtVxQxMSCYqCorDuUSfuIqrdb0VRNRHAAAzpEsrkFDMxrQdEYgnkdinvecKnFZXAkVVA+5dUgQckvsfP9wUwm9I7123YrLYWkbtHsPHmeEruz9spNGC0oGkMwq0OxONwf+DKLSR0KGCYoy7enUVjmZFZSZ0Pab8okdW9DmcPF7pPVb37ZHNCCojEEszoU7O7UxNqhXIjLpViy8TiJnVqTqLMz1pmRXaMRodlvw6IFRWMIZnUoFYIiWlC8WXsokyNZhdwxpmujt90ciAwNoG/7iEZZjyIifURkl9cjT0R+JSJPiMgpr/Krvc55VERSReSgiFxZ37a1oGgMwaybQ7Y0h/LEE0/w7LPP1lhn0aJFPP3hBmLDg/j4hcf46KOP6txObblOAA4cOEBiYiJDhgzh8OHDdbr+mjVr+OGHH+ocV/kOx5VRSjF58mTy8vKqPM+X38c71wu4x1F2nDjH1u07WblyZUX5F198ccFuyZeKUuqgUipRKZUIDAOKgE88h/9RfkwptRJARPoDNwMDgKuAl0WkbgtmPGhB0RiCWTeHLBeUEjHh4KlBDuXKn93OoYJAbhnZGQv1a7u2XCcAK1asYPr06ezcuZMePXrU6fr1FZTqWLlyJYMHDyYiIuKCcqUULpfLp9+nMqO6t6HE7uLLzfsvEJRp06bx2WefNVSul8uBw0qpmna8nA68r5QqVUodBVKBeq1a1VuvaAzB9A6lEbq8/vj5fpLTq/4GXCUOOxSXwNZtYKn6u2D/DhE8fu2AGi/z9NNPs2TJEjp16kRsbCzDhg0D4JVXXmHx4sWUlZXRs2dP3n77bUJDQ/nFP5ZhoS23JnVmt3sbLr777jtefPFFPvnE/cX322+/ZeHChXz88cdVttm1a1e2bdtGQUEBU6dOZezYsfzwww907NiRTz/9lNWrV/P8889jtVpZu3Ytq1ev5p133mHBggWUlZWRlJTEyy+/jNVq5auvvuJ3v/sdTqeTmJgYXnvtNRYtWoTVauWdd97hn//8J3379uXee+/lxIkTgDs512WXXUZ2djazZs0iMzOTkSNHUt1OIUuXLmXevHkAHDt2jKlTpzJp0iQ2btzIihUrmDBhAtu2bSMmJoYnn3ySpUuX0qlTJ2JiYhg2bBi/+c1vAPjwww/5xS9+wfnz53l+4SsAvPrZWoq3f8L69et59NFHmTlzJhMnTuSLL77gpptuqvFvVw9uBrx3vrxfRGYD24BfK6XOAR2BTV510jxldUY7FI0hmN6hWEzoULj0qbrbt2/n/fffZ+fOnXz88ccVGz+Ce8v3rVu3snv3bvr168drr71GQamDVGcMPYPyaBsRXFF38uTJpKSkkJmZCcAbb7zB3LlzfYrh0KFDzJ8/n/3799O6dWuWL1/O1Vdfzb333suDDz7I6tWrSUlJYdmyZWzYsIFdu3ZhtVpZunQpmZmZ3H333Sxfvpzdu3fz4Ycf0rVr14pzd+3axbhx4/jlL3/Jgw8+yNatW1m+fDk///nPAfjjH//I2LFj2blzJ9ddd12F4FRmw4YNFUILcPDgQWbPns3OnTvp0qVLRfm2bdtYvnx5xfvp3cUF4HA42LJlC88//zzP/flJ+rQLp+eYqcycOZNdu3Yxc+ZMwJ0Fc926dT69f4BNRLZ5PeZVVUlEAoHrgA89RQuBHkAikAH8vbxqFafXy4pqh6IxBLM7lBIaXlBqcxIXce4cHD4M/ftDaGi92ly3bh0//elPCfWc7711/b59+3jsscc4f/48BQUFXHnllSzfnoYdKwnBF85OEhFuv/123nnnHebOncvGjRtZsmSJTzF069aNxMRE4MJ8Jd589913bN++vSIXS3FxMW3btmXTpk2MHz+ebt26ARAdHV1lG6tWrbogP0peXh75+fmsXbu2wkVNmzaNqKioKs/PyckhPDy84nWXLl0YNWrURfXWr1/P9OnTKxKPXXvttRccr5yX5Zbu0by3KY++le7hbdu2JT09vcpYqsChlBruQ72pwA6l1BmA8p8AIvIK8IXnZRrQyeu8eMDnYLzRgqJpdJRS7nwo2qHUDT+NoVS3KPGOO+5gxYoVDB48mDfffJM1a9bw+e50oqWQ9gElF9WfO3cu1157LcHBwcyYMaMi+VZtlOcqAXe+Eu80wOUopZgzZw5//vOfLyj/7LPPfFpU6XK52LhxY8WN3htfzrfZbLhcLiyersWwarJk1ra5buW8LEnd2rBk43FyLRcmDSspKaky1ktkFl7dXSISp5TK8Lz8KbDP8/wz4F0ReQ7oAPQCttSnQZ+6vETkKs90slQReaSGeiNExCkiP6tPMJqWgcPlTnjU0h1KnfGDoIwfP55PPvmE4uJi8vPz+fzzzyuO5efnExcXh91uZ+nSpbgQ9qXn0t6SX+W1OnToQIcOHXjqqae444476h1TVVx++eV89NFHnD17FnA7huPHjzN69Gi+//57jh49WlEOF+ZKgYtzn+zypFr2zsfy5Zdfcu5c1Rs29unThyNHjtQa59ixY/n8888pKSmhoKCgIm1xdSR1dzuqdMeFAvXjjz8ycODAWtvzFREJBa4AvAe1nhGRvSKyB5gEPAiglNoPfAAkA18B85VS9foA1CoonuljL+G2T/2BWZ5pZlXV+yvwdX0C0bQc7C73TdusDsUpUIYJU7b6QVCGDh3KzJkzSUxM5MYbb2TcuHEVx5588kmSkpK44oor6Nu3L0W2CErsLmIs1aeZuPXWW+nUqRP9+190S7gk+vfvz1NPPcWUKVNISEjgiiuuICMjg9jYWBYvXswNN9zA4MGDK8Ygrr32Wj755BMSExNZt24dCxYsYNu2bSQkJNC/f38WLVoEwOOPP87atWsZOnQo33zzDZ07V73qf9q0aaxZs6bWOEeMGMF1113H4MGDueGGGxg+fHi1KYsBYloFcePgdmQfSyExMZFly5YBsHr1aqZNm1bHd6l6lFJFSqk2Sqlcr7LblVKDlFIJSqnrvNwKSqmnlVI9lFJ9lFK1ZyuroeEaH8Bo4Guv148Cj1ZR71fAfOBN4Ge1XTc0NFRpWibni88rnkA998NzRodyMQ8/rEpsoq5YckWDXD45Obn+J+flKbV1q1K5uf4LqAaWbTmhujz8hTp8Nr/aOvPnz1evvvpqo8TTmKSnp6uf/OQnPtXNz3e/P4WFhWrYsGFq+/btdWrr9OnTavLkyVUeq+r/BShUtdxfjXr40unZETjp9ToNSPKuICIdcffJTQZGVHchz2yEeQCBgYE+NK1pjpjeodgsFTGaikZeh7I77TzhQTa6tql6/GDYsGGEhYXx97//vcrjTZm4uDjuvvtu8vLyLlqLUpl58+aRnJxMSUkJc+bMYejQoXVq68SJE83mPfRFUHyZUvY88LBSylnTgJdSajGwGNwZG30NUtO8sDs9gmLSMRSHVSpiNBWNLCh7T+UyKD4Si6Xqz/T27dsvKktKSqpIpVvO22+/zaBBgxokxobE1zUh77777iW1Uz6TrTngi6D4MqVsOPC+R0xigKtFxKGUWuGXKDXNCtM7FKt2KKUOJykZedw5tludztu8eXMDRaRpCvgiKFuBXiLSDTiFe+XlLd4VlFIV/3Ui8ibwhRYTTXWY3aE4bZYGdShKqfrlE2lEQTl4Oh+7UzE4Xu8sbBTK4Oyc9aHWWV5KKQdwP+7ZWynAB0qp/SJyr4jc29ABapof5d/+bRYTLoOy23FZGs6hBAcHk52dXb+bRSMKyu409+SgQR2rn7GkaTiUUmRnZxMcHFx7ZRPh0ydauXelXFmpbFE1de+49LA0zZkKh2LWLq8GdCjx8fGkpaVVbFlSJxwOyMpyC4pnfUZDsW5fJhFBFvJPHyPljM7OaATBwcHEx8cbHUadMOFXRE1zp2IMxaRdXq4GnOUVEBBQsW1InTl5EhIS4JVXwLM3VUNx4utMhnRp4/f1JZrmjd4cUtPomN+hWM05yyvA837ZGza24jInP57JJyFed3dp6oYWFE2jY36HYjXnLK9GEpTkjFxcChL0gLymjmhB0TQ6ZncoLmvLdii7T7oH5LVD0dQVLSiaRsfsDkUFtGyHsvdULm3Dg2gX0bRmGGmMRwuKptExvUOx2Vq2Q0k7r7u7NPVCC4qm0TG9QzHrGIrV6v7ZgIKSX2LnSGah7u7S1AstKJpGx+wORdlsOFwO861UFnG7lAYUlL2n9PiJpv5oQdE0OqZ3KAHu5VnlicBMRQMJilKKo1mFfLzjFKBXyGvqh17YqGl0TO9QAtxpW+0uu/li9LOgnCss469fHWDdoSxOnXen4k3qFk2bVkG1nKnRXIwWFE2jY3aHgs2d/8PutIPZQvSzoCz6/jAfbk/jJ/3acu+E7oztFUvXNqF+u76mZaEFRdPomN+huOMy5cC8nwVlfWoWw7tE8a/bh/vtmpqWix5D0TQ6pnco5YJi1qnDfhKU7IJS9qfnMa5XjF+up9FoQdE0OmZ3KLQQh7LhcDYAl/XUgtLcEJFjIrJXRHaJyDZPWbSIfCsihzw/o7zqPyoiqSJyUESurG+7WlA0jY52KJeAHwVl/aFMIoJtehFj82WSUipRKVXen/kI8J1Sqhfwnec1ItIfd+LEAcBVwMsiYq1Pg1pQNI2O2R2KtACHopRi/aEsxvSIwVpNznhNs2M68Jbn+VvA9V7l7yulSpVSR4FUYGR9GtCComl0ym/U1vp9CWpY7HYkMND9tBk7lKNZhaTnljBWj580RWwiss3rMa+KOgr4RkS2ex1vp5TKAPD8bOsp7wic9Do3zVNW98Dqc5JGcynYnXYCLAH1y6vekCjlzooY4BGUZuxQ1qdmATBWj580RRxe3VjVcZlSKl1E2gLfisiBGupW9UGs1zYR2qFoGh1TLhgEt5hAi3Ao6w9lER8VQhe95qRZopRK9/w8C3yCuwvrjIjEAXh+lueRTgM6eZ0eD6TXp10tKJpGp9yhmA7PjdrSzB2Kw+li4+FsxvWKMZ9L1FwyIhImIuHlz4EpwD7gM2COp9oc4FPP88+Am0UkSES6Ab2ALfVpW3d5aRod0zoUz43a9A6lqOiSLrE7LZf8UoeeLtx8aQd84vmyYAPeVUp9JSJbgQ9E5C7gBDADQCm1X0Q+AJIBBzBfKeWsT8NaUDSNjvkdyn/38jIdfnAo6w9lIQKX9dCC0hxRSh0BBldRng1cXs05TwNPX2rbustL0+iY3qEEeQTFrA7lEgVlQ2oWAztEEhUW6KegNBo3WlA0jYpSijLtUOrPJQpKQamDHSfO6e4uTYOgBUXTqDz7zUH2JP/E1A7F2owdyrofM3G4lN6/S9MgaEHRNBqlDidLN5+goDAGmwQbHc7FlDuUwObrUN7dcoK4yGCSukX7MSiNxo0WFE2jsfpAJueL7IAFqzLhN+RyhxLoFrvm5lCOZRWy7lAWN4/ojM2qP/oa/6P/qzSNxsc70ihf9iCOWGODqYoKh+IRlGbmUN7bcgKrRZg5olPtlTWaeqAFRdMo5BSWsfrgWa5J6OAucJrYoQQ1P4dS6nDywbaTXNGvHe0jTdjdqGkWaEHRNApf7EnH7lTcO6E7InZcDhP24Xtu1LZm6FC+2neac0V2bh3VuQGC0mjc+CQoInKVJ/FKqog8UsXxW0Vkj+fxg4hctKhG07JZvuMU/eIiGNAhEovtPE5HVO0nNTbN2KG8s+k4XdqE6sWMmgalVkHxJFp5CZgK9AdmeRKyeHMUmKCUSgCeBBb7O1BN0yX1bAG7T57nxqHuHbEtATk4ykyY1KncoQSFuF+a2aEo3zeDPXg6n63HznFrUmcsOveJpgHxxaGMBFKVUkeUUmXA+7gTslSglPpBKXXO83IT7t0qNRoAPtmZhtUiXJfoHj8RazZlZRGoOtwUG4Wm4lAAnL5vtfTu5uMEWi38bJgejNc0LL7s5VVV8pWkGurfBXxZ1QFPopd5AIGBetuHloDLpfhkxynG94qhbbj7Rq2smbhcgZwrshNtpu0/Kk8bNqtDAXesNhsnc4r4y1cHOJZVSKDNQpDNQqDNitXLiGw+msPVg9qb673WNEt8ERSfk6+IyCTcgjK2quNKqcV4usPCwsJM9vVU0xCsPZRJem4Jj17dr6LMaTkDwImcInPd5Lx2Gw6wBJjaodhLSnltSzrPr/oRiwiju7ehzOmi1OEit9h+gfvr2z6ceyb0MCpiTQvCF0HxKfmKiCQArwJTPbtaalo4P6Rm8T/v7qRj6xCu6N+uotwpZxDcgpLYyURjKeWD3QEBBFgDTOtQUmK78uBrOzmQVcSU/u144roBdGgdYnRkGo1PgrIV6OVJvHIKuBm4xbuCiHQGPgZuV0r96PcoNU2OT3ed4jcf7qZbTBhvzh1JcMB/88fb5TSBwInsQuMCrApvQTGxQ3ly8t1kFpbxyuzhFwi1RmM0tQqKUsohIvcDXwNW4HVPQpZ7PccXAX8A2gAve5K6+JLzWNMMUUqxeO0R/vzlAZK6RbN49nAiQy7cCNKuCgkPKOFEzqUlivI7TcChKFsA+9r34NruEVpMNKbDpwRbSqmVwMpKZYu8nv8c+Ll/Q9M0Rb7ef4Y/f3mAaQlxPHfTYIJs1ovq2J12wkKKzS0oJnUoaQSRF9yKAdEmGnvSaDzojI0av/L6hqN0ig5hwc1DsFaz5sHustMqtJSTOcWNHF0tNAGHst/h3gl5QKT+6GrMh956ReM3UjLy2HI0h9tHdalWTMDtUCJCykjPLabM4WrECGuhskMxo6CUBmJ1OenTSi9Q1JgPLSgav7Fk43GCbBZuGl79AjqlFHaXncgwO0pB2jkTdXtVdigm7PLaX2KlR3Yawcr3hY2aloWIdBKR1SKSIiL7ReSXnvInROSUiOzyPK72OudRz9ZaB0Xkyvq2rX2zxi/kFtlZsfMU1yd2pHVo9f37Ts+NsHWYA3BPHe4e26pRYqyVpuBQiiyMOXMY7HrOi6ZaHMCvlVI7RCQc2C4i33qO/UMp9ax3Zc9WWjcDA4AOwCoR6a1U3b+1aIei8Qsfbj9Jsd3J7aO71Fiv/Ft/dLi7q+ukmQbmywXFajWlQ8kqKOVMGQw4e+SSsjZqmjdKqQyl1A7P83wgBfeOJ9UxHXhfKVWqlDoKpOLecqvOaEHRXDIul2LJxuMM7xLFwI6RNdYt/9YfGSIE2Szmmullt7tXoouY0qHsT88DoP+Zw1pQND4hIl2BIcBmT9H9nl3hXxeR8i2/q9peqyYBqhYtKJpL5vsfMzmRU8ScMV1rrVv+rT/QFkDn6FBzCgqY0qHsT88FYMAZ7VBaODYR2eb1mFdVJRFpBSwHfqWUygMWAj2ARCAD+Ht51SpOr9fWWHoMRXPJvLXxGG3Dg7hyQPta65Z/6w+wlAuKiaYOewuKCR1Kcnoe8aEWIksLtaC0bGpdOC4iAbjFZKlS6mMApdQZr+OvAF94Xvq0vZYvaIeiuSQOnM5jzcFMbknqTKCt9n+n8m/9AdYAOkWHciK70Dzb2JvcoSSn5zEgymu3YY2mCsS9XclrQIpS6jmv8jivaj8F9nmefwbcLCJBni22egFb6tO2diiaeqOU4onP9tM6NIA5o7v6dI63Q+nSJpTCMic5hWW0aRVU67l5JXay8kvpFhOGZ4sf/2Jih1JY6uBodiHTB3nGqLSgaKrnMuB2YK+I7PKU/Q53csRE3N1Zx4B7ADxbaX0AJOOeITa/PjO8QAuK5hL4Yk8Gm47k8NT1A4nycRt6b4fSLjoUcE8d9haUnMIyDmTkceB0PgdP53Mkq4CjWYVkFZQBML53LM/cmED7yGD//kImdigpGXkoBQNiPbsKa0HRVINSaj1Vj4usrKKs/JyngacvtW0tKJp6UVjq4Ol/pzCgQwSzRnb2+bzKYyjgFpQhnaPYdyqXhz7aQ0pGXkX9NmGB9Gjbip/0a0fXmDCcLsWL/0llyj++54/TB3B9Ykf/uZVKDqXMYZ7V6OUzvAa0DXMXaEHRmBAtKJp68eLqVE7nlfDSrdXv2VUV3g4lPsotKMezi1i89jB/+/og0WGB/O7qvvSPi6RP+3Biwy/uCps2KI7ffLibB5ft5t97MvjZsHhGd48hMjTgorp1wktQSorbU3LqNv761QF+e2WfhuliqwP703OJDgukfWvtUDTmRQuKps4cySzg1XVHuGFoR4Z1ia7Tud4OJSTQStvwIF5cnUqZw8WVA9rxlxsSau0+6xoTxrJ7RvPquiP88z+prEo5i0VgUMdIJvRpy03D4yvEqm7B/VdQzmb2BqwsXHOYnIIynv7pQGxW4+aw7E/PY0CHCKQ8dbYWFI0J0YKiqRGlFK+tP0pyRh5BNgtBNivbj58jyGblkal963w9b4cC0Kd9ONuOneOvNw7ipuGdfHYCVotwz4QezL2sG7vTzrP+UBYbUrN48T+HePE/h5jctx23j+7CuJ4xWHx1UB5BKS5zkpkTjzNoA7+6bA7//E8q54rKWDBryAWJwhqLMoeLH8/kc+fYbhfmlNdoTIYWFE21KKX44+fJvPnDMdpFBOFSUGp3ohT84dr+tA2v+6C4t0MBeH5mIi5FlV1bvhBoszCiazQjukbz4BW9STtXxHtbTvD+lpOsSjlDQnwkb9wxwqdZZOWC8vX+0zidAZSFreXXU/4f0WGB/PHzZOa8voVFtw3zeQKCvzh0Nh+7U9E/LkILisbUaEHRVIlSij994RaTOy/rxu+v6eeXcYTKDsWnG30diI8K5aEr+/LA5b34bFc6j63Yx6xXNrH056NqFy2PoCzfkUZocDHnZT8Acy/rRlRoIL/9aA/X/HM9L906lMROrS84tdThpMTuuig7ZX3YdCSbdzYdJzoskI6tQ0g75178OaBDJARY/hurRmMytKBoLkIpxZNfpPDGhmPMvayr38QELnYoDUWQzcqM4Z3o2DqEu97axs2LN/Lu3aNoF1GDq7LbyQiNYn1qFn27nOZ4TlnFoeuHdKR7bBj3vbODGYt+4A/X9Oe2UV3YdyqPD7ad5NNdp8grcdA5OpSE+EgGx7eme2wYMa2CiAkPok1YIOeKyjiaWciRrEIycov5Sb92DOkcdUEI724+wR8+3Ud4sA2HS5Ff4t6VOTIkgG4xYeB0VMSq0ZgNLSiaCyixO/nTF8m8u/kEd4zpyh+u6e/XGU6VHUpDM6ZnDG/dOZK5b2xh5r828s7Pk6ofsLfb+aTdQJSCHh2y+S7zwpt2Qnxr/v3AWH61bBe//3Q/C9ccJj23hCCbhakD29OrXTh703LZcfwcX+zJqDW2l1Yf5qoB7fnNlX3o2iaUp/6dwps/HGNSn1hemDWEiOAA8krsZJwvITTQ6p5NJ7aKWDUas6EFRVPBgdN5/PK9XRw8k889E7rzyFV9/T5dtrEcijcju0Wz5K4k7nh9C5OeXcOU/u2ZMTyecb1iL5jyrOx2lkf3Z0TXKKLDndhddpRSF7wHrUMDeX3OCBatPcwPqdncN6kn1w3ucFFXV2Z+KWnnisgqKCO7oJTswjIiQwLoHhNGt9gwWgXZeGPDMRavPcI3yafpFhPG4cxC7h7XjUem9quIKyI4gIj2XtcWAZtNC4rGlGhB0eByKd784Rh/+eoAEcEBvDF3BJP6tG2QthrboZQzrEsUn/3PWJZsPMaKnaf4994MOkQGc/vorswe3YWwIBu7Q9txODiKu4fGc6jYHZ9TObHJhR8Ti0X4xcSe/GJiz2rbiw0PqnXM5oHLe3FrUmdeXnOYT3ed4pmfJdSY7bKCgAAtKBpTogWlhZN6toDHVuxl05EcftKvLX+5MYEYPw+Ue2OEQymnW0wYj187gEem9mVV8lne23KCv351gFfWHWHe+O4cjh9BkMvB1QlxvLzNHZ/dacdmabiPSZtWQfz+mv78/pr+vp+kBUVjUrSgtFBK7E5eXp3Kou+PEBxg4c83DOLmEb6vA6kvRjkUb4JsVqYlxDEtIY4dJ87x/KpD/OXLA9BlJNcVHiMiOKBC8OwuOyGEGBZrlWhB0ZgULSgtDJdL8U3yGf761QGOZhUyPbEDj03rX+91IHXFSIdSFUM7R7HkzpFsP36Odx/6O/Mi3HtmlQuemTaIrEALisakaEFpISjlFpIXVh0iZWecLQAACdRJREFUOSOP7jFhvH3XSMb1im3UOMzgUKpiWJcohq1/HaZNA7jAoZgOLSgak6IFpQWQkpHHQx/tZt+pPLq2CeW5mwZz3eAOhuxNZTaHcgGVtq8H7VA0mrqgBaWZ8/nudH770R7Cg208O2Mw1ycaIyTlmNWhABdtXw/aoWg0dUELSjPgbF4Jn+1OJyI4gBHdounaJhSXgme+OsC/1h5heJcoXr5taL323vI32qH4AS0oGpOiBaWJoJSi2O4kyOZeMa2UYsvRHJZsOs7X+07jcP03L3tMK/dWHwfP5HP7qC78/pr+PuV7bwzKb9ANORW33miHotFcEib8VGsq88PhLP74WTIHz+QDYLMINqtUbEZ4x5iu3JLUGZdSbDl6jq3Hcjh0Np9nbkzgphE+LJRrROwu97oOoxNWXYRS2qFoNJeIFhQDUEqRVVDG6dwS2kcGE9MqsMobbPr5Yp5emcK/92TQKTqE30zp7d5C3uGk1O6id/twrk3oQEjgf3N09Gwbzi1JvqfkbWzsTrs5u7ucTvdP7VA0mnrjk6CIyFXAC4AVeFUp9ZdKx8Vz/GqgCLhDKbXDz7GaGqdLUVDiIDjQnYTKm/NFZXz/YyZrf8zi0Nl8jmYWkl/qqDgeHmSja0wYHVoH43QpSh0uSh0u9qbl4lKKB3/Sm3smdDckuZO/sbvs5h2QB+1QNM2C2u7ZDUWtgiIiVuAl4AogDdgqIp8ppZK9qk0FenkeScBCz09DUUpRYndR5nTVWOdsfimnzheTfr6YM7klFJU5KXO6KLW7cClFdFigZxvyQKwWC8ezCjma5d6G/ExeCfklDgo8AmER6BgVQreYVnSODuHg6Xy2Hz+HS0F0WCADOkTw06Ed6RYTRlxkMKdzSyqudSSzkACrhaAAC0E2C9ckxPHA5b3oFF2PdLYmxbQOpbKgmN2hlJQYHYXGpPh4z24QfHEoI4FUpdQRABF5H5gOeAc3HViilFLAJhFpLSJxSqna9/CuI9veeJrWv3uy4vWyAdNZmvAzAp12Alx2Ap12XCIUBoRSGBiKs46DvxaXk0BnGQFOB4GuMkRBbnAE9krfqmMLMonPyyChKIuwskJalRUSai+mMCCUtIgOnIqIY0d4e9oVnGVW2nZGpW2nT1YqFtQF1+kJjK0hnlIgtU6/gbm532XnAbHAWwOMDuVCKnd5ef7ec1bMoVVgK6OiqpKX00+Q9GMhaXHGz9rT1I+0GVcyccGnDXV5X+7ZDYIvd9uOwEmv12lc7D6qqtMRuEBQRGQeMA8gMLB+aVSDomLJ7BpT8To2IJsJZ3/AbglwP8SGoAhxFBPqLCbUUUyAy1Ht9QRFZFkeMaU5tCnJIaos96KbvgKKrcHkBkbgEBttSzIJqumba24VZa0gu1WbOv62zZM2ITHQppfRYVzM0KEwdSoAie0TuTPxTvLK8gwO6mI2TWtNYPhxo8PQXAIBcR0v5XSbiGzzer1YKbXY67Uv9+wGwRdBqWo6jqpHHTy/9GKAsLCwi477wqDr58H18ypej67PRTSaWmgV2IrXpr9mdBhVM8PoADQG41D/v727C5GqjOM4/v2hvej2ovaGtYIKYonkaiKWEaUVKmG3CpJEEIGRRhBKEHgfURdRRFlQYZBZLRKpWN10Ub7Xmm4aipnlSlDCFqH57+I8yw6L267b2X2egd8HDjPn7Cz75ezMPDPnnJkTMfc/fj6o5+PhMJgPJ5wEGo89bQVODeE2ZmY2/LI9Hw9mQNkFTJM0RdLlwHKgvc9t2oFHVJkP/DEc+0/MzGxAg3nOHhYDbvKKiPOSngS2UR2CtjEiDkp6Iv38NeBTqkOGj1IdNvzo8CWbmVl/+nvOHom/rerArJHX0tIS3d3dWf62mVmzkvRnRLTk7riYMr7gyczMmp4HFDMzq4UHFDMzq4UHFDMzq0W2nfKSLgB/DfHXRwP9f/w9r1LbSu0Ctw1FqV1QblupXXBpbWMiosg3A9kGlP9D0u4BPimaTaltpXaB24ai1C4ot63ULii77VIUOcqZmVnz8YBiZma1aNYB5fWBb5JNqW2ldoHbhqLULii3rdQuKLtt0JpyH4qZmZWnWd+hmJlZYTygmJlZLZpuQJG0WFKnpKOS1mVu2SipS1JHw7IJknZIOpIux2fomiTpC0mHJB2UtKaENklXSvpG0oHUtaGErj6NoyTtk7S1lDZJxyV9J2l/z5n6SuhKHeMkbZZ0ON3f7iyhTdL0tL56prOS1hbS9nS6/3dI2pQeF9m76tBUA4qkUcArwBJgBrBC0oyMSW8Di/ssWwfsjIhpwM40P9LOA89ExG3AfGB1Wk+52/4GFkbELKANWJzOn5O7q9Ea4FDDfClt90VEW8NnFUrpehn4LCJuBWZRrbvsbRHRmdZXG3AH1Wk1PsrdJukW4ClgbkTMpPp6+eW5u2oTEU0zUZ3xd1vD/HpgfeamyUBHw3wnMDFdnwh0FrDePgEeKKkNGAvspTrXdRFdVGe22wksBLaW8v8EjgPX91lWQtc1wDHSwT0ltfXpeRD4qoQ2es/3PoHq0/FbU19R62yoU1O9Q6H3n9HjZFpWkpsina0yXd6YM0bSZGA28DUFtKVNSvuBLmBHRBTRlbwEPAtcaFhWQlsA2yXtkfR4QV1TgTPAW2kz4RuSWgppa7Qc2JSuZ22LiJ+BF4ATwC9UZ7fdnrurLs02oOgiy3zccz8kXQV8CKyNiLO5ewAi4p+oNkO0AvMkzczdBCDpIaArIvbkbrmIBRExh2pT72pJ9+QOSkYDc4BXI2I20E1hm2rSKXCXAR/kbgFI+0YeBqYANwMtklbmrapPsw0oJ4FJDfOtwKlMLf05LWkiQLrsyhEh6TKqweS9iNhSUhtARPwOfEm1D6qErgXAMknHgfeBhZLeLaEtIk6lyy6q/QDzSuiiejyeTO8yATZTDTAltPVYAuyNiNNpPnfb/cCxiDgTEeeALcBdBXTVotkGlF3ANElT0iuP5UB75qa+2oFV6foqqv0XI0qSgDeBQxHxYiltkm6QNC5dH0P14DqcuwsgItZHRGtETKa6X30eEStzt0lqkXR1z3Wq7e0dubsAIuJX4CdJ09OiRcD3JbQ1WEHv5i7I33YCmC9pbHqcLqI6kCF3Vz1y78S51AlYCvwA/Ag8l7llE9V20HNUr9YeA66j2rF7JF1OyNB1N9WmwG+B/WlamrsNuB3Yl7o6gOfT8uzrrE/nvfTulM+9zqYCB9J0sOc+n7uroa8N2J3+px8D4wtqGwv8BlzbsCx7G7CB6oVUB/AOcEUJXXVM/uoVMzOrRbNt8jIzs0J5QDEzs1p4QDEzs1p4QDEzs1p4QDEzs1p4QDEzs1p4QDEzs1r8C4gATjSuA3m2AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "timeseries = model.results \n",
    "timeseries[\"daily_infected\"] = timeseries[\"total_infected\"].diff(1)\n",
    "timeseries['social_distance'] = timeseries.apply(lambda X: 1 if ( X.time >40) & ( X.time <= 50 ) else 0, axis=1)\n",
    "\n",
    "plt.figure\n",
    "timeseries.lockdown.plot( x = \"time\", kind = \"line\", label = \"lockdown\", legend = True, style = 'g' )\n",
    "timeseries.social_distance.plot( x = \"time\", kind = \"line\", label = \"social distance\", legend = True , style = 'r' )\n",
    "timeseries.daily_infected.plot( x = \"time\", legend = True, secondary_y = True )\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide_input": true
   },
   "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}