{"word_count": 441, "character_count": 4728}
In Statistical Inference in a Stochastic Epidemic SEIR Model with Control Intervention, a more complex model than the one we've seen yesterday was considered (and is called the SEIR model). Consider a population of size [latex]N[/latex], and assume that [latex]S[/latex] is the number of susceptible, [latex]E[/latex] the number of exposed, [latex]I[/latex] the number of infectious, and [latex]R[/latex] for the number recovered (or immune) individuals, [latex display="true"]\displaystyle{\begin{aligned}{\frac {dS}{dt}}&#38;=-\beta {\frac {I}{N}}S\\[8pt]{\frac {dE}{dt}}&#38;=\beta {\frac {I}{N}}S-aE\\[8pt]{\frac {dI}{dt}}&#38;=aE-b I\\[8pt]{\frac {dR}{dt}}&#38;=b I\end{aligned}}[/latex]Between [latex]S[/latex] and [latex]I[/latex], the transition rate is [latex]\beta I[/latex], where [latex]\beta[/latex] is the average number of contacts per person per time, multiplied by the probability of disease transmission in a contact between a susceptible and an infectious subject. Between [latex]I[/latex] and [latex]R[/latex], the transition rate is [latex]b[/latex] (simply the rate of recovered or dead, that is, number of recovered or dead during a period of time divided by the total number of infected on that same period of time). And finally, the incubation period is a random variable with exponential distribution with parameter [latex]a[/latex], so that the average incubation period is [latex]a^{-1}[/latex].Probably more interesting, Understanding the dynamics of ebola epidemics suggested a more complex model, with susceptible people [latex]S[/latex], exposed [latex]E[/latex], Infectious, but either in community [latex]I[/latex], or in hospitals [latex]H[/latex], some people who died [latex]F[/latex] and finally those who either recover or are buried and therefore are no longer susceptible [latex]R[/latex].Thus, the following dynamic model is considered[latex display="true"]\displaystyle{\begin{aligned}{\frac {dS}{dt}}&#38;=-(\beta_II+\beta_HH+\beta_FF)\frac{S}{N}\\[8pt]\frac {dE}{dt}&#38;=(\beta_II+\beta_HH+\beta_FF)\frac{S}{N}-\alpha E\\[8pt]\frac {dI}{dt}&#38;=\alpha E+\theta\gamma_H I-(1-\theta)(1-\delta)\gamma_RI-(1-\theta)\delta\gamma_FI\\[8pt]\frac {dH}{dt}&#38;=\theta\gamma_HI-\delta\lambda_FH-(1-\delta)\lambda_RH\\[8pt]\frac {dF}{dt}&#38;=(1-\theta)(1-\delta)\gamma_RI+\delta\lambda_FH-\nu F\\[8pt]\frac {dR}{dt}&#38;=(1-\theta)(1-\delta)\gamma_RI+(1-\delta)\lambda_FH+\nu F\end{aligned}}[/latex]In that model, parameters are [latex]\alpha^{-1}[/latex] is the (average) incubation period (7 days), [latex]\gamma_H^{-1}[/latex] the onset to hospitalization (5 days), [latex]\gamma_F^{-1}[/latex] the onset to death (9 days), [latex]\gamma_R^{-1}[/latex] the onset to "recovery" (10 days), [latex]\lambda_F^{-1}[/latex] the hospitalisation to death (4 days) while [latex]\lambda_R^{-1}[/latex] is the hospitalisation to recovery (5 days), [latex]\eta^{-1}[/latex] is the death to burial (2 days). Here, numbers are from Understanding the dynamics of ebola epidemics (in the context of ebola). The other parameters are [latex]\beta_I[/latex] the transmission rate in community (0.588), [latex]\beta_H[/latex] the transmission rate in hospital (0.794) and [latex]\beta_F[/latex] the transmission rate at funeral (7.653). Thusepsilon = 0.001 Z = c(S = 1-epsilon, E = epsilon, I=0,H=0,F=0,R=0)p=c(alpha=1/7*7, theta=0.81, delta=0.81, betai=0.588,    betah=0.794, blambdaf=7.653,N=1, gammah=1/5*7,    gammaf=1/9.6*7, gammar=1/10*7, lambdaf=1/4.6*7,    lambdar=1/5*7, nu=1/2*7)If [latex]\boldsymbol{Z}=(S,E,I,H,F,R)[/latex], if we write [latex display="true"]\frac{\partial \boldsymbol{Z}}{\partial t} = SEIHFR(\boldsymbol{Z})[/latex]where [latex]SEIHFR[/latex] isSEIHFR = function(t,Z,p){  S=Z[1]; E=Z[2]; I=Z[3]; H=Z[4]; F=Z[5]; R=Z[6]  alpha=p["alpha"]; theta=p["theta"]; delta=p["delta"]  betai=p["betai"]; betah=p["betah"]; gammah=p["gammah"]  gammaf=p["gammaf"]; gammar=p["gammar"]; lambdaf=p["lambdaf"]  lambdar=p["lambdar"]; nu=p["nu"]; blambdaf=p["blambdaf"]  N=S+E+I+H+F+R  dS=-(betai*I+betah*H+blambdaf*F)*S/N  dE=(betai*I+betah*H+blambdaf*F)*S/N-alpha*E  dI=alpha*E-theta*gammah*I-(1-theta)*(1-delta)*gammar*I-(1-theta)*delta*gammaf*I  dH=theta*gammah*I-delta*lambdaf*H-(1-delta)*lambdaf*H  dF=(1-theta)*(1-delta)*gammar*I+delta*lambdaf*H-nu*F  dR=(1-theta)*(1-delta)*gammar*I+(1-delta)*lambdar*H+nu*F  dZ=c(dS,dE,dI,dH,dF,dR)  list(dZ)}We can solve it, or at least study the dynamics from some starting valueslibrary(deSolve)times = seq(0, 50, by = .1)resol = ode(y=Z, times=times, func=SEIHFR, parms=p)For instance, the proportion of people infected is the followingplot(resol[,"time"],resol[,"I"],type="l",xlab="time",ylab="",col="red")lines(resol[,"time"],resol[,"H"],col="blue")