Case Study 1 — When a Differential Equation Locked Down the World
Field: Epidemiology and public-health policy Calculus used: the SIR system and its early-outbreak linearization (§19.9), $R_0 = \beta N/\gamma$ (§19.9), numerical integration as the successor to Euler's method (§19.6)
In the second week of March 2020 a research team at Imperial College London published a report that, within nine days, helped persuade two governments to confine tens of millions of people to their homes. The report contained no new virus, no new drug, and no new data that the world did not already have. What it contained was a differential equation — the SIR model of §19.9, dressed in extra compartments and calibrated to the early case counts from Wuhan and Lombardy. This is the story of how a rate law became policy, told as a single thread you can follow with the calculus of this chapter.
The question that only a rate law can answer
By early March the raw case counts were doubling every two to three days in several countries. Doubling is the fingerprint of exponential growth, and exponential growth is what you get from the infectious equation of the SIR system at the start of an outbreak. Recall the system from §19.9:
$$\frac{dS}{dt} = -\beta S I, \qquad \frac{dI}{dt} = \beta S I - \gamma I, \qquad \frac{dR}{dt} = \gamma I.$$
At the very beginning almost everyone is susceptible, so $S \approx N$, and the infectious equation linearizes to
$$\frac{dI}{dt} \approx (\beta N - \gamma)\,I,$$
a pure exponential with solution $I(t) \approx I_0\,e^{(\beta N - \gamma)t}$. The sign of the exponent decides everything. Epidemiologists repackage that sign as the dimensionless basic reproduction number
$$R_0 = \frac{\beta N}{\gamma},$$
the average number of secondary infections one case produces in a fully susceptible population. For the ancestral strain of SARS-CoV-2, early estimates put $R_0$ between $2$ and $3$. Take $R_0 = 3$ as the working figure. The crucial point is not that $3 > 1$ — that only tells you the epidemic grows — but how fast, and how far, the growth runs before it stops on its own.
What the model predicted, and why the broad strokes were unavoidable
To turn $R_0 = 3$ into a forecast you must integrate the full nonlinear system, because the linear approximation is valid only while $S \approx N$. Once enough people have been infected, the $\beta S I$ term bends the curve. The Imperial model used a far more elaborate compartment structure than plain SIR — it added an exposed (latent) stage, split the infectious into mild, hospitalized, and intensive-care streams, and stratified by age — but every one of those refinements is a compartment with an inflow rate and an outflow rate, exactly the bookkeeping of §19.9. The whole apparatus, several dozen coupled ODEs, was integrated numerically by the production-grade descendants of the Euler method we hand-coded in §19.6.
Three predictions came out of that integration, and all three follow from the structure of the model rather than from any finely tuned parameter.
First, the curve has a peak, and the peak sits where $S = N/R_0$. This is not an assumption; it is a theorem about the system. The infectious population stops growing exactly when $dI/dt = 0$, and from $dI/dt = (\beta S - \gamma)I$ that happens when $\beta S = \gamma$, i.e.
$$S^{*} = \frac{\gamma}{\beta} = \frac{N}{R_0}.$$
For $R_0 = 3$ the peak arrives when the susceptible pool has fallen to one-third of the population. Before that point infections accelerate; after it, they decline. The entire logic of "flatten the curve" lives in this one equation: every public-health measure that lowers $\beta$ — distancing, masks, closing crowded venues — lowers $R_0 = \beta N/\gamma$, which raises the peak's susceptible threshold $N/R_0$ and so lowers and delays the peak itself.
Second, the epidemic ends with susceptibles to spare. It is tempting to read $R_0 > 1$ as "everyone gets infected," but the SIR model flatly refuses that reading, and the refusal is the model's most counterintuitive gift. Once $S$ drops below $N/R_0$, each case infects fewer than one other person and the chain decays even though susceptibles remain. The leftover susceptible fraction $s_\infty$ satisfies the final-size relation $\ln(s_\infty) = R_0\,(s_\infty - 1)$, which for $R_0 = 3$ gives roughly $s_\infty \approx 0.06$. About $94\%$ are eventually infected if nothing is done — devastating, but tellingly not $100\%$, and the gap between those two numbers is the whole prize for intervening.
Third, the timing and height of the peak set the hospital problem. If the peak infectious fraction exceeds the supply of intensive-care beds, deaths climb not only from the disease but from the unavailability of care. The model's headline scenarios — hundreds of thousands of deaths under no intervention, an order of magnitude fewer under suppression — were qualitative claims of this form. They were robust precisely because they depended on the shape of the SIR solution, not on the third decimal place of any parameter.
From a number to a policy: the herd-immunity inequality
The same $R_0$ that governs whether the epidemic grows also fixes how much immunity stops it. Suppose a fraction $p$ of the population is immune (recovered or vaccinated). Then a typical case meets, on average, only $(1-p)$ as many susceptibles, so the effective reproduction number is $R_0(1-p)$. Transmission can no longer sustain itself once
$$R_0(1-p) < 1 \quad\Longleftrightarrow\quad p > 1 - \frac{1}{R_0}.$$
For $R_0 = 3$ this herd-immunity threshold is $p > 1 - 1/3 \approx 67\%$. For later, more transmissible variants with $R_0$ near $10$–$12$, the threshold climbs above $90\%$ — and because vaccine-derived immunity wanes, that ceiling proved punishing to reach. This single ODE-derived inequality is what turned "model output" into "vaccination target" for entire national programs.
A small worked simulation you can reproduce
Strip the Imperial model back to the bare SIR of §19.9 and you can watch the same logic unfold on a toy population. The code below uses the adaptive Runge–Kutta integrator that succeeds the hand-coded Euler method; the printed outputs are computed by hand from the model's structure, not by running the cell.
# Bare-bones SIR with R0 = 3 on a town of 1,000 people.
import numpy as np
from scipy.integrate import solve_ivp
def sir(t, y, beta, gamma):
S, I, R = y
return [-beta*S*I, beta*S*I - gamma*I, gamma*I]
N, gamma, R0 = 1000, 0.1, 3.0 # 10-day infectious period (1/gamma)
beta = R0 * gamma / N # beta = 3*0.1/1000 = 3e-4
sol = solve_ivp(sir, [0, 160], [N-1, 1, 0], args=(beta, gamma), max_step=1.0)
# Peak occurs when S = N/R0 = 1000/3 ≈ 333; final susceptible s_inf ≈ 0.06*N.
print("Susceptible threshold at peak: N/R0 =", N/R0) # ≈ 333
print("Final susceptible (never infected): ≈", round(0.06*N)) # ≈ 60
The peak lands when the susceptibles fall through $N/R_0 \approx 333$, and the outbreak burns out leaving roughly $60$ people — about $6\%$ — never infected, exactly the final-size prediction. Lower $\beta$ (say, halve it to mimic distancing) and $R_0$ drops to $1.5$: the peak shrinks, arrives later, and the spared fraction grows. That is "flattening the curve," rendered in three lines.
How the forecast aged
The model was qualitatively right about the things that mattered for the March 2020 decision: the exponential early phase, the existence and rough height of a peak, the dominant role of $R_0$, and the order-of-magnitude gap between acting fast and acting slowly. It was quantitatively off in places — it could not anticipate how much voluntary behavior change would bend $\beta$ before any mandate, it did not foresee the variants (each effectively a new disease with its own $\beta$, $\gamma$, and $R_0$), and it underestimated the speed of vaccine development. But those are corrections to a curve whose shape the calculus got right. The lesson of §19.5 applies: even when you cannot pin the numbers exactly, the qualitative behavior of the system — does it grow, where does it peak, where does it settle — is often robust and often enough.
Discussion Questions
- The peak of an outbreak occurs at $S = N/R_0$. Explain, in one sentence each, why a larger $R_0$ pushes the peak earlier and higher. Tie your answer to the linearized growth rate $\beta N - \gamma$.
- The final-size relation predicts that even an uncontrolled $R_0 = 3$ epidemic leaves about $6\%$ of people uninfected. Why does the SIR model forbid a $100\%$ attack rate, and what feature of the dynamics enforces the leftover susceptibles?
- Non-pharmaceutical interventions act by lowering $\beta$. Using $R_0 = \beta N/\gamma$ and the herd-immunity inequality $p > 1 - 1/R_0$, explain why a measure that halves $\beta$ helps on two fronts at once.
- The Imperial model used dozens of compartments; §19.9's SIR uses three. What does a modeller gain by adding an exposed (SEIR) compartment or age structure, and what does she risk by adding too many uncertain parameters?
- For a variant with $R_0 = 10$ and lifelong immunity, compute the herd-immunity threshold. If immunity instead wanes over a year, qualitatively how does the picture change?
A Short Annotated Reading List
- Kermack, W. O., & McKendrick, A. G. (1927). "A contribution to the mathematical theory of epidemics." Proc. R. Soc. A 115, 700–721. The founding paper of the SIR model and the source of the final-size relation used above. Surprisingly readable; the modern $R_0$ falls straight out of its threshold theorem.
- Anderson, R. M., & May, R. M. (1992). Infectious Diseases of Humans. Oxford. The standard reference that turned mathematical epidemiology into a working science; the herd-immunity inequality is developed here in full.
- Ferguson, N. M., et al. (2020). "Impact of non-pharmaceutical interventions..." Imperial College COVID-19 Response Team. The March 2020 report itself — read it to see how the three-compartment skeleton scales into a policy instrument.
- Kucharski, A. (2020). The Rules of Contagion. Basic Books. A lucid popular account of how $R_0$ and the SIR family shape everything from disease to financial panics to viral content; the best non-technical companion to §19.9.
In Chapter 39 you will take the bare SIR model of §19.9, fit $\beta$ and $\gamma$ to a real case-count series, compute $R_0$ from your fit, and simulate an intervention — turning the anchor example you met here into the centerpiece of your modeling portfolio.