Approximate Sampling a Distribution using Counts — sample_approx

Convolves cases by a PMF function. This function will soon be removed or replaced with a more robust stan implementation.

sample_approx_dist(
  cases = NULL,
  dist_fn = NULL,
  max_value = 120,
  earliest_allowed_mapped = NULL,
  direction = "backwards",
  type = "sample",
  truncate_future = TRUE
)

Arguments

cases: A dataframe of cases (in date order) with the following variables: date and cases.
dist_fn: Function that takes two arguments with the first being numeric and the second being logical (and defined as dist). Should return the probability density or a sample from the defined distribution. See the examples for more.
max_value: Numeric, maximum value to allow. Defaults to 120 days
earliest_allowed_mapped: A character string representing a date ("2020-01-01"). Indicates the earliest allowed mapped value.
direction: Character string, defato "backwards". Direction in which to map cases. Supports either "backwards" or "forwards".
type: Character string indicating the method to use to transform counts. Supports either "sample" which approximates sampling or "median" would shift by the median of the distribution.
truncate_future: Logical, should cases be truncated if they occur after the first date reported in the data. Defaults to TRUE.

Value

A data.table of cases by date of onset

Examples

cases <- example_confirmed
cases <- cases[, cases := as.integer(confirm)]
print(cases)
#>            date confirm cases
#>   1: 2020-02-22      14    14
#>   2: 2020-02-23      62    62
#>   3: 2020-02-24      53    53
#>   4: 2020-02-25      97    97
#>   5: 2020-02-26      93    93
#>  ---                         
#> 126: 2020-06-26     296   296
#> 127: 2020-06-27     255   255
#> 128: 2020-06-28     175   175
#> 129: 2020-06-29     174   174
#> 130: 2020-06-30     126   126

# total cases
sum(cases$cases)
#> [1] 240581

delay_fn <- function(n, dist, cum) {
  if (dist) {
    pgamma(n + 0.9999, 2, 1) - pgamma(n - 1e-5, 2, 1)
  } else {
    as.integer(rgamma(n, 2, 1))
  }
}

onsets <- sample_approx_dist(
  cases = cases,
  dist_fn = delay_fn
)

# estimated onset distribution
print(onsets)
#>            date cases
#>   1: 2020-02-12     1
#>   2: 2020-02-13     0
#>   3: 2020-02-14     0
#>   4: 2020-02-15     1
#>   5: 2020-02-16     0
#>  ---                 
#> 136: 2020-06-26   225
#> 137: 2020-06-27   184
#> 138: 2020-06-28   124
#> 139: 2020-06-29    71
#> 140: 2020-06-30    32

# check that sum is equal to reported cases
total_onsets <- median(
  purrr::map_dbl(
    1:100,
    ~ sum(sample_approx_dist(
      cases = cases,
      dist_fn = delay_fn
    )$cases)
  )
)
total_onsets
#> [1] 240620.5


# map from onset cases to reported
reports <- sample_approx_dist(
  cases = cases,
  dist_fn = delay_fn,
  direction = "forwards"
)


# map from onset cases to reported using a mean shift
reports <- sample_approx_dist(
  cases = cases,
  dist_fn = delay_fn,
  direction = "forwards",
  type = "median"
)