Hospital-capacity exceedance probabilities

Moment-based prevalence summaries can inform capacity risk only after a clear model connects infection to healthcare demand.

From infectious prevalence to critical-care demand

Let \(C(t)\) be the number of infectious people requiring critical care. Suppose that, conditional on \(I(t)\),

\[C(t)\mid I(t)\sim\operatorname{Binomial}(I(t),p_c),\]

where \(p_c\) is the critical-care probability. This adds patient-level demand variation to epidemic variation.

Transform the moments

\[\mathbb E[C]=p_c\mathbb E[I],\]

\[ \operatorname{Var}(C) = p_c(1-p_c)\mathbb E[I] + p_c^2\operatorname{Var}(I). \]

The first term is conditional binomial variation; the second transfers epidemic prevalence variation into demand.

Define the capacity event

If capacity is \(L\) beds, exceedance at time \(t\) is

\[\{C(t)>L\}.\]

This is a time-specific probability. The probability of ever exceeding capacity during a period is a pathwise event and cannot be calculated from separate marginal moments alone.

Interactive Python laboratory

Interactive PythonHospital-capacity probability

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
from scipy.stats import norm

N, I0 = 100, 5
beta, gamma = 0.25, 0.15
critical_probability = 0.20
capacity = 6
times = np.linspace(0, 30, 61)

def rhs(t, y):
    m1, m2 = y
    m3 = 3 * m1 * m2 - 2 * m1 ** 3
    return [
        (beta - gamma) * m1 - beta * m2 / N,
        (2 * (beta - gamma) - beta / N) * m2
        - 2 * beta * m3 / N
        + (beta + gamma) * m1
    ]

sol = solve_ivp(
    rhs, (0, 30), [I0, I0 ** 2],
    t_eval=times, rtol=1e-9, atol=1e-11
)
mean_I, second_I = sol.y
variance_I = np.maximum(second_I - mean_I ** 2, 0)

mean_C = critical_probability * mean_I
variance_C = (
    critical_probability * (1 - critical_probability) * mean_I
    + critical_probability ** 2 * variance_I
)
sd_C = np.sqrt(variance_C)

# C > capacity means C >= capacity + 1.
# Continuity boundary is capacity + 0.5.
normal_exceedance = 1 - norm.cdf(
    (capacity + 0.5 - mean_C) / sd_C
)

# CTMC plus conditional binomial benchmark
rng = np.random.default_rng(310)
M = 8000

def ctmc_path():
    I, time = I0, 0.0
    event_times, values = [0.0], [I]
    while I > 0 and time < times[-1]:
        b = beta * (N - I) * I / N
        d = gamma * I
        wait = rng.exponential(1 / (b + d))
        if time + wait > times[-1]:
            break
        time += wait
        I += 1 if rng.random() < b / (b + d) else -1
        event_times.append(time)
        values.append(I)
    pos = np.searchsorted(event_times, times, side="right") - 1
    return np.array(values)[pos]

I_paths = np.array([ctmc_path() for _ in range(M)])
C_samples = rng.binomial(
    I_paths.astype(int), critical_probability
)
empirical_exceedance = (C_samples > capacity).mean(axis=0)

selected = [0, 10, 20, 30]
idx = [int(t * 2) for t in selected]
table = pd.DataFrame({
    "time": selected,
    "mean critical demand": mean_C[idx],
    "SD critical demand": sd_C[idx],
    "normal P(C>6)": normal_exceedance[idx],
    "simulation P(C>6)": empirical_exceedance[idx]
})
print(table.round(4).to_string(index=False))

plt.figure(figsize=(9.5, 5.2))
plt.plot(times, normal_exceedance,
         color="#b42318", linewidth=2,
         label="Moment + normal approximation")
plt.plot(times, empirical_exceedance, "--",
         color="#174c34",
         label="CTMC + binomial simulation")
plt.xlabel("Time")
plt.ylabel("Probability demand exceeds 6 beds")
plt.title("Time-specific hospital-capacity exceedance risk")
plt.ylim(0, 1)
plt.legend()
plt.grid(alpha=0.25)
plt.tight_layout()

Output

Run the code to see the result.

Assumptions that must be stated

Critical-care probability is constant and identical between individuals.
Demand is contemporaneous with infectious prevalence; no admission delay or length of stay is modelled.
Capacity is fixed.
The normal approximation may be poor for small counts or skewed demand.
Parameter and capacity uncertainty are not included.

Decision interpretation

A high exceedance probability identifies periods of capacity risk under the stated model. It is not a guarantee that capacity will be exceeded and should be accompanied by sensitivity analysis for \(p_c\), capacity and epidemic parameters.

What this lesson adds

You can now transform infectious moments into demand moments, include conditional patient-level variation, define a precise capacity event and validate a moment-based exceedance probability against simulation.