Simulate Cases by Date of Infection, Onset and Report

simulate_cases(
  rts,
  initial_cases,
  initial_date,
  generation_interval,
  rdist = rpois,
  delay_def,
  incubation_def,
  reporting_effect,
  reporting_model,
  truncate_future = TRUE,
  type = "sample"
)

Arguments

rts	A dataframe of containing two variables `rt` and `date` with `rt` being numeric and `date` being a date.
initial_cases	Integer, initial number of cases.
initial_date	Date, (i.e `as.Date("2020-02-01")`). Starting date of the simulation.
generation_interval	Numeric vector describing the generation interval probability density
rdist	A function to be used to sample the number of cases. Must take two arguments with the first specfying the number of samples and the second the mean. Defaults to `rpois` if not supplied
delay_def	A single row data.table that defines the delay distribution (model, parameters and maximum delay for each model). See `lognorm_dist_def` for an example of the structure.
incubation_def	A single row data.table that defines the incubation distribution (model, parameters and maximum delay for each model). See `lognorm_dist_def` for an example of the structure.
reporting_effect	A numeric vector of length 7 that allows the scaling of reported cases by the day on which they report (1 = Monday, 7 = Sunday). By default no scaling occurs.
reporting_model	A function that takes a single numeric vector as an argument and returns a single numeric vector. Can be used to apply stochastic reporting effects. See the examples for details.
truncate_future	Logical, should cases be truncted if they occur after the first date reported in the data. Defaults to `TRUE`.
type	Character string indicating the method to use to transfrom counts. Supports either "sample" which approximates sampling or "median" would shift by the median of the distribution.

Value

A dataframe containing three variables: date, cases and reference.

Examples

## Define an initial rt vector 
rts <- c(rep(2, 20), (2 - 1:15 * 0.1), rep(0.5, 10))
rts
#>  [1] 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0
#> [20] 2.0 1.9 1.8 1.7 1.6 1.5 1.4 1.3 1.2 1.1 1.0 0.9 0.8 0.7 0.6 0.5 0.5 0.5 0.5
#> [39] 0.5 0.5 0.5 0.5 0.5 0.5 0.5
## Use the mean default generation interval for covid
generation_interval <- rowMeans(EpiNow::covid_generation_times)

## Sample a report delay as a lognormal
delay_def <- EpiNow::lognorm_dist_def(mean = 5, mean_sd = 1,
                                      sd = 3, sd_sd = 1, max_value = 30,
                                      samples = 1, to_log = TRUE)
                                      

## Sample a incubation period (again using the default for covid)
incubation_def <- EpiNow::lognorm_dist_def(mean = EpiNow::covid_incubation_period[1, ]$mean,
                                          mean_sd = EpiNow::covid_incubation_period[1, ]$mean_sd,
                                          sd = EpiNow::covid_incubation_period[1, ]$sd,
                                          sd_sd = EpiNow::covid_incubation_period[1, ]$sd_sd,
                                          max_value = 30, samples = 1)

## Simulate cases with a decrease in reporting at weekends and an increase on Monday                                     
simulated_cases <- simulate_cases(rts, initial_cases = 100 , initial_date = as.Date("2020-03-01"),
                    generation_interval = generation_interval, delay_def = delay_def, 
                   incubation_def = incubation_def, 
                   reporting_effect = c(1.1, rep(1, 4), 0.95, 0.95))
                   
print(simulated_cases)
#>            date cases reference
#>   1: 2020-03-02    37 infection
#>   2: 2020-03-03    53 infection
#>   3: 2020-03-04    75 infection
#>   4: 2020-03-05    76 infection
#>   5: 2020-03-06    86 infection
#>  ---                           
#> 126: 2020-04-10  2085    report
#> 127: 2020-04-11  1786    report
#> 128: 2020-04-12  1632    report
#> 129: 2020-04-13  1745    report
#> 130: 2020-04-14  1453    report



## Simulate cases with a weekly reporting effect and stochastic noise in reporting (beyond the delay)                                  
simulated_cases <- simulate_cases(rts, initial_cases = 100 , initial_date = as.Date("2020-03-01"),
                    generation_interval = generation_interval, delay_def = delay_def, 
                   incubation_def = incubation_def, 
                   reporting_effect = c(1.1, rep(1, 4), 0.95, 0.95),
                   reporting_model = function(n) {
                      out <- suppressWarnings(rnbinom(length(n), as.integer(n), 0.5))
                      out <- ifelse(is.na(out), 0, out)
                      })
                   
print(simulated_cases)
#>            date cases reference
#>   1: 2020-03-02    51 infection
#>   2: 2020-03-03    57 infection
#>   3: 2020-03-04    70 infection
#>   4: 2020-03-05    96 infection
#>   5: 2020-03-06   148 infection
#>  ---                           
#> 127: 2020-04-10  2613    report
#> 128: 2020-04-11  2385    report
#> 129: 2020-04-12  2083    report
#> 130: 2020-04-13  2297    report
#> 131: 2020-04-14  2002    report