{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Ein einfaches epidemiologisches Modell\n", "_Dr. Thomas Viehmann _\n", "\n", "**Dies ist ein Entwurf, bitte nicht weiter verbreiten!**\n", "\n", "**This is a German draft, please do not copy!**\n", "\n", "\n", "Epidemiologie ist die Wissenschaft von der Ausbreitung von Krankheiten. Wie in vielen wissenschaftlichen Disziplinen gibt es mathematische Modelle, mit denen man versuchen kann, die Ausbreitung zu beschreiben.\n", "Solche Modelle vereinfachen notgedrungen die Realität, dass heißt, es gibt einen _Modellfehler_, und wir kennen die einfließenden Größen nicht genau, sie unterliegen also einem _Schätzfehler_.\n", "\n", "_Hinweis:_ Ich bin kein Epidemiologe, sondern Mathematiker und Modellierer. Ich hoffe, mit diesem Notebook ein wenig zum Verständnis für Ausbreitung von Krankheiten wie CoViD-19 beizutragen. Insbesondere sind meine Annahmen für die Parameter \"aus der Luft gegriffen\" und die Erklärungen allgemeinverständlich gehalten.\n", "\n", "\n", "Eines der einfachsten Modelle ist das sogenannte [SIR-Modell](https://en.wikipedia.org/wiki/Compartmental_models_in_epidemiology) (englisch Susceptible to Infectious to Recovered) dass letzlich auf Kermack und McKendrick zurückgeht (eine Arbeit von 1927). Es trifft auf Krankheiten zu, gegen die der Körper eine Immunität besitzt, wenn man sie einmal hatte oder geimpft ist.\n", "\n", "Wir stellen uns vor, dass man die Bevölkerung in drei Gruppen einteilen kann:\n", "- Ansteckbar (also noch gesund),\n", "- Infiziert (also krank oder ansteckend, auch ohne, dass die Krankheit ausgebrochen ist),\n", "- Genesen bzw. Geimpft (und damit nicht mehr ansteckbar und nicht mehr gesund).\n", "\n", "Jeder Mensch ist zu jeder Zeit in einer der drei Gruppen. Ein Übergang ist möglich von Ansteckbar zu Infiziert sowie von Infiziert zu Genesen. Wir nehmen an, dass keine anderen Übergänge möglich sind. Inbesondere sind alle geimpften direkt zum Start geimpft.\n", "\n", "Wir nehmen an, dass es insgesamt $N$ Personen gibt (und $N$ so groß ist, dass wir in \"Teilen\" modellieren können).\n", "Nach den englischen Begriffen nennen wir $S$ die Anzahl der ansteckbaren Personen, $I$ die der Infizierten, und $R$ die der Genesenen. Diese sind Funktiontionen (also für jede Zeit $t$ definiert), die wir $S(t)$, $I(t)$ und $R(t)$ nennen. Die Gesamtzahl $N$ ist unabhängig von $t$.\n", "\n", "Ein Sachverhalt, den wir hier vernachlässigen, sind Geburten und Todesfälle (oder Umzüge etc.). Meist sind diese gering im Vergleich zur Dynamik zwischen zwischen $S$, $I$, und $R$. Damit wird unser Modell zu einem geschlossenen System, es gilt immer\n", "\n", "$$S(t) + I(t) + R(t) = N.$$\n", "\n", "Wir beschreiben die Geschwindigkeit der Änderungen mit Ableitungen, die wir mit $\\dot S$, $\\dot I$ und $\\dot R$ beschreiben. Aus der Gleichung für das geschlossene System erhalten wir die Gleichgewichtsbedingung\n", "\n", "$$\\dot S(t) + \\dot I(t) + \\dot R(t) = 0.$$\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Ansteckung\n", "\n", "Wir führen eine positive Zahl $\\beta > 0$ ein, die angibt, an wie viele Personen, eine infizierte Person Krankheitserreger in einem relevanten Maß weitergibt. \n", "\n", "Wenn wir davon ausgehen, dass alle Personen aus der Gesamtheit gleich wahrscheinlich der Empfänger ist, treffen die Erreger mit Wahrscheinlichkeit $S(t) / N$ eine ansteckbare Person, die dann in den infizierten Zustand wechselt. Das dies für $I(t)$ Personen gleichzeitig passiert, ist also die Übergangsrate -- mit der die Anzahl der ansteckbaren, gesunden abnimmt -- gerade $\\beta S(t)I(t) / N$, es gilt also \n", "\n", "$$\n", "\\dot S(t) = - \\frac{\\beta S(t)I(t)}{N}.\n", "$$\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Genesung\n", "\n", "Für die Genesung haben verschiedene Möglichkeiten. Wir könnten zum Beispiel eine Annahme darüber treffen, wie lange die Krankheit dauert. Wir machen es hier uns einfach, und führen eine Gesundungsrate $\\gamma$ ein, also der Anteil der infizierten Personen, die zu jedem Zeitpunkt gesunden. Diese erhöht die Zahl Genesenen, es gilt also\n", "\n", "$$\n", "\\dot R(t) = \\gamma I(t).\n", "$$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Gleichung für die Infizierten\n", "\n", "Durch Einsetzen der Gleichungen für Ansteckung und Genesung in die Gleichgewichtsbedingung erhalten wir\n", "\n", "$$\n", "\\dot I(t) = - \\dot S(t) - \\dot R(t) = \\left(\\beta \\frac{S(t)}{N} - \\gamma \\right) I(t).\n", "$$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Nebenbei: Kontinuierliche Verzinsung\n", "\n", "Wie so oft geben wir hier Raten als Ableitungen an, auch wenn wir vielleicht eher daran gewöhnt sind, Raten von einem Zeitschritt als zum nächsten Anzusehen.\n", "Wenn man das mit Verzinsungsprozessen vergleicht, ist dies hier das äquivalent zur kontinuierlichen Verzinsung.\n", "Dennoch haben die Koeffizienten $\\beta$ und $\\gamma$ als Einheit $\\frac{1}{Zeiteinheiten}$, also zum Beispiel Anzahl pro Tag.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Invarianz gegenüber der Gesamtzahl\n", "\n", "Wir können von absoluten Größen $S$, $I$, $R$, zu relativen Größen $s = S/N$, $i = I/N$ und $r = R/N$ übergehen,\n", "dann sind die gleichungen, dieselben wie für $N=1$.\n", "Aufgrund dieser Invarianz können wir auf $N=1$ übergehen und $S$, $I$ und $R$ als Anteile an der Gesamtbevölkerung ansehen. Das machen wir im Folgenden." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Typische Krankheitsdauer\n", "\n", "Wenn wir mit nur Kranken beginnen würden (also auch keine Kranken dazukämen), hätten wir die Gleichung\n", "\n", "$$\n", "\\dot I(t) = - \\gamma I(t),\n", "$$\n", "\n", "also die Gleichung für einen exponentiellen Zerfall. Deren Lösung ist\n", "\n", "$$\n", "I(t) = I(0) exp(- \\gamma t).\n", "$$\n", "\n", "Wir können das startend von einer infizierten Person graphisch darstellen:" ] }, { "cell_type": "code", "execution_count": 75, "metadata": {}, "outputs": [], "source": [ "import torch\n", "%matplotlib inline\n", "from matplotlib import pyplot\n", "import math" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(tensor(2.1212), tensor(0.6543))" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3yV5f3/8dcngySQkBASAiRA2EOWEBRRtA4Ud511jw61aq22tdbavfSnHda667a2ah0F90CqoqCEvSHssJKwE0YIfH5/nIPfFJJwgJzcSc77+Xjcj+Qe5z7vaMLn3Nd93ddl7o6IiMSuuKADiIhIsFQIRERinAqBiEiMUyEQEYlxKgQiIjEuIegABysrK8vz8/ODjiEi0qRMmTKlzN2za9rX5ApBfn4+hYWFQccQEWlSzGx5bfvUNCQiEuNUCEREYpwKgYhIjFMhEBGJcSoEIiIxLmqFwMyeMrMSM5tdy34zswfMrMjMZprZkGhlERGR2kXziuAZYHQd+08HeoaX64BHophFRERqEbVC4O6fABvqOORc4DkPmQRkmFmHaOVZvr6CX78xh12790TrLUREmqQg7xHkAiurrReHt+3HzK4zs0IzKywtLT2kNysqKefpz5bx78LiQ3q9iEhzFWQhsBq21ThLjrs/7u4F7l6QnV3jE9IHdFKfdgzt0oYHxi1ix67dh3QOEZHmKMhCUAx0qraeB6yO1puZGbef1pu1W3bwj0m1PmktIhJzgiwEY4Grwr2HhgOb3X1NNN9weLe2jOyZxUPji9i6Y1c030pEpMmIZvfRfwETgd5mVmxm3zKzG8zshvAhbwNLgCLg78CN0cpS3e2n9Wbjtl08NWFZQ7ydiEijF7XRR9390gPsd+CmaL1/bQbmZTD6iPb8/dMlXHlMFzJbtWjoCCIijUpMPln8w1N7sa2yiofHFwUdRUQkcDFZCHrmpHHh0Dyem7ic4o3bgo4jIhKomCwEALee0gsz+PMHC4OOIiISqJgtBB0zUrjm2Hxen7aKeWu2BB1HRCQwMVsIAG48oQdpSQnc++78oKOIiAQmpgtBestEbjqxB+MXlPL54rKg44iIBCKmCwHA1SPy6ZiezB/enseePTWOcCEi0qzFfCFITozn9tG9mb1qC2NnRG2ECxGRRivmCwHAuYNy6Z/bmvveW6AB6UQk5qgQAHFxxk/P6MuqTdt5+rNlQccREWlQKgRhI7pncXKfdjw8voj15TuDjiMi0mBUCKq584w+bNu1m/s/XBR0FBGRBqNCUE2PdmlcfnRnXvhiOQvWbg06johIg1Ah2Metp/QiNSmB3701l9AAqSIizZsKwT4yW7Xg+6f04tNFZfx3waHNjywi0pSoENTgyuFd6JbVit++NZddu/cEHUdEJKpUCGrQIiGOu87sy5LSCp6fqPmNRaR5UyGoxUl92nF8r2z+8uFCytSdVESaMRWCWpgZvzirH9srd/PH9xYEHUdEJGpUCOrQo10q1x6bz0uFK5lZvCnoOCIiUaFCcAC3nNyTtq2S+NXYORqdVESaJRWCA0hLTuSO0b2ZumITr09bFXQcEZF6p0IQgQuG5DG4UwZ3vzOPzdt3BR1HRKReqRBEIC7O+N3X+7OhopK/aLJ7EWlmVAgi1D83nSuGd+G5icuYs3pz0HFEROqNCsFB+OGo3rRp2YKf/2e2bhyLSLOhQnAQ0lsmcucZfZm6YhOvTCkOOo6ISL1QIThIFwzJZVh+G+5+Zx4bKiqDjiMicthUCA6SmfH78wawdUcVf3h7XtBxREQOmwrBIeiVk8Z3ju/GK1OKmbRkfdBxREQOiwrBIbrlpJ50ykzhrtdnUVmloapFpOlSIThEKS3i+c25/VlcWsFjHy8OOo6IyCFTITgMJ/Zux5kDO/C38UUsLi0POo6IyCGJaiEws9FmtsDMiszsJzXsTzezN8xshpnNMbNro5knGn55dj+SE+K487VZerZARJqkqBUCM4sHHgJOB/oBl5pZv30OuwmY6+6DgK8BfzKzFtHKFA3t0pK568y+fLl0Ay8Vrgw6jojIQYvmFcFRQJG7L3H3SuBF4Nx9jnEgzcwMSAU2AFVRzBQVFxd0Yni3TP7w9jxKtuwIOo6IyEGJZiHIBap/RC4Ob6vuQaAvsBqYBXzf3ffrgmNm15lZoZkVlpaWRivvITMz7j5/IDur9vCLMXOCjiMiclCiWQishm37NqKfBkwHOgKDgQfNrPV+L3J/3N0L3L0gOzu7/pPWg65Zrbj1lJ68O2ctb89aE3QcEZGIRbMQFAOdqq3nEfrkX921wGseUgQsBfpEMVNUXTeyGwNy0/nFmNkafkJEmoxoFoLJQE8z6xq+AXwJMHafY1YAJwOYWQ7QG1gSxUxRlRAfx30XDWTz9l385g01EYlI0xC1QuDuVcDNwHvAPOBld59jZjeY2Q3hw34LjDCzWcA44A53L4tWpobQp31rbjqxB/+ZvpoP564LOo6IyAGZe9Pq+15QUOCFhYVBx6hTZdUeznlwAhsqKvngthNIb5kYdCQRiXFmNsXdC2rapyeLo6BFQhx/vGgQGyoq+ZWaiESkkVMhiJL+uencfFIPXp+2indnrw06johIrVQIouimE3vQP7c1d70+i/XlO4OOIyJSozoLgZnFm9l9DRWmuUmMj+NPFw1my45d/HzMbJra/RgRiQ11FgJ33w0MDQ8BIYegd/s0bhvVi7dnrWXM9H0foxARCV5CBMdMA8aY2b+Bir0b3f21qKVqZq4/vjsfzSvh52NmM6xrJrkZKUFHEhH5SiT3CDKB9cBJwNnh5axohmpu4uOMP188mD17nB+9PEPDVYtIo3LAKwJ3b3JzBDRGndu25JdnH8GPX53JU58t5dsjuwUdSUQEiOCKwMx6mdk4M5sdXh9oZj+LfrTm56KCPEb1y+He9xYwf+2WoOOIiACRNQ39HbgT2AXg7jMJjRskB8nMuOf8AaSnJHLLv6axY9fuoCOJiERUCFq6+5f7bGtyk8c0Fm1Tk/jTRYNYuK6cu9+eF3QcEZGICkGZmXUnPJeAmV0IaMD9w3B8r2y+dVxXnp24nHHzNDCdiAQrkkJwE/AY0MfMVgG3At+NaqoYcPtpvenTPo3bX5mp6S1FJFCRFIJV7n4KkA30cffjAN3pPEzJifH87dIj2VZZxa0vTWe3upSKSEAiKQSvmVmCu1e4+1Yzaw98EO1gsaBnThq/PucIPl+8nkc/Xhx0HBGJUZEUgv8Ar4THHcoH3ifUi0jqwcUFnThnUEf+/MFCCpdtCDqOiMSgAxYCd/87oSuA/wBvADe4+/vRDhYrzIzfn9efvDYp3PKvaWzaprmORaRh1VoIzOwHexcgmdBE9NOB4eFtUk/SkhN58NIhlJVX8kMNQSEiDayuK4K0aksq8DpQVG2b1KMBeencdWZfxs0v4fFPlwQdR0RiSK1jDbn7rxsyiMBVx3Thy6UbuO+9BQzt0oZh+ZlBRxKRGBDpWEOPm9n7ZvbR3qUhwsUaM+OeCwbQqU0KN/9zKmWa1UxEGkAkvYb+TWhOgp8Bt1dbJArSkhN5+PKhbNq2i1v+NY2q3XuCjiQizVwkhaDK3R9x9y/dfcreJerJYli/jq353df78/ni9fzpg4VBxxGRZi6SQvCGmd1oZh3MLHPvEvVkMe6igk5celRnHvnvYt6bszboOCLSjEUyVeXV4a/Vm4Mc0MwqUfarc/oxd/VmfvTyDHrenEq37NSgI4lIMxTJA2Vda1hUBBpAUkI8D18xlMSEOK57fgrlOzX6t4jUv0iahjCz/mZ2sZldtXeJdjAJyc1I4cHLjmRpWQU/eGm6HjYTkXoXSffRXwJ/Cy8nAvcC50Q5l1QzonsWPz2jL+/PXceD44uCjiMizUwkVwQXAicDa8MT2Q8CkqKaSvbzzWPzOe/IXP7y4UI+nKvJbESk/kRSCLa7+x6gysxaAyXoRnGDMzPuPn8A/Tum8/0Xp7Fw3dagI4lIMxFJISg0swxCk9hPAaYC+85hLA0gOTGex68aSsukBL79bCEbKzRSqYgcvjoLgZkZcLe7b3L3R4FRwNXhJiIJQIf0FB67cihrt+zguy9MYZeePBaRw1RnIXB3JzQPwd71Ze4+M+qppE5DOrfhnvMHMGnJBn4xZg6h/00iIocmkqahSWY27FBObmajzWyBmRWZ2U9qOeZrZjbdzOaY2ceH8j6x6PwhedxwQnf+9eUKnpywNOg4ItKERfJk8YnA9Wa2HKgAjNDFwsC6XmRm8cBDhJqTioHJZjbW3edWOyYDeBgY7e4rzKzdIf4cMenHp/VmWVkFv397HvltW3FKv5ygI4lIExTJFcHpQHfgJOBs4CzgoghedxRQ5O5L3L0SeBE4d59jLgNec/cVAO5eEmlwgbg44y/fGEz/junc8uI05qzeHHQkEWmC6pqq8ucA7r68+gJsBB6L4Ny5wMpq68XhbdX1AtqY2X/NbEptTyyb2XVmVmhmhaWlpRG8dexIaRHPE1cXkJ6SyLeeKWTN5u1BRxKRJqauK4KRZvb76hvMrD3wKRDJxDRWw7Z972omAEOBM4HTgJ+bWa/9XuT+uLsXuHtBdnZ2BG8dW3JaJ/PUNcMo31nFtU9PZsuOXUFHEpEmpK5CcA4wyMz+DGBmPYEJwMPu/psIzl1MaML7vfKA1TUc8667V7h7GfAJoSeX5SD17dCaR64YQlFJOTe9MFXdSkUkYrUWAnffAZwHdDGzF4EPgdvdPZJmIYDJQE8z62pmLYBLgLH7HDOG0JVHgpm1BI4G5h3sDyEhI3tm84fzB/DpojJ+8uosdSsVkYjU2mvIzH4Q/vZL4MeEmoS67t3u7n+u68TuXmVmNwPvAfHAU+4+x8xuCO9/1N3nmdm7wExgD/CEu88+3B8qll1c0IlVG7fz13GLaJ+exO2n9Qk6kog0cnV1H02r9v0DNWw7IHd/G3h7n22P7rN+H3DfwZxX6nbrKT0p2bqDh8Yvpl1aMlePyA86kog0YrUWAnf/dUMGkfpjZvz23P6Ubq3kV2/MITstiTMGdAg6log0UhFNTCNNT0J8HH+79EiGdG7DrS9O57OisqAjiUgjpULQjKW0iOfJqwvomtWK654rZGbxpqAjiUgjpELQzGW0bMFz3zqKzNQWXPP0ZIpKyoOOJCKNTCRTVSaZ2WVm9lMz+8XepSHCSf3IaZ3M8988mjgzrnzyC1Zu2BZ0JBFpRCK5IhhDaIygKkKDzu1dpAnJz2rF8986ioqdVVzx5BeUbNkRdCQRaSTsQA8dmdlsd+/fQHkOqKCgwAsLC4OO0WRNW7GRK574go4ZKbx0/TFktmoRdCQRaQBmNsXdC2raF8kVwedmNqCeM0lAjuzchieuHsaKDdu48skv2LxN4xKJxLq6Rh+dZWYzgeOAqeEJZmZW2y5N1DHd2/LolUNZtK6cq57+kq0apE4kptX1ZPFZDZZCGtyJvdvx0OVD+O4/pnDN05N57ptH0SopknmKRKS5qWvQub3zD3QANlRb3wC0b6iAEj2j+uXwt0uPZPrKTVz7zGS2VVYFHUlEAhDJPYJHgOqdzyvC26QZOH1AB+7/xmAKl23gmqcnU7FTxUAk1kRSCMyrdS1y9z1ENtexNBFnD+rIXy85kinLN3KtioFIzImkECwxs1vMLDG8fB9YEu1g0rBCxWAwU1Zs5BrdQBaJKZEUghuAEcCq8HI0cF00Q0kwzhrYkQcuOZJpKzZx5ZNfsnm7ioFILDhgIXD3Ene/xN3bhZfL3L2kIcJJwztzYAcevnwIc1Zv5rK/T2JjRWXQkUQkyiIZayjPzF43sxIzW2dmr5pZXkOEk2CcekR7Hr+qgEUl5Vzy+CRKtmo4CpHmLJKmoacJzTXcEcgF3ghvk2bsxN7tePqa0BPIFz86keKNGqhOpLmKpBBku/vT7l4VXp4BsqOcSxqBY3tk8Y9vH82GikouenQii0s1hLVIcxRJISgzsyvMLD68XAGsj3YwaRyGdmnDi9cdw67de7j40YnMXrU56EgiUs8iKQTfBC4G1gJrgAvD2yRG9OvYmpevP4bkxHgueXwSny/WtJcizUkkvYZWuPs57p4d7jX09fBQExJDumWn8up3R9AxI5lrnprMu7PXBB1JROpJJL2GssOzkz1uZk/tXRoinDQu7dOTefn6Y+if25obX5jK85P0eUCkOYhkqIgxwKfAh8Du6MaRxi6jZQte+PZwbv7nVH7+n9ms3bydH53aGzMLOpqIHKJICkFLd78j6kmkyUhpEc9jVw7l52Nm89D4xazZvIN7zh9Ii4RIbjmJSGMTyV/um2Z2RtSTSJOSEB/HH84bwA9G9eK1qau49hkNSSHSVEVSCL5PqBhsN7MtZrbVzLZEO5g0fmbGLSf35I8XDeKLJRu46NHP9eCZSBMUSa+hNHePc/cUd28dXm/dEOGkabhwaB7PffMo1mzewXkPf86MlZuCjiQiB0GNulIvRvTI4rXvjiApIY6LH5vIWzPVvVSkqVAhkHrTMyeNMTcdy4DcdG7651T+Nm4R1eY0EpFGSoVA6lXb1CRe+M7RnH9kLn/6YCG3vDid7ZXqdSzSmNXafdTMMut6obtvqP840hwkJcTzp4sH0SMnlfveW8DSsnIev7KAjhkpQUcTkRrU9RzBFMCBmp4UcqBbVBJJs2Bm3Pi1HvTOSeP7L07nnAcn8MgVQxmWX+fnCxEJQK1NQ+7e1d27hb/uu0RUBMxstJktMLMiM/tJHccNM7PdZnbhofwQ0nid3DeH/9w0grTkRC59fBLPT1ym+wYijUythcDM+oS/DqlpOdCJzSweeAg4HegHXGpm/Wo57v8B7x3qDyGNW492afznpmM5oVc2Px8zh9tfmcmOXbpvINJY1NU09ANCk9T/qYZ9Dpx0gHMfBRS5+xIAM3sROBeYu89x3wNeBYZFEliapvSURP5+VQF/HbeIv45bxPy1W3jk8qF0ymwZdDSRmFdrIXD368JfTzzEc+cCK6utFwNHVz/AzHKB8wgVFRWCZi4uzrhtVC8G5qVz20vTOfOBT7n/ksGc1Ccn6GgiMS2SYahbmtnPzOzx8HpPMzsrgnPXdpO5uvuBO9y9znYCM7vOzArNrLC0tDSCt5bG7OS+Obz5vZHktWnJN58p5N5351O1e0/QsURiVqST11cCI8LrxcDvInhdMdCp2noesHqfYwqAF81sGaGZzx42s6/veyJ3f9zdC9y9IDtb0yU3B53btuS1G0fwjYJOPPzfxVz2xBes27Ij6FgiMSmSQtDd3e8FdgG4+3Zq/rS/r8lATzPramYtgEuAsdUPCPdAynf3fOAV4EZ3/8/B/ADSdCUnxvP/LhzIny8exKzizZzx10/5eKGu+EQaWiSFoNLMUgg365hZd2DngV7k7lXAzYR6A80DXnb3OWZ2g5ndcBiZpZk5f0geb3zvWLJSk7j6qS+5++15VFapqUikodiB+nSb2anAXYS6gL4PHAtc6+7jox9vfwUFBV5YWBjEW0uU7di1m9++OZcXvljBoLx0Hrj0SLq0bRV0LJFmwcymuHtBjfsiebjHzNoCwwk1CU1y97L6jRg5FYLm751Za7jj1Zns3uP86pwjuHBonqbCFDlMdRWCSHoNjXP39e7+lru/6e5lZjau/mOKhJw+oAPv3Ho8R+Smc/srM7n5n9PYtK0y6FgizVZdTxYnhweeyzKzNmaWGV7ygY4NFVBiU25GCv/6znB+PLo3781Zy+j7P+XTRbqRLBINdV0RXE9o4Lk+4a97lzGEho4Qiar4uNDAda/dOIJWSfFc+eSX/HLMbA1rLVLPIrlZ/D13/1sD5Tkg3SOITTt27ebedxfw1GdL6ZbVivsuGsTQLm2CjiXSZNTHzeIRQD7VhqRw9+fqK+DBUCGIbZ8XlXH7KzNZs3k73xnZjdtG9SI5MT7oWCKN3uHeLH4e+CNwHKHxgIYReiJYpMGN6JHFu7eO5BvDOvPYJ0s484FPmbJccySJHI5ImobmAf28kQwirysC2euThaXc+dosVm/eztXH5PPj0b1p2aKuAXVFYtdhXREAs4H29RtJ5PAd3yub9247niuHd+GZz5dx6l8+4RMNUSFy0CIpBFnAXDN7z8zG7l2iHUwkEqlJCfzm3P68fP0xtIiP46qnvuS2l6azvvyAo6CISFgkTUMn1LTd3T+OSqIDUNOQ1GbHrt08PL6IRz5eTKukBO48vQ8XDe1EXJyeShY5rKYhd/+4+gJUARfXd0iRw5WcGM8PTu3NW7eMpGe7VO54dRYXPzaR+Wu3BB1NpFGLpGkIMxtsZveG5w34HaHRREUapV45abx03THce+FAFpeWc+YDE/jdm3PZumNX0NFEGqVau1iYWS9CcwhcCqwHXiLUlHSoU1eKNJi4OOPigk6c0jeHe9+dz5OfLWXsjNX89Iy+nDu4owaxE6mmriuC+cDJwNnuflz46WI92y9NSmarFtxzwUBev/FY2qcnc+tL0/nGY5OYvWpz0NFEGo26CsEFwFpgvJn93cxOJrKZyUQancGdMnj9xmO5+/wBFJWWc/aDE7jztVmUqXeRSES9hloBXyfURHQS8Czwuru/H/14+1OvITlcm7fv4oFxi3j282WkJMZz00k9uPbYfJISNFSFNF+HPdZQtRNlAhcB33D3k+op30FRIZD6UlRSzh/ensdH80volJnCHaP7cOaADrp/IM1SvRWCxkCFQOrbp4tK+f1b85i/diuDO2Vw15l9GZafGXQskXp1uENMiDRrI3tm89YtI7n3goGs2bydix6dyLefLWThuq1BRxNpELoiEKlme+VunpywhMc+XkJFZRXnD8njtlG9yM1ICTqayGFR05DIQdpQUcnD44t4buJyAC47ujM3ntiddmnJAScTOTQqBCKHaNWm7fxt3CL+PaWYxHjj6hH5XH98dzJbtQg6mshBUSEQOUxLyyq4/8OFjJ2xmpTEeK4Zkc93RnajjQqCNBEqBCL1pKhkK38dV8SbM1fTMjGeK4/J59sju5KVmhR0NJE6qRCI1LOF67by4EdFvDFzNUkJcVx+dBeuO74bOa11D0EaJxUCkSgpKinnofFFjJ2xmngzLhiax3dP6E7nti2DjibyP1QIRKJsxfptPPrJYl4pLKZqzx7OGNCBG07oTv/c9KCjiQAqBCINZt2WHTw1YSkvfLGC8p1VjOyZxXdGdmNkzywNXSGBUiEQaWCbt+/ihS+W88xnyyjZupM+7dP49shunD2ogwa3k0CoEIgEZGfVbsZOX80Tny5lwbqtZKUmceXwLlw+vLN6GkmDUiEQCZi7M6GojKcmLGX8glJaxMdx1qAOXDMin4F5GUHHkxhQVyGodapKEak/ZsbIntmM7JlNUUk5z36+jFenFvPa1FUc2TmDq4/J5/QB7dVsJIGI6hWBmY0G/grEA0+4+z377L8cuCO8Wg58191n1HVOXRFIc7Flxy5eKSzm+UnLWVpWQdtWLbh4WCcuO6oznTLV/VTqVyBNQ2YWDywERgHFwGTgUnefW+2YEcA8d99oZqcDv3L3o+s6rwqBNDd79jifLS7j+YnL+XDeOpzQ0NiXHdWZk/u2IzFeo8XL4QuqaegooMjdl4RDvAicC3xVCNz982rHTwLyophHpFGKi/u/ZqPVm7bz0uSVvDR5JTf8YwrZaUlcODSPbxR0Ij+rVdBRpZmKZiHIBVZWWy8G6vq0/y3gnSjmEWn0OmakcNuoXnzvpB6MX1DKS5NX8NjHi3nkv4s5umsmFxd04vQB7WnZQrf3pP5E87eppqdnamyHMrMTCRWC42rZfx1wHUDnzp3rK59Io5UQH8eofjmM6pfD2s07eHVqMf8uXMkP/z2DX46dwxkD2nPBkDyG5WcSF6cH1eTwRPMewTGE2vxPC6/fCeDud+9z3EDgdeB0d194oPPqHoHEKndn8rKN/LtwJW/PWkNF5W46ZaZw3pF5fH1wR7plpwYdURqxoG4WJxC6WXwysIrQzeLL3H1OtWM6Ax8BV+1zv6BWKgQisK2yivfnrOPVqcV8VlTGHodBeemcOziXswZ10Exqsp/AHigzszOA+wl1H33K3X9vZjcAuPujZvYEcAGwPPySqtqC7qVCIPK/1m3Zwdjpq3l92irmrtlCnMGI7lmcPagDpx3RnoyWmjxH9GSxSMwoKtnK2OmrGTNjNcvXbyMx3jiuRxZnDezIKf1ySE9JDDqiBESFQCTGuDuzVm3mzZlreGvmGlZt2v5VUThjQAdO6ZujaTZjjAqBSAxzd6av3MTbs9bw9qy1rNq0nfg4Y3i3TEYf0Z5R/drTPl33FJo7FQIRAUJFYfaqLbw7Zw3vzF7LktIKAAZ1yuDUfjmc0jeHXjmpmjuhGVIhEJEaFZVs5b0563h/7jpmrNwEQF6bFE7pm8PJfdtxVNdMDYTXTKgQiMgBrduyg3HzShg3bx0TisrYWbWHVi3iOa5nFif2bsfXerdTE1ITpkIgIgdle+VuPl9cxkfzS/hofglrNu8AoE/7NE7onc0JPbMZmt9GVwtNiAqBiBwyd2fhunL+u6CE8QtKmLJ8I7t2Oy1bxHN018zwgHlZ9GinewuNmQqBiNSbip1VTFy8no8XljKhqIylZaEbzu3SkjiuRxYjemQxontbOmakBJxUqtMMZSJSb1olJXBKvxxO6ZcDwMoN25hQVMZnRWX8d2Epr01bBUB+25Yc070tw7uFlpzWur/QWOmKQETqzZ49zvy1W/l8cRmTlqzniyUb2LqzCoBuWa04ulsmR3XNZFh+JnltNAtbQ1LTkIgEomr3Huau2cIXSzYwacl6vly2ga07QoUhNyOFgvw2FORnMiy/DT3bpRGvIbWjRoVARBqF3XucBWu38sXS9RQu28iXyzZQunUnAGlJCRzZpQ1DO7dhSJcMBnfKIC1ZYyPVFxUCEWmU3J3l67cxZflGCpdvZMryDSwqKccdzKBXuzQGd8rgyM4ZDO6coauGw6BCICJNxpYdu5i+YhNTV2xk+spNTFuxic3bdwHQskU8/XPTGZSXzsC8DAbmpdM5s6W6rUZAvYZEpMlonZzI8b2yOb5XNhC6ali2fhvTVmxkZvFmZhRv4tmJy6msWv8fQ7UAAAn7SURBVApAekoi/XNb0z83nQG56fTvGCoOmsIzcioEItKomRlds1rRNasV5w/JA6Cyag8L121lZvFmZq3axKxVm3lqwlJ27Q61cKQmJdCvQ2v6dQwvHVrTMydVT0LXQoVARJqcFglx9M9Np39uOtAZgJ1Vu1m0rpw5qzczZ/UWZq/azMuFK9lWuRuAhDije3YqfTqk0ad96/DXNNq3To75piUVAhFpFpIS4qsVh5Dde5zl6yuYu2YL89ZsYf6arUxeuoEx01d/dUzr5AR6t0+jV05o6ZmTSq+cNLJSk4L4MQKhQiAizVZ8nNEtO5Vu2amcNbDjV9s3b9vFgnVbmb92CwvWbmXhuq2MnbH6q2ccADJbtaBHdio9clJDX9uFlg7pze8KQoVARGJOestEjuoaesp5L3dn3ZadLFwXKgyLS8tZtK6ct2au+arXEoR6LnXNakX37FS6Zbf66vv8rFakJjXNf1KbZmoRkXpmZrRPT6Z9evJXPZYgVCDWV1RSVFJOUUk5i0vLWVJawdQVG3lj5mqq98DPTkuia9tW5Ge1pEvbUJHonNmSLm1bNuqH41QIRETqYGZkpSaRlZrE8G5t/2ffjl27Wb5+G0tKy1lSVsGysgqWra/go/mllJUX/8+xbVu1oEvblnTODC2dqn3NaZ0c6INyKgQiIocoOTGe3u3T6N0+bb995TurWL6+gmVl21i+oYIV67exfP02CpdvZOyM1eypdiWRGG/kZqSQ16YleW1SyGuTQm6b0HpuRkrUC4UKgYhIFKQmJXBEx3SO6Ji+375du/ewauN2Vm7cxsoN21mxYRvFG7dRvHE7H85bR1l55f8cHx9ntG+dzDUj8vnO8d3qPasKgYhIA0uMjyM/qxX5Wa1q3L+9cjerNm0PLRu3szr8fbvW0enSqkIgItLIpLSI/6q7akOIa5B3ERGRRkuFQEQkxqkQiIjEOBUCEZEYp0IgIhLjVAhERGKcCoGISIxTIRARiXFNbvJ6MysFlh/iy7OAsnqME21NKW9TygpNK29TygpNK29TygqHl7eLu2fXtKPJFYLDYWaF7l4QdI5INaW8TSkrNK28TSkrNK28TSkrRC+vmoZERGKcCoGISIyLtULweNABDlJTytuUskLTytuUskLTytuUskKU8sbUPQIREdlfrF0RiIjIPlQIRERiXMwUAjMbbWYLzKzIzH4SdJ7amFknMxtvZvPMbI6ZfT/oTJEws3gzm2ZmbwadpS5mlmFmr5jZ/PB/42OCzlQXM7st/Hsw28z+ZWbJQWeqzsyeMrMSM5tdbVummX1gZovCX9sEmXGvWrLeF/5dmGlmr5tZRpAZq6spb7V9PzIzN7Os+nivmCgEZhYPPAScDvQDLjWzfsGmqlUV8EN37wsMB25qxFmr+z4wL+gQEfgr8K679wEG0Ygzm1kucAtQ4O79gXjgkmBT7ecZYPQ+234CjHP3nsC48Hpj8Az7Z/0A6O/uA4GFwJ0NHaoOz7B/XsysEzAKWFFfbxQThQA4Cihy9yXuXgm8CJwbcKYaufsad58a/n4roX+ocoNNVTczywPOBJ4IOktdzKw1cDzwJIC7V7r7pmBTHVACkGJmCUBLYHXAef6Hu38CbNhn87nAs+HvnwW+3qChalFTVnd/392rwquTgLwGD1aLWv7bAvwF+DFQbz19YqUQ5AIrq60X08j/cQUws3zgSOCLYJMc0P2EfjH3BB3kALoBpcDT4WasJ8ys5tnDGwF3XwX8kdAnvzXAZnd/P9hUEclx9zUQ+mADtAs4T6S+CbwTdIi6mNk5wCp3n1Gf542VQmA1bGvU/WbNLBV4FbjV3bcEnac2ZnYWUOLuU4LOEoEEYAjwiLsfCVTQeJot9hNuWz8X6Ap0BFqZ2RXBpmqezOwuQs2yLwSdpTZm1hK4C/hFfZ87VgpBMdCp2noejewSuzozSyRUBF5w99eCznMAxwLnmNkyQk1uJ5nZP4KNVKtioNjd915hvUKoMDRWpwBL3b3U3XcBrwEjAs4UiXVm1gEg/LUk4Dx1MrOrgbOAy71xP1jVndCHghnhv7c8YKqZtT/cE8dKIZgM9DSzrmbWgtANt7EBZ6qRmRmhNux57v7noPMciLvf6e557p5P6L/rR+7eKD+1uvtaYKWZ9Q5vOhmYG2CkA1kBDDezluHfi5NpxDe3qxkLXB3+/mpgTIBZ6mRmo4E7gHPcfVvQeeri7rPcvZ2754f/3oqBIeHf68MSE4UgfDPoZuA9Qn9IL7v7nGBT1epY4EpCn6ynh5czgg7VjHwPeMHMZgKDgT8EnKdW4SuXV4CpwCxCf6+NakgEM/sXMBHobWbFZvYt4B5glJktItS75Z4gM+5VS9YHgTTgg/Df2qOBhqymlrzRea/GfSUkIiLRFhNXBCIiUjsVAhGRGKdCICIS41QIRERinAqBiEiMUyGQZs3MzqvWDXfvssfMTj/A6z4Pf803s8vqOK6Dmb1pZqdVO395eKTb6Wb2XD3/PPeb2fH1eU4RdR+VmGJm1wGXAye6+wHHRjKzrwE/cvezatl/HzDB3cdU2/bf8GsK6yX0/75fd+BBd6+zkIkcDF0RSMwws16Exmm5cm8RMLPbzWxyeDz6X1c7tjz87T3AyPCn+9tqOO0FwLsHeN/uZvZpeKC7KWZ2dHh7vJk9Gp5v4A0ze9fMvh7eN8zMPg4f/46Z5QC4+2Kgg5llH+Z/DpGvqBBITAiP3/RPQp/UV4S3nQr0JDRM+WBgaA3NLj8BPnX3we7+l33O2RXY6O47D/D2a4BR4YHuLgceCG+/iNAouAOA64FjwudNIjRvwgXuPhT4B/DbauebRtMYc0iaiISgA4g0kN8Cc9z9xWrbTg0v08LrqYQKwycRnrMDoWGtDyQJeNDMBhEa4bJ7ePtxhIY72QOsNrOPw9v7AkcAH4aGGCKe0Lgye5UQGo1UpF6oEEizF27nv4D9Rxo14G53f+wQT70diGTqyB8Smg/jCiAR2NvsVNPw6Hu3z3T3kbXsTw6/t0i9UNOQNGvhMf2fBq4Kz/hW3XvAN8NzP2BmuWa27yQqWwkNSlaThUB+BDHSgTXhIY6v5v8KwATgQgvpQGj2NAiNiJprZkeFc7UwsyOqna8XsN88tiKHSoVAmrsbCM2Q9cg+XUi/EZ7t65/ARDObRWikz33/0Z8JVJnZjH1vFrt7BbDYzHocIMODwLfNbBLQBdh7T+FlQs08swnNqf0FoVnIdgIXAn82sxmEmq723mBOIlR8piFST9R9VOQwmNl5wFB3/9khvj7V3cvDvYC+AI5291rvO5jZRUA/d/91bceIHCzdIxA5DO7+upm1PYxTvGNmrQndO/hlXUUgzAhNXi5Sb3RFICIS43SPQEQkxqkQiIjEOBUCEZEYp0IgIhLjVAhERGLc/wcryxO6sRhzsAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "gamma = 0.2 # komplett ausgedacht(!)\n", "t = torch.linspace(0, 14)\n", "I = torch.exp(-gamma * t)\n", "pyplot.plot(t, I)\n", "pyplot.xlabel(\"Zeit (Tage)\")\n", "pyplot.ylabel(\"Anteil noch Kranker\")\n", "t[15], I[15]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Lesen wir auf der $x$-Achse die Zeit ab (z.B. in Tagen), dann sehen wir auf der y-Achse der Anteil noch kranker Personen: Nach 2 Tagen sind noch ca. 65% der zum Start infizierten Personen krank.\n", "\n", "Wir können aber die Fläche unter dem Graphen gedanklich in horizontale Linien aufteilen. Dann erhalten wir Krankheitsdauern. 65% der zum Start infizierten haben eine Krankheitsdauer von mindestens 2 Tagen.\n", "Damit ist die Fläche unter dem Graphen die durchschnittliche Krankheitsdauer $T_I$. Wir können sie als Integral ausrechnen:\n", "\n", "$$\n", "T_I = \\int_0^\\infty \\exp(-\\gamma t) dt = - \\frac{1}{\\gamma} \\exp(-\\gamma t) \\big|_{t=0}^\\infty = \\frac{1}{\\gamma}.\n", "$$\n", "\n", "Wollen wir also die Annahme für die Genesungsrate $\\gamma$ aus der Krankheitsdauer ableiten, ist also der Kehrwert der durchschnittlichen Krankheitsdauer $\\gamma = 1/T_I$. (Wir machen aber mit dem Modell _auch_ eine strukturelle Annahme.)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Die Basisreproduktionszahl\n", "\n", "Wenn wir davon ausgehen, dass die Weitergabe der Krankheitserreger mit der Rate $\\beta$ erfolgt und die Infektion $T_I$ Tage andauert, dann hat jeder Kranke $R_0 = \\beta T_I = \\frac{\\beta}{\\gamma}$ \"Infektionskontakte\" (wenn die Kontakte alle in $S$ wären, wären es Ansteckungen. Diese Zahl ist immens wichtig und hat daher den Nahmen Basisreproduktionszahl oder Grundvermehrungsrate.\n", "\n", "Man kann sich überlegen, dass wenn man $\\beta$ und $\\gamma$ um denselben Faktor multipliziert, dass eine \"Zeitskalierung\" gleichkommt: Verdoppelt man beide, läuft alles doppelt so schnell ab, d.h. man hat zu $t$ dieselben werte wie zuvor in $2t$. Wir werden sehen, dass die Form des Verlaufs -- insbesondere das Maximum an Infektionen -- wesentlich von der Basisreproduktionszahl $R_0$ abhängen.\n", "\n", "Eine Sache, die man leicht sieht, ist, dass für $R_0 < 1$ die Ableitung von $I$ stets negativ ist (weil $S/N < 1$). Die Anzahl der Infizierten nimmt in diesem Fall sofort ab, es ist dann also -- in diesem Modell und im Sinne der Epidemiologie -- alles unter Kontrolle.\n", "\n", "Auf der englischen Wikipedia gibt der [Artikel zur Basisreproduktionszahl](https://en.wikipedia.org/wiki/Basic_reproduction_number) für einige krankheiten Schätzer an, für CoViD-19 wird sie auf zwischen $1.4$ und $3.9$ geschätzt. Wir kommen darauf aber nochmal zurück." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Das Maximum der Anzahl der Infizierten\n", "\n", "Eine sehr wichtige Frage ist, wieviele infizierte Personen es maximal gleichzeitig gibt. Davon hängt, wenn man annimmt, dass z.B. 5% der Kranken eine intensivmedizinische Behandlung benötigt, zum Beispiel ab, ob das medizinische System überlastet wird. Sobald $5\\% \\max I > Kapazität$, ist die Not groß (andere Faktoren spielen auch eine Rolle, die Kapazität ist zum Beispiel nicht konstant).\n", "\n", "Wir nehmen für die Basisreproduktionszahl $R_0 > 1$ an (siehe oben für $R_0 < 1$).\n", "\n", "Bestimmen wir also das Maximum der gleichzeitig Infizierten $I_{max}$. Man kann zeigen, dass dieses Maximum angenommen wird, wir nennen den Zeitpunkt, zu dem dies geschiet $t^*$.\n", "Dies können wir nicht direkt bestimmen, sondern bestimmen zunächst $S(t^*)$ und $R(t^*)$ und nutzen dann wieder\n", "die Gleichung $S+I+R = 1$ des geschlossenen Systems. (Wir erinnern uns, dass wir $N=1$ annehmen.)\n", "\n", "Dieses ist der vielleicht mathematisch anspruchsvollste Abschnitt, man kann ihn auch überspringen und nur das Ergebnis anschauen.\n", "\n", "### Anzahl $S(t^*)$ der Ansteckbaren am Maximalpunkt für $I$\n", "\n", "Dies ist der Fall, wenn die Veränderung der Zahl der Infizierten gerade ihr Vorzeichen von positiv zu negativ ändert bzw. wenn die Ableitung $\\dot I(t) = 0$. Für $\\dot I(t)$ haben wir eine Gleichung, also gilt an $t^*$\n", "$$\n", "0 = \\dot I(t^*) = \\left(\\beta S(t^*) - \\gamma \\right) I(t^*).\n", "$$\n", "\n", "Da aber $I > 0$ (andernfalls gibt es keine Infizierten) folgt\n", "$$\n", "0 = \\beta S(t^*) - \\gamma\n", "$$\n", "oder\n", "$$\n", "S(t^*) = \\gamma / \\beta.\n", "$$\n", "\n", "Damit haben wir $S(t^*)$ bestimmt.\n", "\n", "### Zusammenhang zwischen den Zahlen der Ansteckbaren $S$ und Genesenen $R$\n", "\n", "Um $R(t^*)$ zu bestimmen, können wir einen Zusammenhang zwischen den \n", "Zahlen der Ansteckbaren $S$ und Genesenen $R$ herstellen.\n", "Dazu eliminieren wir $I$ aus den Gleichungen für $\\dot S$ und $\\dot R$.\n", "\n", "Lösen wir die Gleichungen für $S$ und $R$ nach $I$ auf und setzen sie gleich, so erhalten wir\n", "$$\\dot S / (-\\beta S) = \\dot R / \\gamma$$\n", "\n", "oder, wenn wir den Quotienten $\\dot S/S$ als Ableitung des Logarithmus schreiben,\n", "$$\n", "\\dot {(\\ln S)} = \\dot S / S = -\\frac{\\beta}{\\gamma} \\dot R.\n", "$$\n", "\n", "Die Methode der [Trennung der Veränderlichen](https://de.wikipedia.org/wiki/Trennung_der_Veränderlichen) lässt uns nun auf beiden seiten in der Zeit integrieren:\n", "\n", "$$\n", "\\ln S(t) - \\ln S(0) = -\\frac{\\beta}{\\gamma} (R(t) - R(0))\n", "$$\n", "und das Exponential ist\n", "$$\n", "S(t) = S(0) \\exp \\left( -\\frac{\\beta}{\\gamma}{(R(t) - R(0))}\\right)\n", "$$\n", "\n", "### Anzahlen $R(t^*)$ der Genesenen und $I(t^*)$ am Maximalpunkt für $I$\n", "\n", "Wenn also $S(t^*) = \\frac{\\gamma}{\\beta}$, $R(0) \\approx 0$ und $S(0) \\approx 1$, so ist\n", "$$\n", "R(t^*) = R(0) -\\frac{\\gamma}{\\beta} (\\ln S(t_{crit} - \\ln S(0))\n", " \\approx -\\frac{\\gamma}{\\beta} \\ln \\frac{\\gamma}{\\beta}\n", "$$\n", "und\n", "$$I_{max} = I(t^*) = 1 - S(t^*) - R(t^*) \\approx 1-\\frac{\\gamma}{\\beta} (1-\\ln \\frac{\\gamma}{\\beta} )$$\n", "\n", "\n", "### Das Maximum der gleichzeitig Infizierten in Abhängigkeit der Basisreproduktionszahl\n", "\n", "Wie wir sehen, hängt das (approximative) Maximum unter der Annahme, dass anfangs fast die gesamte Population ansteckbar und (dann ja auch gegeben) nahezu keiner genesen ist, nur von der Basisreproduktionszahl ab und zwar als\n", "$$\n", "I_{max} \\approx 1-\\frac{1}{R_0} (1 + \\ln R_0).\n", "$$\n", "Wir können diesen Zusammenhang graphisch darstellen:" ] }, { "cell_type": "code", "execution_count": 100, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.1515151262283325" ] }, "execution_count": 100, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEKCAYAAAAB0GKPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3wVZdbA8d+hV+k9CUWaSCc0QRexLLYFy65YEQsLim1f99VddW27a3nX3VVXQXRRUAFdFQXFghURkCT0KqEEQuiBQAgh7bx/zMS9hhuYSe7NzU3O9/PJJ3fmmefOmTtwT2bmKaKqGGOMMadSJdIBGGOMiQ6WMIwxxnhiCcMYY4wnljCMMcZ4YgnDGGOMJ9UiHUA4NW3aVNu1axfpMIwxJmokJSXtV9VmwcoqdMJo164diYmJkQ7DGGOihoikFFdmt6SMMcZ4YgnDGGOMJ5YwjDHGeGIJwxhjjCeWMIwxxnhiCcMYY4wnljCMMcZ4UqYJQ0RGiMhGEUkWkQdOsl1/EckXkav81jXGmMoqv0D5cv0eJn2zOSzvX2Yd90SkKvAicAGQCiSIyBxVXRdku6eBz/zWNcaYyujg0RzeSdzBG0tSSD14jDYNazN2SDtqVa8a0v2UZU/vAUCyqm4BEJFZwEig6Jf+ncB7QP8S1DXGmEpjbVoG0xZt48MVaRzPK2BA+8b84aIzuPDMFlSvGvobSGWZMNoAOwKWU4GBgRuISBvgcmA4P08Yp6wb8B7jgHEAcXFxpQ7aGGPKk7z8Aj5ft4fXv9/G0m3p1K5elSv6xjDmrLZ0bXlaWPddlglDgqwrOj/sP4H7VTVf5Gebe6nrrFSdAkwBiI+Pt/lnjTEVwsGjOcxM2M4bi1PYlZFNbOPaPHjxGfwmPpYGdaqXSQxlmTBSgdiA5Rggrcg28cAsN1k0BS4WkTyPdY0xpsLZsPswr3+/jdnLd3I8r4AhHZvwxMjunNu1OVWrBPtbOnzKMmEkAJ1EpD2wExgNXBu4gaq2L3wtIq8DH6nqByJS7VR1jTGmoigoUL7asJfXFm3l++QD1KpehSv6xjB2SDs6t6gfsbjKLGGoap6ITMRp/VQVmKqqa0VkvFs+2W/dsojbGGPKytHjebyblMpr329l24EsWjWoxf0jujK6fyyN6taIdHiIasW9zR8fH682H4YxprzblXGM1xdtY+YP2zmcnUfv2IbcPLQ9F3VvGZbWTicjIkmqGh+srEJPoGSMMeXZ6tQMXl24hY9X7aJAlYu6t+Lmoe3p17ZRpEMLyhKGMcaUocLnE698t4UftqZTr2Y1bjqrHWPOakds4zqRDu+kLGEYY0wZyM7NZ/bynbzy3Ra27DtK6wa1ePDiM7h6QCyn1SqbZrGlZQnDGGPC6FBWDm8sTmHa4m3sz8yhe5vTeG50by7u0arMn0+UliUMY4wJgx3pWfx74VbeTtjBsdx8hnVpxrhzOjC4QxOKdEyOGpYwjDEmhNalHeblBZv5aNUuBBjZuw3jzulAl5aR6z8RKpYwjDGmlFSVH7amM+mbzXz74z7q1qjKzUPacfPQ9rRqUDvS4YVMiRKGiFRV1fxQB2OMMdGkoED5Yv0eJn27meXbD9G0Xg1+/8suXD+wbZmN71SWSnqFMUVE7lTVLBE5R1UXhDQqY4wpx3LzC5i7Mo1J32xm095MYhvX5olR3fl1v5iQz0FRnpQ0YfwJ+Lc7MOAKwBKGMabCy87N5z+JO5j87RZ2HjpG15b1+efVvbm0ZyuqRVmLp5IoacJ4AtgIdADeCV04xhhT/mQez+OtJSm88t1W9mcep29cQx4feSbDuzaP2hZPJVHShPG/qrpfROoCzwG3hjAmY4wpFzKycnl90Tamfr+VjGO5DO3YlDvO7cOgDo0rVaIoVNKEMVFEBuIMNb48hPEYY0zEHcg8zr8XbmX64hQyj+dx/hktmDi8I71jG0Y6tIgqacJoCCwBpgP3hi4cY4yJnL1HsnllwRbeXLKd7Lx8Lu7RionnduSMVuGd+jRalDRhpLt197qvjTEmau05nM3kbzcz44ft5BUoI3u15vZzO9Kxeb1Ih1aulChhqOrjItIaeB5Y47WeiIzAeeZRFXhVVZ8qUj4S54F6AZAH3KOqC92ybcARIB/IK268dmOM8Wp3RjaTvklmZsIO8guUK/u24fZhHWnXtG6kQyuXStpxrxvQFnhEVVM91qkKvAhcgDNHd4KIzFHVdQGbfQnMUVUVkZ44LbC6BpSfq6r7SxKzMcYU2pVxjEnfbGbW0h0UqHJVvxjuOLdjuR9ePNJKekvqMeA/wDgRaauqYzzUGQAkq+oWABGZBYwEfkoYqpoZsH1doOJOB2iMKXO7M7J56ZvknxLFr+NjuH2YJQqvSpow5qvqO/jrg9EG2BGwnAoMLLqRiFwOPAk0By4JKFLgcxFR4GVVnRJsJyIyDhgHEBcX5yM8Y0xFtfdINi99vZkZS7dTUGBXFCVV0oRxlvs84gCwXlX/7qFOsEbLJ1xBqOpsYLaInIPzPON8t2iIqqaJSHNgvohsCDYkiZtIpoAzp7e3wzHGVEQHMo8z+dvNvLEkhdx85xnFncM7WaIooZImjDWq+jcRqQac6bFOKhAbsBwDpBW3saouEJHTRaSpqu5X1TR3/V4RmY1zi8uGJDHGnOBQVg6vfLeF177fRnZuPqP6tOGu4Z3sYXYplTRhXCoix4HPVHWlxzoJQCcRaY/T4W80cG3gBiLSEdjsPvTuC9QADrg9yquo6hH39YXA4yWM3RhTQR3JzmXqwm28+t0WjhzP47Jerbn7vE7WPDZESpowrgb6AFeIyOmqetupKqhqnohMBD7DaVY7VVXXish4t3wycCVwo4jkAseAq93k0QLnNlVhzDNU9dMSxm6MqWCyc/OZvngbk77ZzMGsXC7s1oJ7L+hsHe5CTFRPfptfRF4CVgOrgNWqergsAguF+Ph4TUxMjHQYxpgwyckr4J3EHbzw1Sb2HD7OOZ2bcd+FnekZU7mH8CgNEUkqrp+blyuMFUBPnFtI3UXkCD9PILNCFqkxxniQX6DMWbmTf8zfxPb0LOLbNuL50X0Y2KFJpEOr0E6ZMIo2XxWRGJwE0gOn2aslDGNMmVBVvtqwl//7bCMbdh/hjFan8dpN/RnWpVmlHD22rPl+huH27E4F5oU+HGOMCS5xWzpPf7qBhG0HadukDs+N7s1lPVtTpYolirJS0ofexhhTJjbtOcLTn27ki/V7aFa/Jn8e1Z2r+8dSvRLMcFfeWMIwxpRLuzKO8Y/5P/JuUip1a1Tj97/swtgh7ahTw762IsXzJy/ODcLrgA7uaLVxQEtVXRq26Iwxlc7h7FwmfbOZqQu3ogpjh7Rn4rkdaVS3RqRDq/T8pOqXcIYdH47Tae4I8B7QPwxxGWMqmZy8At5cksILX23iYFYuo3q35n8u7GLDeJQjfhLGQFXtKyLLAVT1oIhYyjfGlIqqMm/1bp7+dAPb07MY0rEJf7joDLq3aRDp0EwRfhJGrjunhQKISDOcKw5jjCmRpJR0/vLxepZtP0SXFvV5fWx/ftHZmsiWV34SxvPAbKC5iPwFuAp4OCxRGWMqtO0Hsnj60w18vHoXzevX5Jkre3JlvxiqWhPZcs1zwlDVt0QkCTgPZ6jyUaq6PmyRGWMqnIxjubz0dTKvfb+NqlWEu8/rxLhzOlC3prV8igZ+Wkk9rar3AxuCrDPGmGLl5Rcwc+l2/vHFJg5m5XBV3xju+2UXWpxWK9KhGR/8pPULgKLJ4aIg64wx5icLftzHEx+tY9PeTAZ1aMxDl3SzB9pR6pQJQ0QmALcDp4vIqoCi+sD34QrMGBPdNu/L5C8fr+erDXtp26QOL9/Qjwu7tbAH2lHMyxXGDJw5LF4FxgasP6Kq6WGJyhgTtTKO5fL8l5uYtmgbtatX5Y8Xd2XMWe2oWa1qpEMzpeRltNoMIENEGqpqShnEZIyJQvkFytsJO/jb5xs5mJXD6P6x/M+FXWhar2akQzMh4mf0rsUiUqpe3SIyQkQ2ikiyiDwQpHykiKwSkRUikigiQ73WNcZEztKt6Vz2wkL+OHs1HZvVY+7EoTx5RU9LFhWMn4fe5wLjRWQbcBSnaa2qak8vld1Ofy/iPDxPBRJEZI6qrgvY7Etgjjsta0/gHaCrx7rGmDK2K+MYT87bwJyVabRuUIsXrunDpT1b2XOKCspPwriolPsaACSr6hYAEZkFjAR++tJX1cyA7evi9ir3UtcYU3aO5+Xz74Vb+ddXyeQVKHcN78j4YafbSLIVnJ+zu50go9UCXp9rtAF2BCynAgOLbiQilwNPAs1xZvTzXNetPw4YBxAXF+cxNGOMV19v3Mvjc9exdf9RLuzWgocv7WYDBFYSfp5hvAQMBq5xl4/g3CbyKtg1qp6wQnW2qnYFRgFP+Knr1p+iqvGqGt+sWTMf4RljTmZHehbjpicy9rUEBJh28wCm3BhvyaISKcvRalOB2IDlGCCtuI1VdYGInC4iTf3WNcaEzvG8fKZ8u4V/fZ1MFRHuH9GVm4daM9nKqCxHq00AOolIe2AnMBq4NnADEekIbHYfevcFagAHgEOnqmuMCb1vf9zHIx+uYduBLC7p0YoHLzmD1g1rRzosEyGlHa32Ia+VVTVPRCbidAKsCkxV1bUiMt4tnwxcCdwoIrnAMeBqVVUgaF0fsRtjfNiVcYwnPlrHvNW7ad+0Lm/cMoCzO9kt3spOnO9jjxuLdOW/o9V+Wd5Hq42Pj9fExMRIh2FM1MjNL+D177fxjy9+JL9AuXN4R247p4PdfqpERCRJVeODlflqA6eqGwgYrdYYU3EkpaTz4Ow1bNh9hOFdm/PYr860B9rmZ7wMPrhQVYeKyBF+3jKpsOPeaWGLzhgTdoeycnj6043MXLqdVg1q2SCBplhexpIa6v6uH/5wjDFlRVX5cEUaT3y0jkPHcrnt7Pbcc35nm8zIFMvPBEpfAs+q6ryAdVNUdVxYIjPGhM22/Ud56IM1LEzeT6/YhrxxeQ+6tbabBebk/Pwp0R64X0T6q+pj7rqgD0aMMeVTbn4BUxZs4fkvN1GjahWeGHkm1w5sa3NpG0/8JIxDOC2knheRucD14QnJGBMOy7Yf5A/vrWbjniNc1L0lj/7qTJsi1fjiJ2GIquYBt4vITcBCoFFYojLGhEzm8Tz+9tlGpi3eRov6tXjlxngu6NYi0mGZKOQnYUwufKGqr4vIauCO0IdkjAmVrzbs4aHZa9h1OJsbB7Xlvl92oX6t6pEOy0QpzwlDVV8uspwE3BzyiIwxpXYg8ziPzV3HnJVpdGpej3fHn0W/tnZDwJSO9cMwpgJRVeasTOPROWvJPJ7Hved3ZsKw06lRzc/A1MYEZ/0wjKkgdmUc48HZa/hqw156xzbkmat60rmF/bc1oePlCuMNVb1BRO5W1efKIihjjHeqyqyEHfz14/XkFhTw0CVnMHZIe2sqa0LOyzOMfiLSFrhZRKZTZDIjVU0PS2TGmFPakZ7FA++v4vvkAwzu0ISnruxB2yZ1Ix2WqaC8JIzJwKdAB2BZkTJ11xtjylBBgfLWDyk8+ckGBPjzqO5cOyCOKnZVYcLIyzOM53E6601S1QllEJMx5iR2pGfx+3dXsmRLOmd3asqTV/QgppGNKmvCz08/jHtF5FqgXWA9VX3c6xuIyAjgOZxJkF5V1aeKlF8H3O8uZgITVHWlW7YNZx7xfCCvuPHajamoCgqUt5Zu58l566kiwpNX9GB0/1gbVdaUGT8J4wMgA0gCjvvdkTu964vABThzdCeIyBxVXRew2VbgF+584RcBU4CBAeXnqup+v/s2JtqlHszi/vecZxVnd2rKU1f2pI1NlWrKmJ+EEaOqI0qxrwFAsqpuARCRWcBI4KeEoaqLArZfAsSUYn/GRD1V5e2EHfz54/WoKn+9vAfXDLCrChMZfhLGIhHpoaqrS7ivNsCOgOVUfn71UNQtwCcBywp8LiIKvKyqU4JVEpFxwDiAuLi4EoZqTOTtOZzN/e+t4puN+xjcoQnPXNXTZsAzEeUnYQwFbhKRrTi3pAp7evf0WD/Yn0RBJxQXkXNxEsbQgNVDVDVNRJoD80Vkg6ouOOENnUQyBZw5vT3GZky5MmdlGg9/sIbjefk8clk3xgxuZy2gTMT5SRgXlXJfqUBswHIMkFZ0IxHpCbwKXKSqBwrXq2qa+3uviMzGucV1QsIwJpodysrhoQ/W8NGqXfSJa8izv+5Fh2b1Ih2WMYC/wQdTSrmvBKCTiLQHdgKjgWsDNxCROOB94AZV/TFgfV2giqoecV9fCHhunWVMNPhm417+991VpB/N4b4LOzP+F6dTraqNAWXKDy9DgxQddPCnInwMPqiqeSIyEfgMp1ntVFVdKyLj3fLJwJ+AJsBL7kO9wuazLYDZ7rpqwAxV/dTLfo0p747l5PPXeet5Y0kKnZrXY+pN/enepkGkwzLmBKJacW/zx8fHa2JiYqTDMKZYq1IPcc+sFWzZf5Rbhrbn97/sQq3qVSMdlqnERCSpuH5ufp5hGGNCJC+/gEnfbOa5LzfRrH5NZtw6kLM6No10WMaclCUMY8rYjvQs7n17BYkpB/lVr9Y8MbI7DerYLHim/LOEYUwZUVVmL9/Jnz5ciwD/vLo3o/q0iXRYxnjmOWGIyO+CrM4AklR1RehCMqbiyTiWy0MfrGHuyjQGtGvMs7/pZZ3wTNTxc4UR7/7MdZcvwWkqO15E/qOqz4Q6OGMqgoRt6dwzawW7D2dz34WdmTCso01uZKKSn4TRBOirqpkAIvII8C5wDs6AhJYwjAmQl1/AC18l88JXm4hpVId3xw+mT1yjSIdlTIn5SRhxQE7Aci7QVlWPiYjv0WuNqchSD2ZxzyznwfYVfdrw+Kju1KtpjwxNdPPzL3gGsEREPnSXLwNmuj2v1xVfzZjK5ZPVu7j/vVUUqD3YNhWLn6FBnhCRT4AhOL28x6tqYa+468IRnDHRJDs3n8c/WseMH7bTK7Yhz4/ubfNrmwrF1zWymyCs67QxRfy45wgTZyzjxz2Z/PYXHbjvwi5Ut3GgTAXjZSyphao6NMiYUr7GkjKmIiqc4OjRuWupV7Ma028ewDmdm0U6LGPC4pQJQ1WHur/rhz8cY6LHkexcHpy9hjkr0xjSsQn/uLo3zevXinRYxoSN52tmEXnayzpjKoM1OzO47IWFfLQqjfsu7Mz0mwdasjAVnp+brBcEWVfaSZWMiSqqyhtLUrjipUVk5xYwa9xgJg7vZB3xTKXg5RnGBOB2oIOIrAooqg98H67AjClvjmTn8sD7q/l41S6GdWnG33/Tm8Z1a0Q6LGPKjJdWUjOAT4AngQcC1h9R1fSwRGVMObMu7TB3zFjG9vQs7h/Rld+e08Hm2DaVzilvSalqhqpuU9VrVDUl4Md3shCRESKyUUSSReSBIOXXicgq92eRiPTyWteYcFBVZi7dzqiXvicrJ49Z4wYxYdjplixMpVSSZrUS+Ntrs1oRqQq8iPMsJBVIEJE5qhrYS3wr8AtVPSgiFwFTgIEe6xoTUlk5eTw0ew3vL9/J2Z2a8o+re9O0Xs1Ih2VMxJRls9oBQLKqbgEQkVnASAKGFVHVRQHbLwFivNY1JpSS92Zy+1tJbNqbyb3nd2bicBth1hg/zWpFRK4XkYfd5VgRGeBjX22AHQHLqe664tyC8+zEV10RGSciiSKSuG/fPh/hGeP4aFUaI/+1kP2ZOUy/eQB3n2+toIwBf81qXwIGA9e6y5k4t4m8CvY/ToOsQ0TOxUkY9/utq6pTVDVeVeObNbMet8a7nLwCHp2zlokzltO11Wl8fNdQzu5k/4aMKeRnLKmBqtpXRJYDuM8Z/LQpTAViA5ZjgLSiG4lIT+BV4CJVPeCnrjEltTsjm9vfSmLZ9kPcPKQ9f7i4q40FZUwRfhJGrvvwWQFEpBlQ4KN+AtBJRNoDO4HR/PdqBfc944D3gRtU9Uc/dY0pqcWbD3DnzGVk5eTzwjV9uKxX60iHZEy55CdhPA/MBpqLyF+Aq4CHvVZW1TwRmQh8BlQFpqrqWhEZ75ZPBv6EM7PfSyICkOfeXgpa10fsxpxAVXnluy08/elG2jWpw8zbBtGphQ2ZZkxxRDXoo4DgG4t0Bc7DeabwpaquD1dgoRAfH6+JiTYauzlR5vE87n93FR+v3sXFPVryzFW9bEY8YwARSVLV+GBlnv+HiMiXwLOq+mLAuimqOi4EMRpTZjbvy2T8G0ls3pfJHy/uym1nd8C9ojXGnISfP6naA/eLSH9VfcxdFzQLGVNezV+3h9+9vYLq1arw5i0DOatj00iHZEzU8NMM5BDO7agWIjJXRBqEKSZjQq6gQPn7/B+5bXoi7ZrWZe6dQy1ZGOOTnysMUdU84HYRuQlYCDQKS1TGhNDh7FzunbWCLzfs5ap+Mfx5VHdqVa8a6bCMiTp+Esbkwheq+rqIrAbuCH1IxoRO8t4jjJuexPb0LB4feSY3DGprzyuMKSHPCUNVXy6ynATcHPKIjAmRz9fu5nfvrKRW9SrMuG0QA9o3jnRIxkQ1P62kagJXAu0C66nq46EPy5iSKyhQnv9qE//8YhM9Yxow+fp+tG5YO9JhGRP1/NyS+hDIAJKA4+EJx5jSyTyex/+8s4LP1u7hir5t+OvlPex5hTEh4idhxKjqiLBFYkwpbT+QxW3TE9m09wgPX9qNm4e0s+cVxoSQn4SxSER6qOrqsEVjTAktSt7P7TOWoQrTbx7I0E7WZNaYUPMy495qnAEHqwFjRWQLzi2pwhn3eoY3RGOKp6q8sSSFx+auo0PTurxyYzztmtaNdFjGVEherjAuDXsUxpRATl4Bj8xZy8yl2zn/jOb84+re1K9VPdJhGVNheZmiNQVARKYBd6vqIXe5EfAs1rTWRED60RwmvJnED1vTuX3Y6dx3YReq2Kx4xoSVn2cYPQuTBfw0gVKfMMRkzElt3H2EW6cnsOfwcZ4b3ZuRvU82068xJlT8JIwqItJIVQ8CiEhjn/WNKbUv1+/hrpnLqVuzGu/8djC9YxtGOiRjKg0/gw8+i9NS6gkReRxYBDzjZ2ciMkJENopIsog8EKS8q4gsFpHjInJfkbJtIrJaRFaIiE1yUcmoKlMWbObW6Ym0b1aXOROHWrIwpoz5GRpkuogkAefitJC6QlXXea3vTu/6InABzhzdCSIyp8h7pAN3AaOKeZtzVXW/132aiiEnr4CHPljNO4mpXNyjJc/+uje1a1hnPGPKmq9bSu60qCWdGnUAkKyqWwBEZBYwEvgpYajqXmCviFxSwn2YCib9aA7j30xi6dZ07hzekXvP72wPt42JEC/9MBaq6lAROYLTH+OnIpx+GKd53FcbYEfAciow0HOkzr4/FxEFXlbVKcXEOw4YBxAXF+fj7U15k7w3k1umJbArI9sebhtTDnhpVjvU/V2/lPsK9meh9wnFYYiqpolIc2C+iGxQ1QUnvKGTSKaAM6d3yUI1kbZw034mvJVEzWpVmHnbIPq1talXjIk0zw+9xXG9iDzsLseKyAAf+0oFYgOWY4A0r5VVNc39vReYjXOLy1RAb/2QwpjXltK6QW1m3z7EkoUx5YSfVlIvAYOBa93lTJyH2F4lAJ1EpL2I1ABGA3O8VBSRuiJSv/A1cCGwxse+TRTIL1D+/NE6Hpy9hqEdm/LuhMHENq4T6bCMMS4/D70HqmpfEVkOP3Xcq+G1sqrmichE4DOgKjBVVdeKyHi3fLKItAQSgdOAAhG5B+gGNAVmuyOPVgNmqOqnPmI35VxWTh53zVzBF+v3MGZwWx6+tBvVqvr5e8YYE25+Ekau2zRWAUSkGVDgZ2eqOg+YV2Rd4NSvu3FuVRV1GOjlZ18meuw5nM0t0xJYl3aYRy/rxk1D2kc6JGNMEH4SxvM4zw6ai8hfgKuAh8ISlak01qUd5pZpCRw+lsurY+IZ3rVFpEMyxhTDT8e9t9yOe+fhtHgaparrwxaZqfC+3riXiW8to36t6vxn/Fl0a+21hbYxJhL8dtzbAGwIUyymEnlzSQqPzFlL15b1+feY/rRsUCvSIRljTsFLx73CDnuF/SgK+zb47bhnDAUFylOfbmDKgi0M79qcF67pQ92aNoalMdHAS8e90nbYMwaA7Nx87n17BZ+s2c2Ng9vyJ2sJZUxU8fynnYj8LsjqDCBJVVeELiRTER3IPM6t0xNZseMQD11yBrcMbY/bTNoYEyX83AuId3/musuX4HTGGy8i/1FVX0Odm8pjy75MbnotgT2Hs5l0XV9GdG8V6ZCMMSXgJ2E0AfqqaiaAiDwCvAucAyThc24MUzkkbEvntumJVBVh5rhB9I2zYT6MiVZ+EkYckBOwnAu0VdVjInI8tGGZiuCjVWn87p2VxDSszWtj+9O2Sd1Ih2SMKQU/CWMGsEREPnSXLwNmumM7eZ5IyVR8qsor323hr/M2EN+2Ea/cGE+jup5HkTHGlFN+Ou49ISLzgKE4TWrHq2rhVKnXhSM4E33yC5TH565l2uIULunRimd/04ta1W12PGMqAr8d95JwnlcYc4JjOfncNWs589ft4dah7fnjxWfY7HjGVCDWY8qExIHM49wyLZGVqYdsAEFjKihLGKbUtu0/yk2vLWVXRjaTruvHiO4tIx2SMSYMLGGYUlm+/SC3TEtEVZlx20D6tW0c6ZCMMWHiZyypE4qwsaQqtfnr9nDnzGU0q1+TaWMH0KFZvUiHZIwJo1MO5KOq9VX1tCA/9f0mCxEZISIbRSRZRB4IUt5VRBaLyHERuc9PXVO23vohhd++kUjnFvV5f8IQSxbGVAK+bkmJSCOgE/DTWNSqusBj3ao4c4BfAKQCCSIyR1UD+3CkA3cBo0pQ15QBVeVvn2/kxa8322izxlQyfgYfvBW4G2cK1RXAIGAxMNzjWwwAklV1i/t+s4CRBHT6U9W9wF4RucRvXRN+OXkFPPD+Kt5ftpNrBsTyxMjuNtqsMZWIn//tdwP9gRRVPRfoA+zzUb8NsCNgOdVdF9K6IjJORBJFJHHfPj/hmZPJPJ7HLdMSeH/ZTn53QWf+enkPSxbGVDJ+/sdnq2o2gIjUdGff6+KjfrAeXMEepjSW4f4AABJ4SURBVJeqrqpOUdV4VY1v1qyZ5+BM8fYezuY3kxezaPMBnrmqJ3ed18mGJjemEvJz8zlVRBoCHwDzReQgkOanPhAbsBzjo35p6ppSSN6byZipSzmYlcO/x8QzrEvzSIdkjIkQP2NJXe6+fFREvgYaAJ/62FcC0ElE2gM7gdHAtWVQ15RQ4rZ0bp2eSLUqwtvjBtMjpkGkQzLGRFCJmreo6rclqJMnIhOBz4CqwFRVXSsi493yySLSEkgETgMKROQeoJuqHg5WtySxG28+XbObu2ctp3XD2kwbO4C4JnUiHZIxJsJE1dtjBBGJBx4E2hKQaFS1Z3hCK734+HhNTEw89YbmZ6Yv3sYjc9bSO7Yh/x7Tn8Y2NLkxlYaIJKlqfLAyP1cYbwG/B1YDBaEIzJQvBQXKM59tZPK3mzn/jBa8cE0fatewocmNMQ4/CWOfqs4JWyQmonLyCrj/vVXMXr6T6wbG8divzrRms8aYn/GTMB4RkVeBL4GfpmRV1fdDHpUpU0eyc5nw5jIWJu/n97/swu3DTrdms8aYE/hJGGOBrkB1/ntLSgFLGFFsz+FsxkxdSvLeTP72615c1S8m0iEZY8opPwmjl6r2CFskpsxt2nOEm15L4FBWDlNv6s85na2jozGmeH5uUi8RkW5hi8SUqR+2HODKSYvIyS/g7d8OtmRhjDklP1cYQ4ExIrIV5xlG4XwY5bZZrQlu7so0/uedlcQ2rs3rYwcQ29j6WBhjTs1PwhgRtihMmVBVXv1uK3+Zt57+7Rrxyo3xNKxjfSyMMd74GRokJZyBmPDKL1Aen7uWaYtTuLhHS/7+m97Uqm59LIwx3nmZonWhqg4NMlWrTdEaJY7l5HPXrOXMX7eHW4e2548Xn0GVKtZs1hjjzykThqoOdX/XD384JtT2Zx7n1mmJrEw9xCOXdWPskPaRDskYE6U8t5IK1kJKRIaFNBoTUpv3ZXLFS4vYsPswk67rZ8nCGFMqfprVviMi94ujtoi8ADwZrsBM6Szdms6VkxZx9HgeM28bxIjuLSMdkjEmyvlJGANxJjFahDM/RRowJBxBmdL5cMVOrn/1BxrXrcHs24fQJ65RpEMyxlQAfprV5gLHgNpALWCrqtqoteWIqvKvr5J5dv6PDGzfmJdv6GfNZo0xIePnCiMBJ2H0x+nEd42IvOtnZyIyQkQ2ikiyiDwQpFxE5Hm3fJWI9A0o2yYiq0VkhYjYJBdF5OQVcN9/VvHs/B+5ok8b3rhloCULY0xI+bnCuEVVC7+odwMjReQGr5VFpCrwInABzhzdCSIyR1XXBWx2EdDJ/RkITHJ/FzpXVff7iLlSOHg0h9++mcTSrencc34n7j6vk402a4wJOT8d9xJFpBHOl3ktd7WfznwDgGRV3QIgIrOAkUBgwhgJTFdnGsAlItJQRFqp6i4f+6lUtuzL5ObXE0jLyOa50b0Z2btNpEMyxlRQnhOGiNwK3A3EACuAQcBiYLjHt2gD7AhYTuXnVw/FbdMG2IXTafBzEVHgZVWdUkyc44BxAHFxcR5Di06Lkvcz4a1lVKsizLxtIP3aNo50SMaYCszPM4y7cZ5fpKjquUAfYJ+P+sHukRSdUPxk2wxR1b44t63uEJFzgu1EVaeoaryqxjdrVnFHYH1zSQo3TF1Ki9Nq8sEdQyxZGGPCzs8zjGxVzRYRRKSmqm4QkS4+6qfiNMstFIPTNNfTNqpa+HuviMzGucW1wMf+K4Tc/AL+/NE6pi1OYXjX5jw3ujf1a1WPdFjGmErAzxVGqog0BD4A5ovIh5z4hX8yCUAnEWkvIjWA0UDROcLnADe6raUGARmquktE6opIfQARqQtcCKzxse8K4eDRHMZMXcq0xSncOrQ9r9wYb8nCGFNm/Dz0vtx9+aiIfA2cBnzqo36eiEwEPgOqAlNVda2IjHfLJwPzgIuBZCALZ1pYgBbAbLflTzVghqp63ndFsGH3YW6bnsiejOM2laoxJiLEaZB0kg1Eil4F/FSEM1rtr0IeVYjEx8drYmL0d9mYt3oX9/1nJXVrVuPlG/rR13puG2PCRESSVDU+WJmXK4zBOC2XZgI/EPzBtAmD/ALlb59vZNI3m+kd25DJ1/ejZYNap65ojDFh4CVhtMTpbHcNcC3wMTBTVdeGM7DK7uDRHO6atZzvNu3nmgFxPPqrbtSsZhMeGWMix8t8GPk4zyo+FZGaOInjGxF5XFVfCHeAldHy7Qe5461l7M/M4ckrenDNgIrdn8QYEx08PfR2E8UlOMmiHfA88H74wqqcVJXpi1P488fraHFaLd6bcBY9YhpEOixjjAG8TdE6DegOfAI8pqqVrjlrWcg4lssf3l/FvNW7Oa9rc579TS8bPNAYU654ucK4ATgKdAbuChjUzub0DpFl2w9y54zl7DmczR8u6sptZ3ewObeNMeWOl2cYfjr3GR/y8guY/O1m/vnFJlo2qMV/xg+2yY6MMeWWn6FBTAhtP5DFve+sICnlIJf2bMVfLu9Bg9rWa9sYU35ZwihjBQXKzITt/PXj9VSpIjYkuTEmaljCKEMpB47ywHurWbzlAEM6NuGZq3rRpmHtSIdljDGeWMIoA7n5BUxduJV/frGJalWEp67owdX9Y21WPGNMVLGEEWZLthzg4Q/WsGlvJuef0ZwnRnWnVQO7qjDGRB9LGGGScuAoz3y2kY9X7SKmUW1evTGe87u1iHRYxhhTYpYwQuxA5nH+9XUyby5JoVqVKtx1Xicm/OJ0atewcaCMMdHNEkaI7Mo4xpQFW5i5dDs5eQVc3T+Oe8/vRPPTbHRZY0zFYAmjFFSVxJSDzPxhO3NXpVGgMKp3GyYMO52OzetFOjxjjAmpMk0YIjICeA5nxr1XVfWpIuXill+MM+PeTaq6zEvdsqKqbN6Xyefr9vD+sp0k782kXs1qXDsgjlvP7kBs4zqRCMsYY8KuzBKGiFQFXsSZWyMVSBCROaq6LmCzi4BO7s9AYBIw0GPdsMjOzWfzvkzW7jzM6p0ZfLdpH9sOZAHQO7Yhz1zZk0t7taJODbtYM8ZUbGX5LTcASFbVLQAiMgsYCQR+6Y8Epqszb+wSEWkoIq1whlQ/Vd2QUFUufWEh6UdzOJSVy7Hc/J/K6tWsRr+2jbjl7A6c17U5ra3TnTGmEinLhNEGZ6rXQqk4VxGn2qaNx7oAiMg4YBxAXJz/iYdEhM4t6lOtitCwTnUa1K5OXJO69GjTgLaN69gossaYSqssE0awb1r1uI2Xus5K1SnAFID4+Pig25zKP67uXZJqxhhToZVlwkgFYgOWY4A0j9vU8FDXGGNMGJXlXBcJQCcRaS8iNYDRwJwi28wBbhTHICBDVXd5rGuMMSaMyuwKQ1XzRGQi8BlO09ipqrpWRMa75ZOBeThNapNxmtWOPVndsordGGMMiNMgqWKKj4/XxMTESIdhjDFRQ0SSVDU+WJlNv2qMMcYTSxjGGGM8sYRhjDHGE0sYxhhjPKnQD71FZB+QUsLqTYH9IQwnkirKsVSU4wA7lvKoohwHlO5Y2qpqs2AFFTphlIaIJBbXUiDaVJRjqSjHAXYs5VFFOQ4I37HYLSljjDGeWMIwxhjjiSWM4k2JdAAhVFGOpaIcB9ixlEcV5TggTMdizzCMMcZ4YlcYxhhjPLGEYYwxxpNKnTBEZKqI7BWRNcWUi4g8LyLJIrJKRPqWdYxeeTiWYSKSISIr3J8/lXWMXohIrIh8LSLrRWStiNwdZJuoOC8ej6XcnxcRqSUiS0VkpXscjwXZJlrOiZdjKffnJJCIVBWR5SLyUZCy0J4XVa20P8A5QF9gTTHlFwOf4Mz4Nwj4IdIxl+JYhgEfRTpOD8fRCujrvq4P/Ah0i8bz4vFYyv15cT/neu7r6sAPwKAoPSdejqXcn5Mi8f4OmBEs5lCfl0p9haGqC4D0k2wyEpiujiVAQxFpVTbR+ePhWKKCqu5S1WXu6yPAepw53QNFxXnxeCzlnvs5Z7qL1d2foq1louWceDmWqCEiMcAlwKvFbBLS81KpE4YHbYAdAcupROF/+ACD3UvxT0TkzEgHcyoi0g7og/NXYKCoOy8nORaIgvPi3vZYAewF5qtq1J4TD8cCUXBOXP8E/hcoKKY8pOfFEsbJSZB10frXyDKcMWJ6AS8AH0Q4npMSkXrAe8A9qnq4aHGQKuX2vJziWKLivKhqvqr2BmKAASLSvcgmUXNOPBxLVJwTEbkU2KuqSSfbLMi6Ep8XSxgnlwrEBizHAGkRiqVUVPVw4aW4qs4DqotI0wiHFZSIVMf5gn1LVd8PsknUnJdTHUs0nRcAVT0EfAOMKFIUNeekUHHHEkXnZAjwKxHZBswChovIm0W2Cel5sYRxcnOAG92WBoOADFXdFemgSkJEWoqIuK8H4Jz7A5GN6kRujP8G1qvq34vZLCrOi5djiYbzIiLNRKSh+7o2cD6wochm0XJOTnks0XBOAFT1D6oao6rtgNHAV6p6fZHNQnpeqpU83OgnIjNxWkQ0FZFU4BGch2Co6mRgHk4rg2QgCxgbmUhPzcOxXAVMEJE84BgwWt1mFOXMEOAGYLV7nxngj0AcRN158XIs0XBeWgHTRKQqzpfnO6r6kYiMh6g7J16OJRrOSbHCeV5saBBjjDGe2C0pY4wxnljCMMYY44klDGOMMZ5YwjDGGOOJJQxjjDGeWMIwxhjjiSUMY4wxnljCMBEhIvnuXAMrRWSZiJxVwvdZFOrYQkFEHhWR+05S3k6CzF0iIg1F5PYi6yJ2jMXF6ZZlFrP+tyKy2z23m0XkxvBGacqKJQwTKcdUtbc7wNsfgCdL8iaq6inRuEMjlPjfe2nr+9AQ+FnC8HqM5UhP4FH33F4DFDfEi4kyljBMeXAacBBARD4QkSRxZkMb566rKyIfu3+xrhGRqwsrikhmceXuX8frReQlnBFIY0XkenFmXFshIi+7Q123E5ENIjJNnFnJ3hWROsXU/527jzUick/gQYjIgyKyUUS+ALoExLAmYJv7ROTRIvU6iDNjWn/gKeB0N77/KzzGgG1P2H9AnK+4n9vnIlI72OciIuPlvzPJbRWRrwPe+4TP3lW16Huf4nz2ADa6r7cCOafY3kSL0sy+ZD/2U9IfIB9YgTPwWwbQz13f2P1dG1gDNAGuBF4JqNsg4HVmceVAO5x5Aga5y2cAc4Hq7vJLwI3udgoMcddPBe4LUr8fsBqoC9QD1gJ9ipTVwUmAyQHvsSYgtvuARwvX4ySW5UDvgJjXFPmsMk+2f7dOXsB7vANcf4rPrTrwHXBZwLpgn33Q9w6MK8i5PQi0xhla+zFgbKT/vdlPaH7sCsNESuEtqa44w0tPFxEB7hKRlcASnGGZO+F8SZ4vIk+LyNmqmlHkvU5WnqLOTGMA5+F86Sa4gwGeB3Rwy3ao6vfu6zeBoUHqDwVmq+pRdYa/fh842y072y3LUmfOizkePoNmwIc4X8ArTrXxKfa/NeA9knC+6E/2uTyHM7rp3IB1wT774t47KBGJxZmOdh7OBEVDgdfdsrruVdwrInKdh+M15YwlDBNxqroYaAr8Gme46cHq3P9eDtRS1R/571/XT4rIn4rUP1n50YDXAkxzE1VvVe2iqo8Wvk3RsIqpf9JDCbIuj5//P6sV8DoDZza0Iad4Xy/7Px7wOh+oVtznIiI3AW1x/vrHXTeMIJ99ce99kjh6AgvUmaCoM9AVGOyWXQG8q6q3Ab86yXuYcsoShok4EekKVMX5YjqoqlnuukFueWsgS1XfBP4G9C1S/6TlAb4ErhKR5m69xiLS1i2LE5HCL7ZrgIVB6i8ARrnPN+oCl+Pc1iksu9x9dlAfuMxdvwdoLiJNRKQmcGnA++UAo3DmK7jWXXcE5y/0YE62/xME+1xEpB/ObbHrVTVwWs8GBPnsS6AHTrJBVQ8CM3DmnAZn8p7C6ULzS/j+JoIq9XwYJqJqy3/niBBgDPAFMF5EVuE8NC28FdQD+D8RKQBygQlF3utU5QCo6joReQj43G3xlAvcAewG1gNjRORlYBMwCWhepP4yEXkdWOquelVVlweUvY3zXCYF94tcVXNF5HGcuby3UmSyHlU9Ks5Um/NF5Kiqfigi37sPyj9R1d+fav/izBceTLDPZSLQGPjauQNIoqreCnxK8M/erx7AJwHLc3Fufz2IM/tbDM5nZH+sRiGbD8NUeu4X7keqWnRuZxNC7lXRv4BsYKGqvhXhkIxPdoVhjCkTqnqU8jsTn/HArjCMMcZ4YvcRjTHGeGIJwxhjjCeWMIwxxnhiCcMYY4wnljCMMcZ4YgnDGGOMJ5YwjDHGePL/oCSX38CnGHEAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "R0 = torch.linspace(1.0, 4.0)\n", "Imax = 1-1/R0*(1+torch.log(R0))\n", "pyplot.plot(R0, Imax)\n", "pyplot.xlabel('Basisreproduktionszahl $R_0$')\n", "pyplot.ylabel('Maximal gleichzeitig Infizierte $I_{max}$')\n", "R0[(Imax-0.01).abs().argmin()].item() # Kritisches R_0, siehe unten" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Schussfolgerungen\n", "\n", "*Unter den Annahmen des Modells* (eine entscheidende Einschränkung) steigt das Maximum der Anzahl der gleichzeitig Infizierten $I(t^*)$ sehr stark mit zunehmender Basisreproduktionszahl $R_0$.\n", "\n", "Wenn man davon ausgeht, dass es in Deutschland ca [25.000 Intensivbetten gibt](https://www.aerzteblatt.de/nachrichten/111029/Ueberlastung-deutscher-Krankenhaeuser-durch-COVID-19-laut-Experten-unwahrscheinlich) das sind ca. 3 pro 10000 Einwohner. Nehmen wir an, die Hälfte davon, also 1.5 pro 100000 Einwohner, stünde für CoViD-Patienten mit kritischen Verläufen zur Verfügung. Ferne seien 1.5% aller Verläufe kritisch (das ist [relativ zu den geschätzten 5% in China (untere Grafik rechts unter Diagrams)](https://en.wikipedia.org/wiki/2019%E2%80%9320_coronavirus_pandemic#Signs_and_symptoms) eine eher zu geringe Annahme, aber ggf. sind kritische Fälle auch nur zeitweise kritisch, ich weiß es nicht). Dann würde eine maximale Infektionsrate von $\\frac{1.5}{10000} / 1.5\\% = 1\\%$ kritisch. Im Modell korrespondiert dies zu $R_0 \\approx 1.15$, weit weniger als die Spanne von 1.4 - 3.9.\n", "\n", "Was also tun?\n", "\n", "### Auswirkung frühzeitiger Erkennung und Isolation der bekannten Fälle\n", "\n", "Es scheint, als sei dies in Hongkong, Singapur, Südkorea und Taiwan funktioniert.\n", "Bayern hat es am 16. März aufgegeben, weil die Infektionsketten vollkommen unklar sind.\n", "\n", "Eine spekulative Anmerkung: Eine wichtige Rolle scheinen (Ski-) Urlaube als Quelle und die Schulen als Übertragungswege gespielt zu haben. Jedenfalls spricht folgendes Zahlenverhältnis aus meiner Sicht dafür: Am 12.3. Tag vor der bekanntgabe der allgemeinen Schulschließung gab es 500 bekannte Fälle und 125 Schulschließungen aufgrund mindestens einem bekannten Fall an den jeweiligen Schulen. Schüler nach den Ferien eine Woche in die Schule gehen zu lassen und dann (erst) unter Quarantäne zu Stellen, wie im Fall von Südtirolurlaubern geschehen, ist vermutlich nicht optimal gewesen.\n", "\n", "Wie würde sich als (grob) eine frühzeitige Erkennung auswirken?\n", "Ist die Kranheitsdauer $T_I$, man diagnostiziert aber schon nach $T_D$ Tagen und erfolgreich isoliert, steckt jeder Kranke nur noch $\\tilde R_0 = \\beta T_D$ statt $R_0 = \\beta T_I$ an.\n", "\n", "Wenn wir $I$ als Anteil der nicht-diagnostiziert Infizierten umwidmen, können wir mit diesem niedrigeren $\\tilde R_0$ bzw. $\\tilde \\gamma = 1/T_D$ rechnen. Es gibt dann $\\frac{T_I}{T_D} I$ Infizierte, von denen $I$ nicht diagnostiziert sind.\n", "\n", "### Allgemeine Isolation\n", "\n", "Nun versuchen wir also, wie andere Länder auch, durch Reduktion der Kontakte insgesamt eine Reduktion zu erzielen.\n", "Im Modell könnte man annehmen, dass eine Reduktion der Kontakte um einen bestimmten Faktor dieselbe Reduktion in $\\beta$ bzw. $R_0$ erzielt. Nur, wie einfach ist es, Kontakte um einen Faktor 3 bis 4 zu reduzieren?\n", "\n", "Wir werden es sehen, hoffen wir das beste.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Grenzen des Modells\n", "\n", "Es gibt einige Fragen, die das Modell nicht beantworten kann:\n", "\n", "Welcher Anteil der Bevölkerung wird erkranken? Im Modell steckt (bei $R_0 > 1$) die implizite Annahme, dass die gesamte Bevölkerung im laufe der Zeit erkrankt. Das bedeutet, dass man fragen muss, ob $R_0$ veränderlich ist (das ist ja auch in den Maßnahmen letztlich das Ziel) und unter $1$ gedrückt werden kann.\n", "\n", "Einige Auslassungen schmälern die Aussagekraft des Modells nicht: Das Modell beinhaltet keine krankheitsbedingte Sterblichkeit. Kennt man aber die Sterblichkeitsrate $q$, gibt es (bis auf zeitliche Ungenauigkeit) $q R$ Todesfälle und $(1-q) R$ Genesene. Ähnlich funktioniert ja auch die obige überlegung zum Nutzen umfassender Verfolgung der Fälle.\n", "\n", "Auch kann man das Modell aufgrund der impliziten Annahmen hinterfragen:\n", "- Dem Modell fehlt die räumliche Dimension. Insbesondere ist die Ansteckung unabhängig von der räumlichen Entfernung von Individuen jenseits der mittelbaren Auswirkung über $\\beta$ bzw. $R_0$.\n", "- Eine Dunkelziffer ist nicht vorgesehen. Jenseits sehr pauschaler Annahmen ist sie auch schwer einzubauen. Vermutlich ist es auch deshalb schwierig, das Modell an Fallzahlen zu kalibrieren.\n", "- Die stochastische Dimension wird nicht beachtet. Das ist für ein so grobes Modell wie dieses vermutlich unerheblich, ist dies aber gegebenenfalls für verfeinerte Modelle (z.B. mit räumliche Dimension, Dunkelziffer) könnte es aber eine größere Rolle spielen.\n", "\n", "Hier kommen wir wieder an einen Punkt, bei dem ich darauf hinweisen muss, dass ich kein Epidemiologe bin." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Anhang: Numerische Lösung der Gleichungen in Python\n", "\n", "Wir können die Gleichungen numerisch lösen\n", "Dafür gehen approximieren wir die Ableitungen mit Differenzenquotienten mit einem Zeitschritt $\\tau$, also zum Beispiel\n", "\n", "$$\n", "\\dot S(t) \\approx \\frac{S(t + \\tau) - S(t)}{\\tau}.\n", "$$\n", "\n", "Damit ergeben sich die Vorschriften (wieder rechnen wir $I(t)$ durch das geschlossene System aus).\n", "\n", "$$\n", "S(t + \\tau) = S(t) - \\tau \\beta \\frac{S(t) I(t)}{N},\n", "$$\n", "$$\n", "R(t + \\tau) = R(t) + \\tau \\gamma I(t)\n", "$$\n", "und\n", "$$\n", "I(t + \\tau) = N - S(t + \\tau) - R(t + \\tau).\n", "$$\n", "\n", "Diese Approximation heißt explizites Euler-Verfahren zur Lösung des Systems von Differentialgleichungen.\n", "Wir können sie programmieren:" ] }, { "cell_type": "code", "execution_count": 80, "metadata": {}, "outputs": [], "source": [ "N = 100000\n", "I_0 = 10\n", "T = 100\n", "beta = 0.45\n", "gamma = 0.1\n", "\n", "tau = 0.01\n", "\n", "t = torch.arange(0, T, tau)\n", "S = torch.full_like(t, N - I_0)\n", "I = torch.full_like(t, I_0)\n", "R = torch.zeros_like(t)\n", "\n", "for i in range(t.size(0) - 1):\n", " S[i + 1] = S[i] - tau * beta * I[i] * S[i] / N\n", " R[i + 1] = R[i] + tau * gamma * I[i]\n", " I[i + 1] = N - S[i + 1] - R[i + 1]" ] }, { "cell_type": "code", "execution_count": 86, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 86, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAD4CAYAAADy46FuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3xV5f3A8c/3jgyyByMkgTCC7BlExa0IioqjKlr3qq11tPbn6NRaW22tq1qFVuuoiIoDlVERFRcyZc+QMMIMIQnZuTf3+f1xTiBAyLzJTW6+79fv/s45z3mec77Hhvu953nOEGMMSiml1PE4Ah2AUkqptk0ThVJKqTppolBKKVUnTRRKKaXqpIlCKaVUnVyBDsDfEhMTTVpaWqDDUEqpdmXZsmX7jTGda1sXdIkiLS2NpUuXBjoMpZRqV0Rk2/HWadeTUkqpOmmiUEopVSdNFEoppeoUdGMUSinVXB6Ph5ycHMrLywMdit+FhYWRkpKC2+1ucBtNFEopdZScnByioqJIS0tDRAIdjt8YY8jLyyMnJ4devXo1uF29XU8i8oqI7BORNTXK4kVknohstqdxNdY9JCKZIrJRRMbXKB8lIqvtdc+J/V9fREJF5G27fJGIpNVoc4O9j80ickODj0oppZqhvLychISEoEoSACJCQkJCo8+UGjJG8Sow4aiyB4H5xph0YL69jIgMBCYDg+w2/xQRp93mReB2IN3+VG/zFiDfGNMXeBp4wt5WPPAHYAxwIvCHmglJKaVaUrAliWpNOa56u56MMV/V/JVvmwScac+/BnwJPGCXTzfGVADZIpIJnCgiW4FoY8xCO9DXgUuAOXabh+1tzQCet882xgPzjDEH7DbzsJLLW40+ygYorfTy0pdbEBEcIjgdHDHvEEFEcAo4HFa59am5DE6H4HQIYS4nYW4nYW6HPT08HxnqIsztrD8opZRqA5o6RtHVGLMbwBizW0S62OXJwPc16uXYZR57/ujy6jY77G15RaQQSKhZXkubI4jI7VhnK/To0aNJB1RaWcU/vsiktV7PEepyENcphNhObuI6hRAfEUL32DBS4zuREhdOalwn0hIjcDv1wjSlOqrHHnuMadOm4XQ6cTgcTJkyhTFjxrR6HP4ezK7tnMbUUd7UNkcWGjMVmAqQkZHRpK/6xMhQsv8yEWMMPgM+Y6yPz5qvMgZTY95nDMZAle/Yea/PUO6potzjs6dVlHmqqPD4KPdWUVzhpaDUQ35JJfmlHgpKK1m/+yCfrd9Lhdd3KKYQp4MTukUxqHs0I3vEcWp6It1jw5tyeEqpdmbhwoV88sknLF++nNDQUPbv309lZeVx61e/hK4lusyamij2ikiSfTaRBOyzy3OA1Br1UoBddnlKLeU12+SIiAuIAQ7Y5Wce1ebLJsbbYNXdS85a81TL8vkM+4sr2JFfxvYDJazfXcTaXYXMXbuH6Uusk6venSOYMKgbl41Mpm+XqFaPUan2zOvzUllVSXlVuTX1llNRVXHMJ84bR355PgZz6Au4ev6IKQbr/44sq61N9XJ1/WpHl1UvL920lLCYMLJLsqEEEDChhsK8wsPtj9pWuCuc3rG9/f7fTRryKlR7jOITY8xge/lvQJ4x5nEReRCIN8bcLyKDgGlYg8/dsQa6040xVSKyBLgLWATMBv5hjJktIncCQ4wxd4jIZOAyY8yV9mD2MmCkHcZyYFT1mMXxZGRkmGB71pMxho17i/hm834WbMrl28z9+AwMS43lJ6f3ZvygbjgdwTnwpjo2Ywyl3lIKKwo5WHnw0LS4sphSbymlnlJKPCWUeq1pmbeMEk/JobJSTyll3jLKvVZi8Bpvg/b7zMBn6NarGwBTP88ja18tv+Rr/HKXGv//mLKj/mkK0LdrOHeP615jU9W1D09LSkq49LxLKSst47SzTmPS5ZM4+bSTj6hzeBfWvNvhJi6s/mt+1q9fz4ABA446HFlmjMmorX69ZxQi8hbWL/tEEcnBuhLpceAdEbkF2A5cAWCMWSsi7wDrAC9wpzGmyt7UT7GuoArHGsSeY5e/DLxhD3wfwLpqCmPMARF5FFhi1/tjfUkiWIkI/btF079bNLee1pt9ReV8vHI3//1+Gz97czm9EyN48Pz+jBvYNWiv1FDBwWd8FFQUsL9sP3lleYemeeV55JXlUVBRcERCOFhxsN4vd5e46OTuRCd3JyJcEUS4Iwh3h5MQlmDNu8IJdYUS6qzlc5xy7x4v6XHpCEJc2AbC3UW1JoOmigqJokd0PeOpUbDqh1V8/fXXfPHFF/zsxp/x+OOPc+ONN/olhsZo0BlFexKMZxTHU+Uz/G/tHp6et4nN+4o5o19n/nLZEB3HUAFhjOFA+QF2l+y2PsW7D8+X7Ca3NJcD5QeoOvTb8bBQZygJYQnEhMYQExpDdEj0sfMhMUSHRhMdEk2E20oIndydCHGE+P0HUm2/uANtxowZvPbaa3z88cfN3pbfzyhU2+V0CBcMSWLcwK68sXAbf/90I+c/+zVPXD6UCYO7BTo8FYSMMeSV57Ht4Da2HdzG1oNb2VZozecU51BRVXFE/XBXON0jutMtshv94/uTEJZAQngCieGJJIRZ08TwRCLcEXo2fJSNGzficDhIT08HYMWKFfTs2TMgsWiiCAJup4ObT+3F2f27cPf0H7jjv8v41Xn9uPOsvvqPTzVZiaeEzfmb2ZS/iY0HNrIpfxOZBZkUe4oP1XE73PSI6kHP6J6cmnwq3SO7kxSRRFJkEkkRSUSHROvfYBMVFxdz1113UVBQgMvlom/fvkydOjUgsWiiCCJpiRHMuOMUHnhvFU9+uoldheX8adJgHDrQrepRUVXB+rz1rMpdxar9q1i7fy05xYdvfYpyR5Eel87E3hPpFdOLtOg0ekb3JCkiCadDbx5tCaNGjeK7774LdBiAJoqgE+Jy8NSVw+gaHcZLC7YQ4nTwh4sG6q86dYTCikKW7FnC0r1LWZW7ivUH1uP1WYPG3SO6MyhxEJf0vYR+cf04If4EkiKS9G+oA9NEEYREhAcmnIC3yse/v8kmrlMI95ybHuiwVACVekpZvm85i3YvYtHuRWw4sAGDIdwVzqCEQVw/8HqGdh7K0MShdO5U62uTVQemiSJIiQi/mTiAA6WVPP3ZJgYkRXHeIB3g7kj2luxlQc4CvtzxJYt2L6LSV4nb4WZ4l+H8bPjPGJM0hsGJg3E7Gv5eAtUxaaIIYiLCny8dQua+Yn75zko++nkkvTtHBjos1YJyinKYkz2Hz7Z/xrq8dQCkRKZw5QlXclrKaYzoMoJwl14+rRpHE0WQC3M7eenaUZz/7Nf86t2VvHvHKXoXd5DJK8vj022fMitrFitzVwIwtPNQ7hl5D2elnkXvmN46vqCaRRNFB9A9Npw/ThrEPdNXMPWrLH56Zp9Ah6SayWd8LNy1kHc3vcuXO76kylTRL64f9468l/N7nU/3yO71b0SpBtJE0UFcPKw7c1bv4enPNnHh0CRS4zsFOiTVBHlleXyQ+QEzNs1gZ/FO4sPiuX7g9VzU5yLS4/SChWASGRlJcXFx/RVbgSaKDkJE+MPFA1nwZC6PzVrPS9eNCnRIqhG2HdzGa2tfY2bmTCp9lYzuNpp7R97L2T3OJsQZEujwVJDTRNGBJMWEc+dZfXjy0018l7mfU/omBjokVY81+9fwyppX+GzbZ7gdbi7qcxHXD7y+RR4lrdTxaKLoYG49rTdvLd7BE//byId9gu/l8cFiU/4mnv/heb7Y8QVRIVHcOuRWrhlwDYnhmtxb3ZwHYc9q/26z2xA4/3H/brMFaaLoYMLcTu46uy8Pvr+aLzfmclb/LvU3Uq1m+8HtPL/ieeZmzyXSHcnPh/+cawdeS4Q7ItChqQ5ME0UHdNnIFJ7/IpNnPtvEmSd01rOKNqC4spipq6byxvo3cDvc3Dz4Zm4afBMxoTGBDk21o1/+LUUTRQcU4nJw51l9eej91XyTuZ/T0vWRDYHiMz4+yfqEp5c9zf6y/VzS9xLuGXmPdjGpNsUR6ABUYFw2MpnEyBD+8+3WQIfSYW0t3MpNc2/iN9/8hu4R3Zl2wTQeHfuoJgnV5ugZRQcV6nJyzZiePDd/M9n7S+iVqH3grcXr8/LGujd4YcULhDhD+OMpf2RS30k4RH+3qcPayj0UoGcUHdq1J/XA7RRe+25roEPpMLILs7lu9nU8tewpxnYfy8xJM7k0/VJNEqpN07/ODqxLVBgXDu3OjGU5lFUe+x5j5T/GGD7Y/AFXfXIVOcU5/O2Mv/HMWc/oI71Vu6CJooO7MiOV4govc9fuDnQoQauosogHvn6A33/3ewYnDmbGRTOYkDZBrzZT7YYmig5uTK94UuPDeXdpTv2VVaNl5mcy+ZPJfLr1U+4acRf/GvcvukZ0DXRYSjWKJooOzuEQfjQyle+25LHjQGmgwwkqn2//nB/P/jGl3lJeGf8Ktw+9Xd8vrdolTRSKy0clIwLvL98Z6FCCgs/4eHHli9zzxT30junN9InTGdl1ZKDDUqrJNFEoUuI6MTotnlmrdwU6lHbPU+Xh19/8mn+u+CcX97mYV89/VbuaVJM4nU6GDx/O4MGDueiiiygoKAhYLJooFAAThySxaW8xmfuKAh1Ku1XiKeFn83/GrKxZ3D3ibv409k+EOkMDHZZqp8LDw1mxYgVr1qwhPj6eF154IWCxaKJQAJw/uBsiMGvVnkCH0i7tL9vPTXNvYsmeJTw69lFuG3qbXtWk/Obkk09m587AdQ3rndkKgC7RYYzuGc/s1bu551x9U1pj7C3Zyy2f3sK+0n384+x/cFrKaYEOSfnRE4ufYMOBDX7dZv/4/jxw4gMNqltVVcX8+fO55ZZb/BpDY+gZhTrkgiHd2Li3iMx9befRAW3dnpI93PS/m9hftp+p46ZqklB+U1ZWxvDhw0lISODAgQOMGzcuYLHoGYU6ZMLgJB7+eB2frd9L3y6RgQ6nzdtdvJub/3czBRUFTBk3hWGdhwU6JNUCGvrL39+qxygKCwu58MILeeGFF7j77rsDEoueUahDusWEMTApms837At0KG3evtJ93PS/myisKGTquKmaJFSLiYmJ4bnnnuPJJ5/E4/EEJIZmJQoR+YWIrBWRNSLyloiEiUi8iMwTkc32NK5G/YdEJFNENorI+Brlo0Rktb3uObFHAUUkVETetssXiUhac+JV9Tu7fxeWbcunsDQwf5DtQWFFIT+Z9xPyy/OZMm4KQzoPCXRIKsiNGDGCYcOGMX369IDsv8mJQkSSgbuBDGPMYMAJTAYeBOYbY9KB+fYyIjLQXj8ImAD8U0Sqb1N9EbgdSLc/E+zyW4B8Y0xf4GngiabGqxrmrP5dqPIZFmzODXQobVKpp5Sfz/852w5u47mzn9MkoVrM0Y8Z//jjj7nuuusCEktzu55cQLiIuIBOwC5gEvCavf414BJ7fhIw3RhTYYzJBjKBE0UkCYg2xiw0xhjg9aPaVG9rBnCO6DWHLWp4aizxESF8od1Px/D4PNy34D5W7V/FE6c/wZikMYEOSalW0eREYYzZCTwJbAd2A4XGmE+BrsaY3Xad3UAXu0kysKPGJnLssmR7/ujyI9oYY7xAIZBwdCwicruILBWRpbm5+ku4OZwO4Yx+nfly4z6qfCbQ4bQZxhj+vOjPfLPzG3570m8Z1zNwV6Ao1dqa0/UUh/WLvxfQHYgQkWvralJLmamjvK42RxYYM9UYk2GMyejcWZ/v31xnntCZ/FIPq3cWBjqUNuPN9W8yY9MMbhl8C1f0uyLQ4SjVqprT9XQukG2MyTXGeID3gVOAvXZ3Eva0ug8jB0it0T4Fq6sqx54/uvyINnb3VgxwoBkxqwYY29d6Z/O3mfsDHEnb8HXO1/xt6d84O/Vs7h4ZmMsTlQqk5iSK7cBJItLJHjc4B1gPfATcYNe5AZhpz38ETLavZOqFNWi92O6eKhKRk+ztXH9Um+pt/Qj43B7HUC0oMTKU/t2i+G6LJootBVv4v6/+j35x/fjLaX/RV5aqDqnJN9wZYxaJyAxgOeAFfgCmApHAOyJyC1YyucKuv1ZE3gHW2fXvNMZUv3/zp8CrQDgwx/4AvAy8ISKZWGcSk5sar2qcU/ok8uaibZR7qghzd8x3KJR4Srj3i3sJc4bxj7P/QSd3p0CHpFRANOvObGPMH4A/HFVcgXV2UVv9x4DHailfCgyupbwcO9Go1jW2bwKvfJvN8m35nGJ3RXUkxhge/u5hthdt59/n/ZtuEd0CHZJSAaPn0apWJ/aKx+kQvu2g3U/TN05n7ta53DXiLkZ3Gx3ocJQKKE0UqlZRYW6Gp8bybWZeoENpdatzV/PXJX/ljJQzuHnwzYEOR3VgU6ZM4Y477gh0GJoo1PGN7ZPAqpwCiso7zuM8Sjwl3P/V/XQJ78Jjpz6mg9cqoFatWsWQIYG/+1//FajjGtM7AZ+B5dsD9wrG1vbE4ifYVbKLv5z2F2JCYwIdjurgVq9ezdChQwMdhj5mXB3f8NRYnA5h6dYDnNEv+G9knL99Ph9kfsCtQ25lZNeRgQ5HtRF7/vxnKtb798VFoQP60+3Xv6633po1a/SMQrVtEaEuBnWPZnF28N/jmFuay8PfPcyA+AH8bNjPAh2OUuzYsYPIyEhiY2MDHYqeUai6jU6L57/fb6PS6yPEFZy/K4wxPLzwYcq8ZTx+2uO4ne5Ah6TakIb88m8J1eMTe/bs4aqrrmLixImsXbuWU045hXnz5vHwww8zePBgHnnkEQ4cOEBsbCyPPPIIubm53H///Tz66KP8/ve/Z8qUKbjdzfubDs5/+cpvRqfFUeH1BfVzn2Znz+arnK+4Z+Q99I7tHehwlAIOj0/88MMPXHbZZdx///0UFhZy2223ccUVV7Bt2zZ27tyJx+MhNjaW77//HoDOnTvTo0cP7rvvPp577rlmJwnQRKHqMapnPABLtwZn91N+eT5PLH6CoYlDuab/NYEOR6lDVq9ezZAhQ1ixYgXjx4/H4/GQkJCAw+E4NHbxu9/9jgceeIAbbriB5GTrodvFxcVkZWXhcrmIjPTPK401Uag6dY4KpXdiBEuCNFH8dclfKfIU8fApD+N0dMxHlai26c033+Saa64hMzOTfv36sXbtWgYMGADA1q1b6dGjB4MGDeLJJ5/kqaeeYsSIEXi9Xu6++27+9Kc/MXz4cL788ku/xCLB9oy9jIwMs3Tp0kCHEVTun7GST9ftZflvx+FwBM97o77Z+Q0//eyn3DHsDu4cfmegw1FtyPr16w99KQej2o5PRJYZYzJqq69nFKpeGWnxFJR62JJbXH/ldqLMW8ajCx+ld0xvbhtyW6DDUapN00Sh6jWyRxwAPwTRjXevrHmFXSW7+N1JvyPEGRLocJRq0zRRqHr1TowgKszFDzuCI1HsKNrBK6tf4fxe55PRrdYzbaVUDZooVL0cDmF4aiwrgiRR/G3J33A6nNw36r5Ah6LasGAbv63WlOPSRKEaZERqLBv3HKS00hvoUJrl253f8sWOL7h96O10jega6HBUGxUWFkZeXl7QJQtjDHl5eYSFhTWqnd6ZrRpkeI9YfAZW5RRyUu+EQIfTJJ4qD48vfpye0T25fuD1gQ5HtWEpKSnk5OSQm5sb6FD8LiwsjJSUlEa10UShGmRYivW8mRU7Ctptonh307tsPbiV589+XgewVZ3cbje9evUKdBhthnY9qQZJiAylR3wnVrTTK5+KK4uZsmoKo7uN5vSU0wMdjlLtip5RqAYbnhrbbp8k++raVzlQfoAXRr2ASPDcNKhUa9AzCtVgw1Nj2XOwnD2F5YEOpVFyS3N5fd3rjE8bz+DEwYEOR6l2RxOFarARParHKfIDHEnjvLjyRTxVHu4ecXegQ1GqXdJEoRpsYPdoQpyOdnXjXXZhNu9vfp8rTriCHtE9Ah2OUu2SJgrVYKEuJwO6R7OyHSWKl1a+RIgzhJ8M/UmgQ1Gq3dJEoRplSHI0a3cexOdr+zciZRVmMSd7Dlf3v5qE8PZ5Sa9SbYEmCtUoQ5JjKKrwsjWvJNCh1OullS8R5grjhkE3BDoUpdo1TRSqUYYkWwPabf3VqFkFWczNnsvV/a8mPiw+0OEo1a5polCNkt41khCXgzVtPFFUn03cOOjGQIeiVLuniUI1itvpYEBSNKty2m6i2FKwhblb53JN/2uIC4sLdDhKtXuaKFSjDU2OYe2utjug/a/V/9KxCaX8qFmJQkRiRWSGiGwQkfUicrKIxIvIPBHZbE/jatR/SEQyRWSjiIyvUT5KRFbb654T+xkLIhIqIm/b5YtEJK058Sr/GJIcQ3EbHdDeWbyTudlzuaLfFXo2oZSfNPeM4llgrjGmPzAMWA88CMw3xqQD8+1lRGQgMBkYBEwA/ikiTns7LwK3A+n2Z4JdfguQb4zpCzwNPNHMeJUfDE6OAdrmgPbra19HEK4beF2gQ1EqaDQ5UYhINHA68DKAMabSGFMATAJes6u9Blxiz08CphtjKowx2UAmcKKIJAHRxpiFxnpLyOtHtane1gzgHNEnugVcetdIQl0OVrexcYr88nze3/w+F/S+gG4R3QIdjlJBozlnFL2BXOA/IvKDiPxbRCKArsaY3QD2tItdPxnYUaN9jl2WbM8fXX5EG2OMFygEjrlzSkRuF5GlIrI0GF800tZUD2i3tTOK6RumU15Vzs2Dbw50KEoFleYkChcwEnjRGDMCKMHuZjqO2s4ETB3ldbU5ssCYqcaYDGNMRufOneuOWvnFkDY2oF3mLWPahmmcmXImfWL7BDocpYJKcxJFDpBjjFlkL8/AShx77e4k7Om+GvVTa7RPAXbZ5Sm1lB/RRkRcQAzQPl+IEGSGpFgD2tltZED7g80fUFBRwE2Dbwp0KEoFnSYnCmPMHmCHiJxgF50DrAM+AqqvS7wBmGnPfwRMtq9k6oU1aL3Y7p4qEpGT7PGH649qU72tHwGfm2B723k7NcQe0G4LN95V+ap4fd3rDO88nJFdRwY6HKWCTnPfcHcX8KaIhABZwE1YyecdEbkF2A5cAWCMWSsi72AlEy9wpzGmyt7OT4FXgXBgjv0Ba6D8DRHJxDqTmNzMeJWfpHc5PKA9aXhy/Q1a0IKcBews3smvMn4V0DiUClbNShTGmBVARi2rzjlO/ceAx2opXwoc8+oxY0w5dqJRbYvL6WBg97YxoP3m+jdJikjizNQzAx2KUkFJ78xWTdYWBrQ35W9i8Z7FTO4/GZdDXwGvVEvQRKGabHBy4Ae0p62fRpgzjMvTLw9YDEoFO00UqskCPaBdUF7ArKxZTOw9kZjQmIDEoFRHoIlCNVn1gHagniT73ub3KK8q55oB1wRk/0p1FJooVJMFckDb6/MyfeN0Tux2Iv3i+rX6/pXqSDRRqGYZkhzDugAMaH+540v2lOzRswmlWoEmCtUsgRrQfnfTu3Tt1JUzUs5o1f0q1RFpolDNMjTFfuR4K45T5BTl8N2u77g8/XK9JFapVqCJQjVL386RhLkdrTpO8d7m93CIg0vTL221fSrVkWmiUM3iauVHjnt8Hj7Y/AGnJ5+u75xQqpVoolDNNiQ5hrU7C1tlQPvLHV+SV57HFSfok12Uai2aKFSzDUmOoaSyiqz9LT+g/e7Gd+kW0Y2x3ce2+L6UUhZNFKrZhqS0zh3aO4p2sHD3Qi5Lvwynw1l/A6WUX2iiUM3WWgPa722yB7H76iC2Uq1JE4VqNpfTwcCk6Ba9RNbj8/Bh5oecnqKD2Eq1Nk0Uyi+sR4633ID2NznfkFeep0+JVSoANFEovxjcwgPaM7fMJCEsgbHJOoitVGvTRKH8onpAe/XOAr9v+0D5ARbsWMCFvS/E7XD7fftKqbppolB+cWhAO+eg37c9J3sOXuPl4r4X+33bSqn6aaJQflE9oN0Sl8h+mPkhAxMG6uPElQoQTRTKb6oHtKv8OKC94cAGNhzYwKQ+k/y2TaVU42iiUH5TPaCdvb/Yb9ucmTkTt8PNBb0u8Ns2lVKNo4lC+c3QlFgAv91456nyMCtrFmemnklsWKxftqmUajxNFMpv+nSO8OuA9tc7vya/Ip9L+l7il+0ppZpGjGndV1i2tIyMDLN06dImtd3z5z9TsX6DnyPqWNbsKkQQBnWPPn4l44PyAqgoBlMFrjAIj7OmNWwuyKTEU8ywzsMQpIUjV6r9Cx3Qn26//nWT2orIMmNMRm3r9PVgyq8iQ13sK6rAQO1f7aX74UA2eCusZXFYiQOgUwLEpYE7HI/PS2FFAV07ddUkoVSAaaKooamZWB22ZFkO//fuSj775en07RJ1eIUxMO938N1MSBoOZz4Efc4CZwgcyILV78L3/wTvdhj3KG9FR/Hnxat57+Ip9NTLYpUKKE0Uyq+GJFffoV14ZKL47GH47h8w+laY8Dg4a9xhndAHznwQRt0EH90Fc/6PWelDSI9N13snlGoDdDBb+VX1gPaqmk+SXTkdvn0GMm6BC548MknUFNUVrp7OjjG3sdJbyAUlZeCrap3AlVLHpYlC+dUxd2jnbYFZ90HPsXDB30DqGW9wOJiTbJ1FXJC1CGb+HHy+Fo5aKVUXTRTK74amxLJ210GqqnzwyS/A4YTLplrTehhjmJU1i5FdRtL91Pth5TT48i+tELVS6nianShExCkiP4jIJ/ZyvIjME5HN9jSuRt2HRCRTRDaKyPga5aNEZLW97jkR62eniISKyNt2+SIRSWtuvKrlDU+NpbSyil0L34HsBXDWbyEmpUFtNxzYQFZhFhN7T4QzHoAR18JXf4VV77Rw1Eqp4/HHGcU9wPoayw8C840x6cB8exkRGQhMBgYBE4B/ikj1T8wXgduBdPszwS6/Bcg3xvQFngae8EO8qoWN6BGLGy9x3z4CXQZBxs0Nbjs7ezYucXFez/OsbqqJT0PPU2HmnbBzeQtGrZQ6nmYlChFJASYC/65RPAl4zZ5/DbikRvl0Y0yFMSYbyAROFJEkINoYs9BYd/+9flSb6m3NAM6pPttQbVeP+E5cH/4dkWW7YNwj4GzYxXVVvipmZ8/m1ORTDz+ywxUCV70BEV3g3RuhzP/vu1BK1a25ZxTPAPcDNUcbuxpjdgPY03solwoAACAASURBVC52eTKwo0a9HLss2Z4/uvyINsYYL1AIJBwdhIjcLiJLRWRpbm5uMw9JNZf4vNzh/JANjnToe26D2y3bu4x9pfusbqeaOsXDFf+BgzutM4sge5qAUm1dkxOFiFwI7DPGLGtok1rKjncDb/U3QV3rDhcYM9UYk2GMyejcuXMDw1EtZvW7dPbu4a/lkygo8zS42ezs2XRydeKM1DOOXZl6Ipz7MGz4BBZN8VuoSqn6NeeMYixwsYhsBaYDZ4vIf4G9dncS9nSfXT8HSK3RPgXYZZen1FJ+RBsRcQExwIFmxKxamjHw/YuUxqTzuW8EP+xoWFdRZVUln277lHN6nEO4K7z2Sif/HPqdD/N+D/vW115HKeV3TU4UxpiHjDEpxpg0rEHqz40x1wIfATfY1W4AZtrzHwGT7SuZemENWi+2u6eKROQke/zh+qPaVG/rR/Y+tN+hLctZAntW4TzpJzhE+GFbfoOafZ3zNUWVRcd2O9UkAhf/A0Kj4P3bwFvpp6CVUnVpifsoHgfGichmYJy9jDFmLfAOsA6YC9xpjKm+7fanWAPimcAWYI5d/jKQICKZwC+xr6BSbdjif0FoNKEjr+aEbtENPqOYlT2L+LB4xiSNqbtiZGe4+DnYsxoWPO6HgJVS9fHLs56MMV8CX9rzecA5x6n3GPBYLeVLgcG1lJcDV/gjRtUKinNh7Qcw+hYIjWRkj1g+WrELn8/gcBz/YrWiyiIW7FjAj/r9CJejAX+S/SfC8Gvhm6chfTz0qCe5KKWaRe/MVv6zajr4PNYznYARPeIoqvCSmVv3q1E/2/YZlb7KurudjjbhL9ZNfB/cbr3XQinVYjRRKP9ZOR2SM6Cz9aymkT2seyF+2F73OMWs7FmkRqUyJHFIw/cVFg2XvAT526zBbaVUi9FEofxjz2rYuwaGTT5U1CsxgthObpbVMaCdW5rL4t2LuaDXBTT6Xsq0sXDynbD0Zcic39TIlVL10ESh/GPldHC4YfDlh4pEhIye8SzOPv4VzXO3zsVguKD3BU3b79m/hcR+1nss9K5tpVqEJgrVfFVe66F9/cZbd1HXMKZXPFvzStl7sLzWprOzZjMgfgC9Y3o3bd/ucKsLqmgPzH2oadtQStVJE4VqvqwvoGQfDLv6mFVjeluJY1EtZxXbD25nTd4aLujVxLOJaimj4LRfWo8k3zCredtSSh1DE4VqvjXvQ2gMpJ93zKqBSdFEhrpYnJ13zLrZ2bMBmNBrwjHrGu30+6HrEPj4Hig5dl9KqabTRKGax1sJG2dB/wusJ70exeV0MKpnHIuyjjyjMMYwJ3sOo7qOoltEt+bH4QqBS1+yxilm/UIfHKiUH2miUM2T/RWUF8LAScetcmKveDbvKyavuOJQ2ab8TWQVZjW/26mmboPhrIdg3UxY857/tqtUB6eJQjXPug8hJAp6n3XcKifZ4xRLth4+q5iVPQuXuBjXc5x/4znlHutejln3WQPcSqlm00Shmq7KYz32+4QJ4A47brUhybGEuR2HBrR9xsfc7Lmc3P1k4sLijtuuSZwuqwvKWw4f3a1dUEr5gSYK1XRbv4Gy/Dq7nQBCXA5G9jg8TrFi3wp2l+zm/F7nt0xcienWuys2/w9++G/L7EOpDkQThWq6dTPBHdGgt9id3DuB9XsOkldcwezs2YQ5wzi7x9ktF9uJP7HetT33ISjY3nL7UaoD0EShmsZXBes/hn7nWTe91ePU9ESMga8y9/Lp1k85I/UMItwRLRefwwGXvAAY6/WpPl+9TZRStdNEoZpm23dQur/ebqdqQ1NiiQ5zMXPDl+RX5Ldct1NNcWkw/jHryqwl/275/SkVpDRRqKZZNxNc4dC3YVctOR3CqemJrMz/kih3FKcln9bCAdpG3mB1jc37PeRtaZ19KhVkNFGoxvP5rKud+p4DoZENbnZSn2g8oSs5sesZhDiPvTmvRVS/PtUVAh/cYXWZKaUaRROFarxdy6FoNwy4qFHN3JEbEGclsVWt/Ea66O5w/t8gZzEsfL51961UENBEoRpv/cfgcFlPi22ERbnzkapotu1KaqHA6jD0Suh/IXz+J9i3vvX3r1Q7polCNY4xVrdT2qkQ3vCb5Yoqi/g652t6h5/C91vyKfe0cheQCFz4DIRGwQc/sW4WVEo1iCYK1Ti5GyEv0/p13gjzts2j0lfJpPQLKfNUsTArAE94jexsJYvdK+GLx1p//0q1U5ooVONs+Nia9p/YqGYfb/mYtOg0Jg8dS6cQJ5+t29sCwTXAwIth5PXwzdOweV5gYlCqndFEoRpn/SfWQ/eiuze4yc7inSzdu5QLe19IeIiL09M789n6vfh8AXoO0/l/ha6D4f3boTAnMDEo1Y5oolANV7ADdq+AAY3rdpqVZb117sI+VrtxA7uy92AFa3YV+j3EBnGHwxWvQVUlzLhZxyuUqocmCtVw1a8Z7d/wy2KNMXy85WMyumaQHJkMwFn9u+AQmBeo7ieAxL5w0bOwYxF8/mjg4lCqHdBEoRpuwyfQub/1JdtAq/evZuvBrVzU53ByiY8IISMtPrCJAmDIjyDjZvj2WdgwO7CxKNWGaaJQDVOSB9u+bfRNdh9v+ZhQZ+gxLyg6b2BXNuwpIiu32J9RNt74v0DScGu8Yt+GwMaiVBuliUI1zKY5YHyNuizWU+VhztY5nJ16NlEhUUesmzg0CRH4eOVuf0faOO4wmDzNGreYfrX1fg2l1BE0UaiGWf8JxPSApGENbvL1zq8prCg8NIhdU1JMOKPT4vlo5U5MoN9CF5MMV71hDdbPuBmqvIGNR6k2RhOFql9FEWz53Lp3QqTBzT7a8hHxYfGc0v2UWtdfPKw7W3JLWLf7oL8ibboeJ8GFT1nH+dkfAh2NUm1KkxOFiKSKyBcisl5E1orIPXZ5vIjME5HN9jSuRpuHRCRTRDaKyPga5aNEZLW97jkR69tIREJF5G27fJGIpDX9UFWTbZwLVRUw6JIGN9lftp8FOxZwYe8LcTlctda5YEgSLofw0cpd/oq0eUZeDyfebj04cOl/Ah2NUm1Gc84ovMB9xpgBwEnAnSIyEHgQmG+MSQfm28vY6yYDg4AJwD9FxGlv60XgdiDd/kywy28B8o0xfYGngSeaEa9qqrXvQ1R3SDmxwU0+3vIxXuPl8vTLj1snPiKEU9MT+WTl7sDdfHe08X+23rEx6z7Y9L9AR6NUm9DkRGGM2W2MWW7PFwHrgWRgEvCaXe01oPpn6CRgujGmwhiTDWQCJ4pIEhBtjFlorM7q149qU72tGcA51WcbqpWUF0LmZ9bZhKNhfy7GGN7f/D4ju4ykd2zvOuteNjKFnQVlfJO53x/RNp/TDVe8Ct0Gw7s3ws7lgY5IqYDzyxiF3SU0AlgEdDXG7AYrmQBd7GrJwI4azXLssmR7/ujyI9oYY7xAIZBQy/5vF5GlIrI0NzfXH4ekqm2cY93BPOiyBjdZtncZWw9u5bL0+tuMH9SVuE5upi/Z3pwo/Ss0Eq55FyISYdqVcCA70BEpFVDNThQiEgm8B9xrjKlrVLK2MwFTR3ldbY4sMGaqMSbDGJPRuXPn+kJWjbHmfYhJhZSMBjd5b/N7RLojOS/tvHrrhrqcXD4yhU/X7iW3qKI5kfpXVFf48XvW4z3euBQOBvgyXqUCqFmJQkTcWEniTWPM+3bxXrs7CXu6zy7PAVJrNE8BdtnlKbWUH9FGRFxADHCgOTGrRijLt64CGnRJg692KqwoZN62eUzsPZFwV3iD2kw+sQden+G95W3sAX2d+8GPZ0BJLrw+CUraSPeYUq2sOVc9CfAysN4Y81SNVR8BN9jzNwAza5RPtq9k6oU1aL3Y7p4qEpGT7G1ef1Sb6m39CPjcBPyi+w5kwyzweWDQpQ1uMitrFhVVFXUOYh+tb5dITkyL563F26lqK4Pa1VJHwzVvQ8E2eOMSvSFPdUjNOaMYC1wHnC0iK+zPBcDjwDgR2QyMs5cxxqwF3gHWAXOBO40x1a85+ynwb6wB7i3AHLv8ZSBBRDKBX2JfQaVayaq3IS4Nuo9sUHVjDG9vfJuBCQMZkDCgUbu6cWwa2/JKA//8p9qknQqT37Re2vTfy6G8Ddz3oVQrqv0C9wYwxnxD7WMIAOccp81jwDGvFjPGLAUG11JeDlzR1BhVMxTsgOyv4cyHGtzt9P3u78kqzOKxUxv/9rjxg7qRGh/Ov77OYsLgbo1u3+L6nmtdDfXO9VY31LXvQaf4QEelVKvQO7NV7VZNBwwMm9zgJtM2TCM+LJ4JaRPqr3wUp0O49dTeLNuWz7JtbXQYqv9EuPIN2LsGXr0QivfV30apIKCJQh3LGFjxFvQ8FeJ6NqjJjqIdLNixgMvTLyfEGdKk3V6RkUJMuJuXFmQ1qX2r6H8BXPMO5GfDf87XN+SpDkEThTpWzhI4sAWGX93gJm9veBuHOLjyhCubvNtOIS5uHtuLeev2siqnoMnbaXF9zoLrPrDOKF6ZAPvWBzoipVqUJgp1rBXTwN0JBk5qUPVSTynvZ77PuT3PpVtE88YXbj41jbhObp78dFOzttPiepwEN35i3Yz48njIWhDoiJRqMZoo1JHKD8Lqd2HgJRAaVX99YMamGRRVFnHtgGubvfuoMDd3nNGHrzblsji7jY5VVEsaBrfOh+ju8N/LrASrVBDSRKGOtOptqCyGE29tUPXKqkpeW/caGV0zGN5luF9CuP7kNLpEhfLY7PVt52GBxxObCrf8z7qE9sOfwrw/gK+q/nZKtSOaKNRhxsCSl6H7CEge1aAmn2R9wr7Sfdw6pGGJpSHCQ5w8eH5/Vu4o4J2lO+pvEGhhMdYd3KNugm+fse61KMkLdFRK+Y0mCnXYtu8gdz2MbtiXfpWviv+s+Q8D4gcc9+VETXXpiGRGp8XxxNwNFJRW+nXbLcLphouegYuft/47Tj0Ddv0Q6KiU8gtNFOqwRS9BWGyDnxQ7b9s8th7cys1DbsbfT38XER65eDCFZR4em9WOrioaeR3cPNeaf3k8fP+idaamVDumiUJZ9mfC+o9h9C0Q0qne6l6fl+dXPE/f2L6M6zGuRUIa2D2aO87ow7vLctrmoz2OJ3kk3L7Auox27oPw5o+gqB3Fr9RRNFEoy3fPgSsUxtzRoOozM2ey7eA27hpxF06Hs/4GTXTvuf0YkBTNQ++vIq+4DT2GvD4RCXD1dJj4d9j6Dbx4svWQRaXaIU0UynrXwsq3YPiPIbJLvdUrqip4ceWLDE0cylmpZ7VoaCEuB89cNZyD5V7umb4Cb5WvRffnVyLWeM9PvrIuoZ1+Dbx7kz76Q7U7migUfPcP8HnhlLsaVP3N9W+yt3Qvd4+82+9jE7U5oVsUf7pkMN9k7udv/9vY4vvzu84nwK2fw1m/hQ2fwPOj4Yf/6tiFajc0UXR0BTtgyb9g+DUQ36ve6vtK9zFl5RROTzmdMUljWiFAy5UZqVx3Uk+mfJXFjGXt8PlKrhA44//gjm+hy0CYeSe8OhF2rwx0ZErVSxNFR/flXwCxHifeAE8tewqPz8ODo1v/1SC/u3AgY/sm8MB7q9rX4HZNnfvBjbPgomchdwNMOQNm/lwHu1WbpomiI9u71hqbOPE2iEmpt/rSPUuZlTWLmwbfRGp0ar31/S3E5WDKdRkMTo7hzmnL+WpTbqvH4BcOB4y6Ee5aDiffCSunwz9GwYK/6kuRVJukiaKjMgZm3WfdN3HaffVWL/WU8ofv/kByZLJf78JurMhQF6/eOJo+nSO55bUlzF69O2CxNFt4LIx/DO5cBL1Ohy8eg2eHwtdPQUVxoKNT6hBNFB3VimmwfSGM+2OD3tT27PJn2V60nUfHPkq4K7wVAjy+uIgQpt9+EsNSYrlz2nJe+Sabdv0q9YQ+cPU0uO0LSBkN8x+BZ4dZCaOsDT9uXXUYmig6oqK9MO93kDrGuiS2Hgt3LWTahmn8eMCPGd1tdCsEWL+YcDdv3DKGcQO68sdP1nHv2ysoq2znD+NLHgk/fhdumQdJQ62E8fQgmPMg5G8LdHSqA5N2/UusFhkZGWbp0qWBDqPtMsa6U3jrN9bdw13611l9T8kervrkKmJDY3lr4lt0ctd/13Zr8vkMLy7YwpOfbqR3YgRPXjGMET3iAh2Wf+xeBQtfgDUzwPhgwEXWgwd7nWGNcyjlRyKyzBiTUds6/WvraBZPhczP4Lw/1ZskPFUefrXgV5R7y3n6rKfbXJIAcDiEO8/qy39vGUNZZRWXv/gdj81aR3GFN9ChNV/SULhsCtyzyrrHJfsreOMS+McI+PrveqWUajV6RtGRZH9tfdH0Pdd6vEQdN8v5jI9ff/NrZmXN4skznmR82vhWDLRpiso9/Hn2et5avIPEyFB+Oa4fV2ak4HIGye8hT7l1w96yV2Hr1yBO6HM2DL4c+k+EsOhAR6jasbrOKDRRdBQHsuBfZ0NEZ7j1M+sdCnX4+9K/8+raV7l7xN3cNvS2VgrSP1bsKOBPn6xj6bZ80hI6cdvpvbl8ZAph7pZ7JlWr258JP7wBa96Hwu3gCoP082DwZdDnHE0aqtE0UXR0BdvhPxOhssh6dWdCn+NWNcbw7PJneXnNy1x1wlX8ZsxvWuUxHf5mjOHTdXt54YtMVuUUkhgZypUZKVyZkUpaYkSgw/MfYyBnCayeAWs/gJJ94HBbb9zrNwH6jW/QHfdKaaLoyPK3wesXQ1k+XP8RdD/+60qrfFU8seQJ3trwFj/q9yN+d9LvcEj77rYxxrBwSx7//iabLzfuw2fgxLR4LhyWxLkDutI9NrCX+vqVrwq2fw+b5lqf/Zus8sQToPcZkHaalUAacDm06ng0UXRUO5bA9KuhqhKu/QBSjv9608KKQh746gG+3fUt1w+8nl9l/KpdnknUZe/Bct5bnsN7y3LYklsCwJDkGM46oTMn9U5gZM+44OqeytsCm/4HmfOsBOIpBQS6DbaunEodAykZ1pNtVYeniaKjMca6uunT31lfAj9+FxLTj1t92d5l/Pab37KndA+/HvNrruh3RSsGGxiZ+4qZt24vn67bw8odBfgMhDgdDEuNYXhqLIOTYxicHEOvhAgcjiBImN5K2LXcunIq+yvYsRiq7Pd7RHW3fkQkZ1iJo9uQesewVPDRRNGR5G2B2b+CLZ9D+ni45EXrJTq1KKwo5MWVLzJt/TS6R3bn8dMeZ3iX43dNBauD5R6Wbc3n+6w8FmUfYN3ug1R6rfdeRIQ46dctit6JkfTuHEGvROvTM6ETnUJcAY68GbwVsGc15CyFnUutaX724fUxqdB1kPXpMhC6DrbGtpzuwMWsWpQmio6geB98+ywsmmK9qW7cI5BxS62XwJZ6Snl749v8a/W/KK4sZnL/ydw78t42eZ9EIHiqfGTuK2bNzkLW7Cxk095isveXsOdg+RH1YsLdJMWEWZ/YcJKiw+gaHUZ8RAhxESEk2NPoMFf76MYrybPOOvausR4YuXetNc7hs+9JESfE9YT4PpDQ10oc8b2taUwqtOCbDlXL00QRrIyBXT9Y19WvnG6NRYz4MZz9e4jqekz1rMIsZmyawQebP6DYU8zY5LH8YuQvOCH+hNaPvR0qqfCSvb+E7P0lbD9Qyp7CcnYXlrG7sJzdheUcKKmstZ3LIcRFhBDfKYSoMBeRYS4iQ13WfKiLiNCay27CQxyEuZyEup2EuR2EuZ2EuZ2Euux5l6P17g3xVsD+zVby2L/JOmM9sAXyssBTcrieOK1uzpgU6xOdbM+nQkwyRCVBeLzeUd6GaaIIJt4Kq3856wtY+6H1j9YVZr146OSfH3Hpq6fKw9q8tSzctZBPt31KZkEmLnExLm0c1/S/pkN2M7Wkck8VuUUV5JdWkldSSX5JJQfsT36pNS2u8FJc7qXInhZXeCltwjOqXA6xE4gDt9OByym4HdbU5XDgdgoupzV1Ox24HIeXXXa9ELudUwSHQ3CI4HQIIlhlh8o5qo71SIcIz35iS7cTW7adyPJdRJbvIaJ8D53K99CpbA8Oc+Td8T5xUhkST2VYAhWhCVSGJVIZZk9DE/CGxuINicETEk1VSAze0BhwhoKAACJiT0EQe4q9/vCySI11VJ9UH7W+xnZqasiJX211GrKdo8sa1KYB+65ZK8ztICWuaT0D7T5RiMgE4FnACfzbGPP48eoGVaIoK4CCbVYXwJ7Vh/uUvWUgDug5FoZeCQMuptQVQvbBbDLzM8ksyGT9gfWsyl1FmbcMQRjRZQTnpZ3HeT3Po3OnzoE+MlVDlc9QXOGlpMJLUbmXck+V9fH6Ds1XeHyUe+1yj48KrzUt91ThqfLhrTJ4fAZvlQ9PlbHKfNa8t8qH12eo9FrT6jrV66t8Bp8x+HwGn4EqYzCmurxpxyT4SKSQ7pJHd8mjixSQKIUkUkiiFNJZCg8th4nnuNspMyEUEEmhiaCQCA6aCIoIp9SEUkIYpSaMYsIoJYwSY08Jo9SEUkw4pSaMctxUEEI5IfiC/KlFw1Nj+fDOsU1qW1eiaPOjcSLiBF4AxgE5wBIR+cgYsy6wkdXBGOshblWV4CkDb7n18ZRbX/LeCqu8oghTegBvWR6VpfupKDtAZekBPEW7qCzaTaWnmEoRih0OilxhHIxNpqj/WA5Gd6cwPIp9FQXsyfmQPRtfoqiy6NDuQxwh9Intw6V9L2V0t9GM7DqS+DC9dr6tcjqEmHA3MeFtb6DYGIOxk4eVTMBnjJVMfDXLjT2PnXBqJhprGwb7n4a9fNBAofEhlcW4yvbjKC/EUVmAq6IAR3khzspCHBWFuCsK6VpRSHJlIc6KAhzeXJyeEhyeEpxV5fUcwZF84sLnDMXnDMHnCMXnDKXKGXbEss8ZSpUjFJ/TjREXRlz4HC6Mw4VPXBiH21qW6rLDy9X1Ds1X18EB4sCIE4PgE6e9fLjMmlpliAOfHJ434sBnrLqIXVcc+KguA3AQ3SnU/38EtINEAZwIZBpjsgBEZDowCfBroigs2Mp1719E9Q8og7H+sAHfoRIOfwwYObJedVtfzXr2aaGRw+Uc2q5QWdull+FAeBQQddSKIigpwlWWRXRINF07dSU5MplRXUbRNaIrPaJ6kB6XTmpUKi5He/ifVrV11d04jmM6QfwpFqj/DYu18lVBZYn9KbY/NZYrig//UPNW4PCW46j+oeatOGLd4fmDUFFu/dDzeaHKY00PzXsOD/C3NckZMGC+3zfbHr5NkoEdNZZzgDE1K4jI7cDtAD169GjSTpzOUNJD462+TuTQVSoOkRplDrDXOeypICCH1ztEDtURcSIO1+GP04U43DWWQwkJjSI0LIYQdwRuh5sQZ4j1cYTgdroJcYQQFRJ1xCfMGdY+rqJRqqU5nNZzrVr72VbGHJk4qrz29Khln9eaNz4wVdbUV3V42edr5Lrq+RrrfFWAsWKKPPYiFn9oD4mitm/EI3pOjTFTgalgjVE0ZSeRUUn8/ccLmtJUKdXRiFj3lHSQ+0raw8hODpBaYzkF2BWgWJRSqsNpD4liCZAuIr1EJASYDHwU4JiUUqrDaPNdT8YYr4j8HPgf1uWxrxhj1gY4LKWU6jDafKIAMMbMBmYHOg6llOqI2kPXk1JKqQDSRKGUUqpOmiiUUkrVSROFUkqpOrWLhwI2hojkAtuasYlEYL+fwmkvOtoxd7TjBT3mjqI5x9zTGFPrE0ODLlE0l4gsPd4TFINVRzvmjna8oMfcUbTUMWvXk1JKqTppolBKKVUnTRTHmhroAAKgox1zRzte0GPuKFrkmHWMQimlVJ30jEIppVSdNFEopZSqkyYKm4hMEJGNIpIpIg8GOp6WICKpIvKFiKwXkbUico9dHi8i80Rksz2NC3Ss/iQiThH5QUQ+sZeD+ngBRCRWRGaIyAb7f++Tg/m4ReQX9t/0GhF5S0TCgvF4ReQVEdknImtqlB33OEXkIfs7baOIjG/qfjVRYH2RAC8A5wMDgatFZGBgo2oRXuA+Y8wA4CTgTvs4HwTmG2PSgfn2cjC5B1hfYznYjxfgWWCuMaY/MAzr+IPyuEUkGbgbyDDGDMZ6HcFkgvN4XwUmHFVW63Ha/7YnA4PsNv+0v+saTROF5UQg0xiTZYypBKYDkwIck98ZY3YbY5bb80VYXx7JWMf6ml3tNeCSwETofyKSAkwE/l2jOGiPF0BEooHTgZcBjDGVxpgCgvu4XUC4iLiATlhvwQy64zXGfAUcOKr4eMc5CZhujKkwxmQDmVjfdY2micKSDOyosZxjlwUtEUkDRgCLgK7GmN1gJROgS+Ai87tngPsBX42yYD5egN5ALvAfu8vt3yISQZAetzFmJ/AksB3YDRQaYz4lSI+3Fsc7Tr99r2misEgtZUF73bCIRALvAfcaYw4GOp6WIiIXAvuMMcsCHUsrcwEjgReNMSOAEoKj26VWdp/8JKAX0B2IEJFrAxtVm+C37zVNFJYcILXGcgrWqWvQERE3VpJ40xjzvl28V0SS7PVJwL5AxednY4GLRWQrVnfi2SLyX4L3eKvlADnGmEX28gysxBGsx30ukG2MyTXGeID3gVMI3uM92vGO02/fa5ooLEuAdBHpJSIhWANAHwU4Jr8TEcHqt15vjHmqxqqPgBvs+RuAma0dW0swxjxkjEkxxqRh/W/6uTHmWoL0eKsZY/YAO0TkBLvoHGAdwXvc24GTRKST/Td+Dtb4W7Ae79GOd5wfAZNFJFREegHpwOKm7EDvzLaJyAVY/dlO4BVjzGMBDsnvRORU4GtgNYf77H+NNU7xDtAD6x/dFcaYowfM2jURORP4lTHmQhFJIPiPdzjWAH4IkAXchPXDMCiPW0QeAa7CurLvB+BWIJIgO14ReQs4E+tx4nuBPwAfcpzjFJHfADdj/Xe51xgzp0n71UShlFKqLtr1pJRSqk6aKJRSStVJE4VSSqk6aaJQSilVJ00USiml6qSJQimlVJ00C6C+IgAAAAxJREFUUSillKrT/wOq/3SbuzdEawAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "pyplot.plot(t, S, label=\"S\")\n", "pyplot.plot(t, I, label=\"I\")\n", "pyplot.plot(t, R, label=\"R\")\n", "pyplot.plot([t[0], t[-1]], 2*[(1-gamma/beta*(1-math.log(gamma/beta)))*N], label=\"$I_{max}$\")\n", "pyplot.legend(loc=\"upper right\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "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.6" } }, "nbformat": 4, "nbformat_minor": 4 }