Remove all objects from workspace.
r 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.
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 = 25000 # susceptible hosts
Y = 7 # 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 SIS 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 Infectious sub-population.
plot (I ~ time, data = output, type='b', col = 'red')
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, 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)
# remain on same frame
par (new = TRUE)
In a Susceptible Infectious Susceptible epidemiological model, a susceptible host can become infectious and can immediately reutrn to becoming susceptible once again. This model computes, analyzes, and illustrates the numbers of susceptible (blue) and infectious hosts over time (red). As time progresses, the two host populations parallel each other.