An Analytically Tractable Framework for Multi-Strain Epidemics: Resolving Algebraic Complexity to Map Oscillatory Dynamics

Abstract: Multi-strain epidemiological systems frequently exhibit self-sustained oscillations, yet severe algebraic complexity has long obstructed a complete analytical characterization of these dynamics. Seeking to bypass these barriers, we have identified a broad, analytically tractable class of two-strain models featuring asymmetric cross-immunity that precludes secondary infections for a single strain. This targeted structural simplification allows us to derive explicit expressions for coexistence equilibria and map their stability. We show that the relative transmission advantage of secondary infections fundamentally governs the onset of robust limit-cycle oscillations, extending previous narrow-borderline results. Crucially, numerical simulations reveal a rich macroscopic landscape where small-amplitude local oscillations coexist with large-amplitude, recurrent outbreak cycles. The structural robustness of these phenomena is confirmed by incorporating effects such as waning immunity and isolation and extending the framework beyond the analytically tractable class. These findings provide a clear mechanistic explanation for multi-strain epidemic cycles and offer a highly tractable baseline for future theoretical developments.

• The paper introduces an analytically tractable asymmetric two-strain SIR framework by setting secondary infection with strain 1 to zero, enabling explicit coexistence equilibria and stability analysis.
• The paper identifies the transmission-advantage parameter q₁₂ᶜ as the key instability threshold, showing that secondary infections with q₁₂ᶜ ≥ 1 can destabilize coexistence and generate oscillations even between symmetric strains.
• Numerical simulations reveal coexisting small- and large-amplitude attractors, with amplitude differences reaching roughly 200-fold in the base model, while waning immunity and isolation reduce but do not eliminate oscillatory behavior.

Motivation and problem statement

Two-strain epidemic models with partial cross-immunity have long presented a puzzle: while competitive exclusion arguments predict dominance of the fitter strain, both numerical studies and empirical observations report persistent coexistence and self-sustained oscillations. In the foundational two-strain SIR framework with lifelong partial cross-immunity, sustained oscillations were shown to occur only under highly restrictive conditions [chung2016dynamics], and a recent analysis identified an oscillatory regime confined to an extremely narrow parametric region [gavish2024newoscillatoryregimetwostrain]. The paper under review addresses two questions raised by that prior work: whether the newly identified oscillations are a mathematical artifact of a specific parameter set, and what biological mechanism drives them.

The central technical obstacle is algebraic: in general asymmetric two-strain models, the coexistence equilibrium is defined only implicitly as the solution of a nonlinear algebraic system, rendering analytical stability criteria inaccessible. The author's strategy is to identify a broad subfamily of models in which this equilibrium admits explicit closed-form expressions, without sacrificing the capacity to sustain oscillatory dynamics.

The analytically tractable family

The base model is the standard eight-compartment two-strain SIR system with partial cross-immunity, extended to allow distinct recovery rates for secondary infections ($\gamma_{ji} \neq \gamma_i$) and distinct relative infectivity ($\eta_i$). The tractable family is obtained by setting $\sigma_{21}=0$, i.e., precluding secondary infections with strain 1. This single structural condition has three consequences that satisfy the stated selection criteria:

• Asymmetry: the model is inherently asymmetric, matching documented strongly asymmetric cross-immunity such as the competitive exclusion between Bordetella pertussis and B. parapertussis [wolfe2007antigen].
• Explicit coexistence equilibrium: since $I_{21}\equiv 0$, the equation for $I_1$ forces $S = 1/\mathcal{R}_1$ at any coexistence state; substituting this removes nonlinear terms and reduces the remaining algebraic system to a linear one.
• Generality within the family: the same reduction applies to substantially enriched models (Section on extensions below).

The author is explicit that $\sigma_{21}=0$ is primarily an analytical device rather than a claim of structural novelty, noting that general two-strain models exhibit Hopf bifurcations and chaos even without forcing.

Explicit equilibria and invasion structure

The paper derives the disease-free equilibrium, two single-strain endemic states $\phi^{EE,i}$ (existing when $\mathcal{R}_i>1$), and a unique coexistence steady state $\phi^{CE}$ in the mutual-invasion region $\eta_i$0. All compartment sizes are given explicitly; for example,

$\eta_i$1

with $\eta_i$2 and $\eta_i$3 obtained from a linear solve. Single-strain endemic states are stable if and only if the competing strain's invasion number is below one, proved via Routh–Hurwitz factorization of the characteristic polynomial.

Stability analysis and the transmission advantage parameter

Direct Routh–Hurwitz analysis of $\eta_i$4 is intractable because the characteristic polynomial is quintic. The author instead works in the singular perturbation regime $\eta_i$5, appropriate for acute infections where infectious periods are negligible relative to host lifespan. Four eigenvalues scale as $\eta_i$6, where the $\eta_i$7 satisfy a biquadratic whose discriminant takes the form

$\eta_i$8

with $\eta_i$9 in $\sigma_{21}=0$0. The pivotal new quantity is

$\sigma_{21}=0$1

which aggregates the relative transmission advantage of secondary infections over primary infections—combining enhanced susceptibility or reduced cross-immunity ($\sigma_{21}=0$2), heightened infectivity ($\sigma_{21}=0$3), and longer infectious period ($\sigma_{21}=0$4). This is the paper's principal conceptual contribution.

The resulting stability picture is threefold:

Regime	Discriminant	Outcome
$\sigma_{21}=0$5	$\sigma_{21}=0$6	Eigenvalues purely imaginary at leading order; numerically stable
$\sigma_{21}=0$7	$\sigma_{21}=0$8 near curve $\sigma_{21}=0$9 (where $I_{21}\equiv 0$0)	Coexistence equilibrium linearly unstable
$I_{21}\equiv 0$1	$I_{21}\equiv 0$2 only on $I_{21}\equiv 0$3	Higher-order expansion shows instability along $I_{21}\equiv 0$4

For the borderline case, the author carries the expansion to $I_{21}\equiv 0$5 and shows that $I_{21}\equiv 0$6, so two eigenvalues acquire positive real parts. A key implication follows immediately: instability does not require weak cross-immunity ($I_{21}\equiv 0$7), contrary to what the earlier narrow-borderline result implied. Rather, oscillations emerge whenever secondary infections confer an overall transmission advantage, $I_{21}\equiv 0$8. In retrospect, the earlier study evaluated precisely the borderline case $I_{21}\equiv 0$9 with $I_1$0, explaining its parametric narrowness; the present analysis generalizes it to a broad, biologically plausible regime.

Notably, oscillations can arise even between completely symmetric strains, since $I_1$1 depends on secondary-infection traits rather than inter-strain asymmetry—a point confirmed numerically via antibody-dependent enhancement (ADE) scenarios.

Numerical study of the base model

Bifurcation diagrams over $I_1$2 with realistic parameters ($I_1$3 per week, i.e., ~76-year lifespan) confirm the theory: for $I_1$4 a substantial instability region appears containing $I_1$5 in its interior; for $I_1$6 the coexistence equilibrium is stable wherever it exists; for $I_1$7 instability is confined to a narrow band around $I_1$8. Appendix results with partial cross-immunity ($I_1$9, $S = 1/\mathcal{R}_1$0) reproduce the qualitative behavior, supporting the claim that $S = 1/\mathcal{R}_1$1 itself—not its constituent parameters—governs stability.

Coexistence of two oscillatory mechanisms

Beyond the local Hopf-like picture, the numerical study uncovers a second, qualitatively distinct oscillatory mechanism. Near the boundary of the instability region, solutions exhibit small-amplitude limit cycles around $S = 1/\mathcal{R}_1$2, consistent with prior work. Deep inside the instability region, however, trajectories display large-amplitude outbreak cycles: short intense surges of strain 1 followed by strain 2 surges separated by low-prevalence intervals.

Most strikingly, at a fixed parameter point ($S = 1/\mathcal{R}_1$3, $S = 1/\mathcal{R}_1$4), the system exhibits bistability: initialization near $S = 1/\mathcal{R}_1$5 converges to a low-amplitude limit cycle, whereas seeding a trace invader population ($S = 1/\mathcal{R}_1$6) near the single-strain equilibrium produces a large-amplitude macroscopic attractor. Numerical continuation confirms this bistability persists under small parameter perturbations. One-parameter bifurcation scans show regular periodic large-amplitude oscillations in sub-intervals of $S = 1/\mathcal{R}_1$7 (e.g., $S = 1/\mathcal{R}_1$8 and $S = 1/\mathcal{R}_1$9) and apparently aperiodic patterns elsewhere. The amplitude ratio between large- and small-amplitude attractors reaches roughly 200 in the base model. The author attributes these global transitions to Shilnikov-like orbits connected to a single-strain equilibrium, deferring full characterization to a companion paper—an honest acknowledgment that the global bifurcation structure remains open.

Structural robustness: waning immunity and isolation

To demonstrate the tractability payoff, the author extends the model with waning immunity (rates $\sigma_{21}=0$0) and isolation compartments ($\sigma_{21}=0$1), yielding a nine-dimensional system. Despite the added complexity, Proposition 5 delivers fully explicit expressions for all coexistence equilibrium components—the direct benefit of the $\sigma_{21}=0$2 assumption. This matters practically: stability can be evaluated exactly from eigenvalues of the full Jacobian rather than heuristically.

Numerically, with biologically plausible parameters (immunity waning over ~two years, ~one-week recovery/isolation periods), the extended model preserves the qualitative structure: an oscillatory region appears for $\sigma_{21}=0$3 using the natural analogue $\sigma_{21}=0$4, and vanishes at $\sigma_{21}=0$5. Two caveats are stated plainly. First, the empirical threshold lies strictly above unity ($\sigma_{21}=0$6): fast susceptible replenishment through waning immunity shifts the true boundary away from the base-model value, and no closed-form transcendental threshold exists in this higher-dimensional setting. Second, the extended model's dynamics differ quantitatively: the amplitude ratio between large- and small-amplitude regimes drops to roughly 3, low-activity periods shorten, and oscillations remain entirely regular (no irregular patterns). The author plausibly attributes this damping to faster susceptible replenishment.

Beyond the tractable family

Initial numerical exploration of the general model with $\sigma_{21}=0$7 (coexistence equilibrium reduced exactly to a cubic in $\sigma_{21}=0$8) tests breadth in three scenarios:

• Symmetric partial cross-immunity ($\sigma_{21}=0$9, influenza-like): $\phi^{EE,i}$0 is stable wherever it exists, consistent with $\phi^{EE,i}$1 on both sides.
• One-directional advantage ($\phi^{EE,i}$2): an instability region emerges with the familiar mix of low-amplitude and regular/irregular high-amplitude oscillations.
• ADE-driven symmetric case (following Ferguson et al.): robust oscillations appear even though strains are identical, confirming that secondary-infection advantage—not asymmetry—is the driver. Here the oscillatory region spans part of the line $\phi^{EE,i}$3, though it lies strictly along the boundary of the existence region of $\phi^{EE,i}$4, unlike all other cases considered.

These results support structural robustness of the instability mechanism beyond the analytically tractable class, though the exploration is admittedly not exhaustive.

Limitations and open questions

Several limitations are conceded explicitly. The analytical instability results hold for sufficiently small $\phi^{EE,i}$5; the claim of robustness at realistic $\phi^{EE,i}$6 rests on numerics. For $\phi^{EE,i}$7, stability of $\phi^{EE,i}$8 is supported only by numerical evidence, as next-order terms were not derived. In the extended model, the threshold quantity $\phi^{EE,i}$9 is a heuristic scaling tool, not an exact bifurcation criterion, and the shift of the empirical threshold above unity is unexplained analytically. The global mechanism behind large-amplitude outbreak cycles—including bistability boundaries and the hypothesized Shilnikov connection—is observed but not characterized, motivating the companion paper. Application to dengue would require temporary cross-immunity structures absent from the present model, and the isolated treatment of ADE and waning immunity leaves composite-model behavior unexplored.

Conclusion

This paper converts a previously narrow, seemingly artifactual oscillatory regime into a broadly characterized phenomenon by identifying $\mathcal{R}_i>1$0—the aggregate relative transmission advantage of secondary infections—as the governing parameter for destabilization of the coexistence equilibrium. The $\mathcal{R}_i>1$1 structural restriction resolves the algebraic obstruction to explicit coexistence equilibria while remaining biologically grounded, and the resulting framework accommodates waning immunity, isolation, ADE, and bidirectional secondary infections. The discovery that small-amplitude Hopf-like cycles coexist with large-amplitude outbreak attractors adds a global-dynamical layer that the local analysis cannot capture, establishing a well-defined open problem for subsequent work.