A basic SIR model for the epidemic.
Total population 327 million, with a single initiating infection. 12 day course of disease, initial R0 of 3.5 (as per the stats from China that Steve linked.) Assume 5% of cases are critical, 2% of those critical cases die with ICU care, 5% with general hospital care, and 50% with no care. 90,000 available ICU beds and 900,000 available hospital beds.
Run this model through to completion and it sweeps through the population in about a year, infecting most everyone (less about 9 million who escape), killing 5.4 million.
Now, suppose we impose infection controls on day 80, right about when there’s 1000 deaths from this thing. And then we vary how strong those controls are: from 0.35 (what the Chinese managed) up to nothing at all.
Here we see how the # of deaths varies with the strength of our controls. If we impose Chinese-style controls, we get away with only 5K deaths, or 1000-fold fewer than without the controls. But the inflection point isn’t just at 1.0. In particular if we can get the R0 below about 1.5 that gets us down under 500K, or a 10-fold reduction. At 1.0 we get down to about 50K.
If we take the slope of this plot, we can see the marginal impact of better controls at each point, in terms of # of deaths saved:
Why the discontinuity at R0 of just below 1.5? Well, here’s the number of critical cases over time for several different scenarios. The vertical line is when we impose controls, the horizontal lines are ICU bed capacity and total hospital bed capacity. Right below an R0 of 1.5, we dip the # of critical cases low enough that we never exceed hospital capacity. The trick here though is that until you get R0 below 1.0 you’re still sweeping through most of the population, just more slowly, slowly enough that the hospitals can handle it. And that means that, for example, if you get R0 to exactly 1.0 you have to keep those infection controls in place for ~5 years or else you’ll be back up close to 5million dead again.
Appendix: code to generate this model.
from __future__ import division import numpy as np import matplotlib.pyplot as plt from collections import namedtuple PopState = namedtuple("PopState", [ "day", "susceptible", "infected", "dead", "recovered", ]) EpiParams = namedtuple("EpiParams", [ "days_infected", "infection_rate", "frac_critical", "frac_dead_icu", "frac_dead_hospital", "frac_dead_home", "start_day", "controlled_infection_rate", ]) ICU_BEDS = 90000 HOSPITAL_BEDS = 900000 def update_state(state, params): day = state.day S = state.susceptible I = state.infected D = state.dead R = state.recovered N = S + I + D + R days_betw_contacts = params.days_infected / params.infection_rate if day >= params.start_day: days_betw_contacts = params.days_infected / params.controlled_infection_rate new_infections = (I * S) / (N * days_betw_contacts) recovered_infections = (I / params.days_infected) * (1 - params.frac_critical) # For critical cases, we split into ICU, hospital, and home care, in descending order C = I * params.frac_critical; C_icu = min(C, ICU_BEDS); C_hospital = 0 if C < ICU_BEDS else min(C - ICU_BEDS, HOSPITAL_BEDS) C_home = 0 if C < (ICU_BEDS + HOSPITAL_BEDS) else C - (ICU_BEDS + HOSPITAL_BEDS) recovered_critical = ( ((C_icu / params.days_infected) * (1 - params.frac_dead_icu)) + ((C_hospital / params.days_infected) * (1 - params.frac_dead_hospital)) + ((C_home / params.days_infected) * (1 - params.frac_dead_home)) ) dead_critical = ( ((C_icu / params.days_infected) * (params.frac_dead_icu)) + ((C_hospital / params.days_infected) * (params.frac_dead_hospital)) + ((C_home / params.days_infected) * (params.frac_dead_home)) ) newS = S - new_infections newI = I + new_infections - recovered_infections - recovered_critical - dead_critical newR = R + recovered_infections + recovered_critical newD = D + dead_critical return PopState(day + 1, newS, newI, newD, newR) def run_model(state, params): infections = [state.infected] deaths = [state.dead] while (state.infected) >= 1: state = update_state(state, params) infections.append(state.infected) deaths.append(state.dead) return infections, deaths, state total_population = 327 * 1e6 init_state = PopState( day=1, susceptible=total_population - 1, infected=1, dead=0, recovered=0) params = EpiParams( days_infected=12, infection_rate=3.5, frac_critical=0.05, frac_dead_icu=0.02, frac_dead_hospital=0.05, frac_dead_home=0.5, start_day=79, # Modify this parameter to vary the strength of infection controls controlled_infection_rate=3.5, ) infections_over_time, deaths_over_time, final_state = run_model(init_state, params)