---
title: "SIR model"
output: rmarkdown::html_vignette
bibliography: references.bib
nocite: "@Hens2012"
link-citations: TRUE
vignette: >
%\VignetteIndexEntry{SIR model}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
---
```{r setup, output=FALSE}
library(serosv)
```
```{r, include = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>"
)
```
## Basic SIR model
**Proposed model**
A transmission model consists of 3 compartments: susceptible (S), infected (I), recovered (R)
With the following assumptions:
- Individuals are born into susceptible group (exposure time is age of the individual) then transfer to infected class and recovered class
- Recovered individuals gained lifelong immunity
- Age homogeneity
And described by a system of 3 differential equations
$$
\begin{cases}
\frac{dS(t)}{dt} = B(t) (1-p) - \lambda(t)S(t) - \mu S(t) \\
\frac{dI(t)}{dt} = \lambda(t)S(t) - \nu I(t) - \mu I(t) - \alpha I(t) \\
\frac{dR(t)}{dt} = B(t) p + \nu I(t) - \mu R(t)
\end{cases}
$$
Where:
- $B(t) = \mu N(t)$
- $\lambda(t) = \beta I(t)$ with $\beta$ is the transmission rate
- $\mu$ is the natural death rate
- $\nu$ is the recovery rate
- $\alpha$ is the disease related death rate
- $p$ is the proportion of newborn vaccinated and moved directly to the recovered compartment
**Fitting data**
To fit a basic SIR model, use `sir_basic_model()` and specify the following parameters
- `state` - initial population of each compartment
- `times` - a time sequence
- `parameters` - parameters for SIR model
```{r}
state <- c(S=4999, I=1, R=0)
parameters <- c(
mu=1/75, # 1 divided by life expectancy (75 years old)
alpha=0, # no disease-related death
beta=0.0005, # transmission rate
nu=1, # 1 year for infected to recover
p=0 # no vaccination at birth
)
times <- seq(0, 250, by=0.1)
model <- sir_basic_model(times, state, parameters)
model$parameters
plot(model)
```
## SIR model with constant Force of Infection at Endemic state
**Proposed model**
A transmission model consists of 3 compartments: susceptible (S), infected (I), recovered (R)
With the following assumptions:
- Time homogeneity
- Age heterogeneity
Described by a system of 3 differential equations
$$
\begin{cases}
\frac{ds(a)}{da} = -\lambda s(a) \\
\frac{di(a)}{da} = \lambda s(a) - \nu i(a) \\
\frac{dr(a)}{da} = \nu i(a)
\end{cases}
$$
Where:
- $s(a), i(a), r(a)$ are proportion of susceptible, infected, recovered population of age group $a$ respectively
- $\lambda$ is the force of infection
- $\nu$ is the recovery rate
**FItting data**
To fit an SIR model with constant FOI, use `sir_static_model()` and specify the following parameters
- `state` - initial proportion of each compartment
- `ages` - an age sequence
- `parameters` - parameters for the model
```{r}
state <- c(s=0.99,i=0.01,r=0)
parameters <- c(
lambda = 0.05,
nu=1/(14/365) # 2 weeks to recover
)
ages<-seq(0, 90, by=0.01)
model <- sir_static_model(ages, state, parameters)
model$parameters
plot(model)
```
## SIR model with sub populations
**Proposed model**
Extends on the SIR model by having interacting sub-populations (different age groups)
With K subpopulations, the WAIFW matrix or mixing matrix is given by
$$
C = \begin{bmatrix}
\beta_{11} & \beta_{12} & ... & \beta_{1K} \\
\beta_{21} & \beta_{22} & ... & \beta_{2K} \\
\vdots & \vdots & ... & \vdots \\
\beta_{K1} & \beta_{K2} & ... & \beta_{KK} \\
\end{bmatrix}
$$
The system of differential equations for the i$th$ subpopulation is given by
$$
\begin{cases}
\frac{dS_i(t)}{dt} = -(\sum^K_{j=1}\beta_{ij}I_j(t)) S_i(t) + N_i\mu_i - \mu_i S_i(t) \\
\frac{dI_i(t)}{dt} = (\sum^K_{j=1}\beta_{ij}I_j(t)) S_i(t) - (\nu_i + \mu_i) I_i(t) \\
\frac{dR_i(t)}{dt} = \nu_i I_i(t) - \mu_i R_i(t)
\end{cases}
$$
**FItting data**
To fit a SIR model with subpopulations, use `sir_subpops_model()` and specify the following parameters
- `state` - initial proportion of each compartment for each population
- `beta_matrix` - the WAIFW matrix
- `times` - a time sequence
- `parameters` - parameters for the model
```{r}
k <- 2 # number of population
state <- c(
# proportion of each compartment for each population
s = c(0.8, 0.6),
i = c(0.2, 0.4),
r = c( 0, 0)
)
beta_matrix <- c(
c(0.05, 0.00),
c(0.00, 0.05)
)
parameters <- list(
beta = matrix(beta_matrix, nrow=k, ncol=k, byrow=TRUE),
nu = c(1/30, 1/30),
mu = 0.001,
k = k
)
times<-seq(0,10000,by=0.5)
model <- sir_subpops_model(times, state, parameters)
model$parameters
plot(model) # returns plot for each population
```
## MSEIR model
**Proposed model**
Extends on SIR model with 2 additional compartments: maternal antibody (M) and exposed period (E)
And described by the following system of ordinary differential equation
$$
\begin{cases}
\frac{dM(a)}{da} = -(\gamma + \mu(a))M(a) \\
\frac{dS(a)}{da} = \gamma M(a) - (\lambda(a) + \mu(a)) S(a) \\
\frac{dE(a)}{da} = \lambda(a) S(a) - (\sigma + \mu(a)) E(a) \\
\frac{dI(a)}{da} = \sigma(a) E(a) - (\nu + \mu(a)) I(a) \\
\frac{dR(a)}{da} = \nu I(a) - \mu(a) R(a)
\end{cases}
$$
Where
- $M(0)$ = B, the number of births in the population
- $\gamma$ is the rate of antibody decaying
- $\lambda(a)$ is the force of infection at age $a$
- $\mu(a)$ is the natural death rate at age $a$
- $\sigma$ is the rate of becoming infected after being exposed
- $\nu$ is the recovery rate
**Fitting data**
To fit a MSEIR, use `mseir_model()` and specify the following parameters
- `a` - age sequence
- And model parameters including `gamma`, `lambda`, `sigma`, `nu`
```{r}
model <- mseir_model(
a=seq(from=1,to=20,length=500), # age range from 0 -> 20 yo
gamma=1/0.5, # 6 months in the maternal antibodies
lambda=0.2, # 5 years in the susceptible class
sigma=26.07, # 14 days in the latent class
nu=36.5 # 10 days in the infected class
)
model$parameters
plot(model)
```