{"cells": [{"cell_type": "markdown", "metadata": {}, "source": ["# 1A - Enonc\u00e9 24 novembre 2020\n", "\n", "Correction de l'examen du 24 novembre 2020."]}, {"cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [{"data": {"text/html": ["
\n", ""], "text/plain": [""]}, "execution_count": 2, "metadata": {}, "output_type": "execute_result"}], "source": ["from jyquickhelper import add_notebook_menu\n", "add_notebook_menu()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["## Exercice 1 : gu\u00e9rison\n", "\n", "On commence par g\u00e9n\u00e9rer des donn\u00e9es artificielles \u00e0 partir de v\u00e9ritables donn\u00e9es."]}, {"cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [{"data": {"text/html": ["\n", "\n", "
\n", " \n", " \n", " \n", " rad \n", " dc \n", " \n", " \n", " jour \n", " \n", " \n", " \n", " \n", " \n", " \n", " 2020-03-18 \n", " 1627 \n", " 435 \n", " \n", " \n", " 2020-03-19 \n", " 2322 \n", " 642 \n", " \n", " \n", " 2020-03-20 \n", " 3128 \n", " 890 \n", " \n", " \n", " 2020-03-21 \n", " 3580 \n", " 1041 \n", " \n", " \n", " 2020-03-22 \n", " 4188 \n", " 1251 \n", " \n", " \n", "
\n", "
"], "text/plain": [" rad dc\n", "jour \n", "2020-03-18 1627 435\n", "2020-03-19 2322 642\n", "2020-03-20 3128 890\n", "2020-03-21 3580 1041\n", "2020-03-22 4188 1251"]}, "execution_count": 3, "metadata": {}, "output_type": "execute_result"}], "source": ["import pandas \n", "df = pandas.read_csv(\"https://www.data.gouv.fr/en/datasets/r/63352e38-d353-4b54-bfd1-f1b3ee1cabd7\", sep=\";\")\n", "gr = df[[\"jour\", \"rad\", \"dc\"]].groupby([\"jour\"]).sum()\n", "gr.head()"]}, {"cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [{"data": {"text/html": ["\n", "\n", "
\n", " \n", " \n", " \n", " jour \n", " rad \n", " dc \n", " \n", " \n", " \n", " \n", " 0 \n", " 2020-03-18 \n", " NaN \n", " NaN \n", " \n", " \n", " 1 \n", " 2020-03-19 \n", " 695.0 \n", " 207.0 \n", " \n", " \n", " 2 \n", " 2020-03-20 \n", " 806.0 \n", " 248.0 \n", " \n", " \n", " 3 \n", " 2020-03-21 \n", " 452.0 \n", " 151.0 \n", " \n", " \n", " 4 \n", " 2020-03-22 \n", " 608.0 \n", " 210.0 \n", " \n", " \n", "
\n", "
"], "text/plain": [" jour rad dc\n", "0 2020-03-18 NaN NaN\n", "1 2020-03-19 695.0 207.0\n", "2 2020-03-20 806.0 248.0\n", "3 2020-03-21 452.0 151.0\n", "4 2020-03-22 608.0 210.0"]}, "execution_count": 4, "metadata": {}, "output_type": "execute_result"}], "source": ["diff = gr.diff().reset_index(drop=False)\n", "diff.head()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["On convertit la date en jour de l'ann\u00e9e puis on simule un loi exponentielle de param\u00e8tre 14 pour avoir la date de sortie."]}, {"cell_type": "code", "execution_count": 4, "metadata": {"scrolled": false}, "outputs": [{"data": {"text/html": ["\n", "\n", "
\n", " \n", " \n", " \n", " entree \n", " sortie \n", " issue \n", " \n", " \n", " \n", " \n", " 0 \n", " 39 \n", " 79 \n", " 1 \n", " \n", " \n", " 1 \n", " 71 \n", " 79 \n", " 1 \n", " \n", " \n", " 2 \n", " 58 \n", " 79 \n", " 1 \n", " \n", " \n", " 3 \n", " 79 \n", " 79 \n", " 1 \n", " \n", " \n", " 4 \n", " 53 \n", " 79 \n", " 1 \n", " \n", " \n", "
\n", "
"], "text/plain": [" entree sortie issue\n", "0 39 79 1\n", "1 71 79 1\n", "2 58 79 1\n", "3 79 79 1\n", "4 53 79 1"]}, "execution_count": 5, "metadata": {}, "output_type": "execute_result"}], "source": ["import numpy.random as rnd\n", "\n", "\n", "def donnees_artificielles(hosp, mu=14, nu=21):\n", " dt = pandas.to_datetime(hosp['jour'])\n", " res = []\n", " for i in range(hosp.shape[0]):\n", " date = dt[i].dayofyear\n", " h = hosp.iloc[i, 1]\n", " delay = rnd.exponential(mu, int(h))\n", " for j in range(delay.shape[0]):\n", " res.append([date - int(delay[j]), date, 1])\n", " h = hosp.iloc[i, 2]\n", " delay = rnd.exponential(nu, int(h))\n", " for j in range(delay.shape[0]):\n", " res.append([date - int(delay[j]), date , 0])\n", " return pandas.DataFrame(res, columns=[\"entree\", \"sortie\", \"issue\"])\n", " \n", " \n", "data = donnees_artificielles(diff[1:].reset_index(drop=True))\n", "data.head()"]}, {"cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": ["data.to_csv(\"examen2021.csv\", index=False)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q1\n", "\n", "On r\u00e9cup\u00e8re les donn\u00e9es."]}, {"cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [{"data": {"text/html": ["\n", "\n", "
\n", " \n", " \n", " \n", " entree \n", " sortie \n", " issue \n", " \n", " \n", " \n", " \n", " 0 \n", " 49 \n", " 79 \n", " 1 \n", " \n", " \n", " 1 \n", " 27 \n", " 79 \n", " 1 \n", " \n", " \n", " 2 \n", " 73 \n", " 79 \n", " 1 \n", " \n", " \n", " 3 \n", " 74 \n", " 79 \n", " 1 \n", " \n", " \n", " 4 \n", " 48 \n", " 79 \n", " 1 \n", " \n", " \n", "
\n", "
"], "text/plain": [" entree sortie issue\n", "0 49 79 1\n", "1 27 79 1\n", "2 73 79 1\n", "3 74 79 1\n", "4 48 79 1"]}, "execution_count": 7, "metadata": {}, "output_type": "execute_result"}], "source": ["import pandas\n", "df = pandas.read_csv(\"http://www.xavierdupre.fr/enseignement/complements/examen2021.zip\")\n", "df.head()"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q2 : dur\u00e9e de gu\u00e9rison"]}, {"cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [{"data": {"text/plain": ["(array([30, 52, 6, 5, 31], dtype=int64), array([1, 1, 1, 1, 1], dtype=int64))"]}, "execution_count": 8, "metadata": {}, "output_type": "execute_result"}], "source": ["duree = df['sortie'] - df['entree']\n", "duree = duree.values # conversion en numpy\n", "issue = df['issue'].values\n", "duree[:5], issue[:5]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q3 : estimateur Kaplan-Meier (1)"]}, {"cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [{"data": {"text/plain": ["0.98965342710248"]}, "execution_count": 9, "metadata": {}, "output_type": "execute_result"}], "source": ["t = 10\n", "nt = duree[(duree >= t)].shape[0]\n", "dt = duree[(duree == t) & (issue == 0)].shape[0]\n", "st = 1. - dt / nt\n", "st"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q4 : courbe de Kaplan-Meier"]}, {"cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": ["T = [0]\n", "St = [1.]\n", "for t in range(0, 150):\n", " nt = duree[(duree >= t)].shape[0]\n", " dt = duree[(duree == t) & (issue == 0)].shape[0]\n", " st = 1. - dt / nt\n", " T.append(t)\n", " St.append(st * St[-1])"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q5 : graphe"]}, {"cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAQUAAAD4CAYAAADl7fPiAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAd1ElEQVR4nO3deVxVdf7H8deHXWRTAUERwV0QF8SldVpmcslcWzSnvRznl+W0TfazX/tYzkyppWXW1NQ0aVZWtqlli5VWouKCK65gqCgiKioi398f9zoDBHLFe++5y+f5ePAQzj3c++ZcfHPuud/zPWKMQSmlTguwOoBSyrNoKSilqtFSUEpVo6WglKpGS0EpVU2QVQ8cGxtrUlJSrHp4pfzWihUr9htj4uq63bJSSElJITs726qHV8pvicjOM92uLx+UUtVoKSilqtFSUEpVo6WglKpGS0EpVU29pSAir4nIPhFZV8ftIiLPi0ieiKwRkUznx1RKuYsjewr/BPqf4fYBQHv7xxjgpXOPpZSySr2lYIxZAhSfYZUhwJvG5kcgRkQSnRFudX4Jzy3axJETFc64O6WUA5xxTKElkF/l6wL7sl8RkTEiki0i2UVFRfXe8drdh3j+qzzKyrUUlHIXZ5SC1LKs1plbjDGzjDFZxpisuLg6R1n+yomTlQ3NppQ6S84ohQKgVZWvk4BfnHC/9GzdhOBA4eEP13GqUmeIUsodnFEK84Eb7e9C9AUOGWMKnXC/dE6M4rHB6Xy7uYgpX2x2xl0qpepR7wlRIjIbuASIFZEC4FEgGMAYMxP4DBgI5AFlwC3ODHh972RW55cw/es8MpKi6Zee4My7V0rVUG8pGGNG1XO7Ae50WqIaRIQnhnRh057D3Dd3NW3vjKBdfISrHk4pv+cVIxrDggN56fc9CQ0K4A//yubw8ZNWR1LKZ3lFKQC0iGnE9Osz2XGgjPvfXU2lHnhUyiW8phQAzmvbjIcGdGJh7l5e+nar1XGU8kleVQoAt12YyqCuiTy7aBPLth6wOo5SPsfrSkFEeGZEV1JiG3PX7FXsKz1udSSlfIrXlQJARGgQL43uyZETJ7lr9ioqTumIR6WcxStLAaBjQiSThmXw0/ZintOBTUo5jdeWAsDwzCRG9W7Fi99sZfGGvVbHUconeHUpADx6VTppiVHc/+5q9urxBaXOmdeXQlhwIC9c34PjJyu5d26Ojl9Q6hx5fSkAtI2L4LHBafyQd4BZ322zOo5SXs0nSgHg2qxWDMxI4O8LN7E6v8TqOEp5LZ8pBRHh6WFdiY8MZfycVTqFm1IN5DOlABAdHszUkT3YVVzGox/lWh1HKa/kU6UA0Du1KeMubcf7Kwv4KGe31XGU8jo+VwoAd1/enszkGB7+YB35xWVWx1HKq/hkKQQFBjBtZA8Axs/RYdBKnQ2fLAWAVk3D+cvwDFbuKuH5xVusjqOU1/DZUgAY3K0FIzKTmP51Hj9t09OslXKET5cCwOND0kluGs497+RwqEyncVOqPj5fChGhQUwb2YN9h0/w0AdrsM0zq5Sqi8+XAkC3VjHc368jn63dwzvL8+v/BqX8mF+UAsCYi9pwQbtmPP7xevL2HbE6jlIey29KISBAeO7a7oQFB3D37FWcqDhldSSlPJLflAJA86gw/nZ1N9YXlvK3BZusjqOUR/KrUgD4bVpzbujbmle/387SrfutjqOUx/G7UgD434GdSY1tzAPvrtGrTSlVg1+WQqOQQJ69thuFh47xxMfrrY6jlEfxy1IAyExuwh8vacu7Kwr4Yr1O+qrUaX5bCgDjL+9A58QoHpq3hgNHTlgdRymP4NelEBIUwJTrulF6rIKJH6zT0Y5K4WApiEh/EdkkInkiMqGW26NF5GMRWS0iuSJyi/OjukanhCju+V0HFuTu4UOdlEWp+ktBRAKBGcAAIA0YJSJpNVa7E1hvjOkGXAI8KyIhTs7qMmMubkPP1k145KNcCg8dszqOUpZyZE+hN5BnjNlmjCkH5gBDaqxjgEgRESACKAa8ZubUwADh2Wu6UXHK8Of39KQp5d8cKYWWQNWziArsy6qaDnQGfgHWAuONMV413VFKbGMmXtmZ77bs560fd1odRynLOFIKUsuymn9K+wE5QAugOzBdRKJ+dUciY0QkW0Syi4qKzjKq643uk8zFHeKY9NlGtu8/anUcpSzhSCkUAK2qfJ2EbY+gqluAecYmD9gOdKp5R8aYWcaYLGNMVlxcXEMzu4yI8NcRXQkOFO6bm8MpvQSd8kOOlMJyoL2IpNoPHo4E5tdYZxdwOYCINAc6Al55/baE6DCeHNqFlbtKeHnJVqvjKOV29ZaCMaYCGAcsBDYAc40xuSIyVkTG2ld7EjhfRNYCi4EHjTFee7bR4G4tGJiRwJQvNrP+l1Kr4yjlVmLVkfasrCyTnZ1tyWM7ovhoOVdMWUJsRAgfjbuA0KBAqyMp5RQissIYk1XX7X49ovFMmjYOYfKIDDbuOczUL3WKeOU/tBTO4PLOzbk2K4mXv92qV7JWfkNLoR4PD0ojPjKMB95brVO4Kb+gpVCPqLBgnh6ewea9R5j+VZ7VcZRyOS0FB1zaKZ4RmUm8+M1W1u0+ZHUcpVxKS8FBjwxKo1njEP783hpO6gVrlQ/TUnBQdHgwTw7twvrCUmYt8cpxWUo5REvhLPRLT+DKjESmLd6iF5RRPktL4Sw9NjidRsGBTHh/DZV6boTyQVoKZykuMpRHBqWRvfMgb/2kp1gr36Ol0ADDM1tycYc4Jn++kYKDZVbHUcqptBQaQESYNKwLBnho3lqdqUn5FC2FBkpqEs5DA20zNc3+WS9vr3yHlsI5GN07mQvaNeMvn64nv1hfRijfoKVwDgIChMkjugLwoL4boXyElsI5SmoSzsOD0li69QD/1ncjlA/QUnCCkb1a/WfC110H9GWE8m5aCk4gIkwekUFQoHD/e6v1ZYTyaloKTpIY3YhHBqXx8/Zi3li2w+o4SjWYloITXd0zics6xTN5gV43QnkvLQUnEhGeHp5BSGAAD7y7Wq8bobySloKTNY8K4/Eh6WTvPMjrP2y3Oo5SZ01LwQWGdm/JZZ3ieXbRZh3UpLyOloILiAhPDu1CgMDED9fpuRHKq2gpuEjLmEY80K8jSzYX8VFOzUtvKuW5tBRc6IbzUuiRHMMTn6yn+Gi51XGUcoiWggsFBgjPDO/K4eMneeqT9VbHUcohWgou1jEhkrG/acu8VbtZsrnI6jhK1UtLwQ3uvLQdbeIaM/HDtZSVV1gdR6kz0lJwg7DgQJ4Z3pX84mNM+WKz1XGUOiMtBTfpndqU6/sk84/vt7O2QK8ypTyXloIbTRjQidiIUB58X68ypTyXQ6UgIv1FZJOI5InIhDrWuUREckQkV0S+dW5M3xAVFswTQ9JZX1jKP77XIdDKM9VbCiISCMwABgBpwCgRSauxTgzwIjDYGJMOXOP8qL6hf5dE+qU3Z8oXm9mhZ1IqD+TInkJvIM8Ys80YUw7MAYbUWOd6YJ4xZheAMWafc2P6lieGdCEkMICJH+r08MrzOFIKLYGqc5gX2JdV1QFoIiLfiMgKEbmxtjsSkTEiki0i2UVF/vueffOoMCYM7MQPeQd4b0WB1XGUqsaRUpBaltX88xYE9ASuBPoB/yciHX71TcbMMsZkGWOy4uLizjqsLxnVK5neKU156tMNFB0+YXUcpf7DkVIoAFpV+ToJqHmGTwGwwBhz1BizH1gCdHNORN8UECBMGp7BsfJTPKFDoJUHcaQUlgPtRSRVREKAkcD8Gut8BFwkIkEiEg70ATY4N6rvaRcfwbjL2vHx6l/4auNeq+MoBThQCsaYCmAcsBDbf/S5xphcERkrImPt62wAFgBrgJ+BV40x61wX23eM/U1bOjSP4OEP1nHkhA6BVtYTq45+Z2VlmezsbEse29Os2HmQq2cu5abzUnhscLrVcZSPE5EVxpisum7XEY0eoGfrJtzYtzVvLNvByl0HrY6j/JyWgod4oH8nEqLCeOj9tZRX6BBoZR0tBQ8RERrEU0O7sGnvYWYt2Wp1HOXHtBQ8yOWdm3Nl10SeX5zH1qIjVsdRfkpLwcM8dlU6jUICeWjeWr0mpbKEloKHiYsMZeLAzvy8vZg5y/Pr/walnExLwQNdk5XEeW2a8fTnG9hXetzqOMrPaCl4IBHbEOjyikoenZ9rdRzlZ7QUPFRqbGPG/7Y9n6/bw8LcPVbHUX5ES8GD3XFRGzonRvHIR+soPX7S6jjKT2gpeLDgwAAmj8ig6PAJJn2q55cp99BS8HBdk2IYc3Fb5izP55tNOqGVcj0tBS/wp9+2p318BBPeX8uhY/oyQrmWloIXCAsO5Llru1N05ASPf6zvRijX0lLwEhlJ0dx5SVvmrdzNx6v10vbKdbQUvMhdl7cnMzmGh+atZZueG6FcREvBiwQHBjD9+kyCAoV73smhQq8ypVxAS8HLtIhpxFNDu7C64BAvL9lmdRzlg7QUvNCgri24smsiU7/czIbCUqvjKB+jpeClnhzShehGwdw3d7XO1KScSkvBSzVtHMKkYRmsLyxl+ldbrI6jfIiWghe7Ij2B4T1aMuObrawpKLE6jvIRWgpe7tGr0omLCOWu2at0tKNyCi0FLxcdHsz063vwS8kx7nknR69irc6ZloIPyEppysSBnflq4z69irU6Z1oKPuLG81LoldJEr2KtzpmWgo8ICBCeHp7B8ZOnGPf2Sk7qaEfVQFoKPqRdfCTPjMjgp+3FTPpMJ2VRDaOl4GOG9Uji5vNTeP2HHXy/Zb/VcZQX0lLwQRMGdKJNbGP+/N5qSsrKrY6jvIyWgg8KCw7kuetsk7KMn5OjV5pSZ8WhUhCR/iKySUTyRGTCGdbrJSKnRORq50VUDdG9VQyPXpXOt5uLmPF1ntVxlBeptxREJBCYAQwA0oBRIpJWx3qTgYXODqkaZnSfZK7q1oKpi7eQk19idRzlJRzZU+gN5BljthljyoE5wJBa1rsLeB/QKYc9hIjw1NAuNI8M5Z53cjh6osLqSMoLOFIKLYGqVzotsC/7DxFpCQwDZp7pjkRkjIhki0h2UVHR2WZVDRDdKJjnruvOjgNHeerT9VbHUV7AkVKQWpbVPHI1FXjQGHPqTHdkjJlljMkyxmTFxcU5GFGdq75tmvGHi9sy++d8FqzTS9CpM3OkFAqAVlW+TgJqTiecBcwRkR3A1cCLIjLUGQGVc9z7uw5ktIxmwrw17DmkV7JWdXOkFJYD7UUkVURCgJHA/KorGGNSjTEpxpgU4D3gf4wxHzo7rGq4kKAApo3szomTldzzTg6n9G1KVYd6S8EYUwGMw/auwgZgrjEmV0TGishYVwdUztMmLoLHB6ezbNsBZumkr6oOQY6sZIz5DPisxrJaDyoaY24+91jKVa7JSuLbzUU8u2gTF7RrRtekGKsjKQ+jIxr9jIgwaVgG8ZGh3D17FUf0bUpVg5aCH4oOD2bKdd3ZVVzGfXN1GLSqTkvBT/Vp04yJV6axMHcvL3ylw6DVf2kp+LFbL0hhRGYSU77czKJcHb+gbLQU/JiI8JdhXeiWFM097+Swee9hqyMpD6Cl4OfCggN5+YYswkODGPNmNofKdJp4f6eloEiIDmPm7zPZXXKMcbNX6tWs/ZyWggKgZ+umPDmkC99t2c9fF26yOo6ykEODl5R/GNk7mfWFpcxaso20xCiG9mhZ/zcpn6OloKr5v0FpbNpzmD+/v4aoRkFc1qm51ZGUm+nLB1VNcGAAM3/fk47NIxnz5gp9q9IPaSmoX2nSOIS37+hDesto7p6zilW7DlodSbmRloKqVWRYMP+4KYv4yDBufyObnQeOWh1JuYmWgqpTbEQor9/Si1PGcMvryzl4VK8h4Q+0FNQZtY2L4JUbsygoOcYdb2Zz/OQZZ9xTPkBLQdWrV0pTplzbneydB3n4w3VWx1EupqWgHHJl10Tuvrw9760o4N3s/Pq/QXktLQXlsPGXt+e8Ns2Y+OE6VuwstjqOchEtBeWwwABhxuhMWsY04vY3stmxX9+R8EVaCuqsNG0cwus390JEuPn1nynWdyR8jpaCOmspsY155cYsCg8d13ckfJCWgmqQnq2bMPW67qzcdZB75+bo6dY+REtBNdiAjEQevjKNz9buYfw7OZzUYvAJepakOie3XZhKZaXhL59toOJUJS+MyiQkSP/WeDN99tQ5u+PiNjwyyDYz9Jh/ZVNWrteS8GZaCsopbr0wlUnDMliyuYhRr/zEoWM616O30lJQTnN9n2Rm/r4n6385xC2v/8xRvfqUV9JSUE51RXoCz4/sQU5+CaNf/UnPrPRCWgrK6QZkJPLi6J6s/6WUa19exp5Dx62OpM6CloJyif5dEvjnrb0oPHScES8tZdeBMqsjKQdpKSiXOb9tLLPv6MvR8gqum7VMz5XwEg6Vgoj0F5FNIpInIhNquX20iKyxfywVkW7Oj6q8UUZSNG/f3pcTFZVcN2sZ24qOWB1J1aPeUhCRQGAGMABIA0aJSFqN1bYDvzHGdAWeBGY5O6jyXmktoph9R18qThmufXkZOfklVkdSZ+DInkJvIM8Ys80YUw7MAYZUXcEYs9QYc3rK3x+BJOfGVN6uY0Ikc8eeR6OQQEbOWsZCnTreYzlSCi2BqlPtFNiX1eU24PPabhCRMSKSLSLZRUVFjqdUPqFtXAQf/M8FdEqIYuxbK3jt++1WR1K1cKQUpJZlptYVRS7FVgoP1na7MWaWMSbLGJMVFxfneErlM2IjQpl9R1+uSGvOE5+s57H5uZyqrPXXSVnEkVIoAFpV+ToJ+KXmSiLSFXgVGGKMOeCceMoXNQoJ5MXRPbn9wlT+uXQHY99awREd/egxHCmF5UB7EUkVkRBgJDC/6goikgzMA24wxmx2fkzlawIDhIcHpfHYVWks3rCX/lOXsHyHzvvoCeotBWNMBTAOWAhsAOYaY3JFZKyIjLWv9gjQDHhRRHJEJNtliZVPufmCVOb+4TwCRBg560de/W4bxujLCSuJVU9AVlaWyc7W7lA2pcdPcv/c1Sxav5crMxKZfHVXIkJ1ug9XEJEVxpisum7XEY3KI0SFBfPyDT15sH8nPl9XyNAZP7Bl72GrY/klLQXlMUSEP17Slrdu68PBo+UMfP47nvl8IycqdGJYd9JSUB7n/HaxLPjTxQzp3pKZ325lxEtL2arDo91GS0F5pLjIUP5+TTdm3dCTXQfKGDD1O6Z9uYVKHdPgcloKyqNdkZ7A4vsuoX+XBKZ8uZnb38ym8NAxq2P5NC0F5fHiIkOZNrI7Tw5J5/st+7n0798w4+s8HQnpIloKyiuICDecl8Li+37DpR3j+dvCTVwzcykbCkutjuZztBSUV2nVNJwXR2cy5bpubN9/lEEvfM/Mb7fqgCcn0lJQXkdEGNYjia/vv4T+6Qk88/lGbnp9Odt1Zien0FJQXismPITp1/fgiSHprNx5kAHTlvDJml+dq6fOkpaC8moiwo32Yw1dWkQz7u1VjJ+zSueDPAdaCsonNI8K49939OHOS9uyMHcPlz77Dbe/sVxnkW4ALQXlM0KDAnmgXyeWPHApd13Wnp+2FdN/2hLeXLZDBz2dBS0F5XPio8K493cdWHjPxWSlNOWRj3K5/tUfda/BQVoKyme1iGnEG7f0YvKIDHJ3l+peg4O0FJRPExGu65Vcba/hiqlLeOvHnZSV6xRwtdFSUH7h9F7DtJHdCQsO4OEP19F30mLeXLZDh0vXoDMvKb9jjGHlroNM/XIL323ZT4fmEfzxkrZc1bUFQYG+/3eyvpmXtBSU3zLG8MmaQl74agub9x6hZUwjBmYk0L9LIj1bN7E6nstoKShVj8pKw+KN+3hj6Q5+3l5M+alKeiTHMOaiNlyRnkBgQG2XPvFeWgpKnYWjJyp4b0UBr36/jfziYyQ3DWd0n2QGdEkkuVm41fGcQktBqQY4VWlYlLuHV77bxspdJQCkJUbRv0sC/bsk0D4+AhHv3IPQUlDqHO06UMbC3D0syN3Dip226yi3iW3MFem2guiWFO1VBaGloJQT7Ss9zqL1e1mYu4dlWw9QUWlIjA6jX3oCI3u3olNClNUR66WloJSLlJSVs3jDPhbk7mHJ5iJOVFTSOTGKzOQYMpOb0Du1KUlNGnncXoSWglJuUFJWzjvL8/k+bz85u0o4bL9gbovoMDJbNyEhKoy0FlFc2C6W+KgwS7NqKSjlZqcqDVv2Hebn7cX8tK2YNbtLKDp8guMnKwkQuKBdLH3bNKN1s3ASoxvRvnkEUWHBbsunpaCUB6isNGzcc5gF6wr5dG0hW4v+OwmMCPRoFcP5bWNp3zyC+MgwurWKJjzENdfS1FJQygOVHj9JYclxCg6WsabgEF9t3Mf6wtL/nIcRFhxAt6QY2sRF0DUpmu6tYmgXH0GwE4Zhayko5SVOVJwiv7iMgoPH+GZTEWt3H2LL3sOUHv/v2ZyRYUHEhAeTGNWIVk3DadW0EclNw8lMbkJKbGOHHqe+UtBrfSvlIUKDAmkXH0m7+Egu6RgP2M7P2HGgjJz8g+w8UEZJ2UkOlpVTWHKcpVv3s2fVcU7/Xe+aFM2bt/YmJjzknHJoKSjlwUSE1NjGpNaxF3B67+LrjUWsyj9IdKNzP2DpUCmISH9gGhAIvGqMeabG7WK/fSBQBtxsjFl5zumUUmdUde/CWeo9aiEigcAMYACQBowSkbQaqw0A2ts/xgAvOS2hUsqtHDmU2RvIM8ZsM8aUA3OAITXWGQK8aWx+BGJEJNHJWZVSbuBIKbQE8qt8XWBfdrbrICJjRCRbRLKLiorONqtSyg0cKYXaBm7XfB/TkXUwxswyxmQZY7Li4uIcyaeUcjNHSqEAaFXl6ySg5gX7HFlHKeUFHCmF5UB7EUkVkRBgJDC/xjrzgRvFpi9wyBhT6OSsSik3qPctSWNMhYiMAxZie0vyNWNMroiMtd8+E/gM29uRedjekrzFdZGVUq7k0DgFY8xn2P7jV102s8rnBrjTudGUUlaw7NwHESkCdjqwaiyw38VxHKE5qvOUHOA5WbwlR2tjTJ1H+i0rBUeJSPaZTt7QHP6dAzwni6/k8P3L4SilzoqWglKqGm8ohVlWB7DTHNV5Sg7wnCw+kcPjjykopdzLG/YUlFJupKWglKrGo0tBRPqLyCYRyRORCW583FYi8rWIbBCRXBEZb1/+mIjsFpEc+8dAN2TZISJr7Y+XbV/WVES+EJEt9n9det10EelY5WfOEZFSEfmTO7aHiLwmIvtEZF2VZXX+/CLykP33ZZOI9HNxjr+JyEYRWSMiH4hIjH15iogcq7JdZtZ5x87JUefz0KDtYYzxyA9sQ6q3Am2AEGA1kOamx04EMu2fRwKbsU0w8xhwv5u3ww4gtsayvwIT7J9PACa7+XnZA7R2x/YALgYygXX1/fz252g1EAqk2n9/Al2Y4wogyP755Co5Uqqu54btUevz0NDt4cl7Co5M7uISxphCY59OzhhzGNhALfNDWGgI8Ib98zeAoW587MuBrcYYR0ajnjNjzBKguMbiun7+IcAcY8wJY8x2bOfi9HZVDmPMImPM6amWf8R2drBL1bE96tKg7eHJpeDQxC2uJiIpQA/gJ/uicfbdxddcvdtuZ4BFIrJCRMbYlzU39rNQ7f/GuyHHaSOB2VW+dvf2gLp/fit/Z24FPq/ydaqIrBKRb0XkIjc8fm3PQ4O2hyeXgkMTt7g0gEgE8D7wJ2NMKba5J9sC3YFC4Fk3xLjAGJOJbR7MO0XkYjc8Zq3sp84PBt61L7Jie5yJJb8zIjIRqAD+bV9UCCQbY3oA9wJvi4grL0dd1/PQoO3hyaVg6cQtIhKMrRD+bYyZB2CM2WuMOWWMqQRewUm7pmdijPnF/u8+4AP7Y+49PQem/d99rs5hNwBYaYzZa8/k9u1hV9fP7/bfGRG5CRgEjDb2F/L23fUD9s9XYHst38FVGc7wPDRoe3hyKTgyuYtL2Kes/wewwRjzXJXlVSejHQasq/m9Ts7RWEQiT3+O7cDWOmzb4Sb7ajcBH7kyRxWjqPLSwd3bo4q6fv75wEgRCRWRVGyzi//sqhBiu/TBg8BgY0xZleVxYpsFHRFpY8+xzYU56noeGrY9XHnk2AlHWgdiO/K/FZjoxse9ENtu1hogx/4xEPgXsNa+fD6Q6OIcbbAdPV4N5J7eBkAzYDGwxf5vUzdsk3DgABBdZZnLtwe2EioETmL7y3fbmX5+YKL992UTMMDFOfKwvWY//Tsy077uCPvztRpYCVzl4hx1Pg8N2R46zFkpVY0nv3xQSllAS0EpVY2WglKqGi0FpVQ1WgpKqWq0FJRS1WgpKKWq+X/5W90+Z0DC+AAAAABJRU5ErkJggg==\n", "text/plain": [""]}, "metadata": {"needs_background": "light"}, "output_type": "display_data"}], "source": ["import matplotlib.pyplot as plt\n", "fig, ax = plt.subplots(1, 1, figsize=(4, 4))\n", "ax.plot(T, St);"]}, {"cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlMAAAEGCAYAAABB6hAxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAA1F0lEQVR4nO3dd3iUVd7/8fc3nRRq6KGJSJdiRIqiPugKNtYOKsqKIrbdZ9vvwa3u6q667trWiiLKWhBRwLbiiihID0iPQKQGqQECJKSf3x9J2JDMTCZtJuXzui4vZuacuecbbgY+nvvc55hzDhERERGpnJBgFyAiIiJSlylMiYiIiFSBwpSIiIhIFShMiYiIiFSBwpSIiIhIFYQF64Pj4+Nd586dg/XxIiIiIn5btWrVIedcS09tQQtTnTt3JikpKVgfLyIiIuI3M9vprU2X+URERESqQGFKREREpAoUpkRERESqIGhzpqrbJ+v28saSHeX2G9W3DT8Z1qXmCxIREfFTbm4uqampZGVlBbuUBi8qKoqEhATCw8P9fk+9CVP7jmWxYsfhcvv1bt84ANWIiIj4LzU1lbi4ODp37oyZBbucBss5R1paGqmpqXTp4v/AS7mX+czsNTM7YGYbvLSbmT1rZilmts7MBlagbhERkQYvKyuLFi1aKEgFmZnRokWLCo8Q+jNn6nVgpI/2UUC3ov8mAi9WqAIRERFRkKolKnMeyr3M55xbaGadfXQZDUx3zjlgmZk1NbO2zrm9Fa4mAKYt3sGDo3oSEaa59yIiEnydJ39S9GhbQD93x2NXBPTz6rPqSBTtgd0lnqcWvVaGmU00syQzSzp48GA1fHTl3PjyUn44ejJony8iIlIf/fWvfw12CUFRHWHK03iY89TROTfFOZfonEts2dLjiuzV6g9h05kR8TAzIh5mbOj8U6+v2X2UK55dxMItwQt0IiIi9Y3CVOWlAh1KPE8AfqiG41ZIi5gIr229bCejQ5ec9tqRzFxun7aCZ77YSkGBx+wnIiLSoEyfPp2zzz6bfv36MW7cOMaPH8+sWbNOtcfGxgKwd+9ehg8fTv/+/enTpw+LFi1i8uTJnDx5kv79+3PLLbcA8OSTT9KnTx/69OnD008/DcCOHTvo0aMHd955J3369OGWW27hiy++YNiwYXTr1o0VK1YE/OeuquoIUx8CtxXd1TcYSA/GfKm+CU3KvPbnvNsYk/N7NrlO9LKdZUaonIOnvtjChDdWkn4yN5DlioiI1CobN27kL3/5C19++SVr167lmWee8dr37bff5rLLLmPNmjWsXbuW/v3789hjj9GoUSPWrFnDW2+9xapVq5g2bRrLly9n2bJlvPLKK3z77bcApKSk8LOf/Yx169bx3Xff8fbbb/PNN9/w97//vU6ObvmzNMI7wFKgu5mlmtkEM5tkZpOKunxK4ay5FOAV4N4aq9aHri1juax3a49tc/OHssl1YnBIMo+GTy0TqhZsPsiPn1/M1v3HA1WuiIhIrfLll19y/fXXEx8fD0Dz5s299j333HOZNm0aDz30EOvXrycuLq5Mn2+++YZrrrmGmJgYYmNjufbaa1m0aBEAXbp0oW/fvoSEhNC7d29GjBiBmdG3b1927NhRIz9fTSo3TDnnxjrn2jrnwp1zCc65qc65l5xzLxW1O+fcfc65rs65vs65pJov27O/XdeP3u3KLsr5Tv4IxuT8ngdzJ7CsoKfHy37bD2Xw4+cX89mGfYEqV0REpNZwzpVZFiAsLIyCgoJT7Tk5OQAMHz6chQsX0r59e8aNG8f06dM9Hs+byMjIU49DQkJOPQ8JCSEvL6/KP0ug1av1AZpEh/Pu3UO44uy2HtuLQ5W3y34ZOflMenMV//h8s+ZRiYhIgzJixAhmzpxJWloaAIcPH6Zz586sWrUKgLlz55KbWzglZufOnbRq1Yq77rqLCRMmsHr1agDCw8NP9Rk+fDhz5swhMzOTjIwMZs+ezQUXXBCEn6zm1ZvtZIrFRobx3NgBJHZqxl8+SSbPQyiamz8UQgsnphNaGLJK+ueXKWzYk87TYwbQpJH/e/OIiIjUVb179+a3v/0tF154IaGhoQwYMIDHH3+c0aNHM2jQIEaMGEFMTAwAX331FU888QTh4eHExsaeGpmaOHEiZ599NgMHDuStt95i/PjxDBo0CIA777yTAQMG1MnLeOUxX8NwNSkxMdElJdXsFcFVO49w31ur2XfM87LwMyIeppftZJPrxNz8oWVCVZf4GKaMO4durcteCxYREakO/120M7C0aKd3ycnJ9OzZ87TXzGyVcy7RU/96HaYA0k5k87MZa/gm5VCZtrGh8xkduoTBIckALCso/I0rGaxiIkL5x439GdmnTY3XKiIiDZOnf7wleCoapurVnClPWsRG8sYdg/jp/5xZpq30xHQouyZV8TyqZ+dv9TmZTkRERBqmeh+mAEJDjF/8qDuv3pZIXGTZaWLFocrX5PQn/7OFX8xcS3ZefiBLFxERkVquQYSpYpf0as2c+4fRtWWM1z7Fa1J5Wj5h9rd7uOWV5aSdyK7pUkVERKSOaFBhCgoX95xz3zAu6el5gc/ylk9I2nmEH7+gBT5FRESkUIMLUwBxUeFMGXcOP7/kLK99fI1Q7T58kmtfXMKirdooWUREpKFrkGEKICTE+Nkl3cqdR+VthOp4Vh7jp63kzWU7A1m2iIiI1DINNkwVK28ela8RqvwCx+/mbODPH20iXyumi4hIHRYbG3vq8aeffkq3bt3YtWsXDz30EGZGSkrKqfannnoKM6O6ljiaM2cOmzZtOvX8D3/4A1988UWVj3v06FFeeOGFKh+nPA0+TIHveVTljVABvLZ4OxOnJ3Eiu+7tJyQiIlLS/PnzeeCBB/jss8/o2LEjAH379mXGjBmn+syaNYtevXpV22eWDlN//vOfueSSS6p83ECFqXq3nUxlFc+jenzed7z89bYy7cVb0AwOSWZwSPKpUariBT7nf3eAG15aytTbE2nXtFGgyxcRkfri35Nh3/rqPWabvjDqsXK7LVq0iLvuuotPP/2Url27nnr9xz/+MXPnzuV3v/sd27Zto0mTJoSH+95u7fPPP+ePf/wj2dnZdO3alWnTphEbG8vkyZP58MMPCQsL40c/+hHXXnstH374IV9//TWPPPII77//Pg8//DBXXnkl119/PZ07d+bmm29mwYIF5ObmMmXKFB588EFSUlL49a9/zaRJkzhx4gSjR4/myJEj5Obm8sgjjzB69GgmT57M999/T//+/bn00kt54okneOKJJ5g5cybZ2dlcc801/OlPf6ryb6/CVAkhIcaDo3pyRnwMv5294bR9/d7JH8E7+SNOrZoOZff2S957jNHPL+bV2xLp16FpMH4EERGRSsnOzmb06NF89dVX9OjR47S2xo0b06FDBzZs2MDcuXO56aabmDZtmtdjHTp0iEceeYQvvviCmJgYHn/8cZ588knuv/9+Zs+ezXfffYeZcfToUZo2bcrVV199Kjx50qFDB5YuXcrPf/5zxo8fz+LFi8nKyqJ3795MmjSJqKgoZs+eTePGjTl06BCDBw/m6quv5rHHHmPDhg2sWbMGKAx4W7duZcWKFTjnuPrqq1m4cCHDhw+v0u+dwpQHN53bkQ7No7nnzdWkn8w9ra04VMF/9/abEfHwqRGqg8ezuWnKUp68sT+X920bjPJFRKQu82MEqSaEh4czdOhQpk6dyjPPPFOmfcyYMcyYMYN58+Yxf/58n2Fq2bJlbNq0iWHDhgGQk5PDkCFDaNy4MVFRUdx5551cccUVXHnllX7VdvXVVwOFlxtPnDhBXFwccXFxREVFcfToUWJiYvjNb37DwoULCQkJYc+ePezfv7/McT7//HM+//xzBgwYAMCJEyfYunVrlcOU5kx5MbRrPLPvHUrnFtFe+3ibnJ6VW8C9b63m+QUp2oJGRETqhJCQEGbOnMnKlSv561//Wqb9qquu4l//+hcdO3akcePGPo/lnOPSSy9lzZo1rFmzhk2bNjF16lTCwsJYsWIF1113HXPmzGHkyJF+1RYZGXmqxuLHxc/z8vJ46623OHjwIKtWrWLNmjW0bt2arKwsj3U9+OCDp+pKSUlhwoQJftXgi8KUD2e0jGX2vcM4r0tzj+3lTU5/Yt5mfvXeOnLyCgJVsoiISKVFR0fz8ccf89ZbbzF16tTT2ho1asTjjz/Ob3/723KPM3jwYBYvXnzqDsDMzEy2bNnCiRMnSE9P5/LLL+fpp58+dfktLi6O48crvxh2eno6rVq1Ijw8nAULFrBz506Px73ssst47bXXOHHiBAB79uzhwIEDlf7cYrrMV45mMRH8a8J5/Gb2ematSvXYp3hyeuk5VADvr05l95FMXr71HJrFRASqbBERkUpp3rw5n332GcOHDyc+Pv60tjFjxvh1jJYtW/L6668zduxYsrMLt2B75JFHiIuLY/To0WRlZeGc46mnnjp13Lvuuotnn32WWbNmVbjmW265hauuuorExET69+9/as5XixYtGDZsGH369GHUqFE88cQTJCcnM2TIEKBwOYg333yTVq1aVfgzS7JgXYZKTEx01bU+RSA453jx6+/522ebvfYpnkO1yXU6NYeqWKcW0bw2/ly6toz1+n4REWmYkpOT6dmzZ7DLkCKezoeZrXLOJXrqr8t8fjIz7r3oTF68ZSBR4Z5/23wt8LkzLZPrXlzCqp2HA1GuiIiIBIjCVAWN6tuWmXcPoWVcZJm28uZQHc3M5eZXljNv475AliwiIlJjzjvvPPr373/af+vXV/M6WbWc5kxVwtkJTZl73zAmvJFE8t5jZdp9zaHKzivgnjdX8efRfbh1cKdAli0iIrWYcw4zC3YZFbZ8+fJgl1CtKjP9SSNTldSuaSNmTRrCJT3LTlorb4SqwMHv5mzg7/M2a+kEEREhKiqKtLQ0/ZsQZM450tLSiIqKqtD7NDJVBTGRYbw8LpFHP03m1W+2l2n3NUIF8NyCFPYdy+LRa/sSHqpcKyLSUCUkJJCamsrBgweDXUqDFxUVRUJCQoXeo7v5qsnby3fx+7kbyC8o+/tZ8i4/oMydfsPPasmLtwwkJlLZVkREpDbS3XwBcPN5HXnjJ4OIiyobiIrv8gM83um3cMtBxkxZxsHj2QGpVURERKqPRqaqWcqB44yftpLUIyc9tvtai6pj82jeuGMQXeJjAlWuiIiI+EEjUwF0Zqs4PrhnKL3aet63qHiUanBIMo+GTz1tcvquw4VrUa3ZfTSAFYuIiEhVKEzVgFaNo3j37sGcf2Z8mbbiO/0ezJ3AsoKeZS77Hc7IYeyUZXz5XdndrkVERKT2UZiqIXFR4bw2/lyuGdDeY7uv5RNO5uZz1/RVzFixK5Ali4iISCUoTNWgiLAQ/nFDPyZd2NVrH29b0OQXOCZ/sJ4n/7NF646IiIjUYn6FKTMbaWabzSzFzCZ7aG9iZh+Z2Voz22hmP6n+UuumkBBj8qgePHRVLzwtbFveAp/Pzt/Kr2etIze/IIBVi4iIiL/KDVNmFgo8D4wCegFjzaxXqW73AZucc/2Ai4B/mFlENddap40f1oXnbx5IhJfFOX1tkjxrVSoT3kgiIzsvEKWKiIhIBfgzMjUISHHObXPO5QAzgNGl+jggzgo3FYoFDgP6l7+Uy/u2ZfoEz2tRlTdCtXDLQW5+ZRlpJ7QWlYiISG3iT5hqD+wu8Ty16LWSngN6Aj8A64GfOefKXJcys4lmlmRmSQ11yfzBZ7Rg1qShtG3ied8fXyNUa1PTueGlpew+nBmIUkVERMQP/oQpT1tYl54RfRmwBmgH9AeeM7MyCy0556Y45xKdc4ktW7asYKn1R/c2cbx/z1DOah1bps3TCFXJUapthzK47sUlJO89FuiyRURExAN/wlQq0KHE8wQKR6BK+gnwgSuUAmwHelRPifVTu6aNeG/SUAaf0dxju68taA4cz+bGl5eyYvvhgNQqIiIi3vkTplYC3cysS9Gk8jHAh6X67AJGAJhZa6A7sK06C62PmjQK5407BnFVv3Zl2opHqLzNozqelce4qcv5fOO+QJctIiIiJZQbppxzecD9wDwgGZjpnNtoZpPMbFJRt4eBoWa2HpgP/J9z7lBNFV2fRIaF8sxN/bn7wjO89vE2jyo7r4BJb67i3ZVa3FNERCRYtNFxLfLGkh089NFGvJ0SX5sk/+6Kntx5gfdAJiIiIpWnjY7riNuHdq70WlSPfJLMk59v1mrpIiIiAaYwVctc3rctr99xLrGRFV+L6tkvU/jTR5soKFCgEhERCRSFqVpoaNd4ZkwcTHxspMd2XyNUry/Zwa9mrSVP28+IiIgEhMJULdWnfRPev2cInVpEl2krb4Tqg9V7uPet1WTn5QeyZBERkQZJYaoW69QihlmThtK7XZn1TwHfI1Sfb9rPhNe1n5+IiEhN0918dcDxrFzu/tcqlnyf5rG95F1+wGl3+g3o2JTXxw+iSXR4wOoVERGpb3Q3Xx0XFxXOa+PP5dJerT22+1ot/dtdR7lpylIOHM8KSK0iIiINjUam6pC8/AL+36x1fPDtHq99vK1F1blFNG/eeR4JzcrOwRIRERHfNDJVT4SFhvD3G/px+5BOXvt4m0e1Iy2TG15aSsqBE4EoVUREpMFQmKpjQkKMh67uzf0Xn+mx3dedfnvTs7jx5aWsSz0awIpFRETqN4WpOsjM+NVl3XlwVA+vfbyNUB3OyGHslGUs+V5bJ4qIiFQHhak67O4Lu/LotX0xK9vma4QqIyef8dNWMm/jvgBXLCIiUv8oTNVxYwd15NkxAwgL8ZCo8D5ClZNXwD1vruK9pN2BKlVERKReUpiqB67q145XbkskMqzs6fQ1QlXg4Nez1vHqom2BLllERKTeUJiqJy7u0YrpdwzyuEEy+F4t/ZFPkvn7vM0Ea5kMERGRukxhqh4574wWvHPXYJp5WO28vP38nluQwu/mbCC/QIFKRESkIhSm6pm+CU14b9IQ2jSO8tjua4TqreW7+NmMb8nJKwhEqSIiIvWCVkCvp1KPZHLb1BVsO5Thsd3Xfn7Dz2rJS7cOJDrC8yVDERGRhkYroDdACc2imTlpCL3bNfbY7ms/v4VbDnLrq8tJz8wNSK0iIiJ1mUam6rljWbnc+UYSK7Yf9trH235+3VvH8a8Jg2jl5ZKhiIhIQ6GRqQascVQ40+8YxCU9W3nt420e1eb9x7nupSXsTPN8qVBEREQUphqEqPBQXrz1HK4d0N5ju687/XYfPsn1Ly0lee+xQJYsIiJSZyhMNRDhoSH8/YZ+/GRYZ699vI1QHTyezU0vL2XVTu+XCkVERBoqhakGJCTE+MOVvfjFpWd5bPc1QnUsK49xU1f4nHslIiLSEClMNTBmxk9HdOPh0b09bpAM3keoMnPyGT9NgUpERKQkhakGatyQzjx9U3+PGyT7GqFSoBIRETmdwlQDNrp/e165PZGocM9/DMoboVq5Q4FKREREYaqBu7h7K96ccB5xUWVXOy85QlVaZk4+t7+mQCUiIqIwJSR2bs7Mu4cQHxvptY+nzZEzc/IZr0AlIiINnMKUANCzbWPev2cIHZo3KtPma3PkjKJAlaRAJSIiDZRfYcrMRprZZjNLMbPJXvpcZGZrzGyjmX1dvWVKIHRqEcN7dw+lS3zMaa97mpBecpQqo+iSnwKViIg0ROWGKTMLBZ4HRgG9gLFm1qtUn6bAC8DVzrnewA3VX6oEQpsmUbxz12A6t4gu0+Zrc+TiQKW7/EREpKHxZ2RqEJDinNvmnMsBZgCjS/W5GfjAObcLwDl3oHrLlEBq0ySKGROHlAlUxSNU3pZNyMjJZ9zU5cxP3h+MskVERILCnzDVHthd4nlq0WslnQU0M7OvzGyVmd3m6UBmNtHMksws6eDBg5WrWALCW6Aq5m0eVXZeARP/tYrZ36YGqlQREZGg8idMeVon25V6HgacA1wBXAb83szK7FninJvinEt0ziW2bNmywsVKYLVpEsU7Ez1f8vO1sGd+geMXM9cqUImISIPgT5hKBTqUeJ4A/OChz2fOuQzn3CFgIdCvekqUYGrbpBHvTBxMpwqOUDkHv5y5lrlr9gSqVBERkaDwJ0ytBLqZWRcziwDGAB+W6jMXuMDMwswsGjgPSK7eUiVY2jZpxAwvgcrXCFWBg/99dw1vL98V6JJFREQCptww5ZzLA+4H5lEYkGY65zaa2SQzm1TUJxn4DFgHrABedc5tqLmyJdB8BSrwPUL1m9nreX5BCs6VvjosIiJS91mw/oFLTEx0SUlJQflsqby96ScZM2UZO9MyPbbPiHiYXraTTa4Tc/OH8k7+iFNtd13Qhd9c3hMzT9PwREREai8zW+WcS/TUphXQpULaNmnEe3cPoWfbxh7bfa2W/sqi7fx61jry8gsCUaqIiEhAaGRKKiX9ZC53vZHECi+rnpccoQJOG6Ua1acNT4/pT2RYaMDqFRERqQqNTEm1a9IonNfvOJdBXZp7bPe1Wvq/N+xj4vRVnMzJD0itIiIiNUkjU1IlGdl5jJ+2gpU7jnjt420e1aDOzZk6PpG4qPBAlSsiIlIpGpmSGhMTGcbrPxnEBd3ivfbxNo9qxY7D3PLqco5k5ASiVBERkRqhMCVVFhMZxqu3J3J53zYe20uuRVXautR0bpqylAPHsmq6TBERkRqhMCXVIjIslH+OHciYczv47Fd6YU+ALftPcOPLS0k94nm5BRERkdpMYUqqTWiI8ei1fZl0YVeP7b6WTdiRlskNLy1l28ETgShVRESk2ihMSbUyMyaP6sH/jexRps3X1jMAe9OzuPHlpSTvPRbIkkVERKpEYUpqxD0XdeXh0b09tvkaoTp0IocbX17Ksm1pgShTRESkyhSmpMaMG9KZf9zQj5BSu8eUN0J1PCuP26au4ON1PwS4YhERkYpTmJIadd05Cbxwy0DCQ8vux+drhConv4D73/6Wp/6zhYICbZAsIiK1l8KU1LiRfdry6u3nEhV++h+38kaoAJ6Zv5VJb2q1dBERqb0UpiQgLjyrJdPvOI/YyLAybb5GqAA+37SfW6cuJz0zNxClioiIVIi2k5GAWpd6lNtfW8ERD8HI1+bIAD3axPHGHYNo3TgqYPWKiIiAtpORWuTshKa8e/cQWsZFlmkruTny4JBkHg2fetqlv+/2HWf0c4vZ+EN6QGsWERHxRSNTEhQ70zIYN3UFuw57XvV8bOh8RocuOTVSNSbn96faoiNCef7mgVzco1WgyhURkQZOI1NS63RqEcP79wzl7IQmHtt9TU7PzMnnrulJfLJubyBLFhER8UhhSoKmZVwk79w1mEt6eh9h8jY5Pa/A8cA7q3l7+a5AlCoiIuKVwpQEVUxkGFPGJfLA/5zpsd3XCFWBg9/MXs//m7WWrFwtnSAiIsGhMCVBFxJi/PJH3Xn6pv6ElV4uvYiv5RNmJqVy05RlHM3MCUS5IiIip1GYklrjxwPa88ptiUSGlf1jWd4Cn2t3H2XMlGUcOpEdyJJFREQUpqR2ubhHK96+azAtYiI8tvsaofpu33FufGkpW/cfD0SpIiIigMKU1ELndGrGnPuG0a1VbJm28kaoth3K4OrnFjP729RAliwiIg2YwpTUSh2aRzNr0lD6d2jqsd3XCNXJ3Hx+/u5anvx8M8FaR01ERBoOLdoptdqJ7DzufGMly7Yd9the3hY0Ywd14OHRfQgL1f83iIhI5WnRTqmzYiPDeP0ng7j+nASP7SW3oPE0SvXOit3c9toKTUwXEZEao5EpqROcc7y7cje/n7uB3HzPf2ZLjlKVHqFq2ySKl8edw9kJTQNUsYiI1CcamZI6z8wYM6gjU7wsnQC+51HtTc9izJRlLNxyMBDliohIA6IwJXXKxd1bMW38uURHhJZpK+9Ov8ycfO54fSUfrNadfiIiUn0UpqTOGXpmPLPvHcYZLWM8thePUA0OSebR8Kmnhaq8AscvZhbe6VdQoDv9RESk6vwKU2Y20sw2m1mKmU320e9cM8s3s+urr0SRsrq3ieOj+8/nirPblmkrHqF6MHcCywp6erzs9+yXKTzwzreczNGefiIiUjXlhikzCwWeB0YBvYCxZtbLS7/HgXnVXaSIJzGRYfxzzADuPL+Lx/byLvt9sn4vY6Ys5cCxrECVLCIi9ZA/I1ODgBTn3DbnXA4wAxjtod8DwPvAgWqsT8SnkBDjd1f24ndX9MTLHsk+J6avTU1n9POL2fhDegCqFRGR+sifMNUe2F3ieWrRa6eYWXvgGuAlXwcys4lmlmRmSQcP6q4qqT53XnAG0+84j+Ye9vTzNEJVcpRqb3oWN7y0lP9s2h/oskVEpB7wJ0x5+v/90jN3nwb+zznncwKKc26Kcy7ROZfYsmVLP0sU8c/53eL56IHzPe7pB74X+MzMyWfiv5KY+s32gNQqIiL1hz9hKhXoUOJ5AvBDqT6JwAwz2wFcD7xgZj+ujgJFKqJ900a8N2kI53RqVqateITK2zwq5+Dhjzfx0IcbydedfiIi4id/wtRKoJuZdTGzCGAM8GHJDs65Ls65zs65zsAs4F7n3JzqLlbEH02jI3hzwnlc3reN1z6+5lG9vmQHt722nIPHtQWNiIiUr9ww5ZzLA+6n8C69ZGCmc26jmU0ys0k1XaBIZTSKCOW5sQO5/+IzPbaXd6ff4pQ0rnh2Ecu3pQWqZBERqaO0N5/Uex+sTmXy++vJyS8o0zY2dD6jQ5ec2tNvTM7vT2sPMfjVZd2ZNLwrId5uFxQRkXpPe/NJg3btwATevsv/O/1KjlAVOPjbZ5uZ8MZKjmTkBLJsERGpIxSmpEFI7NycOfcOK/dOP09b0AAs2HyQK55dxOpdRwJVsoiI1BEKU9JgdGwRzfv3DmX4WWWX5fBnC5of0rO48aWlTP1mO8G6PC4iIrWPwpQ0KI2jwpl6eyJjzu3gsb28y355BY6HP97EpDdXkX4yN1Bli4hILaYwJQ1OeGgIj17bl19f1h2rxBY0APM27ueqf37D+lRtQyMi0tDpbj5p0BZuOcj/vruGw14ml8+IePjUnX5QGLLeyR9xqj3E4NbBnfj5JWfRzMMEdxERqR90N5+IF8PPasmnP72AczuXXTEdTt+CxtPk9AIH05fu5H/+8RWfbdgbsLpFRKT20MiUCJCbX8DfP9/My19v89qnvDWpAG5MTOAPV/UmNjKsJssVEZEA8zUypTAlUsL85P38YuZan5PLS176K33ZD6Bj82ieuqm/x/0BRUSkbtJlPhE/jejZmk9+ej79OjT12qe8Nal2Hc7khpeW8OTnm8n1sOq6iIjULxqZEvEgJ6+AR/+dzLTFO7z28eeyX7+EJjx1U3/OaOl5sVAREakbNDIlUkERYSH88arevHjLQOK8zH8qb00qgLWp6Vz93GLmJ+8PRNkiIhIEClMiPozq25aPf3o+vds19tqnvMt+J7LzuHN6Es8vSKGgQCuni4jUN7rMJ+KH7Lx8XljwPS99/T3ZeZ7nQflz2W9o1xb848Z+tG3SqKZLFhGRaqS7+USqSeqRTB76cCNfJB/w2qe8hT6bNArnr9f05Yqz29Z4vSIiUj00Z0qkmiQ0i+aV2xJ57Nq+NAoP9din5EKfnrajST+Zy31vr+aXM9dyPEv7+4mI1HUamRKppO2HMvjfd9ewdvdRr33KW5OqQ/NGPHVjfxI7N6/hakVEpCo0MiVSA7rExzBr0hB+NqIbIeVsmOxtcvruwye58eWlPPThRp8LhYqISO2lkSmRarBq5xF+/u4adh3O9Njuz+T0+NgIHhzVk2sHtsfMSzoTEZGg0AR0kQBIP5nLT9/5lq+3HPTap7zLfgCJnZrxx6t60zehSU2WKyIiFaDLfCIB0KRROK+NP5e7LzzDa5/yLvsBJO08wlXPfcPPZnxL6hHPI10iIlJ7aGRKpAYs2nqQX85cy4Hj2R7b/bnsBxAdEcqvftSd24d2JtTbxCwREalxuswnEgRHMnJ48IP1fLZxn9c+5a1JVaxfQhP+cFUvzumku/5ERIJBYUokSJxzvJeUykMfbSQzJ79Me/EIFcDgkGQAlhX09Bqqrjy7Lf83sgcdmkfXbOEiInIahSmRINuZVrgm1be7jnrtUxysygtVEWEh3Hl+F+67+ExivGzCLCIi1UthSqQWyMsv4MWvvue5BSle9/cD/0NV2yZR/P7KXozq00ZLKYiI1DCFKZFaZFdaJn/+2Pf+fuD/JPULusXzp6t7c0bL2JooV0REUJgSqZW+2LSfP328kd2HT/rs58/aVBGhIdw+tBP3XnQmzWIiaqpkEZEGS2FKpJbKys1n6jfbefGr7zmRneexj7+X/QDiIsO4+8IzuOP8LkRHaD6ViEh1UZgSqeXSTmTzxLzNzFi522sfT6EKPC+nEB8byU9HnMmYczsSEaa1eUVEqqrKYcrMRgLPAKHAq865x0q13wL8X9HTE8A9zrm1vo6pMCVS1pLvD/Hb2RvYfijDa5+SyymUN5+qQ/NG/PLS7lzdrx0hWvRTRKTSqhSmzCwU2AJcCqQCK4GxzrlNJfoMBZKdc0fMbBTwkHPuPF/HVZgS8Swnr4DpS3fw7PytHMvyfOmvmD/zqQB6tW3ML390Fhd3b6VQJSJSCVUNU0MoDEeXFT1/EMA596iX/s2ADc659r6OqzAl4tvhjByenb+Vfy3bSX6B5+9pReZTAXRuEc24IZ0ZO6iD5lSJiFRAVcPU9cBI59ydRc/HAec55+730v9XQI/i/qXaJgITATp27HjOzp07K/SDiDREKQeO89CHm/gm5ZDXPhUNVc1jIrjrgjO4+byONGkUXiN1i4jUJ1UNUzcAl5UKU4Occw946Hsx8AJwvnMuzddxNTIl4j/nHJ+s38vDH29i/zHPmydDxUNVdEQo1w5sz8QLutKxhbaoERHxJiCX+czsbGA2MMo5t6W8ohSmRCruRHYez87fymvfbCfPy6U/qNidfwBhIcb15yRwx/ldOKt1XM0ULyJSh1U1TIVROAF9BLCHwgnoNzvnNpbo0xH4ErjNObfEn6IUpkQqb/uhDJ78zxY+WvuDz34V3UgZCierjz2vIzeck0BUeGj1Fi4iUkdVx9IIlwNPU7g0wmvOub+Y2SQA59xLZvYqcB1QPAkqz9sHFlOYEqm6DXvS+du8zSzccrDcvhW9BBgfG8Et53Xi+nMS6NBclwBFpGHTop0i9dyS7w/xt882s2b30XL7VjRUAQw7swX3XnQmQ7u20KbKItIgKUyJNADOOeZt3M8T877j+4PeF/0sVplQNbBjU647J4GRvdvQIjay2moXEantFKZEGpC8/ALeX53KP79MIfWI702UoXKhKizEGNmnDT8Z1pmBHZtptEpE6j2FKZEGKL/AseC7A0xbsp3FKT5XKgEqfgdgsbNax3LtwASuG5hAyziNVolI/aQwJdLArdp5hH9+uZWvNvs/UR38vwMQICIshJsSOzBx+BmasC4i9Y7ClIgAhaupv75kBzOTUsnJKyi3f2UuAZrB4C4tuPHcBK48ux3hoSHVVr+ISLAoTInIafalZ/H8ghTeW7WbrNyaCVUAHZo34r6LzuTKfu2IjdRegCJSdylMiYhHx7Jy+Wz9Pt5N2s2qnUfK7e9tXhX4nlsVGRbCBd3iGXZmPBd0i+fMVlplXUTqFoUpESlX0o7DvLJoG/OTD/jcqgZOn1cFVHjSeu92jbl1cCdG929HdIRGrESk9lOYEhG/HTqRzZxv9/Dqou3sO5bl13sqO2k9LjKMawe259bBneimPQFFpBZTmBKRCsvKzWfWqlTeWbGLjT8c8/t9lZ1fNahLc24d3ImRvdsQEaZJ6yJSuyhMiUiV7DiUwfurU3l7+S7SMnL8ek9lQ1V8bATXDGjPhWe1IrFzM222LCK1gsKUiFSLrNx8Ptuwj/dXp/JNyiH8+eujsqEKCteuSuzUjPO7xXPV2e20fpWIBI3ClIhUu11pmby08HtmJaWSk1+55RVK8idgDTuzBVf0bcfFPVrStkmjyhcvIlJBClMiUmMOHs9m7po9zFqVynf7jpfbv/SdgFCxSevFerSJ4+IereiX0JRurWPp0iKGkBDtESgiNUNhSkRqnHOOL787wDPzt7IuNb1C763KpcBiTaPDGXZmPBd2a8kFZ8Vr5EpEqpXClIgEjHOODXuO8cn6vfx7w152pmX6/d7KLgrqSbdWsQw/qyUXdIunf4emNI2O8P+HEBEpRWFKRIJm9+FMFm09xMyk3azZfdSv91R1UVBPEpo1ok+7JvRp35g+7ZvQp30T4mMjK3QMEWm4FKZEpFbYsCedt5bvYu6aPWTm5Pv9Pm+LgharTLgCaNM4ij7tG9O7XRP6FgWs1o0jMdPcKxE5ncKUiNQqx7JymfPtHt5ctpMt+09U6L3ljVpB5cMVFK5z9d9w1ZjEzs01giUiClMiUjs550jaeYQ3l+3k3xv2kZNX/hILpfkKV1UJVSX1bNuYnm3j6NQ8hk4tounQPJou8TE0j9E8LJGGQmFKRGq9jOw8Vuw4zJKUQ3yTkkbyXv+3sCmpvPWsilVH0GoaHc4Z8TGc0TKWtk2iiIkMo3lMBF3iYzgjPoYWGtESqTcUpkSkzilev2pm0u4KXwoEz+tZFavuS4PetG/aiIGdmtGzbRwdmkXTLDqC2KgwWsRE0LpxlPYgFKlDFKZEpE7bdvAECzYf5KvNB1i+7bBfK6774s+8K6iZgFVS68aR9GzbmF5tG9O+WSNax0XRuFE4jRuF0bF5NNERYTX22SJSMQpTIlJvZGTnsWxbGht/OMaW/cdJ2nGEfceyqnTM8lZlL62mQxZAiMEZLWPpEh9DfGwEreKi6Ni8cL5Wx+bRtIqL1IrvIgGkMCUi9ZZzjpQDJ1i49RALtxxk2bY0sisxkb00b5cJAzEXyx8hBo3CQ4mODCM6IpSYiDCaxYTTLDqCFjERNIv5769tGkfRo21jYiM10iVSWQpTItJgZOXms2L7YRZ/f4gNe9LZsOcY6Sdzq+34FZ2LVSxQIcsbs8JV4fslNKVfh6b079CU7m3iCA/VvC0RfyhMiUiD5Zwj9chJNv5QGKzW70lnw5500jJyqv2zKjuaVSzQgSsyLISRfdrwzJgBAftMkbrKV5jSmK+I1GtmRoeiuUYj+7QFCgPW/mPZbNiTzvo96aeCVlXnXr2TP8JjGPI1mlVscEgyg0OSy+1XnYErO6+AUM27EqkyhSkRaXDMjDZNomjTJIpLerU+9frB49ls/CGddanpLNuWRtKOI1W+cxC8h6ySqjNwFfMnePXv0NSvY4mId7rMJyLiRVZuPikHTrAzLZNdhzPZdTiDHYcy2X4oo8qjWJXhT+Aq5s+lxbn5Qxkz6Q/0U6ASKZfmTImIVLOM7Dy2H8pg26EM9qWf5ER2Pkcycth+KIONP6RzJLP6Jr1XRnnBqzhsFbToToWu9IWEQqPmp7/W93pI/EklqhSpO6ocpsxsJPAMEAq86px7rFS7FbVfDmQC451zq30dU2FKROor5xzbDmWwPjWdXYcz2X04k8MZORzPzuNIRg77jmVxPCsvqDWODZ3PzZGL6dvcgSuAMD+3vsnLBitxB+ChzYW/xnevnsI8hbViCm0SRFWagG5mocDzwKVAKrDSzD50zm0q0W0U0K3ov/OAF4t+FRFpcMyMri1j6doy1mufA8ez2PTDMVIOnGBfehb7jmVxLCuPo5k5bN1/gpO5+TVa4zv5I4jodwd9R/eBrf+BrPTKHWjb17B7eeULKR3k8rIh42DZfoc2w85vYNmLlf+sivIV7IJFgbJW8mcC+iAgxTm3DcDMZgCjgZJhajQw3RUOcy0zs6Zm1tY5t7faKxYRqQdaxUXRqnsUF3VvVaYtv8CxIy2DvUezSMvIZl96FruPZLL78El2H85kz9GT1bIw6am5Ut0urfxB+l5ftSL8DXJVDW2VUZDvOdgFSzACZV3R4kwY+3bQPt6fMNUe2F3ieSplR5089WkPnBamzGwiMBGgY8eOFa1VRKRBCA0pf2QrN7+AzJx8Tubkk5GTx7GTuRzOyCEtI4cjGTkcLvpvZ1om6/ekexzpqhUTz/0NclUNbZVRlRG7mhCMQCl+8SdMeZqaWHqilT99cM5NAaZA4ZwpPz5bREQ8CA8NoUmjEJo0Ci+3b15+AVsPnGDt7qOsTT3Kmt3p7Es/SZcWMQGotA6ryohdTQhGoBS/+BOmUoEOJZ4nAD9Uoo+IiARBWGgIPds2pmfbxowZVHhVICevQBsli1QTfzZlWgl0M7MuZhYBjAE+LNXnQ+A2KzQYSNd8KRGR2isiTHvyiVSXckemnHN5ZnY/MI/CpRFec85tNLNJRe0vAZ9SuCxCCoVLI+hWAxEREWkQ/NpOxjn3KYWBqeRrL5V47ID7qrc0ERERkdpP47wiIiIiVaAwJSIiIlIFClMiIiIiVaAwJSIiIlIFfm10XCMfbHYQ2BmAj4oHDgXgc6R66HzVPTpndY/OWd2jcxZ8nZxzLT01BC1MBYqZJXnb5VlqH52vukfnrO7ROat7dM5qN13mExEREakChSkRERGRKmgIYWpKsAuQCtH5qnt0zuoenbO6R+esFqv3c6ZEREREalJDGJkSERERqTEKUyIiIiJVUG/DlJmNNLPNZpZiZpODXY94ZmY7zGy9ma0xs6Si15qb2X/MbGvRr82CXWdDZmavmdkBM9tQ4jWv58jMHiz63m02s8uCU3XD5uWcPWRme4q+a2vM7PISbTpnQWRmHcxsgZklm9lGM/tZ0ev6ntUR9TJMmVko8DwwCugFjDWzXsGtSny42DnXv8QaKpOB+c65bsD8oucSPK8DI0u95vEcFX3PxgC9i97zQtH3UQLrdcqeM4Cnir5r/Z1zn4LOWS2RB/zSOdcTGAzcV3Re9D2rI+plmAIGASnOuW3OuRxgBjA6yDWJ/0YDbxQ9fgP4cfBKEefcQuBwqZe9naPRwAznXLZzbjuQQuH3UQLIyznzRucsyJxze51zq4seHweSgfboe1Zn1Ncw1R7YXeJ5atFrUvs44HMzW2VmE4tea+2c2wuFf8kArYJWnXjj7Rzpu1e73W9m64ouAxZfMtI5q0XMrDMwAFiOvmd1Rn0NU+bhNa0BUTsNc84NpPCS7H1mNjzYBUmV6LtXe70IdAX6A3uBfxS9rnNWS5hZLPA+8L/OuWO+unp4TecsiOprmEoFOpR4ngD8EKRaxAfn3A9Fvx4AZlM4VL3fzNoCFP16IHgVihfezpG+e7WUc26/cy7fOVcAvMJ/LwvpnNUCZhZOYZB6yzn3QdHL+p7VEfU1TK0EuplZFzOLoHCi3odBrklKMbMYM4srfgz8CNhA4bm6vajb7cDc4FQoPng7Rx8CY8ws0sy6AN2AFUGoT0op/ke5yDUUftdA5yzozMyAqUCyc+7JEk36ntURYcEuoCY45/LM7H5gHhAKvOac2xjksqSs1sDswr9HCAPeds59ZmYrgZlmNgHYBdwQxBobPDN7B7gIiDezVOCPwGN4OEfOuY1mNhPYROEdSvc55/KDUngD5uWcXWRm/Sm8HLQDuBt0zmqJYcA4YL2ZrSl67Tfoe1ZnaDsZERERkSqor5f5RERERAJCYUpERESkChSmRERERKpAYUpERESkChSmRERERKpAYUpEgsrMmprZvUWP25nZrGo67kNm9quix382s0uq47giIqVpaQQRCaqivcg+ds71qebjPgSccM79vTqPKyJSmkamRCTYHgO6mtkaM3vPzDYAmNl4M5tjZh+Z2XYzu9/MfmFm35rZMjNrXtSvq5l9VrRZ9iIz61H6A8zsdTO7vujxDjP7k5mtNrP1xf2LVuR/zcxWFn3G6AD+HohIHaYwJSLBNhn43jnXH/h1qbY+wM0U7iP3FyDTOTcAWArcVtRnCvCAc+4c4FfAC3585qGiDbZfLHoPwG+BL51z5wIXA08UbXMkIuJTvdxORkTqjQXOuePAcTNLBz4qen09cLaZxQJDgfeKtiUCiPTjuMUbya4Cri16/CPg6uJ5VkAU0BFIrtqPICL1ncKUiNRm2SUeF5R4XkDh318hwNGiUa3KHDef//49aMB1zrnNlStVRBoqXeYTkWA7DsRV5o3OuWPAdjO7AcAK9atkHfOAB6xoiMvMBlTyOCLSwChMiUhQOefSgMVFE8+fqMQhbgEmmNlaYCNQ2YnjDwPhwLqiWh6u5HFEpIHR0ggiIiIiVaCRKREREZEqUJgSERERqQKFKREREZEqUJgSERERqQKFKREREZEqUJgSERERqQKFKREREZEq+P+mc7wZP0P7mQAAAABJRU5ErkJggg==\n", "text/plain": [""]}, "metadata": {"needs_background": "light"}, "output_type": "display_data"}], "source": ["import numpy\n", "from lifelines import KaplanMeierFitter\n", "fig, ax = plt.subplots(1, 1, figsize=(10, 4))\n", "ax.plot(T, St, label=\"custom\", lw=10)\n", "\n", "kmf = KaplanMeierFitter()\n", "kmf.fit(duree, (issue == 0).astype(numpy.int32))\n", "kmf.plot(ax=ax)\n", "ax.legend();"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q6 : application aux donn\u00e9es publiques\n", "\n", "Les donn\u00e9es accessibles librement sur le portail [data.gouv.fr](https://www.data.gouv.fr/fr/) recensent les entr\u00e9es et les sorties des personnes sans relier une entr\u00e9e et une sortie sp\u00e9cifique. Si N est personnes sont sorties gu\u00e9ries, on ne sait pas quand elles sont entr\u00e9es. Donc le calcul ci-dessus n'est pas possible."]}, {"cell_type": "markdown", "metadata": {}, "source": ["## Exercice 2"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q1 : t + 1"]}, {"cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])"]}, "execution_count": 13, "metadata": {}, "output_type": "execute_result"}], "source": ["import numpy\n", "\n", "N = 10\n", "M = numpy.zeros((N, N))\n", "M[4, 5] = 1\n", "M"]}, {"cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 1., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])"]}, "execution_count": 14, "metadata": {}, "output_type": "execute_result"}], "source": ["def propagation(M):\n", " M2 = M.copy()\n", " M2[1:, :] = numpy.maximum(M2[1:, :], M[:-1, :])\n", " M2[:-1, :] = numpy.maximum(M2[:-1, :], M[1:, :])\n", " M2[:, 1:] = numpy.maximum(M2[:, 1:], M[:, :-1])\n", " M2[:, :-1] = numpy.maximum(M2[:, :-1], M[:, 1:])\n", " return M2\n", " \n", "propagation(M.copy())"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q2 : apr\u00e8s T it\u00e9ration"]}, {"cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 1., 0., 0., 0.],\n", " [0., 0., 0., 1., 1., 1., 1., 1., 0., 0.],\n", " [0., 0., 1., 1., 1., 1., 1., 1., 1., 0.],\n", " [0., 0., 0., 1., 1., 1., 1., 1., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 1., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])"]}, "execution_count": 15, "metadata": {}, "output_type": "execute_result"}], "source": ["def propagation_n(M, t):\n", " for i in range(t):\n", " M = propagation(M)\n", " return M\n", "\n", "propagation_n(M, 3)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q3 : vaccin"]}, {"cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0, 0, 0, 0, 0, 0, 1, 0, 1, 0],\n", " [0, 0, 0, 0, 0, 0, 1, 0, 0, 1],\n", " [1, 0, 0, 0, 0, 0, 1, 0, 1, 1],\n", " [0, 1, 0, 1, 0, 0, 0, 0, 0, 0],\n", " [1, 0, 0, 0, 0, 0, 0, 0, 0, 0],\n", " [0, 0, 1, 1, 0, 0, 1, 0, 1, 0],\n", " [0, 0, 1, 0, 0, 0, 0, 1, 0, 0],\n", " [0, 0, 1, 0, 0, 0, 1, 0, 0, 1],\n", " [0, 0, 1, 0, 0, 0, 0, 1, 0, 0],\n", " [0, 0, 0, 1, 0, 0, 0, 1, 0, 0]])"]}, "execution_count": 16, "metadata": {}, "output_type": "execute_result"}], "source": ["p = 0.3\n", "vaccine = (numpy.random.rand(N, N) <= p).astype(numpy.int32)\n", "vaccine"]}, {"cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 1., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])"]}, "execution_count": 17, "metadata": {}, "output_type": "execute_result"}], "source": ["def propagation_vaccine(M, vaccine):\n", " M2 = M.copy()\n", " M2[1:, :] = numpy.maximum(M2[1:, :], M[:-1, :])\n", " M2[:-1, :] = numpy.maximum(M2[:-1, :], M[1:, :])\n", " M2[:, 1:] = numpy.maximum(M2[:, 1:], M[:, :-1])\n", " M2[:, :-1] = numpy.maximum(M2[:, :-1], M[:, 1:])\n", " M2 = numpy.minimum(M2, 1 - vaccine)\n", " return M2\n", "\n", "vaccine[4, 5] = 0\n", "propagation_vaccine(M, vaccine)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q4 : apr\u00e8s T heures"]}, {"cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [{"data": {"text/plain": ["array([[0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 1., 1., 0., 0.],\n", " [0., 0., 1., 1., 1., 1., 1., 1., 1., 0.],\n", " [0., 0., 0., 0., 1., 1., 0., 1., 0., 0.],\n", " [0., 0., 0., 0., 1., 1., 1., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 1., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.],\n", " [0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]])"]}, "execution_count": 18, "metadata": {}, "output_type": "execute_result"}], "source": ["def propagation_n_vaccine(M, t, vaccine):\n", " for i in range(t):\n", " M = propagation_vaccine(M, vaccine)\n", " return M\n", "\n", "propagation_n_vaccine(M, 3, vaccine)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["### Q5 : variation"]}, {"cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [{"data": {"text/html": ["\n", "\n", "
\n", " \n", " \n", " \n", " p \n", " prop \n", " \n", " \n", " \n", " \n", " 0 \n", " 0.70 \n", " 9.1 \n", " \n", " \n", " 1 \n", " 0.80 \n", " 9.6 \n", " \n", " \n", " 2 \n", " 0.85 \n", " 3.8 \n", " \n", " \n", " 3 \n", " 0.86 \n", " 2.4 \n", " \n", " \n", " 4 \n", " 0.87 \n", " 2.4 \n", " \n", " \n", " 5 \n", " 0.88 \n", " 1.4 \n", " \n", " \n", " 6 \n", " 0.89 \n", " 1.6 \n", " \n", " \n", " 7 \n", " 0.90 \n", " 1.2 \n", " \n", " \n", " 8 \n", " 0.94 \n", " 1.2 \n", " \n", " \n", " 9 \n", " 0.95 \n", " 1.5 \n", " \n", " \n", " 10 \n", " 0.96 \n", " 1.0 \n", " \n", " \n", " 11 \n", " 0.97 \n", " 1.2 \n", " \n", " \n", " 12 \n", " 0.98 \n", " 1.2 \n", " \n", " \n", " 13 \n", " 0.99 \n", " 1.1 \n", " \n", " \n", " 14 \n", " 1.00 \n", " 1.0 \n", " \n", " \n", "
\n", "
"], "text/plain": [" p prop\n", "0 0.70 9.1\n", "1 0.80 9.6\n", "2 0.85 3.8\n", "3 0.86 2.4\n", "4 0.87 2.4\n", "5 0.88 1.4\n", "6 0.89 1.6\n", "7 0.90 1.2\n", "8 0.94 1.2\n", "9 0.95 1.5\n", "10 0.96 1.0\n", "11 0.97 1.2\n", "12 0.98 1.2\n", "13 0.99 1.1\n", "14 1.00 1.0"]}, "execution_count": 19, "metadata": {}, "output_type": "execute_result"}], "source": ["import pandas\n", "\n", "res = []\n", "for p in [0.7, 0.8, 0.85, 0.86, 0.87, 0.88, 0.89, 0.9, 0.94, 0.95, 0.96, 0.97, 0.98, 0.99, 1.0]:\n", " cont = []\n", " for test in range(0, 10):\n", " vaccine = (numpy.random.rand(N, N) <= p).astype(numpy.int32)\n", " M[4, 5] = 1\n", " vaccine[4, 5] = 0\n", " M = propagation_n_vaccine(M, 3, vaccine)\n", " contamine = M.ravel().sum()\n", " cont.append(contamine)\n", " cont = numpy.array(cont)\n", " res.append(dict(p=p, prop=cont.mean()))\n", "\n", "df = pandas.DataFrame(res)\n", "df"]}, {"cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAEKCAYAAAALoA6YAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAg6UlEQVR4nO3de3RV5Z3/8feT68mNJISQ2wkEFAUCguSIIt5RwVuhQNBWR1s73qYXZ1Y7/bW/mV9b2zW1ztRZrbV1xo46ulQsF+u1KpVqLXjBhItcFbmHBBJuAUJCbs/vjxOQhISc5Jxkn73zea2VlZOz906+z9rwyc6zn/08xlqLiIi4T4zTBYiISO8owEVEXEoBLiLiUgpwERGXUoCLiLiUAlxExKW6DXBjzJPGmGpjzLpT3htsjPmzMWZz2+fMvi1TREQ6CuUK/H+BGR3e+wGw1Fo7Clja9rWIiPQjE8qDPMaYIuA1a+24tq8/Ba6w1lYZY/KAd6215/ZppSIi0k5cL4/LsdZWAbSF+NBQDhoyZIgtKirq5Y8UERmYysvL91lrszu+39sAD5kx5m7gboBhw4ZRVlbW1z9SRMRTjDE7Onu/t6NQ9rZ1ndD2ubqrHa21j1trA9baQHb2ab9ARESkl3ob4K8Ad7S9vgN4OTLliIhIqEIZRjgf+AA41xhTYYz5BvAL4BpjzGbgmravRUSkH3XbB26t/UoXm6ZFuBYRkU41NTVRUVFBQ0OD06X0KZ/Ph9/vJz4+PqT9+/wmpohIuCoqKkhLS6OoqAhjjNPl9AlrLfv376eiooIRI0aEdIwepReRqNfQ0EBWVpZnwxvAGENWVlaP/spQgIuIK3g5vE/oaRvVhSIhqzveTFVtPZWHGqiqredAXRO3XFBIZkqC06WJDEgKcAGgoamFPbUNVNbWU9UW0JW1DVQdqqeqtoHKQ/Ucbmg+7bgjDU18f8ZoByoWEQX4ANDU0sreww0ng7iqLZgra4NBXXWogf11jacdNzglgbx0H/7MZCaPGExeehL5GT7y0pPIS/fx/15ex+KVFXz32nOJjfH+n7ci3WlpaSE2Nrbffp4C3OVaWi37jh4/GcwnA/qUro7qI8fpOGdZmi+O/PQk8jJ8jC/IID/dR15G0snPeek+fPFn/od4c6CQ+55byXuba7jy3JCmwxFxre3btzNjxgwuvPBCVq1axTnnnMMzzzzD2LFjufPOO1myZAnf+ta3sNby85//HGstN9xwAw899BAAqamp3HPPPbzzzjtkZmbywgsvEO7T6QrwKGat5UBdY7tgbtfFcaiBvYcbaG5tn85J8bHkZfjIT0/islHZ7YL5xOfUxPBP/bQxOWQmx7OorEIBLv3mgVfXs6HycES/59j8Qfz4puJu9/v000954oknmDp1KnfeeSe/+93vgOD47WXLllFZWclFF11EeXk5mZmZXHvttbz00kvMmjWLuro6Jk2axMMPP8xPf/pTHnjgAR599NGw6laAO8Ray+GG5pNdGJUdPlfVBgP7eHNru+MSYmPITfeRl+5r69Y4JZjbujjSk+L75Y59QlwMs84v4LkPd3KwrlE3M8XzCgsLmTp1KgC33XYbjzzyCAA333wzAB9//DFXXHHFySvrW2+9lffee49Zs2YRExNzcr/bbruN2bNnh12PAryPHGtsPtmFcVpAt/VB1zW2tDsmNsaQk5ZIXkYS4wrSubY4NxjQp/Q9Z6UkEBNF/c2lJYU8tXw7L6/ezdemhvbwgUg4QrlS7isdL4xOfJ2SkgIEL8x6+716QwHeC8eb20ZsnHKl3LEPura+6bTjstMSyU/3cXZ2KpeOGnKyD/pEQGenJhIX666h+WPzBzGuYBALyysU4OJ5O3fu5IMPPmDKlCnMnz+fSy65hFWrVp3cfuGFF3L//fezb98+MjMzmT9/Pt/+9rcBaG1tZdGiRdxyyy08//zzXHLJJWHXowDvoLmllb1Hjn8xSqOTm4P7jp4+YiMzOZ689CT8mUlcUDT4ZB90XrqP/Iwkcgb5SIhzVziHqrSkkB+/sp71lbUU56c7XY5InxkzZgxPP/0099xzD6NGjeK+++7jN7/5zcnteXl5PPjgg1x55ZVYa7n++uuZOXMmELxKX79+PSUlJaSnp/OHP/wh7HpCWlItUgKBgHVyQYfWVktNJyM2Th3/XH2kgdaOIzYS49pdKeedEswnujiSEvpv6FC0OXSskcn/tpSvXjiMn3zJuT9vxbs2btzImDFjHK1h+/bt3Hjjjaxbt677nTuRmprK0aNHu92vs7YaY8qttYGO+3rmCry3IzZ88TEnuzIuGTWk3TC6EwGd5gttZrCBKiM5gWuKc3h59W5+eP1oEuMG7i8zkf7kigDv7YiN+FjTNmIjiQuKMtuN1jjRxZGR3D8jNryutMTP659UsXRjNdePz3O6HJGIKyoq6vXVNxDS1XdPuSLA//WldTz30c5278UYyBkUHE5XXJDONWNz2ndxZPgYkpIYVSM2vOzSUdnkDvKxsGyXAlz6hLXW8xdbPe3SdkWAXz02h+FZye0Cemia+0ZseFlsjGFOSQGPvbuFPbUN5Kb7nC5JPMTn87F//35PTyl7Yj5wny/0/zuuCPArzx2qJ/1cYG5JIb99ZwsvrqrgH6442+lyxEP8fj8VFRXU1NQ4XUqfOrEiT6hcEeDiDiOGpDC5aDCLyiq47/KzPHulJP0vPj4+5FVqBhL1QUhEzQ342bqvjvIdB50uRcTzFOASUTeMzyM5IZaFZRVOlyLieQpwiaiUxDhuGJ/Ha59Ucqzx9AUgRCRyFOAScfMuKKSusYU/rd3jdCkinqYAl4gLDM9kxJAUFpTtcroUEU9TgEvEGWOYW+JnxbYD7Nhf53Q5Ip6lAJc+MXtSATEGFpXrZqZIX1GAS5/IS0/i0lHZLCqvoKXj9I4iEhEKcOkz8wKFVNU2sPzzfU6XIuJJCnDpM1ePHUpGcrxuZor0EQW49JnEuFhmTshnyYa91B47fYk5EQmPAlz6VGmgkMbmVl5Zs9vpUkQ8RwEufWpcQTpj8gaxQI/Wi0ScAlz63LyAn7W7a9lYddjpUkQ8RQEufW7mxALiY40muBKJMAW49LnBKQlcPSaHl1bvprHDuqUi0nsKcOkX8wKFHKhr5C+bqp0uRcQzFODSLy4dNYShaYks1JhwkYgJK8CNMf9kjFlvjFlnjJlvjNFKttKpuNgY5pT4efezGqoPNzhdjogn9DrAjTEFwHeAgLV2HBAL3BKpwsR7Skv8tLRaXlylMeEikRBuF0ockGSMiQOSgcrwSxKvGpmdSsnwTBaW7cJaTXAlEq5eB7i1djfwS2AnUAXUWmuXRKow8aZ5AT9baupYteuQ06WIuF44XSiZwExgBJAPpBhjbutkv7uNMWXGmLKampreVyqecMN5+STFx+pmpkgEhNOFcjWwzVpbY61tAl4ELu64k7X2cWttwFobyM7ODuPHiRekJsZx/fg8Xl1TRX1ji9PliLhaOAG+E7jIGJNsjDHANGBjZMoSLysN+Dl6vJk31lU5XYqIq4XTB/4RsAhYCaxt+16PR6gu8bALRwxmeFayHq0XCVNYo1CstT+21o621o6z1v6dtfZ4pAoT7zLGMHeSnw+27mfn/mNOlyPiWnoSUxwxp8SPMbBopa7CRXpLAS6OyM9I4pKzh7C4vIJWLXos0isKcHFMaaCQ3YfqeX/LfqdLEXElBbg45tqxOQzyxbGwXGPCRXpDAS6O8cXHMnNiAW+u20NtvRY9FukpBbg4qjTg53hzK6+u0TQ6Ij2lABdHjS9IZ3RuGgvLNRpFpKcU4OIoYwxzS/ys2XWIz/YecbocEVdRgIvjvnx+AXExRhNcifSQAlwcl5WayLQxQ/njqt00tWjRY5FQKcAlKpSWFLLvaCPvaNFjkZApwCUqXHFuNtlpibqZKdIDCnCJCnGxMcw+v4C/bKqm5ojmRBMJhQJcokZpILjo8Uta9FgkJApwiRpnD03j/GEZLNCixyIhUYBLVCktKWRz9VHWVNQ6XYpI1FOAS1S5cUIevvgYjQkXCYECXKLKIF88143L45U1lTQ0adFjkTNRgEvUKQ34OdLQzFvr9zhdikhUU4BL1LloRBb+zCQWqBtF5IwU4BJ1YmIMpSWFvL9lPxUHteixSFcU4BKV5pQUALC4XGPCRbqiAJeo5M9M5uKzslhYvkuLHot0QQEuUWteoJCKg/V8uE2LHot0RgEuUWt6cS5pvjgWlmmCK5HOKMAlavniY/nShHzeWFfF4QYteizSkQJcolppoJCGplZeW1PldCkiUUcBLlFtgj+dUUNTWViuMeEiHSnAJaoZY5gXKGTVzkN8Xq1Fj0VOpQCXqDfr/AJiY4xuZop0oACXqJedlshVo4eyeKUWPRY5lQJcXKG0xM++o8f566c1TpciEjUU4OIKV44eypDUBN3MFDmFAlxcIT42hi+fX8DSjdXsP6pFj0VAAS4uUhoopLnV8kcteiwCKMDFRc7JSWNCYQaLyiu06LEICnBxmdISP5v2HGHtbi16LBJWgBtjMowxi4wxm4wxG40xUyJVmEhnbpqQT2JcjMaEixD+FfivgTettaOBCcDG8EsS6Vp6UjwzxuXy8urdWvRYBrxeB7gxZhBwGfAEgLW20Vp7KEJ1iXSptKSQww3NLNmw1+lSRBwVzhX4SKAGeMoYs8oY8z/GmJQI1SXSpYvPyqIgI4mFWvRYBrhwAjwOmAQ8Zq09H6gDftBxJ2PM3caYMmNMWU2NnqKT8MXEGOaU+Fn2+T52H6p3uhwRx4QT4BVAhbX2o7avFxEM9HastY9bawPW2kB2dnYYP07kC6UlfqyFF8t1M1MGrl4HuLV2D7DLGHNu21vTgA0RqUqkG4WDk5kyMouF5RVa9FgGrHBHoXwbeM4Y8wkwEfh52BWJhKg04GfngWOs2H7A6VJEHBFWgFtrV7d1j5xnrZ1lrT0YqcJEunPduDxSE7XosQxcehJTXCspIZabJuTxp7VVHD3e7HQ5Iv1OAS6uVhoopL6phdc/qXS6FJF+pwAXVzu/MIOzslNYoG4UGYAU4OJqxhhKA4WU7zjIlpqjTpcj0q8U4OJ6s9sWPV6kMeEywCjAxfWGDvJxxTnZLC6voFmLHssAogAXTygNFFJ95Dh/27zP6VJE+o0CXDzhqtFDGZySwAJNcCUDiAJcPCEhLoZZEwt4e+NeDtQ1Ol2OSL9QgItnzLvAT1OL5eXVWvRYBgYFuHjG6NxBjC9I15hwGTAU4OIp8wJ+NlYdZp0WPZYBQAEunvKlCQUkxMVotR4ZEBTg4inpyfFcOzaHl9dUcrxZix6LtynAxXPmBQo5dKyJtzdUO12KSJ9SgIvnTD17CHnpPo0JF89TgIvnxMYY5pb4+dvmGqpqteixeJcCXDxpbomfVgsvrtSYcPEuBbh40vCsFCaPGMzCsl1Yq0WPxZsU4OJZ8wKFbN9/jLIdWqpVvEkBLp51/fhcUhJiWfCxbmaKNynAxbOSE+K48bx8Xl9bRZ0WPRYPUoCLp5UG/BxrbOH1tVVOlyIScQpw8bSS4ZmMHJLCIk1wJR6kABdPM8YwN+BnxfYDbNtX53Q5IhGlABfPmzPJT4yBReW6mSneogAXz8sZ5OPyc7JZXL6bllaNCRfvUIDLgFAaKGTP4Qb+trnG6VJEIkYBLgPCtDFDyUyOZ2G5bmaKdyjAZUBIjItl5sQC/rx+L4eOadFj8QYFuAwYpQE/jS2tvLy60ulSRCJCAS4DRnF+OsX5g1io0SjiEQpwGVBKS/ys232YDZWHnS5FJGwKcBlQZk4sICE2Rlfh4gkKcBlQMlMSuGZsDi+t2k1jc6vT5YiERQEuA87cgJ+Dx5pYunGv06WIhEUBLgPOZaOyyR3k05hwcb2wA9wYE2uMWWWMeS0SBYn0tdgYw+xJBbz7aTV7Dzc4XY5Ir0XiCvx+YGMEvo9IvykNFGrRY3G9sALcGOMHbgD+JzLliPSPEUNSuKAoU4sei6uFewX+K+D7gG7ni+uUlhSydV8dK3dq0WNxp14HuDHmRqDaWlvezX53G2PKjDFlNTWaCU6ix/Xn5ZGcEMtCrdYjLhXOFfhU4EvGmO3AC8BVxphnO+5krX3cWhuw1gays7PD+HEikZWaGMf14/N4dU0lxxq16LG4T68D3Fr7Q2ut31pbBNwC/MVae1vEKhPpB/MChdQ1tvDG2j1OlyLSYxoHLgPaBUWZFGUls6BMj9aL+0QkwK2171prb4zE9xLpT8YY5pb4+WjbASoOHnO6HJEe0RW4DHg3nJcPwJL1erRe3EUBLgPeiCEpnJOTypIN6gcXd1GAiwDTi3NZse0AB+q03Jq4hwJcBLh2bC6tFt7WDIXiIgpwEWBcwSAKMpJYsl7dKOIeCnARgqNRrhmbw3ub91F3XA/1iDsowEXaTC/OpbG5lfc+05QP4g4KcJE2FxRlkpkcz1vqRhGXUICLtImLjWHamByWbqrWepniCgpwkVNML87lSEMzH23b73QpIt1SgIuc4tJRQ0iKj1U3iriCAlzkFL74WC4/J5sl6/fS2qqVeiS6KcBFOpg+LofqI8dZU3HI6VJEzkgBLtLBVefmEBdjeEuTW0mUU4CLdJCeHM9FI7NYsn6PFjyWqKYAF+nE9OIctu6r4/Pqo06XItIlBbhIJ64ZmwvAkg3qRpHopQAX6URuuo+JhRkaTihRTQEu0oVri3P4pKKWykP1Tpci0ikFuEgXphcHu1H+rG4UiVIKcJEunJWdytlDU9WNIlFLAS5yBtOLc/ho2wEOaqk1iUIKcJEzuHZsLi2tlqWbqp0uReQ0CnCRMzjPn05euk9LrUlUUoCLnIExhmvH5vDe5hrqG1ucLkekHQW4SDemF+fS0NTKX7XUmkQZBbhINy4YMZj0pHh1o0jUUYCLdCM+NoZpY4aydFM1TS1aak2ihwJcJATTi3OprW9ixbYDTpcicpICXCQEl43Kxhcfo4d6JKoowEVCkJQQy2WjgkutNbe00tJqz/ih5dikP8Q5XYCIW0wvzmXJhr2c/S9vdLtvjIHf3TqJGePy+qEyGagU4CIhunFCHgfqGqlv6n48+MLyXTz6zudML87FGNMP1clApAAXCVFiXCx3XTYypH0HpyTwry+t4+PtB5k8YnAfVyYDlfrARfrAnEl+MpLjeWLZVqdLEQ9TgIv0gaSEWL46eRhLNuxl5/5jTpcjHqUAF+kjt08pItYYnnp/m9OliEf1OsCNMYXGmHeMMRuNMeuNMfdHsjARt8tN93HjeXks+HgXhxuanC5HPCicK/Bm4LvW2jHARcA3jTFjI1OWiDd845KR1DW2sODjXU6XIh7U6wC31lZZa1e2vT4CbAQKIlWYiBeM96czecRgnlq+nWbNoyIRFpFhhMaYIuB84KNIfL9OPXXD6e8Vz4LJd0HjMXiu9PTtE78K598Kdfthwe2nb7/gThg3B2or4MV7Tt9+8bfg3Otg32Z49R9P337Z9+CsK6HqE3jzh6dvn/YjGHYh7PwIlv709O0zHoS882DLO/DeL0/fftOvYMgo+PQNeP/R07fP/m9I98O6xfDxk6dvn/cMpGTBqudg9fOnb791ISQkw4rfw/qXTt/+9deDn5c/Ap+91X5bvA9uWxx8/dd/h61/bb89ORNufjb4+u2fwK6P228flA9zfh98/cYPYM/a9tuzzoIvPRJ8/cp3YP+W9ttzx8N1vwi+XnwXHK5sv73wArj6J8HXf7gNjh1sv33k5XD594Ovn50DTQ3tt58zHaZ+J/g6zH97/9X8Iz47doTax1LJSkkMbte/vYH3b+9EmyIo7JuYxphUYDHwj9baw51sv9sYU2aMKaup0XzKMvBkJseTGBdDVW1D9zuL9ICxtvdzNhhj4oHXgLestf/Z3f6BQMCWlZX1+ueJuNVTy7fxwKsbeOmbU5lYmOF0OeIyxphya22g4/vhjEIxwBPAxlDCW2QgKw0UkpYYxxPLNKRQIiecLpSpwN8BVxljVrd9XB+hukQ8JTUxjlsmF/KntVVUHqp3uhzxiHBGoSyz1hpr7XnW2oltH3+KZHEiXnLHxUVYa3n6g+1OlyIeoScxRfqJPzOZ68blMf+jndQdb3a6HPEABbhIP7rzkhEcbmhm8coKp0sRD1CAi/SjkuGZTCzM4Knl27Vqj4RNAS7Sz75xyQi27avjL5uqQ9q/5shxfrN0M699Utn9zjKgaEEHkX523bhc8tN9PLFsG1ePzelyv8+rj/LEsq0sXrmbxubgY/i7DtRz7+UjtcqPAApwkX4XFxvDHRcX8eAbm1hfWUtxfvrJbdZaPt5+kMff28LbG6tJiIthbomf26cM53fvbOGhNzdRc+Q4/3rDGGJivB3i63bXUpiZTHpyvNOlRC0FuIgDbpk8jF8v3cyTy7bz8LwJNLe08tb6vTz+t62s2XWIzOR4vjNtFLdPGc6Q1OD8Kb+6eSJZqQk8uXwb+44e55elE0iI814v6NHjzfzb6xuZv2Ino3PTWHDvFAb5FOKdUYCLOCA9KZ7SEj/Pr9jJ2UNTeX7FDnYdqKcoK5mfzRrH3El+khJi2x0TE2P40Y1jGZrm46E3N3HwWCOP3VZCaqJ3/ht/uHU/31u4ht2H6pk9qYBXVldy37PlPPW1yZ78ZRWusOZC6SnNhSLyhe376rjy4XexNjg65a5LR3LN2BxiQ+gaWVC2ix++uJbi/EE8+bULTl6lu1VDUwv//uanPLl8G0VZyTw8bwIlwwezuLyC7y5cw5fPL+A/500YsH3/Xc2F4p1f3SIuUzQkhSfuCJCelEDJ8MweHTsvUEhWSgLffH4lpf/1Ac/cOZnCwcl9VGnfWrXzIN9duIatNXXcMWU4/+e60SQnBKNpTomfqtp6frnkM/IzfPzz9NEOVxtd9DeJiIOuGp3T4/A+YdqYHJ77+ws5UNfI7MfeZ0PlabM5R7XjzS38x1ubmPPY+zQ0tvDc31/IAzPHnQzvE7555dl8ZfIwfvvOFp79cIdD1UYnBbiIi5UMH8yie6cQF2O4+b8/4MOt+50uKSQbKg8z89Hl/PadLcwt8fPmP13G1LOHdLqvMYafzSzmqtFD+dHL63h7w95+rjZ6qQ9cxAMqD9Vz+5Mr2HngGN+84mwGJUVv72hVbQNPLd9GRnICv5g9nmljuh4Lf6pjjc3c8viHfLb3CC/cPaVH86pbaynfcZC1u2tD2j8lIY6JwzI4Ozs1KoZrdtUHrgAX8YhDxxq565kyPt5+sPudHXbjeXn8bOY4MlMSenRczZHjzH5sOceOt7D4vospGpJyxv07Ds/sqTRfHJOGZVIyPPgxoTDDkVE/CnCRAcBaS219k9NlnFFsjCEtjHHdW2uOMuex90lPimfxfReT1ckInLrjzSws28UTy7edHJ75jUtHct24XOJCuKLeX9fIqp2HKN9xkJU7DvJZ9RGshRgDo3MHnQz0kuGZ+DOT+nx0jAJcRDyjfMcBvvr7jxiTN4j5d110csx89eEGnv5gO89+uJPa+qYeD8/sSm19E6t3fRHoq3YepK6xBYChaYknw3zS8EyK8weRGBfbzXfsGQW4iHjKm+v2cN9z5UwbncP3pp/Dk8u28dKqSppaW5k+Npe7LhtByfDBffKzW1otn+45QvmOA6xsu1LfeeAYAAlxMZxXkH4y0CcNyyQ7Lbxx+gpwEfGcp9/fzo9fWQ+ALz6GeYFC7pw6otu+8b5QfaSBlTsOsXLnweAN04paGluCk5ANz0rmwdnjufiszkfadEcP8oiI55xYpq6usYWvTh7W45uikTQ0zceMcbnMGJcLBMe5r9t9mJU7goE+NM0X8Z+pK3ARkSjX1RW4HuQREXEpBbiIiEspwEVEXEoBLiLiUgpwERGXUoCLiLiUAlxExKUU4CIiLtWvD/IYY2qA3i6pMQTYF8FynKS2RB+vtAPUlmgVTluGW2uzO77ZrwEeDmNMWWdPIrmR2hJ9vNIOUFuiVV+0RV0oIiIupQAXEXEpNwX4404XEEFqS/TxSjtAbYlWEW+La/rARUSkPTddgYuIyCmiIsCNMTOMMZ8aYz43xvygk+3/bIxZ3faxzhjTYowZHMqx/SnMdmw3xqxt2+b4pOkhtCXdGPOqMWaNMWa9MebroR7b38Jsi9vOS6Yx5o/GmE+MMSuMMeNCPbY/hdmOqDknxpgnjTHVxph1XWw3xphH2tr5iTFm0inbwj8f1lpHP4BYYAswEkgA1gBjz7D/TcBfenNstLaj7evtwBCnz0eobQH+L/BQ2+ts4EDbvlFzTsJti0vPy38AP257PRpY2pt/n9Hajig8J5cBk4B1XWy/HngDMMBFwEeRPB/RcAU+GfjcWrvVWtsIvADMPMP+XwHm9/LYvhROO6JNKG2xQJoxxgCpBEOvOcRj+1M4bYk2obRlLLAUwFq7CSgyxuSEeGx/CacdUcVa+x7Bfy9dmQk8Y4M+BDKMMXlE6HxEQ4AXALtO+bqi7b3TGGOSgRnA4p4e2w/CaQcEQ2SJMabcGHN3n1UZmlDa8igwBqgE1gL3W2tbQzy2P4XTFnDfeVkDzAYwxkwGhgP+EI/tL+G0A6LrnHSnq7ZG5HxEw6LGppP3uhoacxOw3Fp74jdeT47ta+G0A2CqtbbSGDMU+LMxZlPbb3cnhNKW6cBq4CrgLII1/y3EY/tTr9tirT2M+87LL4BfG2NWE/xltIrgXxPRdF7CaQdE1znpTldtjcj5iIYr8Aqg8JSv/QSvhDpzC+27HXpybF8Lpx1YayvbPlcDfyT4J5ZTQmnL14EX2/40/BzYRrCvMprOCYTXFtedF2vtYWvt1621E4HbCfbpbwvl2H4UTjui7Zx0p6u2RuZ8RMFNgDhgKzCCLzrzizvZL51gX1NKT491QTtSgLRTXr8PzIjmcwI8Bvyk7XUOsJvgZD1Rc04i0BY3npcMvrgBexfB/lfX/V85Qzui6py01VFE1zcxb6D9TcwVkTwfjjW6QyOvBz4jeFf2X9reuxe495R9vga8EMqxbmsHwTvRa9o+1jvdjlDaAuQDSwj+ebsOuC0az0k4bXHpeZkCbAY2AS8CmdF4Xnrbjmg7JwT/kq4CmgheVX+jQzsM8Nu2dq4FApE8H3oSU0TEpaKhD1xERHpBAS4i4lIKcBERl1KAi4i4lAJcRMSlFOAiIi6lABcRcSkFuAxoxpgiY8wmY8zTbfM1L2qbbEwk6inAReBc4HFr7XnAYeAfHK5HJCQKcBHYZa1d3vb6WeASJ4sRCZUCXOT0aTw1v4S4ggJcBIYZY6a0vf4KsMzJYkRCpQAXgY3AHcaYT4DBBKeXFYl60bAij4jTWq219zpdhEhP6QpcRMSlNB+4iIhL6QpcRMSlFOAiIi6lABcRcSkFuIiISynARURcSgEuIuJS/x/vKrtmCvBIKgAAAABJRU5ErkJggg==\n", "text/plain": [""]}, "metadata": {"needs_background": "light"}, "output_type": "display_data"}], "source": ["ax = df.plot(x=\"p\", y=\"prop\")\n", "ax.plot([0.7, 1.], [2, 2], \"--\");"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Dans cette configuration, en supposant qu'un seul \u00e9l\u00e8ve est contamin\u00e9, il faudrait vacciner \u00e0 plus de 85% pour avoir une chance de ne pas avoir une nouvelle contamination.\n", "\n", "Ce r\u00e9sultat ne change pas si la taille de la matrice change. Il change en revanche en fonction du nombre de tirages, ici 10. Pour n'avoir aucune propagation, il faut entourer la personne contamin\u00e9 de 4 personnes vaccin\u00e9es. $p$ est la probabilit\u00e9 pour une personn\u00e9 d'\u00eatre vaccin\u00e9e (ou immunis\u00e9e si la vaccination est faite sur 100% des personnes). $p^4$ est la probabilit\u00e9 d'avoir 4 personnes vaccin\u00e9es. $q=1-p^4$ est la probabilit\u00e9 qu'une personne ne soient pas immunis\u00e9es parmi les 4. $1 - (1-q)^{10}$ est la probabilit\u00e9 que 10 cours de 3h se passent bien sans contamination."]}, {"cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [{"data": {"text/plain": ["0.9999"]}, "execution_count": 21, "metadata": {}, "output_type": "execute_result"}], "source": ["p = 0.9\n", "p4 = 1 - (1 - p) **4\n", "p4"]}, {"cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [], "source": ["df2 = df.copy()\n", "df2[\"P4**100\"] = 1 - df2[\"p\"] ** 40"]}, {"cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [{"data": {"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlAAAAEKCAYAAADdH2tJAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAA24ElEQVR4nO3deXyU9bn//9eVScJkI3uALJCwZNhFDKCAmqhVFFvtYrVqrVutrVrb8z2nted8e6w9p9r9Vz1aPVat9WsrtWqtCy6tBRVXQJFFCPsSghDCnpD98/tjAgYIkJCZuSeT9/PxyCMzc98z885NcnPN/fnc123OOURERESk6+K8DiAiIiLS26iAEhEREekmFVAiIiIi3aQCSkRERKSbVECJiIiIdJMKKBEREZFuOm4BZWaPmNk2M1va4bEsM/u7ma1q/54Z3pgiIiIi0aMrR6AeBWYc9thtwGvOuRHAa+33RURERPoE60ojTTMrBl5wzo1tv18JlDvntpjZIGCucy4Q1qQiIiIiUSL+BJ83wDm3BaC9iMrrypNycnJccXHxCb6liPRGCxcu3O6cy/U6R09p/yXS9xxr/3WiBVSXmdkNwA0AgwcPZsGCBeF+SxGJIma2wesMoVBcXKz9l0gfc6z914mehbe1feiO9u/bjraic+5B51yZc64sN7fXfwgVEREROeEC6jnga+23vwb8LTRxRERERKJfV9oYPAG8AwTMrMrMrgN+CnzGzFYBn2m/LyIiItInHHcOlHPuK0dZdHaIs4jEpObmZqqqqmhoaPA6Slj5/X4KCwtJSEjwOoqISNiFfRK5SF9XVVVFWloaxcXFmJnXccLCOUdtbS1VVVWUlJR4msXMHgEuBLYdaL1y2HID7gYuAOqBq51zH0Q2pYj0drqUi0iYNTQ0kJ2dHbPFE4CZkZ2dHS1H2R7lyOa/HZ0PjGj/ugG4PwKZRCTGqIASiYBYLp4OiJaf0Tn3BrDjGKtcBDzmgt4FMg6cVSwi0lUawuvFnHPU7GtkbU0d67bXAXDZpKKo+Y9MJEoVAJs63K9qf2xLqN7ggdfXUN/Y0vUndONvtjt/3d3ZFVg3Xrl7r9t1cXGGWTBLnEGcBe93/B5nwWL908c+vR+8fWD9Dq8Bh94/4jXb37fD/QOvZWbExxkJvjgSfHEk+uJIiG+/Hxe83S/ehy9O+92+RgVUL7CvsYV1NXWs3b6PdduDxdKBomnfYTvpcQXpjC1I9yipSK/Q2f90nV7T6vBGwF318Lx1bN/X2KV1u3A1LekFEnyGP8HX/hWHP95HUmLwfmq/eNL88e3fE0jzf3o/PSmBnNR+ZKcmkpPaD3+Cz+sfRbpIBVSUaG5tY+OOeta1F0Zrt9extiZYMG3b++mO2AwKMpIoyUnhixMLKMlJYWhuKpnJiXz23nnMrdymAkpOSGtrKz5fn9h5VwFFHe4XAtWdreicexB4EKCsrKzLpc78/zinJ/lCoivXOf103W68brgytOdocw7nwOFoO3C/rf07we8H1unse5sLvm/H752v59rfs32dtvbncORrtbQ6WtraaGp1NLe00dwa/GpqdTS1tNHU0sb+5lYamltpbGmlobmN/U2tNLS0sr+pla17GlhT08Lehhb2NjTT3Hr07ZLWL/5gMZWT2o9BGX6GZCUzJDuFwdnJFGYm0S++T/ydRj0VUBHknGPrnsZPjyTVBAulddvr2Lijnta2T/+oslISKclJ4YzSXEpyUhiWm0JJTipDspOP+gllfGE6cytruPmsEZH6kaSXWL9+PTNmzGDKlCl8+OGHlJaW8thjjzF69GiuvfZaXn31VW6++Wacc9x5550455g5cyY/+9nPAEhNTeUb3/gGc+bMITMzk1mzZtGLryzwHHCzmc0CpgC7D1zbM5Z0Zyg/fKP+GtbqTENzK/sagwXVrvomavc1sX1fI7V1TdTsbQze3tfE6pp9vL6yhv3NrQefawb56UkMzkqmJDeFiycUMKk4U1M3PKACKgz2NDR/eiSpZt/BImnd9jrqmz79Q/AnxFGcncKoQWlcMG4gQ3NSKclNYWhOChnJid1+3/JAHvf+cxW76ptO6PkSfnc8v4yPq/eE9DVH5/fn9s+OOe56lZWVPPzww0ybNo1rr72W3/72t0Cwf9O8efOorq7m1FNPZeHChWRmZnLuuefy7LPPcvHFF1NXV8fEiRP51a9+xY9//GPuuOMO7r333pD+HKHS3vy3HMgxsyrgdiABwDn3ADCbYAuD1QTbGFzjTVLpqw4M9eWk9gNSjrnugbmuG2vr2VBbz4Yd9WysrWPDjnqeX1TNn97byPjCdK6bXsIF4waR4NO5YZGiAuoENba0smlHPWsPHEXqMPTWce5DnEFhZjIlOSlMKs46eCSpJDeFQf39xIVw4mFFIJd7XlvFG6u287mT8kP2uhIbioqKmDZtGgBXXnkl99xzDwCXXnopAPPnz6e8vPzgkaUrrriCN954g4svvpi4uLiD61155ZV84Qtf8OAn6JpjNP89sNwBN0UojkiPmBl5aX7y0vyUFWcdsmx/UytPf1DFI/PWceusRfzspRVcPa2YSycNJj1JDW3DTQXUMbS1OT7Z03DEnKR12+vYtKOeDiNu5KQGh9zOGpnL0NzU4NyknOCYdaTGq8cXZpCZnMDcFdtUQEWprhwpCpfDD/EfuJ+SEvwE3J35KhouEPFeUqKPK08dwuWTBzOnchsPvbmOO2ev4O5/rOLLk4q4dloJRVnJXseMWSqggN31zQfnJa3tcCRp/fa6Q8aekxJ8lOSkMLYgnYtOyqfkwNGknJSoqPZ9ccaZpbm8vrKGtjYX0qNb0vtt3LiRd955h9NOO40nnniC6dOn8+GHHx5cPmXKFG699Va2b99OZmYmTzzxBLfccgsAbW1tPPXUU1x22WX86U9/Yvr06V79GCJymLg44+xRAzh71ACWbt7NI/PW8f/e2cAf3l7PjLEDuW76UE4Zkul1zJjTZwqohuZWNu6o/3ROUodCaUdd08H1fHFGUWbwLLepw7IPHkkampvKgP79ov6Td8XIPJ5dVM2Szbs5qSjD6zgSRUaNGsUf/vAHvvGNbzBixAi++c1v8j//8z8Hlw8aNIi77rqLiooKnHNccMEFXHTRRUDwKNWyZcs45ZRTSE9P589//rNXP4aIHMPYgnR+fekEvjdjJH94Zz1/fHcDs5d8wsmDM7h++lDOGzOAeM2TCgnrzmH7niorK3MLFiwI2+u3tTmqd+8/4kjS2pp9bN61/5BTdXPT+nU4uy14JGlobgpFmckkxvfeX66ddU1M/O+/c+vZI/jOOaVexxFg+fLljBo1ytMM69ev58ILL2Tp0qUn9PzU1FT27dt33PU6+1nNbKFzruyE3jiKhHv/JRIOdY0tPP1BFQ/PW8eG2noKMpK4Zloxl04qIs3v/chJtDvW/qtXHoHaWdd0xJyktTV1rK+to7Gl7eB6KYk+SnJTmDg4ky9OLGRobgpDc1IpzkmO2V+czJRETi7KYE5ljQooEZE+LqVfPFedVswVU4bw2vKtPPTmOv77xeXc/Y9VXDa5iKunlVCQkeR1zF4paguohuZW1td2OJLUoRP3rvrmg+vFxxmDs5MZmpPCGaU5B48kDc1JITct+ofcwqE8kMf/94+VbN/X2H6arPR1xcXFJ3z0CejS0ScRiV6+OOPcMQM5d8xAPtq0i4fnreORt9bzyFvruWDcIK6fXqJpH90UlQXUnBXbuObR+Yc8NrC/n5KcFC4YN6h9TlJw2K0wM0l9Lw5TEcjj139fyRsra/jCxEKv4wjBM9xivZiP5HQAETlxJxVlcM9XTua280fyh7fX86f3N/L8R9VMKs7kuulD+czoAbq2XxdEZQFVOjCN755T2l4kBb9S+kVl1Kg0Jr8/Oan9mFupAioa+P1+amtryc7OjtkiyjlHbW0tfr/f6ygi0kX5GUn84IJR3HL2CJ6cv4lH3lrHjY8vZHBWMtdOK+aSsiL933sMUbllCjKSuPUcXY7kRMXFGeWBXP7+8VZa25w+SXissLCQqqoqampqvI4SVn6/n8JCFewivU1qv3iunV7C16YW8+qyT3ho3jp+9PzH/PrvK7l8yhC+NnUIg9I1T+pwUVlASc9VBPJ4amEVizbt5JQhWcd/goRNQkICJSUlXscQETkmX5xx/rhBnD9uEB9s3MnDb67jwTfW8NCba7lw/CCuP32oLlbfgQqoGDV9RA6+OGPOihoVUCIi0i0TB2cy8YpMNu2o59G31/Pn+Zt4dlE1pw7N4vrpQzlrZF6fb9as2dcxKj0pgVMGZzKncpvXUUREpJcqykrmhxeO5u0fnMX/nTmKTTv2c/1jCzj/7jfZtqfB63ieUgEVw8pH5rKsek+f/yUXEZGe6e9P4PrTh/L6v5Vz92UT2LSznmv/MJ+6xhavo3lGBVQMqwjkATB3ZWxPXhYRkciI98Vx0YQC7rtiIsu37OWmP31AS2vb8Z8Yg1RAxbCRA9MY2N/PXA3jiYhICFUE8vivi8Yyt7KGH/5taZ/sA6dJ5DHMLNjO4MXFW2hubVPDURERCZnLpwxm86567puzhsLMZG6qGO51pIjS/6gxrjyQx97GFhZu2Ol1FBERiTH/em6Az59cwC9eqeSvH1Z5HSeiVEDFuGnDs0nwGXMrNQ9KRERCy8z42RfHc9rQbL731GLeXr3d60gRowIqxqX5E5hUnKV5UCIiEhaJ8XE88NVTKMlJ4RuPL6Tyk71eR4oIFVB9QHkglxWf7KV6136vo4iISAxKT0rg99dMJinBxzW/f5+tfaB9jgqoPuBgOwMN44mISJgUZCTxyNWT2L2/mZv++EHMn5mnAqoPGJ6XSkFGkrqSi4hIWI0tSOf2z45hwYad/G1RtddxwkoFVB9gZlSMzOXt1dtpbGn1Oo6IiMSwL51SyPjCdO56aXlMdypXAdVHVATyqGtqZcF6tTMQEZHwiYszbv/sGLbuaeS3c1d7HSdsVED1EacNyybRF8ecFRrGExGR8DplSCYXT8jnd2+uY2NtvddxwkIFVB+RnBjPlKFZmgclIiIRcdv5o4iPM34y+2Ovo4SFCqg+pCKQx5qaOjbtiM1PAyIiEj0Gpvu5qWI4ryzbylsx2GBTBVQfUjHyQDsDHYUSEZHwu256CYOzkrnj+WW0tLZ5HSekelRAmdl3zWyZmS01syfMzB+qYBJ6JTkpFGcnM0f9oEREJAL8CT7+Y+YoVm7dx+PvbvA6TkidcAFlZgXAt4Ey59xYwAdcFqpgEh7lgTzeXrOdhma1MxARkfA7d/QApg/P4dd/X8mOuiav44RMT4fw4oEkM4sHkoHY7poVA8oDuTQ0t/Hu2lqvo4iISB9gZtz+2dHUNbXyq1crvY4TMidcQDnnNgO/BDYCW4DdzrlXQxVMwuPUodn4E+J0WReJaWY2w8wqzWy1md3WyfJ0M3vezD5qn4ZwjRc5RfqKEQPS+OqpQ3ji/Y18XL3H6zgh0ZMhvEzgIqAEyAdSzOzKTta7wcwWmNmCmhr9p+01f4KPqcNyNJFcYpaZ+YD7gPOB0cBXzGz0YavdBHzsnDsJKAd+ZWaJEQ0q0sd895xS0pMSuOP5ZTFxnbyeDOGdA6xzztU455qBZ4Cph6/knHvQOVfmnCvLzc3twdtJqFQEcllfW8+67XVeRxEJh8nAaufcWudcEzCL4Ie9jhyQZmYGpAI7gNi95oRIFEhPTuBfzwvw3rodvPrxVq/j9FhPCqiNwKlmlty+EzobWB6aWBJO5YFgOwN1JZcYVQBs6nC/qv2xju4FRhGct7kEuNU5F1vnWItEocsmDaYgIykmzsjryRyo94CngA8I7oDigAdDlEvCqCgrmWG5KepKLrHKOnns8PGC84BFBKcfTADuNbP+R7yQpiCIhJQvzrikrJB5q7f3+qbOPToLzzl3u3NupHNurHPuq865xlAFk/CqCOTx3rod1Ddp1EJiThVQ1OF+IUeeIXwN8IwLWg2sA0Ye/kKagiASepeUBf88/7KwyuMkPaNO5H1Uxcg8mlraeGeN2hlIzJkPjDCzkvaJ4ZcBzx22zkaC0w4wswFAAFgb0ZQifVRBRhKnj8jlqQWbaG3rvZPJVUD1UWXFmaQk+jSMJzHHOdcC3Ay8QnBe5pPOuWVmdqOZ3di+2n8BU81sCfAa8H3nXOxdrEskSl1aVkT17gbm9eJr5MV7HUC80S/ex9ThOcxZUYNzjuB5ACKxwTk3G5h92GMPdLhdDZwb6VwiEnTO6DwykxN4cv4mziztncPjOgLVh1UE8ti8az+rt+3zOoqIiPQh/eJ9fP7kQl79+JNee3kXFVB9WHkgWPWrK7mIiETapZOKaG51PPNB75xMrgKqD8vPSGLkwDTNgxIRkYgLDExjQlEGTy7Y1Cs7k6uA6uPKA3nMX7+DvQ3NXkcREZE+5tJJRazcuo9Fm3Z5HaXbVED1ceWBXJpbHW+tVjsDERGJrAvHDyIpwceTCzYdf+UoowKqjztlSCZp/eJ1cWEREYm4NH8CM8cP4rlF1dQ19q7Gziqg+rgEXxynl+Ywt7KmV45Bi4hI73bppCLqmlp5cckWr6N0iwoooTyQxyd7GljxyV6vo4iISB9TNiSTobkpPDm/dw3jqYASytubmOlsPBERiTQz49KyIhZs2Nmr+hKqgBLy+vsZk9+fuSvUD0pERCLvCxMLiY8z/tKLJpOrgBIg2JV84cad7N6vdgYiIhJZuWn9OGtkHk9/UEVza5vXcbpEBZQAUDEyl9Y2x7xVvffCjiIi0ntdOqmI7fuaeG1575hOogJKAJhQlElGcoLmQYmIiCfOLM0lL61fr+kJpQJKAPDFGaePyGVuZQ1tbWpnICIikRXvi+NLpxQyt3IbNXsbvY5zXCqg5KCKQC7b9zWyrHqP11FERKQPOm/MQNocvLM2+q+OoQJKDjqjNBcz1JVcREQ8MbYgnTR/PG+vjv75uCqg5KCc1H6ML8zQPCgREfGEL844dWg2b6/RESjpZSoCuXy4aRc76pq8jiIiIn3Q1GHZbNxRz6Yd9V5HOSYVUHKI8kAezsGbq9RUU0REIm/a8BwA3onyo1AqoOQQ4wvSyU5JZM4KDeOJiEjkjchLJSc1kbfXRPc8KBVQcoi4OOPM0lzeWLWdVrUzEBGRCDMzThuWw1tranEuev8fUgElRygfmceOuiYWV+3yOoqIiPRB04ZlU7O3kTU10XtxYRVQcoQzRuQQZzCnUvOgREQk8qYOC86Diuaz8VRAyREykhM5eXCm+kGJiIgnBmcnU5iZxFtR3A9KBZR0qiKQy+Kq3b2inb6IiMSeqcOyeXftjqidj6sCSjpVHsgD4I2VGsYTEZHImzosh937m1m+JTovL6YCSjo1Jr8/eWn91JVcREQ8MXVYNkDUDuOpgJJOmRnlgVzeWFlDS2ub13FERKSPyevvZ3heatROJFcBJUdVHshjT0MLH27a5XUUERHpg6YOy2b++h00tUTfB3kVUHJU00fk4IszdSUXERFPTB2WQ31TKx9FYV9CFVByVP39CZQNyWSu+kFJL2NmM8ys0sxWm9ltR1mn3MwWmdkyM3s90hlF5PhOHZqFGby9OvqG8VRAyTFVjMzj4y17+GR3g9dRRLrEzHzAfcD5wGjgK2Y2+rB1MoDfAp9zzo0BLol0ThE5vozkRMbmp0fldfFUQMkxVbS3M3h9pYbxpNeYDKx2zq11zjUBs4CLDlvncuAZ59xGAOecfsFFotTUYdl8uHEX+5tavY5yCBVQckylA1IZlO5nzgoN40mvUQBs6nC/qv2xjkqBTDOba2YLzeyqiKUTkW45bVg2Ta1tLNiww+soh+hRAWVmGWb2lJmtMLPlZnZaqIJJdAi2M8hj3urtUXkWhEgnrJPHDm9lHA+cAswEzgN+aGalR7yQ2Q1mtsDMFtTU6EOEiBcml2QRH2dR186gp0eg7gZeds6NBE4Clvc8kkSbikAu+xpbWLhhp9dRRLqiCijqcL8QqO5knZedc3XOue3AGwT3YYdwzj3onCtzzpXl5uaGLbCIHF1yYjwnD87g7ShrqHnCBZSZ9QfOAB4GcM41Oed2hSiXRJFpw3NI8JkuLiy9xXxghJmVmFkicBnw3GHr/A043czizSwZmII+AIpErdOG5bBk825272/2OspBPTkCNRSoAX5vZh+a2UNmlhKiXBJFUvrFM7kkS5d1kV7BOdcC3Ay8QrAoetI5t8zMbjSzG9vXWQ68DCwG3gcecs4t9SqziBzbtGHZtDl4f130zIPqSQEVD0wE7nfOnQzUAUf0W9EcgthQEchj5dZ9VO2s9zqKyHE552Y750qdc8Occz9pf+wB59wDHdb5hXNutHNurHPuN56FFZHjmjA4A39CXFRdF68nBVQVUOWce6/9/lMEC6pDaA5BbChvb2egppoiIhJp/eJ9TCrO4p0omkh+wgWUc+4TYJOZBdofOhv4OCSpJOoMy02hKCtJBZSIiHhi6rAcKrfupWZvo9dRgJ6fhXcL8EczWwxMAO7scSKJSmZGRSCPt1Zvp7ElupqZiYhI7Js6LBuAd9ZGx1GoHhVQzrlF7cNz451zFzvndJ57DCsP5LK/uTWqJvGJiEjfMLYgnTR/PO9EyWVd1Ilcuuy0oTkkxsepK7mIiEScL86YXJzF/PXRcaxGBZR0WVKij9OGZqsflIiIeGJ0fn/Wba+LiqkkKqCkWyoCuazdXseG2jqvo4iISB9TOiCN1jbHmm3e/x+kAkq6Re0MRETEKyMHpgGwcutej5OogJJuKs5JoSQnRV3JRUQk4opzUkjwGSs+UQElvVB5IJd31tTS0Oz9GLSIiPQdCb44huWmUvnJHq+jqICS7qsI5NHY0hY1vThERKTvCAxMY+XWfV7HUAEl3Te5JIukBB9zV2gYT0REIiswMI3Nu/azp6HZ0xwqoKTb/Ak+pg3PZk5lDc45r+OIiEgfcmAi+SqPJ5KrgJITcmYgj4076lm73ftTSUVEpO8oHRAsoLyeSK4CSk5IeWkuAHM0jCciIhFUkJFEar94KlVASW9UlJXMiLxUXl+pflAiIhI5ZkbpgFQVUNJ7VYzM4721O6hrbPE6ioiI9CGBgf2p3LrX03m4KqDkhJUHcmlqbePtNWpnICIikRMYkMqu+mZq9jZ6lkEFlJywsiFZpCT61JVcREQiKjCwP+DtRHIVUHLCEuPjmD4ih7krtqmdgYiIREygvZWBl/OgVEBJj1QE8qje3cCqbd53hRURkb4hKyWR3LR+VHrYC0oFlPRIeSAPUDsDERGJrJED03QESnqvgel+Rg3qr3lQIiISUaUD0li1bS+tbd5MIVEBJT1WHshlwfqdnl+XSERE+o7AwDQamtvYuKPek/dXASU9VhHIo6XN8daq7V5HERGRPiIwwNuJ5CqgpMcmDs4gzR/P3Ep1JRcRkcgoHZCGmQoo6cXifXGcUZrLnEq1MxARkchISvQxJCuZyq17PHl/FVASEuWluWzb28jHW7z5RRYRkb6ndIB3Z+KpgJKQODOQC6BhPBERiZiRA9NYX1tPQ3NrxN9bBZSERF6an3EF6eoHJSIiEVM6MI3WNseamsg3c1YBJSFTEcjlg4072V2vdgYiIhJ+Iz28pIsKKAmZ8pF5tDl4Y5WG8cRbZjbDzCrNbLWZ3XaM9SaZWauZfSmS+UQkNIZkp5Doi1MBJb3bSYUZZCYnqCu5eMrMfMB9wPnAaOArZjb6KOv9DHglsglFJFQSfHEMy0v15Jp4KqAkZHxxxhmlubxeWUObR631RYDJwGrn3FrnXBMwC7iok/VuAZ4GVPGL9GJeXRNPBZSEVEUgj9q6JpZs3u11FOm7CoBNHe5XtT92kJkVAJ8HHohgLhEJg9IBaWzZ3cDu/ZGdf6sCSkLqjNJczNTOQDxlnTx2+CHR3wDfd84d89xnM7vBzBaY2YKaGv1Oi0SjAxPJV0Z4GE8FlIRUVkoiE4oyNA9KvFQFFHW4XwhUH7ZOGTDLzNYDXwJ+a2YXH/5CzrkHnXNlzrmy3NzcMMUVkZ4obS+gVkR4GE8FlIRceWkeH1XtonZfo9dRpG+aD4wwsxIzSwQuA57ruIJzrsQ5V+ycKwaeAr7lnHs24klFpMfy0/2k9YtnpQoo6e0qRubi1M5APOKcawFuJnh23XLgSefcMjO70cxu9DadiISamVHqwUTy+Ii+m/QJY/PTyUlNZG5lDZ8/udDrONIHOedmA7MPe6zTCePOuasjkUlEwicwMI0XF2/BOYdZZ9MgQ09HoCTk4uKMM0vzeH1lDa1qZyAiImE2cmAau/c3s3VP5KaO9LiAMjOfmX1oZi+EIpDEhoqRueyqb2bRpl1eRxERkRhXOqD9ki4RPBMvFEegbiU4z0DkoNOH5xJnMFdn44mISJgFDhRQn+yJ2Hv2qIAys0JgJvBQaOJIrEhPTuCUIZlqZyAiImGXmZJIXlq/iLYy6OkRqN8A3wPaeh5FYk15II+lm/ewbW+D11FERCTGBQamRbSZ5gkXUGZ2IbDNObfwOOupk28fVRHIA+B1dSUXEZEwCwxIY9XWfRE7eaknR6CmAZ9r7+Q7CzjLzB4/fCV18u27Rg1KY0D/frqsi4iIhF1gYBqNLW1sqK2LyPudcAHlnPuBc66wvZPvZcA/nXNXhiyZ9HpmRnlpHm+sqqG5VaO8IiISPiMH9geIWENN9YGSsKoYmcvehhY+2LDT6ygiIhLDhuelArBq276IvF9ICijn3Fzn3IWheC2JLdOG5xAfZ7qsi4iIhFVSoo/ctH5U79ofkffTESgJqzR/AmML0nl/3Q6vo4iISIzLz0hiswooiRWTS7L4aNNuGppbvY4iIiIxrCDDrwJKYsek4iyaWtv4SJd1ERGRMMpPT6J6136cC38rAxVQEnZlQzIBmL9ew3giIhI++RlJNDS3sbO+OezvpQJKwi4zJZHAgDTeX68z8UREJHzyM5IAIjKRXAWURMSkkkw+2LAzYh1iRUSk7yloL6AiMQ9KBZRExKTiLPY1trB8S+SulC0iIn1LfoYf0BEoiSGTS7IA1M5ARETCJislEX9CnAooiR2D0pMozEzSRHIREQkbMyM/I4nqXQ1hfy8VUBIxk4uzeH/djoicXioiIn1TQYSaaaqAkoiZXJJFbV0Ta7dH5krZIiLS9xzoBRVuKqAkYia1z4Oar3lQIiISJvkZSWzb20hjS3ivfqECSiJmaE4KOamJvK95UCIiEiYHzsTbursxrO+jAkoixswoG5KlieQiIhI2B3pBVe2qD+v7qICSiJpUksWmHfvZsjsyF3sUEZG+5dNu5OE9E08FlETU5GL1gxIRkfAZmB6ZZpoqoCSiRg1KI7VfvIbxREQkLPwJPnJS+6mAktgS74tj4pBM5q/ThYVFRCQ8CjL8Ye8FpQJKIm5ycSaVW/eyq77J6ygiIhKDgt3IVUBJjJnUPg9qwXodhRIRkdA7cDmXcF75QgWURNxJRRkk+uI0D0pERMIiPyOJ/c2t7KpvDtt7qICSiPMn+BhfmM57OhNPwsTMZphZpZmtNrPbOll+hZktbv9628xO8iKniITHgV5Q4ZwHpQJKPDG5JIulm3dT39TidRSJMWbmA+4DzgdGA18xs9GHrbYOONM5Nx74L+DByKYUkXAqONgLSgWUxJhJJVm0tDkWbdzldRSJPZOB1c65tc65JmAWcFHHFZxzbzvnDkzCexcojHBGEQmjA5dzUQElMeeUIZmYoeviSTgUAJs63K9qf+xorgNe6myBmd1gZgvMbEFNTU0II4pIOGWlJNIvPo7q3eHrRq4CSjzR35/AqIH9NZFcwsE6eazTU3HMrIJgAfX9zpY75x50zpU558pyc3NDGFFEwsnMKMhI0hwoiU2TS7L4YMMumlvbvI4isaUKKOpwvxCoPnwlMxsPPARc5JyrjVA2EYmQcPeCUgElnplUnMX+5laWbt7tdRSJLfOBEWZWYmaJwGXAcx1XMLPBwDPAV51zKz3IKCJhlp/hZ/NOFVASgyaVZAJoGE9CyjnXAtwMvAIsB550zi0zsxvN7Mb21f4TyAZ+a2aLzGyBR3FFJEzyM5LYtreRxpbWsLx+fFheVaQL8tL8lOSk8P66ndxwhtdpJJY452YDsw977IEOt68Hro90LhGJnPz2VgZbdzcyODs55K+vI1DiqUnFmSzYsIPm1jZa21ynX+FsxS8iIrEp3M00dQRKPDW5JJsnF1Qx4j86PYscgAlFGfz1W1Mx6+zkKhERkSPlh7mZpgoo8dTMcYPYUddIQ3PnZ+Kt2raP5z+qZvmWvYzO7x/hdCIi0lsNSg9vM00VUOKppEQfN5wx7KjLa/c1MnvJFl5cUq0CSkREusyf4CMntR/Vu8NTQGkOlES17NR+nDY0mxcXb9FcKBER6ZaCDD+bd4WnG7kKKIl6M8cPYn1tPR9v2eN1FBER6UXC2UxTBZREvfPGDMQXZ7y4eIvXUUREpBc5UECFYwTjhAsoMysyszlmttzMlpnZraEMJnJAVkoiU4dlM3uJhvFERKTr8jOSqG9qZff+5pC/dk+OQLUA/8c5Nwo4FbjJzEaHJpbIoWaOCw7jLavWMJ6IiHRNQUbwTLxw9II64QLKObfFOfdB++29BC+ZUBCqYCIdnXtgGG+JhvFERKRrPu0FFfqJ5CFpY2BmxcDJwHuheD0Afj/zyMfGXAyTvw5N9fDHS45cPuFyOPkKqKuFJ686cvmka2HsF2F3FTzzjSOXT70ZAufD9lXw/HeOXH7Gv8KwCtiyGF7+wZHLz/5PGDwFNr4Hr/34yOUz7oJB42HNHHjjl0cu/+xvIGcEVL4Eb9975PIv/C+kF8LSp2H+I0cu//JjkJINH/4RFv3pyOVX/AUSk+H938GyZ49cfs2Lwe9v3QMrXzl0WYIfrnw6ePv1n8Pa1w9dnpwJlz4evP2PH8Gm+Ycu758PX/xd8PZLt8EnSw5dnj0MPndP8PZz34baNYcszho4jqnDPs/sJVv4Xt2vsD3Vhz6/aBKc86Pg7T9fCfU7D10+9Ew483vB249/EZoP+2MqPQ+mfTt4W797Ry7v+Ls39otHLhcRiUIHCqjNO+tD/to9nkRuZqnA08B3nHNHjK+Y2Q1mtsDMFtTU1PT07aQPmzluEBtq68Myli0iIrEnOyWRxPg4qneH/giU9WRSrpklAC8Arzjnfn289cvKytyCBbrouZyYnXVNlP3kH9xwxlC+P2Ok13Gki8xsoXOuzOscPaX9l0jvVPHLuYzO7899l0/s9nOPtf/qyVl4BjwMLO9K8STSU5kpiUwbnqOmmiIi0mX5Gf6w9ILqyRDeNOCrwFlmtqj964IQ5RLp1MxxA9m4o56lm3U2noiIHF9+eniaafbkLLx5zjlzzo13zk1o/5odynAihzt39EDidTaeiIh0UX5GEtv2NtLU0vlF60+UOpFLr5KZksjU4Tm8uKRaw3giInJcBZlJOAdb94R2IrkKKOl1Lhw3iE079msYT0REjqvgQCuDEA/jqYCSXufcMQOIjzNeWFJ9/JVFRKRP+7SZpgoo6eMykrt2Nt7KrXu566XlbN/XGMF0IiISTQalBy/nogJKBJg5fhBVO/ezZPPuI5btbWjmv1/4mAvufpP/fX0tNzy2gIbmVg9SioiI1/wJPnJSE9kc4su5hORSLiKRdu7oAfx7nPHi4i2ML8wAwDnHcx9V85MXl1Ozr5HLJhVxUmEGtz2zhB88s4Rff/kkgu3Loptzjo+qdjMmvz8JPn3GiWbNzc1UVVXR0BD6Lsexyu/3U1hYSEJCgtdRpA/Jzwh9KwMVUNIrZSQnMn1EDi8u2cJt549k5dZ9/OfflvLeuh2ML0znwavKmFCUAUDN3kZ+9feVDM9L5aaK4d4GP47WNscdzy/jsXc2cNVpQ/jxRWO9jiTHUFVVRVpaGsXFxb2iOPeac47a2lqqqqooKSnxOo70IfnpSayp2RfS19THW+m1LhgXHMa75YkPueCeN6ncupc7Pz+Ov35r2sHiCeDms4bzuZPy+cUrlby8NHr7RzU0t3Lznz7gsXc2UDoglcfe2cDba7Z7HUuOoaGhgezsbBVPXWRmZGdn64idRNyBI1ChbH+jAkp6rfNGDyTBF2yqeemkIub8n3IunzIYX9yh/5mZGT//0ngmFGXw3T9/xNJO5k15bff+Zq565H1eWvoJ/3fmKP5203RKclL43lOLqWts8TqeHIOKp+7R9hIv5Gf4qWtqDenF6FVASa+VnpzA76+ezPM3T+fOz48jMyXxqOv6E3w8eNUpZCYn8PXHFrAtxA3VeuKT3Q18+YF3+HDjTu75yslcf/pQkhJ9/OJL49m8az93vbTc64gSxXw+HxMmTGDs2LFccskl1NfXH1zW2trKySefzIUXXnjIcx599FHWr19/yKfxwx/7y1/+wpgxY4iLi+PwiyjfddddDB8+nEAgwCuvvHLw8YULFzJu3DiGDx/Ot7/9bTW7lagRjl5QKqCkV5s+IoexBeldWjcvzc/vvlbGrvpmbvh/C6PizLzV2/byhd++xeZd+3n0msl87qT8g8vKirO4bloJj7+7kbdXayhPOpeUlMSiRYtYunQpiYmJPPDAAweX3X333YwaNerg/c2bN3PdddexceNG5s2bx4033tjpYwBjx47lmWee4Ywzzjjk/T7++GNmzZrFsmXLePnll/nWt75Fa2vwb+mb3/wmDz74IKtWrWLVqlW8/PLLEdgCIsf3aS+o0H14VgElfcqY/HR+c9kEFm3axfeeWuzpJ+QF63fwxfvfobnN8edvnMq04TlHrPOv5wUYmpPCvz21mH0aypPjOP3001m9ejUQnOD+4osvcv311x9cXlBQwJ133skjjzzCrFmzuP/++zt9DGDUqFEEAoEj3uNvf/sbl112Gf369aOkpIThw4fz/vvvs2XLFvbs2cNpp52GmXHVVVfx7LPPRuTnFjmecDTT1Fl40uecN2Yg35sR4OcvV5Kc6CMwMC3iGfY1tHDvnNXkZyTx2LWTKcpK7nQ9f4KPX1wyni898A53zV7OTz4/7rivPX/9jmPO8zJgXGEGE4oyjpgvJifujueX8XF1aC8vNDq/P7d/dkyX1m1paeGll15ixowZAHznO9/h5z//OXv37j24TnV1NbfffjvXXnstJSUl3HTTTfzwhz884rEDRVRnNm/ezKmnnnrwfmFhIZs3byYhIYHCwsIjHheJBtkpiSTGx6mAEumpb545jI219cyav8mzDBMHZ/C7q8rITu13zPVOGZLF108fyoNvrOX8sYOYPuLII1UAW3bv579fWM6LS7p2pmF2SiLlgTzOGZXH9BE5pPnVl6c32r9/PxMmTACCR6Cuu+46XnjhBfLy8jjllFOYO3fuwXXz8/P53e9+x6OPPsrpp5/OlVdeiZkd8dixdHbU1syO+rhINIiLMwoykkI6B0oFlPRJZsZPvzieH1wwyrNhvPSkhC7/B/MvnynlH8u38v2nF/Pyd04/pNhpamnj4Xnr+J9/rqK1zfHdc0q54tTBxB/l6FJTSxvvrK3lnyu28Y/lW3n6gyoSfMapQ7M5a2Qe54wacNQjYr2Fmc0A7gZ8wEPOuZ8ettzal18A1ANXO+c+6Ml7dvVIUagdmAPV0VtvvcVzzz3H7NmzaWhoYM+ePVx55ZU8/vjjAFx99dVHvE5nj3WmsLCQTZs+/eBRVVVFfn4+hYWFVFVVHfG4SLTIz/DrCJRIqKQn9Y6jLv4EH7+85CS+dP/b3Dl7BXd9ITiU98bKGn703DLWbq/jM6MH8J8Xju5S8XPRhAIumlBAS2sbCzfs5LUV23ht+VbueP5j7nj+Y0bkpXL2qAGcMyqPkwdn9qqhPjPzAfcBnwGqgPlm9pxz7uMOq50PjGj/mgLc3/49Jtx1113cddddAMydO5df/vKXB4unnvrc5z7H5Zdfzr/8y79QXV3NqlWrmDx5Mj6fj7S0NN59912mTJnCY489xi233BKS9xQJhfz0JN5cFboTclRAifQSEwdn8vXTh/K/b6xl4uAM/rliGy8t/YTi7GR+f80kKgJ53X7NeF8cU4ZmM2VoNv9+wSjWb6/jH8u38s8V23jozbU88PoaMpMTqAjkcdaoPM4ozaV/9A/1TQZWO+fWApjZLOAioGMBdRHwmAsefnzXzDLMbJBzLno7rUbYX//6V2655RZqamqYOXMmEyZM4JVXXmHMmDF8+ctfZvTo0cTHx3Pffffh8/kAuP/++7n66qvZv38/559/Pueff77HP4XIp/Izkti6t4Hm1raQXCbLIjl8UVZW5g7vJyIiXdfQ3MrMe95kTU0d/oQ4bq4YzvWnD8Wf4Av5e+1paOaNlTW8tnwbcyq3sau+mfg44+xRefzvV8u6/DpmttA51/Un9JCZfQmY4Zy7vv3+V4EpzrmbO6zzAvBT59y89vuvAd93zh11B9XZ/mv58uWHtAmQrtF2Ey88OX8T33t6MW9+r6LL0xSOtf/SESiRXsSf4OO+Kyby5/mbuG56CYWZ4Zur1N+fwIXj87lwfD6tbY4PNu7kteXbekNzxM7GGw8P3ZV1MLMbgBsABg8e3PNkIuKZwMA0Zo4fRKh2YSqgRHqZkQO7fmp7qPjijEnFWUwqzoro+56gKqCow/1CoPoE1sE59yDwIASPQIU2pohE0klFGdx3+cSQvZ4aaYpIrJkPjDCzEjNLBC4DnjtsneeAqyzoVGC35j+JSHfoCJSIxBTnXIuZ3Qy8QrCNwSPOuWVmdmP78geA2QRbGKwm2Mbgmh68n/oddUMvGAIW6RIVUCISc5xzswkWSR0fe6DDbQfc1NP38fv91NbWkp2drSKqC5xz1NbW4vf7vY4i0mMqoERETtCB5pE1NTVeR+k1/H7/IZd8EemtVECJiJyghIQESkpKvI4hIh7QJHIRERGRblIBJSIiItJNKqBEREREuimil3IxsxpgQxdXzwFCd9W/0IrWbNGaC6I3W7TmgujN1t1cQ5xzueEKEynd3H9B7Pz7RVK0ZovWXBC92WIl11H3XxEtoLrDzBZE8vpZ3RGt2aI1F0RvtmjNBdGbLVpzRZto3U7RmguiN1u05oLozdYXcmkIT0RERKSbVECJiIiIdFM0F1APeh3gGKI1W7TmgujNFq25IHqzRWuuaBOt2ylac0H0ZovWXBC92WI+V9TOgRIRERGJVtF8BEpEREQkKnlSQJnZDDOrNLPVZnZbJ8v/zcwWtX8tNbNWM8vqynM9zLXezJa0L1sQylxdzJZuZs+b2UdmtszMrunqcz3M5fU2yzSzv5rZYjN738zGdvW5HuYK2zYzs0fMbJuZLT3KcjOze9pzLzaziV39mWJJtO6/QpAtnL9bUbn/CkE2L7eZ9l+Hvm/k91/OuYh+AT5gDTAUSAQ+AkYfY/3PAv88kedGKlf7/fVAjlfbDPh34Gftt3OBHe3rerrNjpYrSrbZL4Db22+PBF6Lht+zo+WKwDY7A5gILD3K8guAlwADTgXeC/f2iravaN1/9TRbOH+3onX/1dNsUbDNtP869H0jvv/y4gjUZGC1c26tc64JmAVcdIz1vwI8cYLPjVSucOtKNgekmZkBqQT/yFu6+FwvcoVbV7KNBl4DcM6tAIrNbEAXn+tFrrByzr1B8N/naC4CHnNB7wIZZjaI8G6vaBOt+6+eZgunaN1/9TRbOGn/1U1e7L+8KKAKgE0d7le1P3YEM0sGZgBPd/e5Ec4FwT+yV81soZndEKJM3cl2LzAKqAaWALc659q6+FwvcoH32+wj4AsAZjYZGAIUdvG5XuSC8G6z4zla9nBur2gTrfuvnmaD8P1uRev+q6fZwNttpv1X94R8/xUfsmhdZ508drRTAT8LvOWcO1BVdue53dWTXADTnHPVZpYH/N3MVrRXxJHKdh6wCDgLGNae4c0uPjfiuZxze/B+m/0UuNvMFhHcMX5I8JOl19vsaLkgvNvseI6WPZzbK9pE6/6ru68fyX1YtO6/epQtzPsw7b9CL+T7Ly+OQFUBRR3uFxKs7DtzGYceYu7OcyOZC+dcdfv3bcBfCR4WDJWuZLsGeKb98ORqYB3B8Wevt9nRcnm+zZxze5xz1zjnJgBXEZzfsK4rz/UoV7i32fEcLXs4t1e0idb9V0+zhfN3K1r3Xz3N5uk20/6r20K//3IhmLzVnS+CR73WAiV8OmFrTCfrpRMcz0zp7nM9yJUCpHW4/TYwI5LbDLgf+FH77QHAZoIXTfR0mx0jVzRssww+nQz6dYLj457/nh0jV1i3WfvrFnP0SZgzOXQS5vvh3l7R9tXD/URYt1MPs4Xtd6uH+wnPt9kxsnm9zY62n9D+q/NlId9/hSx4N3/IC4CVBGe+/0f7YzcCN3ZY52pgVlee63UugrP3P2r/WhbqXF3JBuQDrxI8ZLoUuDIattnRckXJNjsNWAWsAJ4BMqNkm3WaK9zbjOARiS1AM8FPZdcdlsuA+9pzLwHKIrG9ou3rRPcTkdhOJ5otAr9bUbn/6km2KNhm2n8dmivi+y91IhcRERHpJnUiFxEREekmFVAiIiIi3aQCSkRERKSbVECJiIiIdJMKKBEREZFuUgElIiIi0k0qoERERES6SQWUhJWZFZvZCjP7g5ktNrOn2i9kKiIS1bT/kmNRASWREAAedM6NB/YA3/I4j4hIV2n/JZ1SASWRsMk591b77ceB6V6GERHpBu2/pFMqoCQSDr9ekK4fJCK9hfZf0ikVUBIJg83stPbbXwHmeRlGRKQbtP+STqmAkkhYDnzNzBYDWcD9HucREekq7b+kU/FeB5A+oc05d6PXIUREToD2X9IpHYESERER6SZzTvPhRERERLpDR6BEREREukkFlIiIiEg3qYASERER6SYVUCIiIiLdpAJKREREpJtUQImIiIh00/8P0ot4jWXZ8MEAAAAASUVORK5CYII=\n", "text/plain": [""]}, "metadata": {"needs_background": "light"}, "output_type": "display_data"}], "source": ["import matplotlib.pyplot as plt\n", "fig, ax = plt.subplots(1, 2, figsize=(10, 4))\n", "df2.plot(x=\"p\", y=\"prop\", ax=ax[0])\n", "ax[0].plot([0.7, 1.], [2, 2], \"--\")\n", "df2.plot(x=\"p\", y=\"P4**100\", ax=ax[1]);"]}, {"cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [], "source": []}, {"cell_type": "code", "execution_count": 24, "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.7.2"}}, "nbformat": 4, "nbformat_minor": 2}