Objective

To analyze the dynamics of Susceptibles-Infectious-Recovered (SIR) Epidemiological Models

SIR Model

Remove all objects from workspace

remove (list = objects() )

Load add-on packages - deSolve - contains lsoda function - differential equation solver

library (deSolve) 
## 
## Attaching package: 'deSolve'
## 
## The following object is masked from 'package:graphics':
## 
##     matplot

Function to compute derivatives of the differential equations

sir_model = function (current_timepoint, state_values, parameters)
{
  # create state variables (local variables)
  S = state_values [1]        # susceptibles
  I = state_values [2]        # infectious
  R = state_values [3]        # recovered
  
  with ( 
    as.list (parameters),     # variable names within parameters can be used 
         {
           # compute derivatives
           dS = (-beta * S * I)
           dI = ( beta * S * I) - (gamma * I)
           dR = (gamma * I)
           
           # combine results
           results = c (dS, dI, dR)
           list (results)
         }
    )
}

Parameters

contact_rate = 10                     # number of contacts per day
transmission_probability = 0.07       # transmission probability
infectious_period = 5                 # infectious period

Compute values of beta (tranmission rate) and gamma (recovery rate)

beta_value = contact_rate * transmission_probability
gamma_value = 1 / infectious_period

Compute Ro - Reproductive number

Ro = beta_value / gamma_value

Disease dynamics parameters

parameter_list = c (beta = beta_value, gamma = gamma_value)

Initial values for sub-populations

X = 9999        # susceptible hosts
Y = 1           # infectious hosts
Z = 0           # recovered hosts

Compute total population

N = X + Y + Z

Initial state values for the differential equations

initial_values = c (S = X/N, I = Y/N, R = Z/N)

Output timepoints

timepoints = seq (0, 50, by=1)

Simulate the SIR epidemic

output = lsoda (initial_values, timepoints, sir_model, parameter_list)

Plot dynamics of Susceptibles sub-population

plot (S ~ time, data = output, type='b', col = 'blue')

Plot dynamics of Infectious sub-population

plot (I ~ time, data = output, type='b', col = 'red')

Plot dynamics of Recovered sub-population

plot (R ~ time, data = output, type='b', col = 'green')

Plot dynamics of Susceptibles, Infectious and Recovered sub-populations in the same plot

# susceptible hosts over time
plot (S ~ time, data = output, type='b', ylim = c(0,1), col = 'blue', ylab = 'S, I, R', main = 'SIR epidemic') 

# remain on same frame
par (new = TRUE)    

# infectious hosts over time
plot (I ~ time, data = output, type='b', ylim = c(0,1), col = 'red', ylab = '', axes = FALSE) 

# remain on same frame
par (new = TRUE)  

# recovered hosts over time
plot (R ~ time, data = output, type='b', ylim = c(0,1), col = 'green', ylab = '', axes = FALSE)

Discussion and Results

This is a SIR model, or Susceptible-Infectious-Recovered epidemiological model. A host is susceptible to infection, then becomes infectious, and finally is recovered, typically with lifelong immunity. An example of a disease that follows a SIR model is chickenpox.

SEIR Model

Remove all objects from workspace

remove (list = objects() )

Load add-on packages - deSolve - contains lsoda function - differential equation solver

library (deSolve)

Function to compute derivatives of the differential equations

seir_model = function (current_timepoint, state_values, parameters)
{
  # create state variables (local variables)
  S = state_values [1]        # susceptibles
  E = state_values [2]        # exposed
  I = state_values [3]        # infectious
  R = state_values [4]        # recovered

  with ( 
    as.list (parameters),     # variable names within parameters can be used 
         {
           # compute derivatives
           dS = (-beta * S * I)
           dE = (beta * S * I) - (delta * E)
           dI = (delta * E) - (gamma * I)
           dR = (gamma * I)
           
           
           # combine results
           results = c (dS, dE, dI, dR)
           list (results)
         }
    )
}

Parameters

contact_rate = 10                     # number of contacts per day
transmission_probability = 0.07       # transmission probability
infectious_period = 5                 # infectious period
latent_period = 3                     # latent period

Compute values of beta (tranmission rate), gamma (recovery rate), and delta (exposion rate)

beta_value = contact_rate * transmission_probability
gamma_value = 1 / infectious_period
delta_value = 1 / latent_period

Compute Ro - Reproductive number

Ro = beta_value / gamma_value

Disease dynamics parameters

parameter_list = c (beta = beta_value, gamma = gamma_value, delta = delta_value)

Initial values for sub-populations

X = 9998        # susceptible hosts
A = 1           # exposed hosts
Y = 1           # infectious hosts
Z = 0           # recovered hosts

Compute total population

N = X + A + Y + Z

Initial state values for the differential equations

initial_values = c (S = X/N, E = A/N, I = Y/N, R = Z/N)

Output timepoints

timepoints = seq (0, 50, by=1)

Simulate the SEIR epidemic

output = lsoda (initial_values, timepoints, seir_model, parameter_list)

Plot dynamics of Susceptibles sub-population

plot (S ~ time, data = output, type='b', col = 'blue')

Plot dynamics of Exposed sub-population

plot (E ~ time, data = output, type='b', col = 'purple')

Plot dynamics of Infectious sub-population

plot (I ~ time, data = output, type='b', col = 'red')

Plot dynamics of Recovered sub-population

plot (R ~ time, data = output, type='b', col = 'green')

Plot dynamics of Susceptibles, Exposed, Infectious and Recovered sub-populations in the same plot

# susceptible hosts over time
plot (S ~ time, data = output, type='b', ylim = c(0,1), col = 'blue', ylab = 'S, E, I, R', main = 'SEIR epidemic') 

# remain on same frame
par (new = TRUE)    

# exposed hosts over time
plot (E ~ time, data = output, type='b', ylim = c(0,1), col = 'purple', ylab = 'S, E, I, R', main = 'SEIR epidemic') 

# remain on same frame
par (new = TRUE)

# infectious hosts over time
plot (I ~ time, data = output, type='b', ylim = c(0,1), col = 'red', ylab = '', axes = FALSE) 

# remain on same frame
par (new = TRUE)  

# recovered hosts over time
plot (R ~ time, data = output, type='b', ylim = c(0,1), col = 'green', ylab = '', axes = FALSE)

Discussion and Results

This is a SEIR model, or Susceptible-Exposed-Infectious-Recovered epidemiological model. A host is susceptible to infection, exposed to the pathogen, then becomes infectious, and finally is recovered. An example of a disease that follows a SIR model is influenza.

SIS Model

Remove all objects from workspace

remove (list = objects() )

Load add-on packages - deSolve - contains lsoda function - differential equation solver

library (deSolve)

Function to compute derivatives of the differential equations

sis_model = function (current_timepoint, state_values, parameters)
{
  # create state variables (local variables)
  S = state_values [1]        # susceptibles
  I = state_values [2]        # infectious

  with ( 
    as.list (parameters),     # variable names within parameters can be used 
         {
           # compute derivatives
           dS = (-beta * S * I) + (gamma * I)
           dI = ( beta * S * I) - (gamma * I)
           
           
           # combine results
           results = c (dS, dI)
           list (results)
         }
    )
}

Parameters

contact_rate = 10                     # number of contacts per day
transmission_probability = 0.07       # transmission probability
infectious_period = 5                 # infectious period

Compute values of beta (tranmission rate) and gamma (recovery rate)

beta_value = contact_rate * transmission_probability
gamma_value = 1 / infectious_period

Compute Ro - Reproductive number

Ro = beta_value / gamma_value

Disease dynamics parameters

parameter_list = c (beta = beta_value, gamma = gamma_value)

Initial values for sub-populations

X = 9999        # susceptible hosts
Y = 1           # infectious hosts

Compute total population

N = X + Y

Initial state values for the differential equations

initial_values = c (S = X/N, I = Y/N)

Output timepoints

timepoints = seq (0, 50, by=1)

Simulate the SI epidemic

output = lsoda (initial_values, timepoints, sis_model, parameter_list)

Plot dynamics of Susceptibles sub-population

plot (S ~ time, data = output, type='b', col = 'blue')

Plot dynamics of Infected sub-population

plot (I ~ time, data = output, type='b', col = 'red')

Plot dynamics of Susceptibles and Infectious sub-populations in the same plot

# susceptible hosts over time
plot (S ~ time, data = output, type='b', ylim = c(0,1), col = 'blue', ylab = 'S, I, S', main = 'SIS epidemic') 

# remain on same frame
par (new = TRUE)

# infectious hosts over time
plot (I ~ time, data = output, type='b', ylim = c(0,1), col = 'red', ylab = '', axes = FALSE)

Discussion and Results

This is a SIS model, or Susceptible-Infectious-Susceptible epidemiological model. A host is susceptible to infection, is infectious, and is then susceptible once again. An example of a pathogen that follows a SEIR model is rotaviruses, which can result in gastrointestinal illness.

SEIS Model

Remove all objects from workspace

remove (list = objects() )

Load add-on packages - deSolve - contains lsoda function - differential equation solver

library (deSolve)

Function to compute derivatives of the differential equations

seis_model = function (current_timepoint, state_values, parameters)
{
  # create state variables (local variables)
  S = state_values [1]        # susceptibles
  E = state_values [2]        # exposed
  I = state_values [3]        # infectious

  with ( 
    as.list (parameters),     # variable names within parameters can be used 
         {
           # compute derivatives
           dS = (-beta * S * I) + (gamma * I)
           dE = (beta * S * I) - (delta * E)
           dI = (delta * E) - (gamma * I)
           
           
           # combine results
           results = c (dS, dE, dI)
           list (results)
         }
    )
}

Parameters

contact_rate = 10                     # number of contacts per day
transmission_probability = 0.07       # transmission probability
infectious_period = 5                 # infectious period
latent_period = 3                     # latent period

Compute values of beta (tranmission rate) and gamma (recovery rate)

beta_value = contact_rate * transmission_probability
gamma_value = 1 / infectious_period
delta_value = 1 / latent_period

Compute Ro - Reproductive number

Ro = beta_value / gamma_value

Disease dynamics parameters

parameter_list = c (beta = beta_value, gamma = gamma_value, delta = delta_value)

Initial values for sub-populations

X = 9998        # susceptible hosts
Y = 1           # infectious hosts
Z = 1           # exposed hosts

Compute total population

N = X + Y + Z

Initial state values for the differential equations

initial_values = c (S = X/N, E = Z/N, I = Y/N)

Output timepoints

timepoints = seq (0, 50, by=1)

Simulate the SEIS epidemic

output = lsoda (initial_values, timepoints, seis_model, parameter_list)

Plot dynamics of Susceptibles sub-population

plot (S ~ time, data = output, type='b', col = 'blue')

Plot dynamics of Exposed sub-population

plot (E ~ time, data = output, type='b', col = 'purple')

Plot dynamics of Infectious sub-population

plot (I ~ time, data = output, type='b', col = 'red')

Plot dynamics of Susceptibles, Exposed, and Infectious sub-populations in the same plot

# susceptible hosts over time
plot (S ~ time, data = output, type='b', ylim = c(0,1), col = 'blue', ylab = 'S, E, I, S', main = 'SEIS epidemic') 

# remain on same frame
par (new = TRUE)    

# exposed hosts over time
plot (E ~ time, data = output, type='b', ylim = c(0,1), col = 'purple', ylab = 'S, E, I, S', main = 'SEIS epidemic') 

# remain on same frame
par (new = TRUE)

# infectious hosts over time
plot (I ~ time, data = output, type='b', ylim = c(0,1), col = 'red', ylab = '', axes = FALSE)

Discussion and Results

This is a SEIS model, or Susceptible-Exposed-Infectious-Susceptible epidemiological model. A host is susceptible to infection, is exposed and then infectious, and is then susceptible once again. An example of a disease that follows a SEIR model is tuberculosis.

SI Model

Remove all objects from workspace

remove (list = objects() )

Load add-on packages - deSolve - contains lsoda function - differential equation solver

library (deSolve)

Function to compute derivatives of the differential equations

si_model = function (current_timepoint, state_values, parameters)
{
  # create state variables (local variables)
  S = state_values [1]        # susceptibles
  I = state_values [2]        # infectious

  with ( 
    as.list (parameters),     # variable names within parameters can be used 
         {
           # compute derivatives
           dS = (-beta * S * I)
           dI = ( beta * S * I) 
           
           
           # combine results
           results = c (dS, dI)
           list (results)
         }
    )
}

Parameters

contact_rate = 10                     # number of contacts per day
transmission_probability = 0.07       # transmission probability
infectious_period = 5                 # infectious period

Compute values of beta (tranmission rate) and gamma (recovery rate)

beta_value = contact_rate * transmission_probability
gamma_value = 1 / infectious_period

Compute Ro - Reproductive number

Ro = beta_value / gamma_value

Disease dynamics parameters

parameter_list = c (beta = beta_value, gamma = gamma_value)

Initial values for sub-populations

X = 9999        # susceptible hosts
Y = 1           # infectious hosts

Compute total population

N = X + Y

Initial state values for the differential equations

initial_values = c (S = X/N, I = Y/N)

Output timepoints

timepoints = seq (0, 50, by=1)

Simulate the SI epidemic

output = lsoda (initial_values, timepoints, si_model, parameter_list)

Plot dynamics of Susceptibles sub-population

plot (S ~ time, data = output, type='b', col = 'blue')

Plot dynamics of Infected sub-population

plot (I ~ time, data = output, type='b', col = 'red')

Plot dynamics of Susceptibles and Infectious sub-populations in the same plot

# susceptible hosts over time
plot (S ~ time, data = output, type='b', ylim = c(0,1), col = 'blue', ylab = 'S, I, S', main = 'SI epidemic') 

# remain on same frame
par (new = TRUE)

# infectious hosts over time
plot (I ~ time, data = output, type='b', ylim = c(0,1), col = 'red', ylab = '', axes = FALSE)

Discussion and Results

This is a SI model, or Susceptible-Infectious epidemiological model. A host is susceptible to infection, and is then infectious with no recovery possible. An example of a disease that follows a SI model is HIV.

SIR/C Model

Remove all objects from workspace

remove (list = objects() )

Load add-on packages - deSolve - contains lsoda function - differential equation solver

library (deSolve)

Function to compute derivatives of the differential equations

sirc_model = function (current_timepoint, state_values, parameters)
{
  # create state variables (local variables)
  S = state_values [1]        # susceptibles
  I = state_values [2]        # infectious
  R = state_values [3]        # recovered
  C = state_values [4]        # carrier
  
  with ( 
    as.list (parameters),     # variable names within parameters can be used 
         {
           # compute derivatives
           dS = (-beta * S * I) - (epsilon * beta * S * C)
           dI = ( beta * S * I) + (epsilon * beta * S * C) - (gamma * I)
           dC = (rho * gamma * I) - (tau * C)
           dR = (1-rho) * gamma * I + (tau * C)
           
           # combine results
           results = c (dS, dI, dC, dR)
           list (results)
         }
    )
}

Parameters

contact_rate = 10                     # number of contacts per day
transmission_probability = 0.07       # transmission probability
infectious_period = 5                 # infectious period
reduced_transmission_rate = 0.2       # reduced transmission rate
proportion_acute_infections = 0.4     # proportion of acute infections that become carriers
time_in_carrier_state = 5             # length of time spent in carrier state

Compute values of beta (tranmission rate), gamma (recovery rate), epsilon (reduced transmission rate), rho (proportion of acute infections that become carriers), and tau (rate at which individuals leave the carrier class)

beta_value = contact_rate * transmission_probability
gamma_value = 1 / infectious_period
epsilon_value = 0.2
rho_value = 0.4
tau_value = 1/time_in_carrier_state

Compute Ro - Reproductive number

Ro = beta_value / gamma_value

Disease dynamics parameters

parameter_list = c (beta = beta_value, gamma = gamma_value, epsilon = epsilon_value, rho = rho_value, tau = tau_value)

Initial values for sub-populations

X = 9987        # susceptible hosts
Y = 1           # infectious hosts
Z = 10           # recovered hosts
A = 2           # carriers

Compute total population

N = X + Y + Z + A

Initial state values for the differential equations

initial_values = c (S = X/N, I = Y/N, R = Z/N, C = A/N)

Output timepoints

timepoints = seq (0, 50, by=1)

Simulate the SIR/C epidemic

output = lsoda (initial_values, timepoints, sirc_model, parameter_list)

Plot dynamics of Susceptibles sub-population

plot (S ~ time, data = output, type='b', col = 'blue')

Plot dynamics of Infectious sub-population

plot (I ~ time, data = output, type='b', col = 'red')

Plot dynamics of Recovered sub-population

plot (R ~ time, data = output, type='b', col = 'green')

Plot dynamics of Carrier sub-population

plot (C ~ time, data = output, type='b', col = 'purple')

Plot dynamics of Susceptibles, Infectious, Recovered and Carrier sub-populations in the same plot

# susceptible hosts over time
plot (S ~ time, data = output, type='b', ylim = c(0,1), col = 'blue', ylab = 'S, I, R, C', main = 'SIR/C epidemic') 

# remain on same frame
par (new = TRUE)    

# infectious hosts over time
plot (I ~ time, data = output, type='b', ylim = c(0,1), col = 'red', ylab = '', axes = FALSE) 

# remain on same frame
par (new = TRUE)  

# recovered hosts over time
plot (R ~ time, data = output, type='b', ylim = c(0,1), col = 'green', ylab = '', axes = FALSE)  

#remain on same frame
par (new = TRUE)

#carrier hosts over time
plot (C ~ time, data = output, type='b', ylim = c(0,1), col = 'purple', ylab = '', axes = FALSE)

Discussion and Results

This is a SIR/C model, or Susceptible-Infectious-Recovered/Carrier epidemiological model. A host is susceptible to infection, then becomes infectious, and finally is recovered and this model takes into account a carrier population, meaning those that carry a disease-causing pathogen but are not sick themselves. An example of a disease that follows a SIR/C model is hepatitis B.