Vignette on parameter uncertainty

CITATION.cff

@@ -361,6 +361,8 @@ references:
     email: sebastian.funk@lshtm.ac.uk
   - family-names: Willem
     given-names: Lander
+  - family-names: Gruson
+    given-names: Hugo
   year: '2023'
 - type: software
   title: spelling
@@ -413,6 +415,42 @@ references:
   - family-names: Chang
     given-names: Winston
   year: '2023'
+- type: software
+  title: dplyr
+  abstract: 'dplyr: A Grammar of Data Manipulation'
+  notes: Suggests
+  url: https://dplyr.tidyverse.org
+  repository: https://CRAN.R-project.org/package=dplyr
+  authors:
+  - family-names: Wickham
+    given-names: Hadley
+    email: hadley@posit.co
+    orcid: https://orcid.org/0000-0003-4757-117X
+  - family-names: François
+    given-names: Romain
+    orcid: https://orcid.org/0000-0002-2444-4226
+  - family-names: Henry
+    given-names: Lionel
+  - family-names: Müller
+    given-names: Kirill
+    orcid: https://orcid.org/0000-0002-1416-3412
+  - family-names: Vaughan
+    given-names: Davis
+    email: davis@posit.co
+    orcid: https://orcid.org/0000-0003-4777-038X
+  year: '2023'
+- type: software
+  title: EpiEstim
+  abstract: 'EpiEstim: Estimate Time Varying Reproduction Numbers from Epidemic Curves'
+  notes: Suggests
+  url: https://github.com/mrc-ide/EpiEstim
+  repository: https://CRAN.R-project.org/package=EpiEstim
+  authors:
+  - family-names: Cori
+    given-names: Anne
+    email: a.cori@imperial.ac.uk
+    orcid: https://orcid.org/0000-0002-8443-9162
+  year: '2023'
 - type: software
   title: BH
   abstract: 'BH: Boost C++ Header Files'

DESCRIPTION

@@ -24,6 +24,13 @@ Authors@R:
             email = "edwin.vanleeuwen@ukhsa.gov.uk",
             comment = c(ORCID = "https://orcid.org/0000-0002-2383-5305")
         ),
+        person(
+            given = "Adam",
+            family = "Kucharski",
+            role = c("ctb", "rev"),
+            email = "adam.kucharski@lshtm.ac.uk",
+            comment = c(ORCID = "https://orcid.org/0000-0001-8814-9421")
+        ),
         person(
             given = "Tim",
             family = "Taylor",
@@ -74,7 +81,9 @@ Suggests:
     socialmixr,
     spelling,
     testthat (>= 3.0.0),
-    withr
+    withr,
+    dplyr,
+    EpiEstim
 Config/testthat/edition: 3
 VignetteBuilder: knitr
 LinkingTo: 

_pkgdown.yml

@@ -12,6 +12,7 @@ articles:
   - multiple_interventions
   - rate_interventions
   - time_dependence
+  - parameter_uncertainty
 - title: Guides to library models
   navbar: Guides to library models
   contents:

inst/WORDLIST

@@ -14,6 +14,7 @@ ETUs
 EVD
 Ebolavirus
 Eigen
+EpiEstim
 Epirecipes
 Epiverse
 Eurosurveillance
@@ -75,6 +76,7 @@ doi
 ebola
 ebolavirus
 epicookbook
+epiparameter
 epirecipes
 etc
 ggplot
@@ -88,7 +90,6 @@ org
 packagename
 pkg
 prepper
-rspb
 stochasticity
 susceptibles
 svg

vignettes/parameter_uncertainty.Rmd

@@ -0,0 +1,245 @@
+---
+title: "Generating and modelling R uncertainty"
+output:
+  bookdown::html_vignette2:
+    fig_caption: yes
+    code_folding: show
+pkgdown:
+  as_is: true
+bibliography: references.json
+link-citations: true
+vignette: >
+  %\VignetteIndexEntry{Generating and modelling R uncertainty}
+  %\VignetteEngine{knitr::rmarkdown}
+  %\VignetteEncoding{UTF-8}
+---
+
+```{r, include = FALSE}
+knitr::opts_chunk$set(
+  collapse = TRUE,
+  comment = "#>",
+  message = FALSE,
+  warning = FALSE,
+  fig.width = 5,
+  fig.height = 4,
+  dpi = 300
+)
+```
+
+Uncertainty in the characteristics of an epidemic is a key element and roadblock in epidemic response [@shea2020].
+This vignette aims to show how to model parameter uncertainty using _epidemics_, as well as how it might affect conclusions regarding the implementation of response measures such as non-pharmaceutical interventions.
+
+::: {.alert .alert-warning}
+**New to _epidemics_?** It may help to read the ["Get started"](epidemics.html) vignette first!
+:::
+
+```{r}
+# load epidemics
+library(epidemics)
+
+# Load CRAN packages
+library(EpiEstim)
+library(dplyr)
+library(ggplot2)
+```
+
+## Prepare $R$ estimates
+
+For this example, we consider influenza with pandemic potential [@ghani2010], and prepare multiple samples of the estimated $R$.
+
+We obtain the probability distribution function (PDF) from the distribution of the serial intervals; this is a Gamma distribution with shape $k$ = 2.622 and scale $\theta$ = 0.957 [@ghani2010].
+
+::: {.alert .alert-info}
+The [forthcoming Epiverse package _epiparameter_](https://epiverse-trace.github.io/epiparameter/) is expected to make it substantially easier to access and use epidemiological parameters, such as the serial interval, reported in the literature, making it easier to model scenarios differing in the intrinsic characteristics of the pathogen causing the outbreak.
+:::
+
+We use this PDF to estimate the $R$ of the 2009 influenza pandemic in the U.K., using [the _EpiEstim_ package](https://cran.r-project.org/package=EpiEstim).
+We use the $R$ estimate (the mean and standard deviation) from _EpiEstim_ to generate 100 samples of $R$, assuming that it is normally distributed.
+
+```{r fig.width=6, fig.height=3}
+# Get 2009 influenza data for school in Pennsylvania
+data(Flu2009)
+flu_early_data <- filter(Flu2009$incidence, dates < "2009-05-10")
+
+# get the PDF of the distribution of serial intervals
+serial_pdf <- dgamma(seq(0, 25), shape = 2.622, scale = 0.957)
+serial_pdf <- serial_pdf / sum(serial_pdf) # check sum to 1
+
+# Use EpiEstim to estimate R with uncertainty
+# Uses Gamma distribution by default
+output_R <- estimate_R(
+  incid = flu_early_data,
+  method = "non_parametric_si",
+  config = make_config(list(si_distr = serial_pdf))
+)
+
+# Plot output to visualise
+plot(output_R, "R")
+```
+
+Generate 100 samples of $R$.
+
+```{r}
+# Generate 100 R samples
+set.seed(1)
+r_samples <- rnorm(
+  n = 100, mean = output_R$R$`Mean(R)`, sd = output_R$R$`Std(R)`
+)
+```
+
+## Prepare epidemic model inputs
+
+We define inputs for the default, SEIR-V model in _epidemics_ in the form of `<population>` and `<infection>` objects.
+We omit non-pharmaceutical interventions and vaccinations for this example.
+
+Since this step is identical to that in previous vignettes, the code is folded, but can be expanded for clarity.
+
+```{r class.source = 'fold-hide'}
+# load contact and population data from socialmixr::polymod
+polymod <- socialmixr::polymod
+contact_data <- socialmixr::contact_matrix(
+  polymod,
+  countries = "United Kingdom",
+  age.limits = c(0, 20, 40),
+  symmetric = TRUE
+)
+
+# prepare contact matrix
+contact_matrix <- t(contact_data$matrix)
+
+# prepare the demography vector
+demography_vector <- contact_data$demography$population
+names(demography_vector) <- rownames(contact_matrix)
+
+# initial conditions: one in every 1 million is infected
+initial_i <- 1e-6
+initial_conditions <- c(
+  S = 1 - initial_i, E = 0, I = initial_i, R = 0, V = 0
+)
+
+# build for all age groups
+initial_conditions <- rbind(
+  initial_conditions,
+  initial_conditions,
+  initial_conditions
+)
+rownames(initial_conditions) <- rownames(contact_matrix)
+
+# prepare the population to model as affected by the epidemic
+uk_population <- population(
+  name = "UK",
+  contact_matrix = contact_matrix,
+  demography_vector = demography_vector,
+  initial_conditions = initial_conditions
+)
+
+# an intervention to close schools
+close_schools <- intervention(
+  type = "contacts",
+  time_begin = 200,
+  time_end = 260,
+  reduction = matrix(c(0.5, 0.01, 0.01))
+)
+```
+
+## Run with $R$ uncertainty
+
+Loop over sampled $R$ values to generate epidemic model outputs with uncertainty.
+
+```{r}
+# iterate over the R samples, running the epidemic model for each value of R
+output_samples <- Map(
+  r_samples,
+  seq_along(r_samples),
+  f = function(x, i) {
+    # Prepare an <infection> class object to store the parameters of
+    # the infection which is causing the epidemic which is being modelled.
+
+    # simulate a pandemic, with an R0,
+    # an infectious period, and an pre-infectious period
+    pandemic_influenza <- infection(
+      r0 = x,
+      preinfectious_period = 3,
+      infectious_period = 7
+    )
+
+    # run an epidemic model using `epidemic()`
+    output <- epidemic_default_cpp(
+      population = uk_population,
+      infection = pandemic_influenza,
+      intervention = list(contacts = close_schools),
+      time_end = 600, increment = 1.0
+    )
+
+    # calculate new infections
+    output <- new_infections(output)
+
+    # assign scenario number
+    output[, c("scenario", "R") := list(i, x)]
+
+    output
+  }
+)
+```
+
+We view the output of a few scenario runs.
+
+```{r}
+# View the firt few rows of the first two runs
+lapply(output_samples[1:2], head)
+
+# combine to prepare for plotting
+output_samples <- bind_rows(output_samples)
+```
+
+
+We obtain the epidemic runs and plot the number of new infections in each scenario; the plotting code is folded for brevity.
+
+```{r class.source = 'fold-hide', fig.cap="Modelling scenarios of the spread of influenza with pandemic potential in the U.K., with uncertainty in the reproduction number $R$. Each scenario represents one of 1,000 runs with different values of $R$ estimated from the case data for the 2009 influenza pandemic in the U.K. Scenarios are coloured by $R$, which is divided into four bins for visual clarity. Values of $R <$ 1.0 do not lead to large epidemics. The duration of a non-pharmaceutical intervention to close schools for two months (affecting the social contacts of all individuals under age 20) is also shown. For example, in high $R$ scenarios (green), school closures have little effect, as the epidemic has peaked before the intervention can take effect; this demonstrates how uncertainty in epidemic characterists may influence the effectiveness of response measures.", fig.width=6}
+ggplot(output_samples) +
+  geom_vline(
+    xintercept = c(200, 260),
+    colour = alpha("black", 0.5),
+    linetype = "dashed"
+  ) +
+  annotate(
+    geom = "text",
+    x = 230, y = 600e3,
+    label = "Schools closed",
+    colour = alpha("black", 0.5),
+    angle = 90,
+    hjust = "inward"
+  ) +
+  geom_line(
+    aes(
+      time, new_infections,
+      group = scenario,
+      colour = R
+    ),
+    alpha = 0.3
+  ) +
+  scale_colour_fermenter(
+    palette = "Dark2",
+    breaks = c(0, 1, 1.5, 2.0, 3.0),
+    limits = c(0, 3)
+  ) +
+  scale_y_continuous(
+    name = "New infections",
+    labels = scales::label_comma(scale = 1e-3, suffix = "K")
+  ) +
+  labs(
+    x = "Time (days since start of epidemic)"
+  ) +
+  facet_grid(
+    cols = vars(demography_group)
+  ) +
+  theme_bw() +
+  theme(
+    legend.position = "top",
+    legend.key.height = unit(2, "mm")
+  )
+```
+
+The [forthcoming Epiverse package _scenarios_](https://epiverse-trace.github.io/scenarios/) is expected to include functionality that wraps the iteration step shown in this vignette, to make it easier to model parameter uncertainty.
+
+## References

vignettes/references.json

@@ -4,5 +4,7 @@
   {"id":"bjornstad2020","accessed":{"date-parts":[[2023,3,14]]},"author":[{"family":"Bjørnstad","given":"Ottar N."},{"family":"Shea","given":"Katriona"},{"family":"Krzywinski","given":"Martin"},{"family":"Altman","given":"Naomi"}],"citation-key":"bjornstad2020","container-title":"Nature Methods","DOI":"10.1038/s41592-020-0856-2","ISSN":"1548-7105","issue":"6","issued":{"date-parts":[[2020,6,1]]},"language":"en","license":"2020 Springer Nature America, Inc.","number":"6","page":"557-558","publisher":"Nature Publishing Group","source":"www.nature.com","title":"The SEIRS model for infectious disease dynamics","type":"article-journal","URL":"https://www.nature.com/articles/s41592-020-0856-2","volume":"17"},
   {"id":"bjornstad2020a","accessed":{"date-parts":[[2023,3,14]]},"author":[{"family":"Bjørnstad","given":"Ottar N."},{"family":"Shea","given":"Katriona"},{"family":"Krzywinski","given":"Martin"},{"family":"Altman","given":"Naomi"}],"citation-key":"bjornstad2020a","container-title":"Nature Methods","DOI":"10.1038/s41592-020-0822-z","ISSN":"1548-7105","issue":"5","issued":{"date-parts":[[2020,5,1]]},"language":"en","license":"2020 Springer Nature America, Inc.","number":"5","page":"455-456","publisher":"Nature Publishing Group","source":"www.nature.com","title":"Modeling infectious epidemics","type":"article-journal","URL":"https://www.nature.com/articles/s41592-020-0822-z","volume":"17"},
   {"id":"getz2018","accessed":{"date-parts":[[2023,7,13]]},"author":[{"family":"Getz","given":"Wayne M."},{"family":"Dougherty","given":"Eric R."}],"citation-key":"getz2018","container-title":"Journal of Biological Dynamics","DOI":"10.1080/17513758.2017.1401677","ISSN":"1751-3758","issue":"1","issued":{"date-parts":[[2018,1,1]]},"page":"16-38","PMID":"29157162","publisher":"Taylor & Francis","source":"Taylor and Francis+NEJM","title":"Discrete stochastic analogs of Erlang epidemic models","type":"article-journal","URL":"https://doi.org/10.1080/17513758.2017.1401677","volume":"12"},
-  {"id":"li2019","accessed":{"date-parts":[[2023,6,14]]},"author":[{"family":"Li","given":"Shou-Li"},{"family":"Ferrari","given":"Matthew J."},{"family":"Bjørnstad","given":"Ottar N."},{"family":"Runge","given":"Michael C."},{"family":"Fonnesbeck","given":"Christopher J."},{"family":"Tildesley","given":"Michael J."},{"family":"Pannell","given":"David"},{"family":"Shea","given":"Katriona"}],"citation-key":"li2019","container-title":"Proceedings of the Royal Society B: Biological Sciences","DOI":"10.1098/rspb.2019.0774","issue":"1905","issued":{"date-parts":[[2019,6,19]]},"page":"20190774","publisher":"Royal Society","source":"royalsocietypublishing.org (Atypon)","title":"Concurrent assessment of epidemiological and operational uncertainties for optimal outbreak control: Ebola as a case study","title-short":"Concurrent assessment of epidemiological and operational uncertainties for optimal outbreak control","type":"article-journal","URL":"https://royalsocietypublishing.org/doi/full/10.1098/rspb.2019.0774","volume":"286"}
+  {"id":"ghani2010","accessed":{"date-parts":[[2023,10,20]]},"author":[{"family":"Ghani","given":"Azra"},{"family":"Baguelin","given":"Marc"},{"family":"Griffin","given":"Jamie"},{"family":"Flasche","given":"Stefan"},{"family":"Van Hoek","given":"Albert Jan"},{"family":"Cauchemez","given":"Simon"},{"family":"Donnelly","given":"Christl"},{"family":"Robertson","given":"Chris"},{"family":"White","given":"Michael"},{"family":"Truscott","given":"James"},{"family":"Fraser","given":"Christophe"},{"family":"Garske","given":"Tini"},{"family":"White","given":"Peter"},{"family":"Leach","given":"Steve"},{"family":"Hall","given":"Ian"},{"family":"Jenkins","given":"Helen"},{"family":"Ferguson","given":"Neil"},{"family":"Cooper","given":"Ben"}],"citation-key":"ghani2010","container-title":"PLoS Currents","container-title-short":"PLoS Curr","DOI":"10.1371/currents.RRN1130","ISSN":"2157-3999","issued":{"date-parts":[[2010,6,13]]},"page":"RRN1130","source":"DOI.org (Crossref)","title":"The Early Transmission Dynamics of H1N1pdm Influenza in the United Kingdom","type":"article-journal","URL":"https://currents.plos.org/influenza/index.html%3Fp=1653.html","volume":"1"},
+  {"id":"li2019","accessed":{"date-parts":[[2023,6,14]]},"author":[{"family":"Li","given":"Shou-Li"},{"family":"Ferrari","given":"Matthew J."},{"family":"Bjørnstad","given":"Ottar N."},{"family":"Runge","given":"Michael C."},{"family":"Fonnesbeck","given":"Christopher J."},{"family":"Tildesley","given":"Michael J."},{"family":"Pannell","given":"David"},{"family":"Shea","given":"Katriona"}],"citation-key":"li2019","container-title":"Proceedings of the Royal Society B: Biological Sciences","DOI":"10.1098/rspb.2019.0774","issue":"1905","issued":{"date-parts":[[2019,6,19]]},"page":"20190774","publisher":"Royal Society","source":"royalsocietypublishing.org (Atypon)","title":"Concurrent assessment of epidemiological and operational uncertainties for optimal outbreak control: Ebola as a case study","title-short":"Concurrent assessment of epidemiological and operational uncertainties for optimal outbreak control","type":"article-journal","URL":"https://royalsocietypublishing.org/doi/full/10.1098/rspb.2019.0774","volume":"286"},
+  {"id":"shea2020","accessed":{"date-parts":[[2023,3,14]]},"author":[{"family":"Shea","given":"Katriona"},{"family":"Bjørnstad","given":"Ottar N."},{"family":"Krzywinski","given":"Martin"},{"family":"Altman","given":"Naomi"}],"citation-key":"shea2020","container-title":"Nature Methods","DOI":"10.1038/s41592-020-0943-4","ISSN":"1548-7105","issue":"9","issued":{"date-parts":[[2020,9,1]]},"language":"en","license":"2020 Springer Nature America, Inc.","number":"9","page":"867-868","publisher":"Nature Publishing Group","source":"www.nature.com","title":"Uncertainty and the management of epidemics","type":"article-journal","URL":"https://www.nature.com/articles/s41592-020-0943-4","volume":"17"}
 ]