[CovsirPhy] COVID-19 Python-Paket für die Datenanalyse: S-R-Trendanalyse

Introduction

Wir erstellen ein Python-Paket CovsirPhy, mit dem Sie COVID-19-Daten (z. B. die Anzahl der PCR-Positiven) einfach herunterladen und analysieren können.

Einführungsartikel:

• SIR model
• SIR-F model
• Data loading

** Dieses Mal möchten wir die S-R-Trendanalyse einführen (eine Analysemethode für Trends bei der Ausbreitung von Infektionen). ** ** **

Die englische Version des Dokuments lautet CovsirPhy: COVID-19-Analyse mit phasenabhängigen SIRs, [Kaggle: COVID-19-Daten mit SIR-Modell]( Weitere Informationen finden Sie unter https://www.kaggle.com/lisphilar/covid-19-data-with-sir-model.

1. Ausführungsumgebung

CovsirPhy kann mit der folgenden Methode installiert werden! Bitte verwenden Sie Python 3.7 oder höher oder Google Colaboratory.

--Stabile Version: pip install covsirphy --upgrade

• Entwicklungsversion: pip install" git + https://github.com/lisphilar/covid19-sir.git#egg=covsirphy "

import covsirphy as cs
cs.__version__
# '2.8.2'

Die Tabellen und Grafiken in diesem Artikel wurden mit Daten vom 11.09.2020 erstellt. Klicken Sie hier für den Code [^ 2], um die tatsächlichen Daten von COVID-19 Data Hub [^ 1] herunterzuladen:

data_loader = cs.DataLoader("input")
jhu_data = data_loader.jhu()
population_data = data_loader.population()

Die Anzahl der Infizierten in Japan unterscheidet sich vom Wert [^ 3] des Ministeriums für Gesundheit, Arbeit und Soziales HP. Es scheint, dass. Wenn Sie den Wert des Ministeriums für Gesundheit, Arbeit und Soziales verwenden möchten, laden Sie die Daten aus COVID-19-Datensatz in Japan mit dem folgenden Code herunter. Bitte ersetzen. Wenn zu diesem Zeitpunkt ImportError / ModuleNotFoundError: Anfragen und aiohttp auftreten, installieren Sie bitte by pip install request aiohttp (während der Ursachenuntersuchung sollten Anfragen in das abhängige Paket aufgenommen werden ...).

[^ 3]: Wenn die vom Ministerium für Gesundheit, Arbeit und Soziales festgelegte "Anzahl der Positiven", "Anzahl der entlassenen oder annullierten Personen" und "Anzahl der Todesfälle" als bestätigte Fälle behandelt, wiederhergestellt, wiederhergestellt und tödlich sind. ..

japan_data = data_loader.japan()
jhu_data.replace(japan_data)
# Citation -> str
print(japan_data.citation)

2. Überprüfen Sie die Daten

Bevor wir uns mit der Erläuterung der S-R-Trendanalyse befassen, betrachten wir die tatsächlichen Daten (japanische Daten als Beispiel) anhand von "Scenario.records ()".

snl = cs.Scenario(jhu_data, population_data, country="Japan")
#Wenn das Diagramm nicht angezeigt wird
snl.records(show_figure=False)

#Grafik anzeigen
# filename:Standardwert=None (Nicht als Bilddatei speichern)
_ = snl.records(filename=None)

3. Zweck der S-R-Trendanalyse

Beim Umgang mit Daten zu Infektionskrankheiten unter Verwendung simultaner normaler Differentialgleichungen wie dem SIR-Modell ist es üblich, eine Analyse unter einer der folgenden Annahmen durchzuführen.

• Die Modellparameter sind vom Beginn bis zum Ende der Epidemie konstant
• Die Modellparameter schwanken täglich und folgen der Poisson-Verteilung usw.

Bei der Analyse von COVID-19 dachte ich jedoch, dass es notwendig sei, beide anzupassen, anstatt sie zu übernehmen.

Erstens hat bei diesem Ausbruch jedes Land und jeder Einzelne einige Maßnahmen ergriffen, um den Ausbruch ab etwa Februar und März 2020 zu kontrollieren, als sich die Infektion auf die Welt ausbreitete. Parameter wie die Infektionsrate schwanken von Zeit zu Zeit und können bis zum Ende nicht als konstant angenommen werden (klar aus den oben genannten tatsächlichen Daten und der Wellenform des SIR-F-Modells [^ 4]).

[^ 4]: [CovsirPhy] COVID-19 Python-Paket für die Datenanalyse: SIR-F-Modell

Wenn andererseits die Parameter jeden Tag schwanken, hängt die Analyse von jedem Wert der erfassten Daten ab. Insbesondere im Fall von COVID-19 ist es notwendig, die Umrisse der Epidemie anstelle der täglichen Daten zu erfassen, da das Inspektionssystem im Anfangsstadium nicht eingerichtet wurde und die Möglichkeit einer Infektiosität besteht, auch wenn kein Symptom vorliegt. Ich dachte.

Daher entschied sich CovsirPhy, beide Vorteile zu nutzen, indem die Parameter in Perioden unterteilt wurden, die konstant zu sein scheinen. Der "Zeitraum, in dem die Parameter konstant sind" wird als "Phase" bezeichnet. Vom Beginn bis zum Ende der Epidemie gibt es mehrere Phasen mit unterschiedlichen Parametern. Ich werde den Code in einem anderen Artikel schreiben, aber ich denke zum Beispiel, dass sich die Ausführungsreproduktionsnummer $ R_ {\ mathrm {t}} $ des SIR-F-Modells schrittweise ändert, wie in der Grafik gezeigt.

Wie stellen Sie die Phase ein? Da die Berechnung der Parameter einige Zeit in Anspruch nimmt [^ 5], möchte ich den Verzweigungspunkt der Phase (das Datum, an dem sich der Wert des Parameters ändert) ermitteln, ohne die Parameter zu berechnen.

[^ 5]: Die parallele Berechnung von 10 Phasen und 8 CPUs dauert ca. 2 Minuten. Sie können auf die Analyse eines Landes warten, dies ist jedoch zeitaufwändig und schwierig, da die Daten täglich zunehmen und die Daten vieler Länder möglicherweise analysiert werden.

4. Wie es funktioniert

Mithilfe der S-R-Trendanalyse ist es möglich, den Wendepunkt der Phase aus dem Übergang der Anzahl der Anfälligkeitsinhaber $ S $ und der Anzahl der Wiederhersteller $ R $ zu ermitteln. Ich werde den Prozess auslassen, bis ich ihn gefunden habe, weil er lang sein wird, aber [Kaggle: COVID-19-Daten mit SIR-Modell, SR-Trendanalyse](https://www.kaggle.com/lisphilar/covid-19-data-with- Sir-Modell # SR-Trend-Analyse), beziehen Sie sich also bitte darauf.

In den simultanen normalen Differentialgleichungen des SIR-Modells [^ 6] und des SIR-F-Modells [^ 4] ist dies für $ S $ und $ R $ wie folgt.

[^ 6]: [CovsirPhy] COVID-19 Python-Paket für die Datenanalyse: SIR-Modell

\begin{align*}
& \frac{\mathrm{d}S}{\mathrm{d}T}= - N^{-1}\beta S I  \\
& \frac{\mathrm{d}R}{\mathrm{d}T}= \gamma I  \\
\end{align*}

Teilen wir $ \ frac {\ mathrm {d} S} {\ mathrm {d} T} $ durch $ \ frac {\ mathrm {d} R} {\ mathrm {d} T} $.

\begin{align*}
\cfrac{\mathrm{d}S}{\mathrm{d}R} &= - \cfrac{\beta}{N \gamma} S  \\
\end{align*}

Wenn Sie dies integrieren und als $ a = \ frac {\ beta} {N \ gamma} $ [^ 7] angeben,

[^ 7]: Pionier bei der Analyse von SR-Ebenen im SIR-Modell: Balkew, Teshome Mogessie, "Das SIR-Modell, wenn S (t) eine multiexponentielle Funktion ist." (2010) .Elektronische Thesen und Dissertationen.Papier 1747 .

S_{(R)} = N e^{-a R}

Teilen Sie beide Seiten logarithmisch

\log S_{(R)} = - a R + \log N

** Wenn die Parameter des Modells konstant sind, führt das Zeichnen eines Ein-Log-Diagramms mit der x-Achse als Anzahl der Wiederhersteller $ R $ und der y-Achse als Anzahl der Empfindlichkeitsinhaber $ S $ zu einer geraden Linie! ** ** **

(Grafik: Aus dem oben erwähnten Kaggle-Notizbuch der Unterschied in der Neigung der geraden Linie, wenn die Parameter geändert werden)

5. Codebeispiel

Es kann in einer Zeile mit der Methode Scenario.trend () ausgeführt werden. Das einseitige S-R-Protokolldiagramm wird angezeigt. Die Phase heißt Anfangsphase (0. Phase), 1. Phase, 2. Phase, ... [^ 8].

[^ 8]: Verwandte Technologie: [Python] Natürliche Zahlen in Bestellungen konvertieren

# show_figure:Grafik anzeigen
# filename:Standardwert=None (Nicht als Bilddatei speichern)
snl.trend(show_figure=True, filename=None)

Eine Liste der Start- und Enddaten für jede Phase kann im Datenrahmenformat mit der Methode "Scenario.summary ()" abgerufen werden.

snl.summary()

Da sich der Wert der Anzahl der Anfälligkeitsinhaber $ S $ in Abhängigkeit vom Wert der Gesamtbevölkerung ändert, wird er zusammen als wichtige Analysebedingung angezeigt. Obwohl in diesem Artikel nicht gezeigt, können die geschätzten Werte von Parametern mit der Methode "Scenario.summary ()" aufgelistet werden. Wenn es eine Phase gibt, für die die tatsächlichen Daten nicht registriert sind, lautet die Spalte "Typ" der entsprechenden Zeile "Zukunft".

Obwohl diesmal weggelassen, wird die obige Phaseneinstellung von "Scenario.add ()", "Scenario.combine ()", "Scenario.delete ()", "Scenario.separate ()", "Scenario.clear ()" durchgeführt. Es ist auch möglich zu bearbeiten.

6. Implementierung

Ich verwende das Ruptures-Paket, um die Anzahl der Verzweigungspunkte in Phase zu ermitteln und die Anzahl der Wiederhersteller zu berechnen, die die Verzweigungspunkte sein werden. Ich habe ursprünglich fbprophet verwendet, das als Zeitreihenanalysepaket bekannt ist, aber es gab ein Problem, dass die Anzahl der Verzweigungspunkte manuell festgelegt werden musste. Ilyass Tabial of Colaborator sagte mir [^ 9], und ich ermöglichte es, die Anzahl der Verzweigungspunkte mithilfe von Brüchen automatisch zu schätzen.

Ab covsirphy.ChangeFinder-Implementierungscode lauten die Hauptschritte wie folgt.

1. Holen Sie sich die S- und R-Werte für jedes Datum
2. Berechnen Sie log10 (S)
3. Setzen Sie den Wert von log10 (S) für R = 1, 2, 3 .: Wenn er nicht in den Daten enthalten ist, geben Sie die Werte vorher und nachher ein.
4. Stichproben zur Reduzierung des Rechenaufwands: Anzahl der Stichproben = Anzahl der Originaldaten
5. Verwenden Sie das Ruptures-Paket, um die Menge der Punkte (R, log10 (S)) zu ermitteln, an denen sich log10 (S) in Bezug auf das Ausmaß der Änderung der Anzahl der Wiederhersteller erheblich ändert.
6. Sortieren Sie log10 (S) des Satzes von Änderungspunkten mit den Originaldaten und ermitteln Sie das Datum jedes Änderungspunkts (ungefährer Wert, falls nicht in den Originaldaten. Log10 (S), da die Anzahl der bestätigten Fälle eine größere Abweichung als R aufweist. Es wurde benutzt)
7. Berechnen Sie das Start- und Enddatum jeder Phase ab dem Änderungsdatum

Nachwort

Das nächste Mal werde ich erklären, wie die Parameter für jede Phase geschätzt werden.

Danke für deine harte Arbeit!