Infection, recovery and no-change events

We now connect biological event mechanisms to transition probabilities. This lesson uses a one-event-per-step SIR DTMC: during one sufficiently short step, either one infection occurs, one recovery occurs, or the state does not change.

Biological scenario

A closed population contains 1,000 people and is currently in state

\[(S,I,R)=(990,10,0).\]

We use \(\beta=0.3\) per day, \(\gamma=0.1\) per day and a step of \(\Delta t=0.1\) day. We want the probabilities of infection, recovery and no change during the next step.

The model permits three outcomes per step

One infection

\[(S,I,R)\to(S-1,I+1,R).\]

One recovery

\[(S,I,R)\to(S,I-1,R+1).\]

No change

\[(S,I,R)\to(S,I,R).\]

The outcomes are mutually exclusive in this particular DTMC construction: exactly one is assigned to each step.

Rates and probabilities are different quantities

Event rate

Expected event frequency per unit time at the current state.

× \(\Delta t\) →

Short-step probability

Approximate chance of one event during the next interval.

For a sufficiently short step,

\[\Pr(\text{one event in }\Delta t)\approx\text{event rate}\times\Delta t.\]

A rate can exceed one per day. A probability cannot exceed one. Multiplication by the time length converts the rate into a dimensionless short-interval probability.

Infection probability

The deterministic infection-event rate at state \((S,I,R)\) is

\[a_{\mathrm{inf}}(S,I)=\beta\frac{SI}{N}.\]

The corresponding one-step probability is approximated by

\[p_{\mathrm{inf}}(S,I)=\beta\frac{SI}{N}\Delta t.\]

At the scenario state,

\[p_{\mathrm{inf}}=0.3\times\frac{990\times10}{1000}\times0.1=0.297.\]

Recovery probability

The recovery-event rate is

\[a_{\mathrm{rec}}(I)=\gamma I.\]

Therefore,

\[p_{\mathrm{rec}}(I)=\gamma I\Delta t.\]

For the scenario,

\[p_{\mathrm{rec}}=0.1\times10\times0.1=0.100.\]

No-change probability

The three probabilities must sum to one, so no change receives the remaining probability:

\[p_{\mathrm{stay}}=1-p_{\mathrm{inf}}-p_{\mathrm{rec}}.\]

Hence,

\[p_{\mathrm{stay}}=1-0.297-0.100=0.603.\]

Event	Next state	Probability
Infection	\((989,11,0)\)	0.297
Recovery	\((990,9,1)\)	0.100
No change	\((990,10,0)\)	0.603
Total	—	1.000

The one-event approximation

The formulas above approximate the chance of one infection or one recovery during a short interval. They do not represent multiple infections or recoveries within the same step.

Making \(\Delta t\) smaller reduces the neglected probability of multiple events, but increases the number of simulation steps.

This is one way to construct a DTMC from continuous event rates. Other discrete-time epidemic models may use binomial numbers of infections and recoveries per step; those are different model constructions.

Time-step validity at the current state

The combined event probability is

\[p_{\mathrm{event}}=\left(\beta\frac{SI}{N}+\gamma I\right)\Delta t.\]

We require

\[p_{\mathrm{event}}\leq1.\]

If the total event rate is positive, the largest step satisfying this probability condition at the current state is

\[\Delta t_{\max}=\frac{1}{\beta SI/N+\gamma I}.\]

Necessaryall probabilities lie in \([0,1]\)

Necessaryprobabilities sum to one

Accuracythe step is sufficiently short for the one-event approximation

Satisfying \(p_{\mathrm{event}}\leq1\) is necessary but does not by itself guarantee that the approximation is biologically accurate.

Current-state validity versus whole-model validity

A step can be valid at one state but invalid at another because event rates change with \(S\) and \(I\).

To use one fixed \(\Delta t\) for an entire DTMC, it should satisfy the probability conditions for every state that the model may visit—not only the initial state.

Boundary states

No susceptible people

If \(S=0\), then \(p_{\mathrm{inf}}=0\) because infection cannot move anyone out of \(S\).

No infectious people

If \(I=0\), then both event probabilities are zero and \(p_{\mathrm{stay}}=1\). The epidemic is absorbing in this closed model.

Program the three probabilities

infection_rate = beta * S * I / N
recovery_rate = gamma * I

p_infection = infection_rate * dt
p_recovery = recovery_rate * dt
p_stay = 1 - p_infection - p_recovery

The rate variables retain units of events per day. The probability variables are dimensionless.

Run the event-probability calculator

Interactive PythonSIR event probabilities

import matplotlib.pyplot as plt

# Current SIR state
S = 990
I = 10
R = 0
N = S + I + R

# Event-rate parameters
beta = 0.3
gamma = 0.1

# Fixed DTMC time step
dt = 0.1

# Validate inputs
if min(S, I, R, beta, gamma) < 0:
    raise ValueError("Counts and parameters cannot be negative.")

if N <= 0:
    raise ValueError("Population must be positive.")

if dt <= 0:
    raise ValueError("dt must be positive.")

# Calculate event rates
infection_rate = beta * S * I / N
recovery_rate = gamma * I
total_event_rate = infection_rate + recovery_rate

# Convert rates into short-step probabilities
p_infection = infection_rate * dt
p_recovery = recovery_rate * dt
p_stay = 1 - p_infection - p_recovery

# Validate probabilities
probabilities = {
    "Infection": p_infection,
    "Recovery": p_recovery,
    "No change": p_stay
}

for event, probability in probabilities.items():
    if probability < -1e-12 or probability > 1 + 1e-12:
        raise ValueError(
            f"Invalid probability for {event}: "
            f"{probability:.3f}. Reduce dt."
        )

# Maximum step satisfying the current-state sum condition
if total_event_rate > 0:
    maximum_current_dt = 1 / total_event_rate
else:
    maximum_current_dt = float("inf")

# Possible next states
next_states = {
    "Infection": (S - 1, I + 1, R),
    "Recovery": (S, I - 1, R + 1),
    "No change": (S, I, R)
}

# Display rates and probabilities
print("Current state:", (S, I, R))
print(f"Time step: {dt:.3f} day")
print("\nEvent       Rate per day   Probability   Next state")
print("------------------------------------------------------")

for event in probabilities:
    if event == "Infection":
        rate = infection_rate
    elif event == "Recovery":
        rate = recovery_rate
    else:
        rate = 0.0
    
    print(
        f"{event:10s} {rate:12.3f} "
        f"{probabilities[event]:12.3f} "
        f"{next_states[event]}"
    )

print("------------------------------------------------------")
print(f"Probability total: {sum(probabilities.values()):.3f}")
print(
    "Largest dt satisfying the current-state "
    f"sum condition: {maximum_current_dt:.3f}"
)

# Plot the event distribution
plt.figure(figsize=(8, 4.6))
plt.bar(
    probabilities.keys(),
    probabilities.values(),
    color=["indianred", "steelblue", "seagreen"]
)
plt.xlabel("Next-step outcome")
plt.ylabel("Probability")
plt.title("One-event SIR DTMC probabilities")
plt.ylim(0, 1)
plt.grid(axis="y", alpha=0.25)
plt.tight_layout()

Output

Run the code to see the result.

Understand the new code

Code	Meaning
`float("inf")`	Represents positive infinity when the total event rate is zero and no upper bound is required by the probability sum at that state.
`if ... elif ... else`	Selects the matching displayed rate for each event label.
`probabilities.keys()`	Returns event labels for the horizontal graph axis.
`probabilities.values()`	Returns the corresponding probability heights.
`1 + 1e-12`	Allows only a tiny floating-point tolerance above one; it does not permit a meaningfully invalid probability.

Try changing the time step

Run the program with dt = 0.01, 0.1, 0.25 and 1.0.

Compare the three probabilities.
Identify when no-change probability becomes negative.
Compare each value with maximum_current_dt.
Explain why merely obtaining non-negative probabilities does not prove the one-event approximation is sufficiently accurate.

Boundary of this lesson

The code calculates all possible event probabilities and next states but does not choose one. Lesson 4 will map these probabilities to uniform-random-number intervals.