Themen aus der Numerik mit aktuellem Bezug

Zusammenstellung von Parametern von SARS-CoV-2 (Achtung, diese Parameter sind keine Naturkonstanten):

Serial interval/Serielles Interval: ca. 4 Tage (z.B. hier aus chinesischen Daten);

Das RKI nimmt zur Berechnung von R₀ eine Generationenzeit von 4 Tagen an.

Als „Serielles Interval“ wird die Zeit zwischen Symptombeginn bei Person A und Symptombeginn bei einer angesteckten Person B bezeichnet. Als Generationenzeit die Zeit zwischen der Infektion von A und der Infektion von B. Diese beide Zahlen werden ähnlich oder fast gleich sein.

Basisreproduktionszahl R₀: Laut RKI zwischen 2,4 und 3,3

Rückwärtsschau zum zeitlichen Verlauf von R von Anfang März bis Anfang April 2020 in Deutschland (Berechnungen RKI, die angenommene Generationenzeit ist 4 Tage). Hier die Antwort des RKI zur einer oft gestellten Frage.

Inkubationszeit: im Median 5-6 Tage bei starker Spreizung

• Letalität: Laut RKI am 25.4.2020 in Deutschland 3,6 Prozent. Weitere Zahlen hat das RKI hier zusammengestellt. Die tatsächliche Letalität ist aufgrund der Dunkelziffer schwer zu bestimmen. Sie ist zudem stark altersabhängig.
  Frühe chinesische Zahlen lauteten (finde auf die Schnelle nur diese Sekundärquelle):
• 10-19 Jahre: 0,2 Prozent
• 20-29 Jahre: 0,2 Prozent
• 30-39 Jahre: 0,2 Prozent
• 40-49 Jahre: 0,4 Prozent
• 50-50 Jahre: 1,3 Prozent
• 60-69 Jahre: 3,6 Prozent
• 70-79 Jahre: 8 Prozent
• 80 Jahre: 14,8 Prozent

Bemerkung: Es gibt bisher wenig Daten zu Situationen ohne Dunkelziffer. Auf einem Flugzeugträger wurde die gesamte Besatzung von 4800 Besatzungsmitgliedern getestet und 580 Tests fielen positiv auf SARS-CoV-2 aus. Ein Besatzungsmitglied ist inzwischen an Covid-19 gestorben.

 

Um Entscheidungen zu Epidemien zu treffen (oder zu verstehen) muss man zuerst exponentielles Wachstum verstehen: aus 100 Infizierten werden schon nach 20 Verdopplungen 100 Millionen Infizierte.

Dieses explosive Wachstum führt dazu, dass Krankenhäuser überlastet werden. Um die Geschwindigkeit des Wachstums zu reduzieren, kam es in fast allen Ländern der Welt mit exponentiell schnellen Ausbrüchen zu radikalen Maßnahmen, wie etwa die Kontaktbeschränkungen (+weiteren Maßnahmen) in Deutschland. Einige Länder kamen ohne solche Maßnahme aus, weil Sie den jeden Ausbruch früh genug unter Kontrolle bekommen haben. Italien, Frankreich und Spanien hatten dagegen deutlich strengere Maßnahmen als Deutschland (d.h. echte Ausgangssperren), mit ähnlichem Erfolg wie Deutschland.

Exponentielles Wachstum tritt am Anfang einer Epidemie auf. Durch erworbene Immunität schwächt sich die Epidemie ab, sobald ein substantieller Anteil der Bevölkerung die Infektion durchgemacht hat.

Die einfachsten Modelle, die dieses Verhalten beschreiben, sind die schon erwähnten SIR-Modelle und verfeinerte Varianten davon. Solche Modelle sind nötig, um Fragen zu beantworten wie: Wie hoch wird die Zahl der gleichzeitig Infizierten maximal sein und wann wird das sein? Wie hoch wird der Anteil der Bevölkerung sein, der sich im schlimmsten Fall ansteckt? Auch das RKI verwendet Varianten von SIR-Modellen für Grippesimulationen und jetzt auch für SARS-CoV-2.

Es gibt ein paar Online-Rechner zu SIR-Modellen, zum Teil mit sehr vielen Parametern. Das einfachste Modell kommt jedoch mit zwei Parametern aus: alpha und beta. Die sogenannte Basisreproduktionszahl R₀, über die zur Zeit so viel diskutiert wird, ist dann der Quotient beta/alpha.

Das Modell hat nur drei Variablen:

• S (susceptibles; Ansteckbare)
• I (infected; Infizierte)
• R (recovered/removed; Geheilte oder Verstorbene)

Die Werte der Variablen liegen immer zwischen 0 und 1, was sich dann als 0 Prozent bis 100 Prozent der Bevölkerung interpretieren lässt.

Das SIR-Modell wird durch ein System gewöhnlicher Differentialgleichungen beschrieben (vielleicht folgt in einem weiteren Beitrag die diskretisierte Variante, die dann ohne sichtbare Ableitung auskommt), in der einfachsten Variante:

• S'(t)=-beta*S(t)*I(t)
• I'(t)=beta*S(t)*I(t) – alpha*I(t)
• R'(t)=alpha*I(t)

Hier gibt die Konstante beta an, wieviele andere ein Infizierter pro Zeiteinheit (z.B. pro Tag) ansteckt. Aus I(t) Infizierten (z.B. 10 Prozent der Bevölkerung) werden also nach einem Tag beta*I(t) Infizierte (also z.B. 15 Prozen, wenn beta=3/2), sofern jeder infizierbar ist. Da aber nur S(t) infizierbar ist (90 Prozent sind noch infizierbar, weil 10 Prozent ja schon infiziert sind), wird noch mit S(t) multipliziert.

Die mittlere Gleichung sagt also: Pro Tag kommen beta*S(t)*I(t) Infizierte dazu (immer in Prozent gerechnet, da alle Variablen zwischen 0 und 1 sind).

Die erste Gleichung sagt, dass die Anzahl der Infizierbaren entsprechend sinkt.

Die letzte Gleichung sagt, dass pro Zeiteinheit (z.B. pro Tag) ein Anteil alpha der Infizierten gesundet (oder stirbt, das wird hier nicht unterschieden). Daher ist 1/alpha ist die durchschnittliche Krankheitsdauer (z.B. 7 Tage).

Wenn beta die Infektionen pro Tag sind und 1/alpha die Krankheitsdauer in Tagen, dann sind beta/alpha die Infektionen, die ein Kranker während seiner gesamten Krankheit verursacht (wenn niemand immun oder bereits infiziert ist). Das ist die Basisreproduktionszahl R₀.

Das Modell weicht natürlich von der Realität der SARS-CoV-2-Epidemie ab. Beispiel: Ein Infizierter ist im Modell während seiner Krankheit jeden Tag gleich ansteckend. Bei SARS-CoV-2 ist es so, dass ein Infizierter in den ersten Tagen nach seiner eigenen Infektion besonders ansteckend ist (mit einem Höhepunkt am Tag unmittelbar vor dem Symptombeginn), obwohl er noch keine Symptome hat (im Median 5 Tage lang). In der zweiten Woche ist ein Infizierter kaum mehr ansteckend (außer bei engem Kontakt z.B. mit Krankenhauspersonal) – obwohl er sich dann erst richtig krank fühlt.

Aus diesem Grund wollen wir für unser Modell 1/alpha = 7 Tage Krankheitsdauer annehmen. Das heißt, wir nehmen an, dass alle Ansteckungen während der ersten 7 Tage passiert (und ignorieren die zweite Krankheitswoche des durchschnittlichen Covid-19-Patienten).

Wenn wir nun versuchen, beta so zu wählen, dass es auf die SARS-CoV-2-Infektionsdaten in Deutschland vom Anfang März 2020 passt, bekommen wir den folgenden Verlauf (S blau, I rot, R grün). Tag 0 ist wieder der 1. März 2020. Die genauen Parameter sind beta=2.85/7 und alpha=1/7. Dementsprechend ist R₀=2.85, also knapp unter 3. Das passt grob zur (rückwärts berechneten) Schätzung des RKI von R₀ für Anfang März 2020 und zu anderen Quellen.

Dieses Szenario beschreibt, was passiert wäre, hätte es gar keine eindämmenden Maßnahmen gegeben – für diesen Fall ist dieses einfache Modell sogar relativ gut geeignet. Man sieht, dass kurz nach Tag 50 (genauer: am Tag 54, das wäre der 23. April 2020 gewesen, also vorgestern) das Maximum bei den Infizierten erreicht wird (rote Kurve): Es sind an diesem Tag fast 30 Prozent (genauer: 28 Prozent) der Bevölkerung gleichzeitig infiziert. Eine Katastrophe für die Krankenhäuser, denn selbst wenn nur ein Prozent schwere Verläufe angenommen werden, bedeutet das über 200000 Schwerkranke (bei aktuell 30000 Intensivbetten in Deutschland, von denen maximal 2/3 für Covid-19 zur Verfügung stehen werden.)

Ingesamt werden innerhalb der 100 Tage Berechnungszeit über 90 Prozent der Bevölkerung infiziert (grüne Kurve).

Ich habe oben geschrieben, dass beta=2.85/7 (bei angenommenen alpha=1/7, also 7 Tage Erkrankungsdauer) so gewählt wurde, dass es zu den deutschen Daten von Anfang März passt. Hier ist dargestellt, wie gut es passt:

Der Vergleich von I(t) (rote Kurve) mit den Infektionszahlen (Daten von NTV, siehe den allerersten Blog-Beitrag) geteilt durch 80 Millionen (schwarze Punkte). Die Kurve passt gut.

Hier die Messwerte (schwarze Punkte) noch einmal, aber eingetragen ein einen semilogarithmischen Plot von S, I und R (in der Standarddarstellung wie oben, könnte man nichts zu erkennen).

Man sieht, dass die rote Kurve sehr gut zu den Infektionszahlen in Deutschland (1. März bis 21. März 2020) passt. Da die y-Achse logarithmsch aufgetragen ist (semilog) ist, ist die Phase des exponentiellen Wachstums als Geradenabschnitt zu sehen. Die Steigung passt gut zusammen.

Im ersten Beitrag haben wir gesehen, dass das exponentielle Modell für den Tag 40 etwa 4,4 Millionen Infizierte prognostiziert hat. In der semilogarithmischen Plot (y-Achse logarithmisch) sieht man, dass das exponentielle Modell ungefähr bis Tag 45 gültig ist, hier der Ausschnitt bis Tag 45:

Das bedeutet aber, dass auch dieses Modell (ohne Eindämmungsmaßnahmen) vorhersagt, dass innerhalb von 40 Tagen aus 100 Infizierten mehr als 4 Millionen werden (etwa 5 Prozent der Bevölkerung bzw. 5e-2).

Nach Tag 45 ist die Phase des exponentiellem Wachstum (mit konstantem Wachstumsfaktor) beendet. Allerdings sind nach Tag 46 fast 15 Prozent der Bevölkerung infiziert (entsprechend 12 Mio Infizierten) – erst dann kommt es zu einer sichtbaren Abschwächung des Wachstums.

Bemerkung zu R₀: Wenn man aus echten Daten den Wert von R0 schätzen will, muss man immer eine Generationenzeit annehmen. Das RKI berechnet R₀ so (S. 14):

Schlussbemerkungen: Wir haben nur zwei Parameter (beta und alpha), dennoch gibt es noch andere Kombinationen von beta und alpha, die ähnlich gut zu den Daten passen. Sie liefern aber auch Ergebnisse, die ähnlich sind.

Eigentlich wollte ich Ihnen Octave-Implementierung (entweder mit Euler oder mit dem eingebautem ODE-Löser) an die Hand geben, damit Sie die Rechnungen nachvollziehen können – das verschieben wir aber.

Korrektur (06.07.2020): In der zweiten Gleichung des (einfachsten) SIR-Modells fehlte „- alpha*I(t)“. Danke für den Hinweis.

 

 

 

Nach Medienberichten hat der Regierungssprecher Helge Braun die Strategie der Bundesregierung mit den Worten beschrieben: „Die Zahl der täglichen Neuinfektionen sollte auf einem gleichbleibenden Niveau sein. Das entspricht einer Basisreproduktionszahl R von 1,0 und bedeutet, dass ein Infizierter nur einen weiteren Menschen ansteckt.“ (Zitiert nach T-online).

Eine Herdenimmunität durch Durchseuchung wird demnach nicht angestrebt (wegen sicherer Überlastung des Gesundheitssystems): „Um nur die Hälfte der deutschen Bevölkerung in 18 Monaten zu immunisieren, müssten sich jeden Tag 73.000 Menschen mit Corona infizieren. […] So hohe Zahlen würde unser Gesundheitssystem nicht verkraften und könnten auch von den Gesundheitsämtern nicht nachverfolgt werden.“ (Zitiert nach T-online).

Und „Daher lautet die Strategie, Ansteckungen zu vermeiden und bezüglich der Immunität auf die Einsatzfähigkeit eines Impfstoffs zu warten.“ (Zitiert nach T-online).

Damit ist klar auch, dass eine (kurzfristig teure) Eindämmung (mit R deutlich kleiner als 1, etwa nach dem Vorbild asiatischer Länder) von der Bundesregierung nicht angestrebt wird. Als Grund wird angegeben: „Politisch müssen wir bedenken, dass Deutschland sich mitten in Europa aufgrund der Pendlerströme und Wirtschaftsverkehre nicht so gut abschotten kann und will. […] Selbst wenn wir das Virus stark zurückdrängen, kommt es dann aus dem Ausland immer wieder zurück.“ (Zitiert nach T-online).

Dazu ist zu bemerken, dass die bislang geltenden Maßnahmen einschneidend und (kurzfristig) sehr teuer sind. Allerdings bedeutet ein R nahe 1 eine instabile Lage, die jederzeit außer Kontrolle geraten kann – und selbst im besten Fall ein langes Andauern der Epidemie-Lage bedeutet, bis ein Impfstoff zur Verfügung steht (Simulationen dazu folgen).

Es ist auch festzustellen, dass die Strategie „R nahe 1“ (+Spekulieren auf einen Impfstoff) als auch die Umsetzung der Stategie den wissenschaftlichen Empfehlungen widerspricht:

• Diese Strategie widerspricht der Empfehlung der Helmholtz-Iniative.
• Die Bedingungen, die die Leopoldina für eine Lockerung der Maßnahmen (Maskenpflicht, Verfolgungsapp) sind nicht erfüllt worden.

Bemerkung (20.4.2020): Inzwischen hat nach Sachsen und Mecklenburg-Vorpommern auch Bayern eine Maskenpflicht angekündigt. Jena hat schon länger eine Maskenpflicht und ist damit erfolgreich.

Bemerkung (22.4.2020): Weiter Bundesländer haben inzwischen nachgezogen und werden eine (teilweise) Maskenplicht einführen. Bislang galt nur eine dringende Aufforderung zum Maskentragen.

Bemerkung zur App (22.4.2020): Die Daten einer Infektionsverfolgungsapp lassen sich am besten schützen, wenn sie gar nicht erst zentral gespeichert werden.

Schlussbemerkung: Eine Optimierung von R mit der Anzahl der Intensivbetten (30000, davon aktuell 13000 frei) oder der Anzahl der Beatmungsgeräte als obere Schranke (Optimierung unter Nebenbedingungen) vernachlässigt, dass auch unter diesen Bedingungen viele Patienten sterben oder bleibende Schäden behalten. So sterben offenbar aktuell in New York über 80 Prozent der Covid-19-Patienten, die beatmet werden, in Großbritannien zwei Drittel.

Bemerkung (zu R): Heute wird in den Medien die Schätzung von R durch das RKI viel diskutiert. Gemeint ist dabei Abbildung 5 im schon oben verlinkten Bulletin des RKI von 15.4.2020 (hier der direkte Link zum PDF, Abb. 5 ist auf S. 14). Wichtig dabei ist, dass die Daten bearbeitet wurden (Nowcasting, eine Form der Extrapolation), d.h. sie wurden rechnerisch „um Diagnose-, Melde-, und Übermittlungsverzug“ korrigiert. Das heißt, die Korrektur soll etwa kompensieren, dass Gesundheitsämter Fälle noch nicht in das Meldesystem des RKI eingetragen haben. Außerdem dauert es einige Zeit, bis nach einer Erkrankung eine sichere Diagnose gestellt wird.
Wichtig in unserem Zusammenhang ist, wie das RKI begründet, dass R nach Mitte März nicht weiter gesunken ist (S. 17): 1. Die Testkapazitäten wurden erweitert, so dass mehr Fälle sichtbar wurden. 2. Mehr ältere Menschen haben sich angesteckt, auch durch Ausbrüche in Heimen und Krankenhäusern.
Während der erste Effekt rechnerisch ist, ist der zweite Effekt real: Auf dem Weg durch die Bevölkerungsgruppen trifft das Virus auch auf neue, schnellere Verbreitungsmöglichkeiten, wie etwa Pflegeheime. Es ist zu erwarten, dass dieser Effekt andauert.

Bemerkung (25.4.2020): Inzwischen hat die Bundeskanzlerin in einer Rede im Bundestag die Lockerungen als „zu forsch“ bezeichnet und auch das Robert-Koch-Institut warnt nun offiziell in einer Pressekonferenz vor weiteren Lockerungen. Weitere Lockerungen seien erst bei „wenigen Hundert“ Infektionsfällen pro Tag möglich, weil erst dann die Gesundheitsämter wieder zur Kontaktverfolgung (+Quarantäne) fähig sind. Aber: „Derzeit melden die Ämter im Durchschnitt täglich noch rund 2.000 Neuinfektionen in Deutschland.“. – Allerdings scheint es, als ob ein großer Teil der Medien, Wirtschaft und ein großer Teil der Bundesländer nun andere Prioritäten im Vordergrund sehen wollen.

Bemerkung (26.4.2020): Schon am 22.4. hat sich das RKI vorsichtig kritisch zu Lockerungen geäußert: „Eine unkontrollierte Lockerung der Maßnahmen und eine Rückkehr zu „prä-Pandemie-Verhalten“ würde somit zu einem erneuten Anstieg der täglichen Fallzahlen und einer Annäherung der effektiven Reproduktionszahl an die Basisreproduktionszahl führen.“

Bemerkung (26.4.2020): Christian Drosten äußert sich in der SZ zu R und äußert sich kritisch über die Lockerungen.

 

Anfang März 2020 verdoppelte sich die Zahl der SARS-CoV-2-Infizierten alle 2 bis 3 Tage (siehe auch die Abbildung hier).

Die Verdopplungszeit ist eine griffige Art, exponentielles Wachstum und seine Geschwindigkeit greifbarer zu machen. Relativ schnell haben fast alle Medien Verdopplungszeiten zusätzlich zu den Infektionszahlen angegeben, offenbar (auch?) inspiriert von einem Blog. (Der kleine Bruder der Verdopplungszeit ist die Halbwertszeit – Ältere erinnern sich noch an Tschernobyl, aber das ist eine andere Geschichte.)

Bei exponentiellem Wachstum f(t)=e^at für ein a>0 ist die Verdopplungszeit t_V=ln(2)/a, wie man schnell nachrechnet. Eine schnelle Probe: Für a=ln(2), also f(t)=e^ln(2)t=2^t ist die Verdopplungszeit dann natürlich t_V=1.

Folgt eine Größe (etwa die Infektionszahlen) einem exponentiellem Wachstumsgesetz f(t)=e^at, ist die Verdopplungszeit konstant (s.o.). Ändert sich die Verdopplungszeit, ist es keine Exponentialfunktion e^at mehr.

Übrigens, fast alle BWLer und Bänker kennen die sogenannte 70er-Regel: Wächst die Zahl der Infizierten jeden Tag um x Prozent, dann ist die Verdopplungszeit (ungefähr) 70/x. Bei 5 Prozent Zinsen verdoppelt sich die Einlage also in 14 Jahren. Die Regel kommt daher, dass ln(2) ungefähr 0,69315 ist, also fast 0,7.

Nochmal Übrigens: Jeder Informatiker hat im Kopf, dass 2¹⁰=1024, also ungefähr 1000. Demnach ist 2²⁰ unfähr eine Million. Das wiederum heißt, dass aus einem einzigen Infiziertem nach 26 Verdopplungen mehr als 64 Millionen Infizierte werden. Übrigens, 1 KiB=1024 Byte und 1 kB=1000 Byte.

Bei einer Verdopplungszeit von 2,6 Tagen, wie sie SARS-CoV-2 sie in Deutschland Anfang März hatte (siehe hier), hätten also etwas mehr als 2 Monate genügt, bis die gesamte deutsche Bevölkerung infiziert gewesen wäre. (Allerdings folgt eine Durchseuchung nur in der Anfangsphase einem exponentiellen Gesetz.)

Die Verdopplungszeit machte Karriere, als ihre Verlängerung auf 10 oder mehr Tage erklärtes Ziel der Bundesregierung wurde. Spätestens jetzt fingen Diskussionen an, denn wie die Verdopplungszeit berechnet werden soll, wenn gar keine Exponentialfunktion e^at vorliegt. Die Zeitung „Die Zeit“ etwa gibt in ihren Graphiken an, wieviele Tage es her ist, dass die Infiziertenzahl halb so hoch war wie heute. Andere Medien scheinen die Verdopplungszeit anhand der Infiziertenzahlen der letzten (beiden?) Tage zu berechnen.

Ein Analysis-Professor aus Bielefeld hat in den Medien (etwa hier und hier) und auf Youtube die Verdopplungszeit(en) und ihre Berechnung diskutiert und ein mathematisch relativ naheliegendes Verfahren zur Berechnung vorgeschlagen: Da die Verdopplungszeit sich ändert, sollte man Daten aus den letzten Tagen nehmen, aber mehr als nur die letzten zwei Tage, um Schwankungen wegen Wochenendeffekten etc. zu verringern. Zum Beispiel könnte man daher die letzten 7 Tage nehmen.

Diese Zahlen werden dann logarithmiert. Anders als hier, zu Anfang der Epidemie und ohne Maßnahmen, werden die Datenpunkte nun eher nicht auf einer Geraden liegen. Daher ist es sinnvoll, mit der Methode der Kleinsten Quadrate eine (möglichst einfache) Kurve die Daten anzupassen. Einfach ist etwa ein Polynom zweiten Grades. Letztlich ist man an dann der Steigung der Kurve am heutigen Tag interessiert, also an der Ableitung des Polynoms zum heutigen Zeitpunkt.

Wenn man die Verdopplungszeiten als Leistungsparameter definieren will, ist das sicher eine gute Methode. Das Verfahren ist z.B. exakt, wenn (nach Logarithmierung) tatsächlich ein Polynom zweiten Grades vorliegt (anders als wenn man nur die letzten zwei Tage auswertet).

Letztlich ist aber die Verdopplungszeit vielleicht gar kein gutes Maß, denn wir wollen die Epidemie ja eindämmen, also bestenfalls überhaupt kein Wachstum mehr haben. Tatsächlich weisen eigene Medien für asiatische Länder, die erfolgreiche Eindämmung betrieben haben (erfolgreicher als Europa?) , Verdopplungszeiten von unendlich aus.

Stattdessen wird jetzt über die Reproduktionszahl R bzw. R₀ gesprochen. Dieser Begriff hat aber auch seine Schwierigkeit, weil dazu eigentlich eine Zeitskala gehört, etwa die Generationenzeit. Hier ein Artikel zur Schätzung von R₀ und der Generationenzeit für SARS-CoV-2 zum Zeit des exponentiellen Wachstums in China.