SEIS and SEISd

Compartment Model as an ODE

The SEIS model is a human model modified from the SIS model to include the Exposed group of individuals ($E$). It is incorporated within the ‘ramp.xds’ with the fulfillment of the generic interface of the human component.

The model has three parameters:

• $b$ is the fraction of infective bites that cause an infection;

• $\nu$ is the transition rate from exposed to infectious: the duration of the latent period is $1/\nu$

• $r$ is the clearance rate for infections: the average duration of infection in this model is $1/r$

These are coupled systems of ordinary differential equations forced by the force of infection, denoted $h(t)$, where $h = F_h(E);$ here, we assume that exposure is Poisson, so we let:

$$h = bE$$

$$\frac{dS}{dt} = - h S + r I \\ \frac{dE}{dt} = h S - \nu E \\ \frac{dI}{dt} = \nu E - r I\\ $$

HPop = 1000
MYZo = list(Z = 2000/365, f=1, q=1)
Xo = list(b= 0.5)

seis = xds_setup(Xname = "SEIS", MYZname = "trivial", MYZopts = MYZo, Xopts=Xo, HPop=HPop)

seis = xds_solve(seis, 3650, 15)

unlist(list_Xvars(seis$outputs$last_y, seis, 1)) -> seis_inf
seis_inf
#>          S          E          I          H 
#>  652.95170   25.04472  322.00358 1000.00000

This model has the steady state…

$$\bar E = H \frac {hr}{h(r+\nu) + r \nu} \\ \bar I = H \frac {h\nu}{h(r+\nu) + r \nu} \\ \bar S = H - \bar E - \bar I$$

xde_steady_state_X(1/365, 1000, seis$Xpar[[1]]) -> seis_ss
seis_ss
#>          S          E          I          H 
#>  652.95170   25.04472  322.00358 1000.00000

sum((seis_inf-seis_ss)^2) < 1e-9
#> [1] TRUE

Compartment Model as an DDE

In the delay differential equation model, we let $\nu$ denote the duration of the incubation period, and we let $h_\nu = h(t-\nu).$

$$\frac{dS}{dt} = - h S + r I \\ \frac{dE}{dt} = h S - h_\nu S_\nu \\ \frac{dI}{dt} = h_\nu S_\nu - r I\\ $$ and $H = S+E+I.$ At steady state, $h=h_{v}$, so

$$\bar S = \frac{H}{1+h \nu + h/r} \\ \bar E = h \bar S \nu \\ \bar I = \frac hr \bar S \\ $$

SEISd = xds_setup(Xname = "SEISd", MYZname = "trivial", MYZopts = MYZo, Xopts=Xo, HPop=HPop)
SEISd = xds_solve(SEISd, Tmax = 3650, dt=15)

unlist(list_Xvars(SEISd$outputs$last_y, SEISd, 1))[1:4] -> seisd_inf
seisd_inf
#>          S          E          I          H 
#>  652.95173   25.04467  322.00359 1000.00000

xde_steady_state_X(1/365, 1000, SEISd$Xpar[[1]]) -> seisd_ss
seisd_ss
#>          S          E          I          H 
#>  652.95170   25.04472  322.00358 1000.00000

sum((seisd_inf-seisd_ss)^2) < 1e-5
#> [1] TRUE