A primer on stochastic epidemic models: Formulation, numerical simulation, and analysis

Abstract

Some mathematical methods for formulation and numerical simulation of stochastic epidemic models are presented. Specifically, models are formulated for continuous-time Markov chains and stochastic differential equations. Some well-known examples are used for illustration such as an SIR epidemic model and a host-vector malaria model. Analytical methods for approximating the probability of a disease outbreak are also discussed.

Keywords: 60H10; 60J28; 92D30; Branching process; Continuous-time Markov chain; Minor outbreak; Stochastic differential equation.

Fig. 1

The dashed curve is the ODE solution of $I\left(t\right)$ in the SIR model and the other curves are four sample paths of the SIR CTMC model. Parameter values are $\beta =0.3$, $\gamma =0.15$, and either $N=100$ (top graphs) or $N=750$ (bottom graphs). The initial conditions are $S\left(0\right)=98$ and $I\left(0\right)=2$. The graphs on the right are a close-up view on the time interval $\left[\mathrm{0,15}\right]$ of the graphs on the left. The value $\mathrm{\mathbb{P}}=0.25$ is the estimate for the probability of a minor outbreak. On average, one sample path out of four hits zero, given ${\mathcal{R}}_0=2$ and $I\left(0\right)=2$.

Fig. 2

The dashed curve is the ODE solution of $I\left(t\right)$ in the SIR model and the other solid curves are four sample paths of the SIR SDE model. Parameter values and initial conditions are the same as for the CTMC model in Fig. 1 for $N=100$ (left) and $N=750$ (right).

Fig. 3

The smooth curve in each figure is the solution of the host-vector ODE model. The other curves are one sample path of the CTMC model (top two figures) and one sample path of the SDE model (bottom two figures). The parameters are $H=100,M=1000$, $\alpha =0.2$, $\gamma =0.1$, $\beta =0.2$, and $\mu =0.5$. Initial conditions are $I\left(0\right)=0$, $V\left(0\right)=10$, $S\left(0\right)=100$ and $U\left(0\right)=990.$ The value of ${\mathcal{R}}_0=8$ with $q_1=1/6$ and $q_2=3/4$. The probability of extinction in the CTMC is approximately $\mathrm{\mathbb{P}}\approx 0.06$. The stable endemic equilibrium in the ODE model is $\left(\overline{S},\overline{I},\overline{U},\overline{V}\right)$ = (16.7, 83.3, 750, 250).

Fig. 4

The two histograms (computed from 10,000 sample paths of the host-vector SDE model) approximate the probability density for I and V at $t=50$. Parameter values are the same as in Fig. 3 but with initial conditions sufficiently large such that the probability of hitting zero is close to zero. The computed mean values and standard deviations are ${\mu }_I=83.2$, ${\mu }_V=249.6$, ${\sigma }_I=3.94$, and ${\sigma }_V=17.8$. The mean values are close to the endemic equilibrium values of the host-vector ODE model: $\left(\overline{I},\overline{V}\right)=\left(83.3,250\right)$.

Fig. 5

The smooth curve is the ODE solution of $K\left(t\right)$ (Wiener process neglected in (17)) and the random curve is one sample path of the SDE model (17). Parameter values are $\overline{K}=100$, $r=0.3$, $\hat{r}=4$, $\omega =2\pi /365$, and $\sigma =5$. The initial condition is $K\left(0\right)=110.$