Materials licensed under CC-BY 4.0 by the authors
Template licensed under CC-BY 4.0 by The Carpentries
Built with sandpaper (0.14.0), - pegboard (0.7.1), +
Built with sandpaper (0.14.0), + pegboard (0.7.1.9000), and varnish (0.3.1).
Last updated on 2023-10-26 | +
Last updated on 2023-10-31 | Edit this page
Content from Lesson title
Content from Simulating transmission
Last updated on 2023-10-31 | - Edit this page
Using the R package epidemics, learn how to:
epidemics
List (and hyperlink) the lessons/packages which need to be covered -before this lesson
This tutorial has the following concept dependencies:
Modelling : Components of +infectious disease models e.g. state variables, parameters, initial +conditions, Ordinary differential equations.
Epidemic theory : Transmission, Reproduction +number.
Introductory text
Mathematical models are useful tools for generating future +trajectories of disease spread. Models can be used to evaluate the +implementation of non-pharmaceutical and pharmaceutical interventions +while accounting for factors such as age.
In this tutorial, we will use the R package {epidemics} +to generate trajectories of influenza spread. By the end of this +tutorial, you will be able to generate the trajectory below showing the +number of infectious individuals in different age categories through +time.
{epidemics}
The first step is to install the R packages +epidemics.
+if (!require("pak")) install.packages("pak") +pak::pak("epiverse-trace/epidemics")
if (!require("pak")) install.packages("pak") +pak::pak("epiverse-trace/epidemics")
Lesson content
To generate predictions of infectious disease trajectories, we must +first select a mathematical model to use.
There is a library of models to choose from in +epidemics. Models are prefixed with epidemic and suffixed +by the infection name. In this tutorial, we will use the default +epidemic model, epidemic_default() which is described in +the next section.
epidemic_default()
When using existing model structures always check the model +assumptions. Ask questions such as:
There can be subtle differences in model structures for the same +infection or outbreak type which can be missed without studying the +equations.
We want to generate disease trajectories of an influenza strain with +pandemic potential. We will use the default epidemic model in +epidemics, an age-structured SEIR model described by a +system of ordinary differential equations. For each age group \(i\), individuals are classed as either +susceptible \(S\), infected but not yet +infectious \(E\), infectious \(I\) or recovered \(R\).
The model parameters and equations are as follows :
\[ +\begin{aligned} +\frac{dS_i}{dt} & = - \beta S_i \sum_j C_{i,j} I_j \\ +\frac{dE_i}{dt} &= \beta S_i\sum_j C_{i,j} I_j - \alpha E_i \\ +\frac{dI_i}{dt} &= \alpha E_i - \gamma I_i \\ +\frac{dR_i}{dt} &=\gamma I_i \\ +\end{aligned} +\]
The contact matrix is a square matrix consisting of +rows/columns equal to the number age groups. Each element represents the +frequency of contacts between age groups. If we believe that +transmission of an infection is driven by contact, and that contact +rates are very different for different age groups, then specifying a +contact matrix allows us to account for age specific rates of +transmission.
From the model structure we see that :
This model also has the functionality to include vaccination and +tracks the number of vaccinated individuals through time. We will cover +the use of interventions in future tutorials.
Confusion sometimes arises when referring to the terms ‘exposed’, +‘infected’ and ‘infectious’ in mathematical modelling. Infection occurs +after a person has been exposed, but in modelling terms individuals that +are ‘exposed’ are treated as already infected.
We will use the following definitions for our state variables:
To generate trajectories using our model, we need the following :
To run our model we need to specify the model parameters:
We will learn how to specify the contact matrix \(C\) in the next section.
We will simulate a strain of influenza with pandemic potential with +\(R_0=1.5\), a preinfectious period of +3 days and infectious period of 7 days.
In epidemics, we use the function +infection() to create an infection object containing the +values of, \(R_0\), the preinfectious +period (\(1/\alpha\)) and the +infectious period (\(1/\gamma\)) as +follows.
infection()
+influenza <- infection( + name = "influenza", + r0 = 1.5, + preinfectious_period = 3, + infectious_period = 7 +)
influenza <- infection( + name = "influenza", + r0 = 1.5, + preinfectious_period = 3, + infectious_period = 7 +)
The basic reproduction number, \(R_0\), for the SEIR model is:
\[ R_0 = +\frac{\beta}{\gamma}.\]
Therefore, we can rewrite the transmission rate, \(\beta\), as:
\[ \beta = R_0 \gamma.\]
Contact matrices can be estimated from surveys or contact data, or +synthetic ones can be used. We will use the R package +socialmixr to load in a contact matrix estimated from +POLYMOD survey data (Mossong et +al. 2008).
Load contact and population data from socialmixr::polymod
Using the R package socialmixr, run the following lines +of R code to obtain the contact matrix for the United Kingdom for the +year age bins:
socialmixr
polymod <- socialmixr::polymod contact_data <- socialmixr::contact_matrix( - polymod, + survey = polymod, countries = "United Kingdom", age.limits = c(0, 20, 40), symmetric = TRUE ) -contact_data
Using POLYMOD social contact data. To cite this in a publication, use the 'get_citation()' function
+contact.age.group [,1] [,2] [,3] + [0,20) 7.883663 2.794154 1.565665 + [20,40) 3.120220 4.854839 2.624868 + 40+ 3.063895 4.599893 5.005571
The result is a square matrix with rows and columns for each age +group. Contact matrices can be loaded from other sources, but they must +be in the correct format to be used in epidemics.
One of the arguments of the function contact_matrix() is +symmetric=TRUE. This means that the total number of +contacts of age group 1 with age group 2, should be the same as the +total number of contacts of age group 2 and age group 1 (see the +socialmixr vignette +for more detail). However, when contact matrices are estimated from +surveys or other sources, the reported number of contacts may +differ by age group resulting in a non-symmetric contact matrix (Prem et al +2021).
contact_matrix()
symmetric=TRUE
We have prepared our parameter values, contact matrix and demography +vector. Now we must set the initial conditions, prepare the population +and run the model.
The initial conditions are the proportion of individuals in each +disease state \(S\), \(E\), \(I\) +and \(R\) for each age group at time 0. +In this example, we have three age groups age between 0 and 20 years, +age between 20 and 40 years and over. Let’s assume that in the youngest +age category, one in a million individuals are infectious, and the +remaining age categories are infection free.
The initial conditions in the first age category are \(S(0)=1-\frac{1}{1,000,000}\), \(E(0) =0\), \(I(0)=\frac{1}{1,000,000}\), \(R(0)=0\). This is specified as a vector as +follows:
+initial_i <- 1e-6 +initial_conditions_inf <- c( + S = 1 - initial_i, E = 0, I = initial_i, R = 0, V = 0 +)
initial_i <- 1e-6 +initial_conditions_inf <- c( + S = 1 - initial_i, E = 0, I = initial_i, R = 0, V = 0 +)
For the age categories that are free from infection, the initial +conditions are \(S(0)=1\), \(E(0) =0\), \(I(0)=0\), \(R(0)=0\). We specify this as follows,
+initial_conditions_free <- c( + S = 1, E = 0, I = 0, R = 0, V = 0 +)
initial_conditions_free <- c( + S = 1, E = 0, I = 0, R = 0, V = 0 +)
We combine the three initial conditions vectors into one matrix,
+# build for all age groups +initial_conditions <- rbind( + initial_conditions_inf, + initial_conditions_free, + initial_conditions_free +) +rownames(initial_conditions) <- rownames(contact_matrix) +initial_conditions
# build for all age groups +initial_conditions <- rbind( + initial_conditions_inf, + initial_conditions_free, + initial_conditions_free +) +rownames(initial_conditions) <- rownames(contact_matrix) +initial_conditions
S E I R V +[0,20) 0.999999 0 1e-06 0 0 +[20,40) 1.000000 0 0e+00 0 0 +40+ 1.000000 0 0e+00 0 0
To run the model we need the following inputs:
We have already created our infection object influenza. +The population object requires a vector containing the demographic +structure of the population. The demographic vector must be a named +vector containing the number of individuals in each age group of our +given population. In this example, we can extract the demographic +information from the contact_data object that we obtained +using the socialmixr package.
influenza
contact_data
+demography_vector <- contact_data$demography$population +names(demography_vector) <- rownames(contact_matrix) +demography_vector
demography_vector <- contact_data$demography$population +names(demography_vector) <- rownames(contact_matrix) +demography_vector
Removing participants that have contacts without age information. To change this behaviour, set the 'missing.contact.age' option
[0,20) [20,40) 40+ +14799290 16526302 28961159
To create our population object, we call the function +population() specifying a name, the contact matrix, the +demography vector and the initial conditions.
population()
+uk_population <- population( + name = "UK", + contact_matrix = contact_matrix, + demography_vector = demography_vector, + initial_conditions = initial_conditions +)
uk_population <- population( + name = "UK", + contact_matrix = contact_matrix, + demography_vector = demography_vector, + initial_conditions = initial_conditions +)
No we are ready to run our model. We will specify +time_end=600 to run the model for 600 days.
time_end=600
+output <- epidemic_default_cpp( + population = uk_population, + infection = influenza, + time_end = 600 +) +head(output)
output <- epidemic_default_cpp( + population = uk_population, + infection = influenza, + time_end = 600 +) +head(output)
$matrix - contact.age.group - [0,20) [20,40) 40+ - [1,] 7.883663 3.120220 3.063895 - [2,] 2.794154 4.854839 4.599893 - [3,] 1.565665 2.624868 5.005571 - -$demography - age.group population proportion year -1: [0,20) 14799290 0.2454816 2005 -2: [20,40) 16526302 0.2741283 2005 -3: 40+ 28961159 0.4803901 2005 - -$participants - age.group participants proportion -1: [0,20) 404 0.3996044 -2: [20,40) 248 0.2453017 -3: 40+ 359 0.3550940
time demography_group compartment value +1: 0 [0,20) susceptible 14799275 +2: 0 [20,40) susceptible 16526302 +3: 0 40+ susceptible 28961159 +4: 0 [0,20) exposed 0 +5: 0 [20,40) exposed 0 +6: 0 40+ exposed 0
Our model output consists of the number of individuals in each +compartment in each age group through time. We can visualise the +infectious individuals only (those in the \(I\) class) through time.
+ggplot(output[compartment == "infectious", ]) + + geom_line( + aes(time, value, colour = demography_group), + linewidth = 1 + ) + + scale_colour_brewer( + palette = "Dark2", + labels = rownames(contact_matrix), + name = "Age group" + ) + + scale_y_continuous( + labels = scales::comma, + name = "Infectious indivduals" + ) + + labs( + x = "Model time (days)" + ) + + theme_classic() + + theme( + legend.position = "top" + ) + + theme_grey( + base_size = 15 + )
ggplot(output[compartment == "infectious", ]) + + geom_line( + aes(time, value, colour = demography_group), + linewidth = 1 + ) + + scale_colour_brewer( + palette = "Dark2", + labels = rownames(contact_matrix), + name = "Age group" + ) + + scale_y_continuous( + labels = scales::comma, + name = "Infectious indivduals" + ) + + labs( + x = "Model time (days)" + ) + + theme_classic() + + theme( + legend.position = "top" + ) + + theme_grey( + base_size = 15 + )
Note that there is a default argument of increment = 1. +This relates to the time step of the ODE solver. When the parameters and +maximum number of time steps is days, the default increment is one +day.
increment = 1
The choice of increment will depend on the time scale of the +parameters, and the rate at which events can occur. In general, the +increment should smaller than the fastest event. For example, if +parameters are on a monthly time scale, but some events will occur +within a month, then the increment should be less than one month.
As the epidemic model is deterministic, we have one trajectory for +our given parameter values. In practice, we have uncertainty in the +value of our parameters. To account for this, we must run our model for +different parameter combinations.
We ran our model with \(R_0= 1.5\). +However, we believe that \(R_0\) +follows a normal distribution with mean 1.5 and standard deviation 0.05. +To account for uncertainty we will run the model for different values of +\(R_0\). The steps we will follow to do +this are:
+R0_vec <- rnorm(100, 1.5, 0.05)
R0_vec <- rnorm(100, 1.5, 0.05)
+output_samples <- Map( + R0_vec, + seq_along(R0_vec), + f = function(x, i) { + # create infection object for R0 value + influenza <- infection( + name = "influenza", + r0 = x, + preinfectious_period = 3, + infectious_period = 7 + ) + + # run an epidemic model using `epidemic()` + output <- epidemic_default_cpp( + population = uk_population, + infection = influenza, + time_end = 600, increment = 1.0 + ) + + # extract infectious individuals + output <- output[compartment == "infectious"] + + # assign scenario number + output[, c("scenario", "R") := list(i, x)] + + output + } +) + +# combine to prepare for plotting +output_samples <- bind_rows(output_samples)
output_samples <- Map( + R0_vec, + seq_along(R0_vec), + f = function(x, i) { + # create infection object for R0 value + influenza <- infection( + name = "influenza", + r0 = x, + preinfectious_period = 3, + infectious_period = 7 + ) + + # run an epidemic model using `epidemic()` + output <- epidemic_default_cpp( + population = uk_population, + infection = influenza, + time_end = 600, increment = 1.0 + ) + + # extract infectious individuals + output <- output[compartment == "infectious"] + + # assign scenario number + output[, c("scenario", "R") := list(i, x)] + + output + } +) + +# combine to prepare for plotting +output_samples <- bind_rows(output_samples)
+ggplot(output_samples ,aes(time, value)) + + stat_summary(geom = "line", fun = mean) + + stat_summary(geom = "ribbon", + fun.min = function(z) { quantile(z, 0.025) }, + fun.max = function(z) { quantile(z, 0.975) }, + alpha = 0.3) + + facet_grid( + cols = vars(demography_group) + ) + + theme_grey( + base_size = 15 + )
ggplot(output_samples ,aes(time, value)) + + stat_summary(geom = "line", fun = mean) + + stat_summary(geom = "ribbon", + fun.min = function(z) { quantile(z, 0.025) }, + fun.max = function(z) { quantile(z, 0.975) }, + alpha = 0.3) + + facet_grid( + cols = vars(demography_group) + ) + + theme_grey( + base_size = 15 + )
Deciding which parameters to include uncertainty in depends on a few +factors: how well informed a parameter value is e.g. consistency of +estimates from the literature; how sensitive model outputs are to +parameter value changes; and the purpose of the modelling task.
In this tutorial, we have learnt how to generate disease trajectories +using a mathematical model. Once a model has been chosen, the parameters +and other inputs must be specified in the correct way to perform model +simulations. In the next tutorial, we will consider how to choose the +right model for different tasks.
Content from Choosing an appropriate model
Last updated on 2023-10-31 | + + Edit this page
Using mathematical models in outbreak analysis does not necessarily +require developing a new model. There are existing models for different +infections, interventions and transmission patterns which can be used to +answer new questions. In this tutorial, we will learn how to choose an +existing model to generate predictions for a given scenario.
When deciding whether an existing model can be used, we must consider +:
A model may already exist for your study disease, or there may be a +model for an infection that has the same transmission pathways and +epidemiological features that can be used.
Model structures differ for whether the disease has pandemic +potential or not. When predicted numbers of infection are small, +stochastic variation in output can have an effect on whether an outbreak +takes off or not. Outbreaks are usually smaller in magnitude than +epidemics, so its often appropriate to use a stochastic model to +characterise the uncertainty in the early stages of the outbreak. +Epidemics are larger in magnitude than outbreaks and so a deterministic +model is suitable as we have less interest in the stochastic variation +in output.
The outcome of interest can be a feature of a mathematical model. It +may be that you are interested in the predicted numbers of infection +through time, or in a specific outcome such as hospitalisations or cases +of severe disease.
Finally, interventions such as vaccination may be of interest. A +model may or may not have the capability to include the impact of +different interventions on different time scales (continuous time or at +discrete time points). We will discuss interventions in detail in the +next tutorial.
The R package epidemics contains functions to run +existing models. For details on the models that are available, see the +package vignettes. +To learn how to run the models in R, read the documentation using +?epidemics::epidemic_ebola. Remember to use the ‘Check +model equation’ questions to help your understanding of an existing +model.
?epidemics::epidemic_ebola
Add additional maths (or epi) content for novice learners
Load contact and population data for Poland from socialmixr::polymod -using the following age bins:
You have been asked to explore the variation in numbers of infected +individuals in the early stages of an Ebola outbreak.
Which of the following models would be an appropriate choice for this +task:
epidemic_default
epidemic_ebola
Consider the following questions:
-polymod <- socialmixr::polymod -contact_data <- socialmixr::contact_matrix( - polymod, - countries = "Poland", - age.limits = c(0, 15, 50), - symmetric = TRUE -)
polymod <- socialmixr::polymod -contact_data <- socialmixr::contact_matrix( - polymod, - countries = "Poland", - age.limits = c(0, 15, 50), - symmetric = TRUE -)
A deterministic SEIR model with age specific direct transmission.
The model is capable of predicting an Ebola type outbreak, but as the +model is deterministic, we are not able to explore stochastic variation +in the early stages of the outbreak.
A stochastic SEIHFR (Susceptible, Exposed, Infectious, Hospitalised, +Funeral, Removed) model that was developed specifically for infection +with Ebola.
As this model is stochastic, it is the most appropriate choice to +explore how variation in numbers of infected individuals in the early +stages of an Ebola outbreak.
Content from Modelling interventions
Outbreak response : Intervention +types.
Mathematical models can be used to generate predictions for the +implementation of non-pharmaceutical and pharmaceutical interventions at +different stages of an outbreak. In this tutorial, we will introduce how +to include different interventions in models.
Non-pharmaceutical interventions (NPIs) are measures put in place to +reduce transmission that do not include taking medicine or vaccines. +NPIs aim reduce contact between infectious and susceptible individuals. +For example, washing hands, wearing masks and closures of school and +workplaces.
In mathematical modelling, we must make assumptions about how NPIs +will affect transmission. This may include adding additional disease +states or reducing the value of relevant parameters.
We want to investigate the effect of school closures on reducing the +number of individuals infectious with COVID-19 through time. We assume +that a school closure will reduce the frequency of contacts within and +between different age groups.
Using an SEIR model (epidemic_default() in the R package +{epidemics}) we set \(R_0 = +2.7\), preinfectious period \(= +4\) and the infectious_period \(= +5.5\) (parameters adapted from Davies et +al. (2020)). We load a contact matrix with age bins 0-18, 18-65, 65 +years and older using socialmixr and assume that one in +every 1 million in each age group is infectious at the start of the +epidemic.
We will assume that school closures will reduce the contacts between +school aged children (aged 0-15) by 0.5, and will cause a small +reduction (0.01) in the contacts between adults (aged 15 and over).
The contact matrix is scaled down by proportions for the period in +which the intervention is in place. To explain the reduction, consider a +contact matrix for two age groups with equal number of contacts:
[,1] [,2] +[1,] 1 1 +[2,] 1 1
If the reduction is 50% in group 1 and 10% in group 2, the contact +matrix during the intervention will be:
[,1] [,2] +[1,] 0.25 0.45 +[2,] 0.45 0.81
The contacts within group 1 are reduced by 50% twice to accommodate +for a 50% reduction in outgoing and incoming contacts (\(1\times 0.5 \times 0.5 = 0.25\)). +Similarly, the contacts within group 2 are reduced by 10% twice. The +contacts between group 1 and group 2 are reduced by 50% and then by 10% +(\(1 \times 0.5 \times 0.9= +0.45\)).
To include an intervention in our model we must create an +intervention object. The inputs are the name of the +intervention (name), the type of intervention +(contacts or rate), the start time +(time_begin), the end time (time_end) and the +reduction (reduction). The values of the reduction matrix +are specified in the same order as the age groups in the contact +matrix.
intervention
name
contacts
rate
time_begin
time_end
reduction
-contact_data
rownames(contact_matrix)
$matrix - contact.age.group - [0,15) [15,50) 50+ - [1,] 8.4882943 5.047958 1.765826 - [2,] 1.5998125 14.049041 3.624419 - [3,] 0.9468591 6.132289 4.314050 - -$demography - age.group population proportion year -1: [0,15) 6372248 0.1661054 2005 -2: [15,50) 20106632 0.5241197 2005 -3: 50+ 11883794 0.3097749 2005 - -$participants - age.group participants proportion -1: [0,15) 299 0.2960396 -2: [15,50) 469 0.4643564 -3: 50+ 242 0.2396040
[1] "[0,15)" "[15,65)" "65+"
Therefore, we specify +reduction = matrix(c(0.5, 0.01, 0.01)). We assume that the +school closures start on day 50 and are in place for a further 100 days. +Therefore our intervention object is :
reduction = matrix(c(0.5, 0.01, 0.01))
+close_schools <- intervention( + name = "School closure", + type = "contacts", + time_begin = 50, + time_end = 50 + 100, + reduction = matrix(c(0.5, 0.01, 0.01)) +)
close_schools <- intervention( + name = "School closure", + type = "contacts", + time_begin = 50, + time_end = 50 + 100, + reduction = matrix(c(0.5, 0.01, 0.01)) +)
To run the model with an intervention we set +intervention = list(contacts = close_schools) as +follows:
intervention = list(contacts = close_schools)
+output_school <- epidemic_default_cpp( + population = uk_population, + infection = covid, + intervention = list(contacts = close_schools), + time_end = 300, increment = 1.0 +)
output_school <- epidemic_default_cpp( + population = uk_population, + infection = covid, + intervention = list(contacts = close_schools), + time_end = 300, increment = 1.0 +)
We see that with the intervention (solid line) in place, the +infection still spreads through the population, though the epidemic peak +is smaller than the baseline with no intervention in place (dashed +line).
We can model the effect of other NPIs as reducing the value of +relevant parameters. For example, we want to investigate the effect of +mask wearing on the number of individuals infectious with COVID-19 +through time.
We expect that mask wearing will reduce an individual’s +infectiousness. As we are using a population based model, we cannot make +changes to individual behaviour and so assume that the transmission rate +\(\beta\) is reduced by a proportion +due to mask wearing in the population. We specify this proportion, \(\theta\) as product of the proportion +wearing masks multiplied by the proportion reduction in transmissibility +(adapted from Li +et al. 2020)
We create an intervention object with type = rate and +reduction = 0.161. Using parameters adapted from Li et al. 2020 +we have proportion wearing masks = coverage \(\times\) availability = \(0.54 \times 0.525 = 0.2835\), proportion +reduction in transmissibility = \(0.575\). Therefore, \(\theta = 0.2835 \times 0.575 = 0.163\). We +assume that the mask wearing mandate starts at day 40 and is in place +for 200 days.
type = rate
reduction = 0.161
+mask_mandate <- intervention( + name = "mask mandate", + type = "rate", + time_begin = 40, + time_end = 40 + 200, + reduction = 0.163 +)
mask_mandate <- intervention( + name = "mask mandate", + type = "rate", + time_begin = 40, + time_end = 40 + 200, + reduction = 0.163 +)
To implement this intervention on the parameter \(\beta\), we specify +intervention = list(beta = mask_mandate).
intervention = list(beta = mask_mandate)
+output_masks <- epidemic_default_cpp( + population = uk_population, + infection = covid, + intervention = list(beta = mask_mandate), + time_end = 300, increment = 1.0 +)
output_masks <- epidemic_default_cpp( + population = uk_population, + infection = covid, + intervention = list(beta = mask_mandate), + time_end = 300, increment = 1.0 +)
Models can be used to investigate the effect of pharmaceutical +interventions, such as vaccination. In this case, it is useful to add +another disease state to track the number of vaccinated individuals +through time. The diagram below shows an SEIRV model where susceptible +individuals are vaccinated and then move to the \(V\) class.
The equations describing this model are as follows:
\[ +\begin{aligned} +\frac{dS_i}{dt} & = - \beta S_i \sum_j C_{i,j} I_j -\nu_{t} S_i \\ +\frac{dE_i}{dt} &= \beta S_i\sum_j C_{i,j} I_j - \alpha E_i \\ +\frac{dI_i}{dt} &= \alpha E_i - \gamma I_i \\ +\frac{dR_i}{dt} &=\gamma I_i \\ +\frac{dV_i}{dt} & =\nu_{i,t} S_i\\ +\end{aligned} +\] Individuals are vaccinated at an age group (\(i\)) specific time dependent (\(t\)) vaccination rate (\(\nu\)). The SEIR components of these +equations are described in the tutorial Simulating transmission.
To explore the effect of vaccination we need to create a vaccination +object. As vaccination is age group specific, we must pass an age groups +specific vaccination rate \(\nu\) and +age group specific start and end times of the vaccination program. Here +we will assume all age groups are vaccinated at the same rate and that +the vaccination program starts on day 40 and is in place for 150 +days.
+# prepare a vaccination object +vaccinate <- vaccination( + name = "vaccinate all", + time_begin = matrix(40, nrow(contact_matrix)), + time_end = matrix(40 + 150, nrow(contact_matrix)), + nu = matrix(c(0.01, 0.01, 0.01)) +)
# prepare a vaccination object +vaccinate <- vaccination( + name = "vaccinate all", + time_begin = matrix(40, nrow(contact_matrix)), + time_end = matrix(40 + 150, nrow(contact_matrix)), + nu = matrix(c(0.01, 0.01, 0.01)) +)
We pass our vaccination object using +vaccination = vaccinate:
vaccination = vaccinate
+output_vaccinate <- epidemic_default_cpp( + population = uk_population, + infection = covid, + vaccination = vaccinate, + time_end = 300, increment = 1.0 +)
output_vaccinate <- epidemic_default_cpp( + population = uk_population, + infection = covid, + vaccination = vaccinate, + time_end = 300, increment = 1.0 +)
Here we see that the total number of infectious individuals when +vaccination is in place is much lower compared to school closures and +mask wearing interventions.
Modelling interventions requires assumptions of how interventions +affect model parameters such as contact matrices or parameter values. +Next we want quantify the effect of an interventions. In the next +tutorial, we will learn how to compare intervention scenarios against +each other.
Content from Comparing public health outcomes of interventions
In this tutorial we will compare intervention scenarios against each +other. To quantify the effect of the intervention we need to compare our +intervention scenario to a counter factual scenario. The counter +factual is the scenario in which nothing changes, often referred to +as the ‘do nothing’ scenario. The counter factual scenario may include +no interventions, or if we are investigating the potential impact of an +additional intervention in the later stages of an outbreak there may be +existing interventions in place.
We must also decide what our outcome of interest is to make +comparisons between intervention and counter factual scenarios. The +outcome of interest can be:
The Vacamole model is a deterministic model based on a system of ODEs +in Ainslie +et al. 2022. The model consists of 11 compartments, individuals are +classed as one of the following:
The diagram below describes the flow of individuals through the +different compartments.
See ?epidemics::epidemic_vacamole for detail on how to +run the model.
?epidemics::epidemic_vacamole
Coming soon
Estimated time 15 minutes
Estimated time 75 minutes
Inline instructor notes can help inform instructors of timing -challenges associated with the lessons. They appear in the “Instructor -View”.
READ THESE LINES AND ERASE:
The Workbench-related sections that the developer must keep are:
By the end of this tutorial, learners should be able to replicate the +above image on their own computers.
The Epiverse-TRACE sections that we encourage to keep are:
Take a look to the Contributing.md file for more writing -guidelines.
Estimated time 30 minutes
The focus of this tutorial is understanding existing models to decide +if they are appropriate for a defined question.
In this tutorial different types of intervention and how they can be +modelled are introduced. Learners should be able to understand the +underlying mechanism of these interventions (e.g. reduce contact rate) +as well as how to implement the code to include such interventions.
In this tutorial we introduce the concept of the counter factual and +how to compare scenarios (counter factual versus intervention) against +each other.
+rownames(contact_matrix)
+polymod <- socialmixr::polymod +contact_data <- socialmixr::contact_matrix( + survey = polymod, + countries = "United Kingdom", + age.limits = c(0, 20, 40), + symmetric = TRUE +) +# prepare contact matrix +contact_matrix <- t(contact_data$matrix) +contact_matrix
polymod <- socialmixr::polymod +contact_data <- socialmixr::contact_matrix( + survey = polymod, + countries = "United Kingdom", + age.limits = c(0, 20, 40), + symmetric = TRUE +) +# prepare contact matrix +contact_matrix <- t(contact_data$matrix) +contact_matrix