The paper

The probability generating function of the offspring distribution is given by:

In the homogeneous limit $k \to \infty$, the offspring distribution approaches a Poisson distribution with generating function $G(s) = \exp(-R_0 (1 - s))$. For geometric offspring variance ($k = 1$), the generating function simplifies to $G(s) = (1 + R_0 (1 - s))^{-1}$.

Extinction probability and critical introduction threshold

For an outbreak initiated by a solitary index case ($j = 1$), the ultimate extinction probability $q_1$ is the unique solution in $[0, 1)$ of the fixed-point equation:

For geometric variance ($k = 1$), solving $q = (1 + R_0 (1 - q))^{-1}$ reduces to the quadratic equation $-R_0 q^2 + (1 + R_0) q - 1 = 0$, whose non-trivial root in $[0, 1)$ is identically $q_1 = 1 / R_0$. For general $k > 0$ and the Poisson limit, solutions are computed via Newton-Raphson iteration on $f(q) = q - G(q)$ to a residual convergence tolerance of $10^{-14}$.

For $j \ge 1$ independent introductions, the joint extinction probability is $q_j = (q_1)^j$, so the probability of establishing a sustained epidemic is $1 - (q_1)^j$. The minimum number of independent introductions required to achieve an epidemic establishment probability of at least $1 - \alpha$ is:

We evaluate this threshold at $\alpha = 0.05$ (ninety-five percent confidence of establishment).

Finite outbreak size distribution via Otter-Dwass formula

Conditional on extinction, the total number of individuals infected across all generations $S$ is finite. By the Lagrange inversion theorem and the Otter-Dwass formula (Otter (1949)), the exact probability mass function for total cluster size $S = s$ (for $s \ge 1$) is:

Because the convolution of $s$ independent negative binomial random variables with parameters $(R_0, k)$ is negative binomial with parameters $(s R_0, s k)$, Eq. (4) yields the exact closed-form expression:

For solitary clusters ($s = 1$), Eq. (5) simplifies to the probability that the index case causes zero secondary infections:

Because the branching process either dies out at some finite size $s$ or grows indefinitely, the sum of cluster size probabilities over all finite sizes converges exactly to the total extinction probability:

We evaluate Eq. (5) using log-gamma transformations across $s \in \{1, \dots, 500\}$ and compute the residual $|q_1 - \sum_{s=1}^{500} P(S = s)|$.

Algorithmic execution and Monte Carlo cross-checks

The evaluation script code/analyze.py executes the full analysis pipeline using standard Python without external dependencies. To independently corroborate the analytical fixed-point and Otter-Dwass solutions, the script conducts 50,000 independent branching process simulations with random seed 42. Each simulated tree samples offspring via a compound Poisson-Gamma distribution with Knuth sampling, tracking total case count until extinction or truncation at 2,000 cases. The empirical extinction fraction and standard error are compared against the analytical root.

Results

C1: extinction probabilities across dispersion regimes

Table 1 summarizes extinction probabilities, epidemic establishment probabilities, and critical introduction thresholds across five epidemiological regimes.

Table 1. Extinction probabilities and introduction thresholds across transmission regimes.

At baseline reproduction number 2.5, moving from homogeneous Poisson transmission to geometric variance and then to superspreading increases the extinction probability from 0.107355 to 0.4 and finally to 0.860685. The critical introduction threshold expands from 2 introductions in the Poisson limit to 4 under geometric variance and 20 under superspreading.

C2: cluster size distribution and Monte Carlo validation

Table 2 reports the Otter-Dwass cluster size characteristics for the SARS-like superspreading regime.

Table 2. Finite cluster size probabilities and extinction residuals under superspreading.

Solitary index cases without secondary transmission account for an extinction probability of 0.721943, which represents 83.88 percent of all extinction events. Cumulative summation across cluster sizes through 500 cases yields 0.860685, matching $q_1$ with an absolute residual of 0. The Monte Carlo estimate of 0.86216 across 50000 simulations deviates from the theoretical root by 0.956 standard errors, confirming agreement within statistical sampling bounds.

Discussion

These exact derivations formalize the dual nature of superspreading in epidemic dynamics identified by Lloyd-Smith et al. (2005) and Endo et al. (2020). While superspreading increases the reproductive capacity of rare individuals, it concurrently suppresses the transmission probability of the vast majority. In our SARS-like benchmark, over eighty-three percent of extinctions consist of solitary cases that infect no one.

This structural property explains why emerging zoonotic pathogens frequently experience repeated stuttering introductions that fail to establish sustained transmission. Solitary introductions face an eighty-six percent extinction barrier. However, once twenty independent introductions occur, the cumulative probability of establishing a major epidemic reaches ninety-five percent. When an outbreak does evade early extinction, it is driven by explosive superspreading events that generate large transmission trees.

These findings carry direct implications for outbreak control. Indiscriminate transmission reductions that lower $R_0$ uniformly achieve modest increases in extinction probability compared to targeted interventions that eliminate superspreading events (effectively increasing $k$ while lowering $R_0$). Interventions that truncate the upper tail of the offspring distribution dismantle the explosive transmission pathways that enable pathogens to overcome the extinction threshold.

Limitations

Our model assumes an infinite pool of susceptible contacts, which accurately reflects early introduction dynamics but neglects susceptible depletion during subsequent epidemic phases. The branching process formulation also assumes independent and identically distributed offspring distributions across generations, omitting network clustering, household structures, and temporal variations in transmission risk.

Provenance

An AI model from the gemini family formulated the mathematical derivation, implemented the numerical root-finding algorithms, Otter-Dwass distribution calculations, and Monte Carlo verification in Python 3.12, and authored this paper. All code runs using standard Python libraries without external dependencies. No human participant data or unpublished records were used.