[PYTHON] (Erster Beitrag) Eine Geschichte über die numerische Berechnung von Influenza und neuem Lungenentzündungs-Coronavirus mit Tensorflow
Einführung
Abstrakt
- Influenza
- Corona-Virus
Ich habe versucht, mit einem mathematischen Modell zu analysieren.
Hintergrund
Dieses Mal habe ich versucht, das Modell für Infektionskrankheiten numerisch zu berechnen, aber bitte haben Sie Verständnis dafür, dass der Autor kein Experte für mathematische Modelle wie das Modell für Infektionskrankheiten ist. Seit meinem Hauptfach Physik an der Universität und der Graduiertenschule habe ich mich nur teilweise mit Informatik befasst. Der Hintergrund dieses Beitrags war, dass ich neugierig war, ob ich Coronavirus quantitativ analysieren könnte, was derzeit ein heißes Thema ist. In den vergangenen Fällen von Infektionskrankheiten gab es SARS und MERS, aber ich dachte, es lohnt sich zu untersuchen, welche Art von mathematischem Modell zur Vorhersage der Konvergenzzeit von Infektionskrankheiten verfügbar ist.
Modell für Infektionskrankheiten
SIR-Modell
Es gibt ein SIR-Modell als Modellgleichung, das den epidemischen Prozess von Infektionskrankheiten beschreibt. Es gibt drei Variablen: Anfällig (vermutet infiziert), Infiziert (infiziert) und Wiederhergestellt (Immunitätsträger). Dieses Modell wurde 1927 in einem Artikel veröffentlicht („Ein Beitrag zur mathematischen Theorie der Epidemien“). Insbesondere wird es häufig in Bereichen eingesetzt, die sich mit bekannten Phänomenen wie der Sozialphysik befassen.
S+I \rightarrow 2I \\
I \rightarrow R \\
Der elementare Prozess sind die beiden oben genannten Prozesse. Die Oberseite ist der Infektionsprozess und die Unterseite ist der Wiederherstellungsprozess. Betrachten Sie es als eine chemische Reaktionsformel. Wenn es tatsächlich in eine Differentialgleichung aufgenommen wird, sollte dieselbe Berechnungslogik wie die Reaktionsgeschwindigkeitstheorie in der Chemie verwendet werden.
Formel
Das SIR-Modell ist wie folgt formuliert.
\frac{dS(t)}{dt} = -\beta S(t)I(t) \\
\frac{dI(t)}{dt} = \beta S(t)I(t)-\gamma I(t)\\
\frac{dR(t)}{dt} = \gamma I(t) \\
Für $ t $ zu einem bestimmten Zeitpunkt ist $ S (t) $ die Anzahl der vermuteten Infektionen, $ I (t) $ die Anzahl der infizierten Personen und $ R (t) $ die Anzahl der Immunträger.
- Obere Formel ➤ Die Anzahl der Personen, bei denen der Verdacht auf eine Infektion besteht, nimmt proportional zur Anzahl der Kontakte mit der infizierten Person ab.
- Mittlere Formel ➤ Die Anzahl der infizierten Personen ist proportional zur Anzahl der Kontakte, nimmt jedoch mit zunehmender Anzahl der heilenden Personen ab.
- Niedrigere Formel ➤ Wenn die Anzahl der Personen, die sich erholen, zunimmt, nimmt die Anzahl der infizierten Personen ab.
Hier steht $ \ beta $ für die Infektionsrate und $ \ gamma $ für die Wiederherstellungsrate und die Isolationsrate. Ebenfalls
-\beta I(t) \\
Entspricht der Infektiosität.
Insbesondere, warum es wie die obige Differentialgleichung beschrieben werden kann http://www2.meijo-u.ac.jp/~tnagata/education/react/2019/react_03_slides.pdf Ich denke, Sie sollten sich darauf beziehen. Es ist die gleiche Logik wie bei der Formulierung der Reaktionsgeschwindigkeitsformel für die Chemie der High School.
Hier löse ich die obigen simultanen Differentialgleichungen, und als ich das sah, erinnerte ich mich an die Lorenz-Gleichungen, die ich im College gelernt hatte. Die Lorenz-Gleichung ist eine gewöhnliche Differentialgleichung (allgemein bekannt als ODE), über die im Bereich des Chaos in der Physik gesprochen wird. Ich denke, ich bin ein ziemlicher Wahnsinniger, also werde ich nicht auf Details eingehen, sondern nur die Formel schreiben.
\frac{dx}{dt} = -px+py \\
\frac{dy}{dt} = -xz+rx-y\\
\frac{dz}{dt} = xy-bz \\
Ich habe diese numerische Simulation für eine Universitätsaufgabe durchgeführt, daher dachte ich, ich sollte dieselbe Simulation wie damals durchführen. Ich erinnere mich, dass ich es damals mit Rungekutta gelöst habe, aber dieses Mal werde ich versuchen, TensorFlow zum Lernen zu verwenden.
Was ist TensorFlow?
TensorFlow ist eine Open Source-Softwarebibliothek für maschinelles Lernen, die von Google entwickelt wurde. Kurz gesagt, es ist eine Python-Bibliothek für die numerische Hochgeschwindigkeitsanalyse. Dieses Mal werde ich kein maschinelles Lernen durchführen, aber ich werde es für die numerische Berechnung verwenden. Ich glaube nicht, dass Sie den Tensor hören können, aber ich verstehe, dass es eine Erweiterung der Matrix ist. Da es außerhalb des Rahmens dieses Beitrags liegt, werde ich es irgendwo veröffentlichen. Das Merkmal von TensorFlow ist, dass der Quellcode einfach geschrieben werden kann. Andererseits ähnelt der größte Teil der Berechnungsverarbeitung einer Black Box, sodass die Berechnungsverarbeitung nicht klar ist.
Aktuelle numerische Simulation (1)
Quellcode
Parameter etc. werden im Ergebnisteil erklärt.
sir_model.py
import numpy as np
import tensorflow.compat.v1 as tf
tf.disable_v2_behavior()
#Abhängig von der Version von TensorFlow funktioniert dies möglicherweise nicht(Wahrscheinlich aus Kompatibilitätsgründen)
#import tensorflow as tf
import matplotlib.pyplot as plt
import time
from scipy.integrate import solve_ivp
%matplotlib inline
def main():
x=tf.placeholder(dtype=tf.float64,shape=(3),name='x')
class SIR(object):
def __init__(self, beta=0.00218, gamma=0.452, **kwargs):
self.beta = float(beta)
self.gamma = float(gamma)
def __call__(self, t, x):
dx_dt = -self.beta * x[0]*x[1]
dy_dt = self.beta*x[0] *x[1] - self.gamma*x[1]
dz_dt = self.gamma*x[1]
dX_dt = tf.stack([dx_dt, dy_dt, dz_dt])
return dX_dt
f_sir = SIR()
fx =f_sir(None, x) # ODE
x0 = np.array([762, 1, 0], dtype=np.float64) # initial value
sess = tf.Session()
sess.run(tf.initializers.global_variables())
print('###value of f(x)###')
print(sess.run(fx, feed_dict={x:x0}))
dt=0.1
tstart=0.0
tend=30.0
ts=np.arange(tstart, tend+dt, dt) # 0.Ausgabe 30 Tage in einem Schritt
def f_sir_tf(t,xt):
return sess.run(fx, feed_dict={x: xt})
start_time =time.time() # measurement of time
sol_sir = solve_ivp(fun=f_sir_tf,
t_span=[tstart, tend], y0=x0, t_eval=ts)
integration_time_tf = time.time() - start_time
t_lo = sol_sir['t'] #Holen Sie sich die Zeit für jeden Schritt
x_lo = sol_sir['y'] #X jedes Schritts(t)Holen Sie sich den Wert von
#Verarbeitungszeit berechnen
print("processing time (tf) %.5f" % integration_time_tf)
fig=plt.figure(figsize=(14,8))
ax=fig.add_subplot(111)
ax.plot(t_lo,x_lo[0],'*',color='blue',alpha=0.2,label="Susceptible",markersize=10)
ax.plot(t_lo,x_lo[1],'^',color='red',alpha=0.2,label="Infected",markersize=10)
ax.plot(t_lo,x_lo[2],'.',color='green',alpha=0.2,label="Recovered",markersize=10)
ax.grid()
ax.legend()
ax.set_xlabel('day')
ax.set_ylabel('number')
if __name__ == "__main__":
main()
Simulationsergebnis
Für diese Simulation habe ich auf der folgenden Folie auf Seite 14 verwiesen.
http://www.actuaries.jp/lib/meeting/reikai20-7-siryo.pdf
Dies ist ein Fall von Influenza in einem 1978 untersuchten Internat. Von den insgesamt 763 Personen sind $ S = 762 $ Personen und $ I = 1 $ Personen die Ausgangsbedingungen. Auch für die Infektionsrate $ \ beta $ und die Wiederherstellungsrate $ \ gamma $
Übrigens ist der Graph der Beziehung zwischen $ S $ und $ I $ zu diesem Zeitpunkt wie folgt. Es konvergiert zum Ursprung, wobei das Bild von rechts unten nach links verläuft.
Interpretation
Die Anzahl der wiedergewonnenen grünen $ R $ ist hier nicht sehr wichtig. Am wichtigsten ist die Anzahl der mit rotem $ I $ infizierten Personen. In diesem Fall wird ab dem Anfangszustand von $ I = 1 $ $ S $ in Kontakt mit der infizierten Person zu $ I $, steigt schnell an, erreicht eine Woche später einen Spitzenwert und konvergiert nach etwa 14 bis 15 Tagen gegen 0. Sie können sehen, wie es läuft. Sie können sehen, dass es tatsächlich gut passt wie das weiße Kreisdiagramm auf der Referenzfolie S. 14.
Aktuelle numerische Simulation (2)
Nachdem wir die numerischen Werte für Influenza berechnet haben, ist der nächste "Coronavirus". Ich habe nichts über Influenza erklärt, aber ich werde es anhand der "Grundreproduktionsnummer $ R_0 $" erklären, die ein wichtiger Parameter beim Umgang mit Infektionskrankheiten ist.
Grundreproduktionsnummer
Die Definition lautet "die erwartete Anzahl von Sekundärinfektionen, die von einer infizierten Person während des gesamten Infektionszeitraums verursacht wurden". Wenn beispielsweise $ R_0 = 2 $ ist, infiziert eine Person durchschnittlich zwei Personen. Ich werde nicht auf Details eingehen, aber es ist allgemein bekannt, dass der Ausbruch endet, wenn
Ich werde die folgende URL für dieses $ R_0 $ zitieren. https://wired.jp/2020/01/26/scientists-predict-wuhans-outbreak-will-get-much-worse/
Ein wesentlicher unsicherer Faktor ist die Infektiosität von 2019-nCoV. Wie ansteckend ist es tatsächlich? In Reeds Modell wird die Anzahl der Personen, die eine einzelne infizierte Person infizieren kann (Grundreproduktion des Virus), auf 3,6-4,0 geschätzt. Im Vergleich dazu beträgt SARS (Severe Acute Respiratory Syndrome) 2-5 und die infektiösesten Masern beim Menschen 12-18.
Berechnen wir diesen Koronavirus mit der Grundreproduktionsnummer $ R_0 = 3.4 $. Obwohl die Berechnung weggelassen wird, kann die Infektionsrate $ \ beta $ nach der folgenden Formel berechnet werden.
Betrachten Sie nun die Anfangsbedingungen. Wuhan City, wo das Corona-Virus weit verbreitet ist, ist eine große Stadt in China mit 11,08 Millionen Einwohnern.
Dann erhalten Sie $ \ beta = 2.19 \ times10 ^ {-8} $. Hier setze ich $ \ gamma = 0.071428 $ (1/14 entspricht 2 Wochen). $ S (0) $ ist 11.079.999.
Ergebnis
Der Quellcode ist grundsätzlich der gleiche wie für Influenza. Nur die Parameter wurden geändert.
Interpretation
Die Anfangsbedingungen setzen den Zeitpunkt voraus, zu dem die erste infizierte Person gefunden wird. Realistisch gesehen ist es der 1. Dezember 2019. Es ist erst Anfang Februar. Wenn Sie also glauben, dass ungefähr zwei Monate vergangen sind, sind es 60 Tage. Es ist wirklich beängstigend, wenn ich darüber nachdenke. .. .. Dieses Mal habe ich eine ziemlich große Anzahl von R = 3,4 verwendet, also nehme ich den schlimmsten Fall an. Tatsächlich ist es laut WHO (Weltgesundheitsorganisation) $ R_0 $ = 2,2 ~ 2,8, also sollte es früher konvergieren. Der entscheidende Punkt ist, wann die erste infizierte Person zu finden ist. Unter der Annahme, dass dies der Höhepunkt infizierter Personen ist (ungefähr 100 Tage in der Grafik), wird erwartet, dass es ungefähr 40 Tage dauern wird, bis die Anzahl infizierter Personen von hier aus auf 0 konvergiert. Wie viele von Ihnen vielleicht bemerkt haben, berücksichtigt dieses Modell die Toten nicht. Daher bleibt das S + I + R-System erhalten. Wenn Sie sich ein realistischeres Modell vorstellen möchten, können Sie das SEIR-Modell verwenden.
Zusammenfassung
Dieser Beitrag
Ich habe zum ersten Mal bei Qiita gepostet. Dieses Mal haben wir das SIR-Modell numerisch berechnet, das das grundlegendste der "Infektionsmodelle" ist. Die tatsächlich verwendeten Parameter wurden aus früheren Fällen gelernt, aber eine fortgeschrittenere Mathematik ist erforderlich, um die Parameterschätzung im Detail durchzuführen. Wenn ich Zeit habe, werde ich in Zukunft einen solchen Artikel schreiben. Als Einzelperson interessiere ich mich für maschinelles Lernen und Wettbewerbsprogrammierung, daher denke ich darüber nach, ein Repository zu erstellen, das Memoranden davon enthält.
Ausblick
- Numerische Berechnung des SEIR-Modells (Modell einschließlich Latenzzeit in SIR)
- SIR im Netzwerkmodell
- Spielen Sie mit TensorFlow
- Wasche ernsthaft deine Hände
Referenzdokument
- https://www.ism.ac.jp/editsec/toukei/pdf/54-2-461.pdf
- http://www.orsj.or.jp/archive2/or56-12/or56-12_728.pdf
Recommended Posts