Lessons from Epidemiological Models

Richard Wiener, Research Corporation for Science Advancement, rwiener@rescorp.org

Two or greater is bad

In early February I attended a presentation by Myron Cohen. Cohen is the chief architect of a clinical trial recognized as Science journal’s 2011 Breakthrough of the Year. The trial showed remarkable success for a new protocol to prevent transmission from HIV infected people to non-infected sexual partners. This breakthrough was the culmination of three decades of Cohen’s inspiring work to reduce the horrific toll HIV exacts on humanity.

At the end of his talk, Cohen offered to take questions on his HIV research or on the novel coronavirus outbreak. Not surprisingly most of the questions pertained to the latter, given its currency in the news. Cohen said several things that stuck with me. He said not to travel to Asia until we have a better understanding of the virus (even though at that time there were few confirmed cases outside China). He also said prophetically that epidemiologists are greatly concerned for the safety of health care workers and the potential for community spread, since estimates of the reproduction number $R_0$ were between 2 and 5.

That last statement got me to sit up and take notice. It’s been 15 years since I taught nonlinear dynamics, but I remember enough to know $R_0\ge 2$ is bad. Really bad. I dusted off Strogatz’ Nonlinear Dynamics and Chaos (1) and looked up the SIR model, which describes the dynamics of an unconstrained epidemic. I also read the Wikipedia page on Compartmental Models in Epidemiology (2). These resources are excellent background to understand how epidemiological models can be applied to the COVID-19 pandemic. A couple examples which stand out: Goldenfeld and Maslov using an SEIR model to inform policymakers in Illinois (3), and Maier and Brockmann (4) introducing the SIR-X model, which provides important insights into the effects of containment measures.

In this article I’ll first provide a brief summary of the SIR model, then discuss lessons from slightly more involved models applied to COVID-19, and end with some open questions.

The SIR Model

The SIR model posits a total population comprised of three compartments (i.e. types of individuals) with respect to an infectious disease: susceptible to infection, infectious, and removed from the transmission process, either through recovery with conferred immunity or by death. The model can be written as a third-order system of nonlinear ordinary differential equations with time-dependent variables $S(t)$, $I(t)$, and $R(t)$ representing the fraction of the population in each compartment:

${\partial }_{t}S= -\alpha SI$ (1)

${\partial }_{t}I= \alpha SI- \beta I$ (2)

${\partial }_{t}R= \beta I$ (3)

where ${\partial }_t$ means differentiation with respect to time, and positive constants $\alpha$ and $\beta$ represent the transmission and removal rates. The ratio of these rates defines the reproduction number

$R_0=\raisebox{1ex}{$\alpha $}\!\left/ \!\raisebox{-1ex}{$\beta $}\right.$ (4)

(By perhaps unfortunate convention $R(t)$ is used to represent the dynamic variable for the removed fraction, and $R_0$ to represent a parameter.) The inverse of these constants, ${\alpha }^{-1}$ and ${\beta }^{-1}$, are the time between transmissions, when most of the population is susceptible, and the duration of infectiousness, respectively. The ratio of the latter to the former time equals $R_0$ and represents the mean number of new infections which result from an infectious person’s interactions with susceptible people during the early part of an outbreak, which is a good way to conceptualize the meaning of $R_0$.

The model assumes people can only change compartments irreversibly from $S$ to $I$ to $R$ and the population is well-mixed so that susceptible and infectious people interact at a rate proportional to the size of these groups and each interaction has a constant probability of transmission. The model also assumes the rate at which infectious people are removed (i.e. recover or die) is constant, and the time a person is infectious coincides with the time he or she is infected. For a real population these may be poor assumptions. For example, voluntary changes in population behavior and government interventions such as requiring people to shelter in place reduce interactions between susceptible and infectious people. Also, people may vary in infectiousness, or those who recover may not acquire full immunity and be completely removed from the transmission process. The SIR model illustrates the fundamental dynamics of an unconstrained outbreak without the complications of a real-world epidemic. Other epidemiological models contain various more realistic assumptions.

The initial condition of most interest is when almost the entire population is susceptible with only a few people infectious, as when a viral infection is introduced from a wild animal or a handful of infectious people travel to an uninfected region. Setting the initial susceptible fraction to its approximate value $S_{in}=1$ and integrating Eq. (2) yields

$I\left(t\right)\cong I_{in}{e}^{(R_0-1)\beta t}$ (5)

as the approximate behavior of $I$ when the system first starts evolving. Eq. (5) shows $R_0=1$ is the threshold for an outbreak, with initial exponential growth for $R_0>1$ and exponential decay for $R_0<1$. This makes sense, since a reproduction number greater than unity means on average each infectious person transmits more than one new infection when the population is mostly susceptible.

Since $R$ plays no role in the dynamics of $S$ and $I$, the SIR model is equivalent to a second-order system consisting of Eqs. (1) and (2). (Solving this system for $S$ and $I$, the fraction removed can be determined from $R=1-S-I$.) Dividing Eq. (2) by Eq. (1), separating variables and integrating, and then substituting approximate initial values of $S_{in}=1$ and $I_{in}=0$ gives the trajectories (i.e. solution curves) in the phase plane which intersect the point $\left(S,I\right)=(\mathrm{1,} 0)$:

$I=1-S+R_0^{-1}lnS$ (6)

The system has fixed points $(S^{*},I^{*})$ when ${\partial }_{t}S={\partial }_{t}I=0$, which occurs when $I^{*}=0$, and can be found by setting Eq. (6) equal to zero and solving for $S^{*}$ as a function of $R_0$. For $R_0=2$, the minimum estimate for COVID-19 given by Cohen, $\left(S^{*}, I^{*}\right)=(\mathrm{1,} 0)$ and $(0.2, 0)$. Trajectories in the phase plane flow from the former to the latter fixed point, since the dynamics of the system can only irreversibly reduce $S$, while $I$ initially grows nearly exponentially, so long as $R_0>1$, before peaking and then decreasing to zero, at which point the epidemic is over. A nonzero fraction of infectious people is enough to spark an epidemic. With $R_0=2$ only 20% of those susceptible at the start of the outbreak would escape infection. This simple model provides a sobering back of the envelope estimate of the devastating effect an unconstrained epidemic running its course could have on a population.

Taking the derivative of $I$ with respect to $S$ in Eq. (6) and setting it equal to zero yields ${S=R}_0^{-1}$ when $I$ peaks. For $R_0=2$, at the peak of infection 50% of the population is still susceptible, 15% is infectious, and 35% is removed from the transmission process. Such a peak could massively overwhelm healthcare systems for a disease such as COVID-19 in which even a small percentage of those infected become critically ill.

Numerical solutions of the SIR model depend on the size of the initial infectious fraction. A solution with $I_{in}{=10}^{-6}$ and $R_0=2$ shows $I$ peaking after a time of about ${13\beta }^{-1}$. For ${\beta }^{-1}$ of roughly 10 days, as for COVID-19, this would mean peaking four months after an outbreak begins, if the disease is left to spread unconstrained.

Stay at Home

According to APS News (3), just before spring break in early March at the University of Illinois, Urbana-Champaign, UIUC physicists Goldenfeld and Maslov collaborated on an SEIR model of COVID-19. An SEIR model, which is only slightly more involved than the SIR model, uses an additional compartment $E$ (for exposed) when a disease has a nonnegligible latency period during which an infected individual is not yet infectious, and assumes people move irreversibly from $S$ to $E$ to $I$ to $R$. According to the model, which took just a few hours to analyze, if students returned to campus after spring break, there would be a huge wave of infections. Alerting university administrators to the model’s results led to a rapid decision to move classes online. Eleven days later, citing Goldenfeld and Maslov’s modeling as part of the rationale, the governor of Illinois issued a statewide stay-at-home order.

An obvious question is, based on epidemiological modeling to what extent do containment measures such as a stay-at-home order flatten the curve?

The SIR-X Model

According to the SIR model, if a population could be constrained to avoid all interaction between susceptible and infectious people (i.e. effectively $\alpha =0$), then no new infections would occur, and all infectious people would be removed and the epidemic completely suppressed after a time of ${\beta }^{-1}$. Unfortunately, the danger posed by an outbreak typically goes unrecognized until a substantial number of people have been infected, at which point it may not be practical to identify all infectious individuals, including those who are asymptomatic, and prevent their interaction with those who are susceptible. Instead, as happened with COVID-19, it becomes necessary to implement blunter measures to reduce interaction, such as population-wide isolation through stay-at-home orders as well as quarantine of symptomatic infectious individuals. Maier and Brockmann include the effects of containment policies that deplete the susceptible and infectious fractions of a population, thereby reducing their interaction, by introducing the SIR-X model (4):

${\partial }_{t}S= -\alpha SI-{\kappa }_{0}S$ (7)

${\partial }_{t}I= \alpha SI- \beta I-{\kappa }_{0}I-\kappa I$ (8)

${\partial }_{t}R= \beta I+{\kappa }_{0}S$ (9)

${\partial }_{t}X=(\kappa +{\kappa }_0)I$ (10)

For this model, physical distancing measures applied to the whole population deplete individuals from both the $S$ and $I$ compartments at the same rate ${\kappa }_0$. Additionally, quarantining of those who test positive or are symptomatic depletes infectious individuals at a rate $\kappa$. Both ${\kappa }_0$ and $\kappa$ are positive constants. There is an additional compartment $X$ which quantifies infectious individuals who have been separated from the transmission process by containment or quarantine. The model assumes this new fraction of the population $X(t)$ is proportional to confirmed cases of a disease.

The effective reproduction number is then

$R_{\mathrm{0,}eff}=\alpha /(\beta +{\kappa }_0+\kappa )$ (11)

The key result of the SIR-X model is that the explosiveness of an outbreak is damped because ${R_{\mathrm{0,}eff}<R}_0$. Protection of the susceptible fraction of the population by containment and quarantine leads to initial subexponential growth for $X(t)$, with a power law scaling $t^{\mu }$, for a wide range of model parameters. After a period of algebraic growth, saturation sets in, primarily as a consequence of the separation of those who are susceptible from unidentified infectious individuals, leading to a lower peak that occurs earlier than the peak for an unconstrained epidemic. The curve is flattened. This contrasts with initial exponential growth of confirmed cases for an outbreak with little containment, which is expected from the SIR model and observed for some epidemics such as the Ebola outbreak in West Africa.

Remarkably, using the parsimonious SIR-X model Maier and Brockmann are able to reproduce quantitative growth behavior observed in data from the COVID-19 epidemic in nine Chinese provinces including Hubei Province, the epicenter of the outbreak. The parameter choices that best fit the data are a reproduction number $R_0=6.2$ for an unconstrained epidemic and a mean duration of infectiousness of ${\beta }^{-1}=8d$. The model’s growth curves follow the observed scaling of data for the provinces with exponents $\mu \approx 2$. However, wide variations in parameter choices produce similar scaling and thus the model doesn’t permit inference of specific parameter values. More sophisticated epidemiological models coupled with serological studies are needed to yield reliable estimates for epidemiological parameters.

Nonetheless, the SIR-X model’s mathematical form for the growth of confirmed cases implies that the observed subexponential growth is a result of basic epidemiological processes, caused by a balance between transmission and containment. It offers a guide to judge the expected effectiveness of various containment measures from voluntary stay-at-home advisories up to mandatory curfews and hard lockdowns, as was used in China and some European nations in response to COVID-19. Containment efforts can be evaluated by whether the resulting growth of cases exhibits power law scaling and, if so, by the size of the exponent.

Open questions

The epidemiological models discussed above offer important lessons that can be used to inform public health policy on a basic level. The SIR and SEIR models show that doing nothing allows for explosive growth with potentially devastating consequences. The SIR-X model shows that containment measures and quarantine procedures produce a predictable flattening of the curve with algebraic instead of exponential growth.

But these simple models leave many questions unanswered for a real-world disease such as COVID-19 that need to be addressed by more sophisticated models. How many cases of infection go unidentified through lack of testing, particularly of asymptomatic individuals who are nonetheless infectious? Does the transmission rate for those with documented infections differ from the rate for those with unidentified infections or for symptomatic versus asymptomatic infectious individuals? What is the average latency period and average duration of infection and by how much do these vary? Is there a seasonality effect on the reproduction number? Is immunity conferred through infection and how variable and long lasting is the strength of acquired immunity? How does susceptibility to and severity of infection vary based on medical history, geographical location, socioeconomic status, and genetics? How do real populations of susceptible and infectious individuals interact to make transmission more or less likely? Determining which epidemiological models best answer these and other questions when compared to data from serological studies will better prepare us to respond to the next potential pandemic.

Science provides the exit strategy from epidemics. Testing, therapeutics, vaccine development, and, importantly, epidemiological modeling are essential. Hopefully, one lesson from COVID-19 is the unequivocal need for modeling to inform public health policy.

References

1. S. H. Strogatz, Nonlinear Dynamics and Chaos. Westview Press, 2nd Ed. (2015).
2. https://en.wikipedia.org/wiki/Compartmental_models_in_epidemiology.
3. https://www.aps.org/publications/apsnews/202005/challenge.cfm.
4. B. F. Maier and D. Brockmann, Science 10.1126/science.abb4547 (2020).

---