systems.builtin.stochastic.continuous.StochasticSIR
systems.builtin.stochastic.continuous.StochasticSIR(*args, **kwargs)Stochastic SIR epidemic model with demographic noise.
Models disease spread in a finite population with inherent randomness from discrete infection and recovery events. Essential for understanding outbreak extinction, variability in epidemic outcomes, and uncertainty in disease forecasting.
Stochastic Differential Equations
Diffusion approximation of discrete SIR:
dS = -β·S·I/N·dt + σ_S·dW_S
dI = (β·S·I/N - γ·I)·dt + σ_I·dW_I
dR = γ·I·dt + σ_R·dW_R
where: - S(t): Number susceptible (can be infected) - I(t): Number infected (infectious) - R(t): Number recovered (immune) - N = S + I + R: Total population (constant) - β: Transmission rate [1/time] - γ: Recovery rate [1/time] - σ_S, σ_I, σ_R: Noise intensities - W_S, W_I, W_R: Wiener processes
Physical Meaning:
Transmission: β·S·I/N - Rate of new infections - Mass action: Proportional to S·I - β: Contacts per time × transmission probability
Recovery: γ·I - Rate of recoveries - First-order: Proportional to I - γ = 1/(infectious period)
Demographic Noise:
True stochastic SIR has √ diffusion: Diffusion ~ √(rate)
From Poisson statistics of discrete events.
This implementation: Simplified additive noise for tractability.
Physical Interpretation
Susceptibles S:
Decrease only (monotonic): - Start: S(0) ≈ N (nearly all susceptible) - End: S(∞) > 0 (some escape infection) - Never increases (no loss of immunity modeled)
Infected I:
Non-monotonic (rise then fall): - Start: I(0) = small (index cases) - Peak: I_max at t_peak (outbreak peak) - End: I(∞) = 0 (disease dies out)
Recovered R:
Increase only (monotonic): - Start: R(0) = 0 (no immunity initially) - End: R(∞) = final outbreak size - Measure of epidemic impact
Conservation: S + I + R = N (total population constant)
Key Features
Nonlinearity: S·I term creates threshold behavior and epidemic curve.
Positivity: Must have S, I, R ≥ 0 (counts of people).
Conservation: S + I + R = N always (no births/deaths).
Extinction: I → 0 eventually (disease dies out). Time random, final size random.
Threshold (R₀): Probabilistic in stochastic model.
Finite Time: Epidemic is transient (not steady state).
Mathematical Properties
Basic Reproduction Number: R₀ = β/γ
Critical threshold (deterministic): - R₀ < 1: Dies out - R₀ > 1: Epidemic
Stochastic Threshold: Even R₀ > 1: Can die out with probability (1/R₀)^{I₀}
Final Size Relation:
Deterministic: R_∞ satisfies R_∞ = N - S₀·exp(-R₀·R_∞/N)
Stochastic: Distribution around this value.
Epidemic Peak:
Deterministic: I_max ≈ I₀ + S₀ - N/R₀ - (N/R₀)·ln(S₀·R₀/N) t_peak ≈ (1/γ)·ln(R₀·S₀/N)
Stochastic: Random variables (compute via simulation).
Physical Interpretation
Transmission Rate β: - Units: [1/time] - β = contact rate × transmission probability - Typical: 0.2-2.0 per day
Examples: - Influenza: β ≈ 0.5 per day - Measles: β ≈ 1.5 per day - COVID-19: β ≈ 0.3-0.6 per day (varies)
Recovery Rate γ: - Units: [1/time] - γ = 1/(infectious period) - Typical: 0.1-1.0 per day
Examples: - Influenza: γ ≈ 0.5 per day (2 days infectious) - COVID-19: γ ≈ 0.1 per day (10 days) - Measles: γ ≈ 0.1 per day (10 days)
Basic Reproduction Number: R₀ = β/γ
Examples: - Influenza: R₀ ≈ 1-2 - COVID-19: R₀ ≈ 2-5 (variant dependent) - Measles: R₀ ≈ 12-18
Noise Intensity:
From demographic stochasticity (Poisson): σ ≈ √(rate/N)
Relative noise: 1/√N
State Space
State: X = [S, I, R] ∈ ℝ₊³ - S, I, R ≥ 0 (non-negative counts) - S + I + R = N (conservation) - Bounded: 0 ≤ S, I, R ≤ N
Control: u (optional, interventions) - Reduce β (social distancing) - Increase recovery (treatment) - Vaccination (move S → R)
Noise: w = [w_S, w_I, w_R] ∈ ℝ³ - Demographic stochasticity - Should be correlated (conservation) - This implementation: Simplified independent
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| beta | float | Transmission rate [1/day] - Typical: 0.2-2.0 - Disease and behavior dependent | 0.5 |
| gamma | float | Recovery rate [1/day] - 1/γ = infectious period - Typical: 0.1-1.0 | 0.1 |
| N | float | Total population - Larger N: More deterministic - Smaller N: More stochastic | 1000.0 |
| sigma_S | float | Susceptible noise intensity [1/√day] - From √(rate/N) for demographic noise - Typical: √(β·S·I/N²) | 0.1 |
| sigma_I | float | Infected noise intensity [1/√day] | 0.1 |
| sigma_R | float | Recovered noise intensity [1/√day] | 0.1 |
Stochastic Properties
- System Type: NONLINEAR
- Noise Type: ADDITIVE (simplified)
- SDE Type: Itô
- Noise Dimension: nw = 3
- Stationary: No (epidemic is transient)
- Positive: Should be (may need projection)
- Conserved: S + I + R = N (approximately)
Applications
1. Epidemiology: - COVID-19, influenza, measles modeling - Outbreak prediction with uncertainty - Intervention timing and effectiveness
2. Parameter Estimation: - Estimate β, γ from outbreak data - Bayesian inference with uncertainty - Real-time estimation (particle filter)
3. Public Health Policy: - Vaccination strategies - Social distancing timing - Resource allocation (ICU beds)
4. Extinction Analysis: - Probability small outbreak dies out - Early intervention effectiveness - Import risk assessment
5. Rare Events: - Superspreading events - Large outbreaks in small populations - Timing of peak (ICU planning)
Numerical Integration
Recommended: - Euler-Maruyama: dt = 0.01-0.1 days - Project to non-negative: max(X, 0) - Check conservation: S+I+R ≈ N
Event Detection: - Extinction: I < 0.5 (declare extinct) - Peak: Max I(t) - Duration: Time I > threshold
Monte Carlo Guidelines
Ensemble Analysis: - N_runs = 100-1,000 - Compute: Extinction probability, mean final size - Histogram: Final R (bimodal if near threshold)
Comparison with Deterministic
Deterministic: - Smooth epidemic curve - Single final size - R₀ sharp threshold
Stochastic: - Variable epidemic curves - Distribution of final sizes - R₀ probabilistic threshold - Extinction possible
Limitations
- Additive noise (not √ diffusion)
- Independent noise (should be correlated)
- Homogeneous mixing (no network structure)
- Constant β, γ (no seasonality)
- Closed population (no births/deaths)
Extensions
- SEIR: Add Exposed class
- Age structure: Multiple age groups
- Spatial: Geographic spread
- Network: Contact structure
- Time-varying: Seasonal β(t)
- Vaccination: Control via S → R
See Also
CoxIngersollRoss : Similar √ diffusion structure
Methods
| Name | Description |
|---|---|
| check_conservation | Check conservation constraint S + I + R = N. |
| define_system | Define stochastic SIR epidemic dynamics. |
| estimate_extinction_probability | Estimate probability of stochastic extinction (approximation). |
| get_basic_reproduction_number | Get basic reproduction number R₀ = β/γ. |
| get_herd_immunity_threshold | Get herd immunity threshold H = 1 - 1/R₀. |
check_conservation
systems.builtin.stochastic.continuous.StochasticSIR.check_conservation(x)Check conservation constraint S + I + R = N.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| x | np.ndarray | State [S, I, R] | required |
Returns
| Name | Type | Description |
|---|---|---|
| float | Absolute error |S+I+R - N| |
Examples
>>> sir = StochasticSIR(N=1000)
>>> x = np.array([500, 200, 300])
>>> error = sir.check_conservation(x)
>>> print(f"Conservation error: {error:.2f}")define_system
systems.builtin.stochastic.continuous.StochasticSIR.define_system(
beta=0.5,
gamma=0.1,
N=1000.0,
sigma_S=0.1,
sigma_I=0.1,
sigma_R=0.1,
)Define stochastic SIR epidemic dynamics.
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| beta | float | Transmission rate [1/day] - β = contact rate × transmission probability - Typical: 0.2-2.0 per day - Higher → faster spread | 0.5 |
| gamma | float | Recovery rate [1/day] - γ = 1/(infectious period) - Typical: 0.1-1.0 per day - Higher → faster recovery | 0.1 |
| N | float | Total population (constant) - Larger N: More deterministic - Smaller N: More stochastic - Typical: 100-1,000,000 | 1000.0 |
| sigma_S | float | Susceptible noise intensity [1/√day] - From demographic noise: ~ √(β·S·I/N²) - Typical: 0.01-1.0 | 0.1 |
| sigma_I | float | Infected noise intensity [1/√day] - From demographic noise: ~ √(β·S·I/N² + γ·I/N) - Most critical (affects outbreak dynamics) | 0.1 |
| sigma_R | float | Recovered noise intensity [1/√day] - From demographic noise: ~ √(γ·I/N) | 0.1 |
Notes
Basic Reproduction Number: R₀ = β/γ
Determines epidemic threshold: - R₀ < 1: Dies out (subcritical) - R₀ = 1: Critical (boundary) - R₀ > 1: Epidemic (supercritical)
Stochastic Extinction:
Even if R₀ > 1, can randomly die out: P(extinction | I₀) ≈ (1/R₀)^{I₀}
Example: R₀ = 2, I₀ = 1 → P_ext = 50%
Noise Scaling:
Physical demographic noise: σ_physical ~ √(rate/N)
For N = 1000, typical rates ~ 0.1-1: σ ~ 0.01-0.1
Larger population → smaller relative noise.
Infectious Period: T_inf = 1/γ
Examples: - Influenza: 2-3 days (γ ≈ 0.4) - COVID-19: 7-14 days (γ ≈ 0.1) - Measles: 7-10 days (γ ≈ 0.12)
Herd Immunity: H = 1 - 1/R₀
Fraction immune needed to prevent epidemic.
Validation:
Check conservation: S + I + R should equal N. If drifts: Renormalize periodically.
estimate_extinction_probability
systems.builtin.stochastic.continuous.StochasticSIR.estimate_extinction_probability(
I0,
)Estimate probability of stochastic extinction (approximation).
For R₀ > 1: P_ext ≈ (1/R₀)^{I₀}
Parameters
| Name | Type | Description | Default |
|---|---|---|---|
| I0 | float | Initial number infected | required |
Returns
| Name | Type | Description |
|---|---|---|
| float | Extinction probability |
Examples
>>> sir = StochasticSIR(beta=0.6, gamma=0.2) # R₀=3
>>> P_ext_1 = sir.estimate_extinction_probability(I0=1)
>>> P_ext_10 = sir.estimate_extinction_probability(I0=10)
>>> print(f"1 infected: P_ext = {P_ext_1:.2%}")
>>> print(f"10 infected: P_ext = {P_ext_10:.2%}")get_basic_reproduction_number
systems.builtin.stochastic.continuous.StochasticSIR.get_basic_reproduction_number(
)Get basic reproduction number R₀ = β/γ.
Returns
| Name | Type | Description |
|---|---|---|
| float | R₀ |
Examples
>>> sir = StochasticSIR(beta=0.5, gamma=0.1)
>>> R0 = sir.get_basic_reproduction_number()
>>> print(f"R₀ = {R0:.2f}")get_herd_immunity_threshold
systems.builtin.stochastic.continuous.StochasticSIR.get_herd_immunity_threshold(
)Get herd immunity threshold H = 1 - 1/R₀.
Returns
| Name | Type | Description |
|---|---|---|
| float | Herd immunity threshold (fraction) |
Examples
>>> sir = StochasticSIR(beta=0.5, gamma=0.1)
>>> H = sir.get_herd_immunity_threshold()
>>> print(f"Herd immunity: {H:.1%}")