The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.

What corncrake() does

Corncrakes are famously detected far more often by ear than by eye — the overwhelming majority of field records are calls, not sightings. Surveyors have long applied call-based correction factors to convert a count of detections into an estimate of the true, largely-unheard population behind it.

corncrake() does the same thing to a surveillance count: every notification system under-ascertains true disease burden to some degree, and that degree is rarely constant — it typically varies by age group, by region, and over time as testing behaviour, clinical thresholds, or case definitions change. corncrake() takes an observed count, together with a factor that may vary by stratum and by time, and returns an estimate of the true total sitting behind it — with uncertainty bounds wherever they can be derived.

It is designed to run directly on roost() output (count_col defaults to "n", and time_col auto-detects the aggregation column from a roost_tbl), but works on any tidy data frame with a count column, per the ecosystem’s “no function requires another to have run first” design.


Three ways to get a factor

corncrake() supports three sources for the ascertainment factor, controlled by method:

  1. "user_supplied" (default): you supply factor_table, a lookup table of factors that can vary by group_by stratum and by a date_start/date_end validity window. This is the right choice when a factor comes from a published estimate, an external evaluation study (e.g. a capture-recapture study), or expert judgement.
  2. "ratio_estimate": the factor is derived internally as secondary_count / count at each stratum/time point, from a second, more-complete data stream. This is the classic surveillance “multiplier method” — for example, dividing all positive laboratory tests by notified cases to estimate under-notification.
  3. "severity_anchor": the factor is derived by comparing an observed severity ratio already in your data (e.g. a case-fatality or case-hospitalisation rate) against a reference_rate representing the believed-true rate from a well-ascertained source. This is the case-fatality/infection-fatality-rate anchor inversion described below, and it’s the method to reach for early in a novel outbreak, before any seroprevalence survey exists.

Synthetic data

set.seed(42)
n <- 400
diag <- data.frame(
  onset_date = as.Date("2024-01-01") + sample(0:59, n, replace = TRUE),
  age        = sample(0:90, n, replace = TRUE),
  stringsAsFactors = FALSE
)
diag <- preening(diag, age_col = "age", scheme = "flucan_sentinel")

cases_monthly <- roost(
  diag,
  date_col   = "onset_date",
  time_unit  = "month",
  group_cols = "age_group"
)

cases_monthly
#> # A tibble: 10 × 3
#>    age_group month          n
#>    <ord>     <date>     <int>
#>  1 0-4       2024-01-01     7
#>  2 0-4       2024-02-01    13
#>  3 5-15      2024-01-01    20
#>  4 5-15      2024-02-01    28
#>  5 16-49     2024-01-01    87
#>  6 16-49     2024-02-01    80
#>  7 50-64     2024-01-01    37
#>  8 50-64     2024-02-01    28
#>  9 65+       2024-01-01    48
#> 10 65+       2024-02-01    52
#> 
#> -- roost_meta --------------------------------------
#>  time_unit  : month 
#>  date_range : 2024-01-01 to 2024-02-29 
#>  group_cols : age_group 
#>  hemisphere : southern 
#>  n_rows_in  : 400

User-supplied factors

Suppose an evaluation study estimated ascertainment separately for two broad age bands, with a slightly higher (and less certain) multiplier for younger ages, tightening over time as testing improved:

factors <- data.frame(
  age_group    = rep(c("0-4", "5-15", "16-49", "50-64", "65+"), each = 2),
  date_start   = rep(as.Date(c("2024-01-01", "2024-02-01")), 5),
  date_end     = rep(as.Date(c("2024-01-31", "2024-02-29")), 5),
  factor       = c(3.2, 2.8,  2.6, 2.3,  1.8, 1.6,  1.5, 1.4,  1.3, 1.2),
  factor_lower = c(2.4, 2.1,  2.0, 1.8,  1.4, 1.3,  1.2, 1.1,  1.1, 1.0),
  factor_upper = c(4.2, 3.6,  3.3, 2.9,  2.3, 2.0,  1.9, 1.7,  1.6, 1.5),
  source       = "Illustrative multiplier, SCPHU surveillance evaluation 2025"
)

knitr::kable(factors)
age_group date_start date_end factor factor_lower factor_upper source
0-4 2024-01-01 2024-01-31 3.2 2.4 4.2 Illustrative multiplier, SCPHU surveillance evaluation 2025
0-4 2024-02-01 2024-02-29 2.8 2.1 3.6 Illustrative multiplier, SCPHU surveillance evaluation 2025
5-15 2024-01-01 2024-01-31 2.6 2.0 3.3 Illustrative multiplier, SCPHU surveillance evaluation 2025
5-15 2024-02-01 2024-02-29 2.3 1.8 2.9 Illustrative multiplier, SCPHU surveillance evaluation 2025
16-49 2024-01-01 2024-01-31 1.8 1.4 2.3 Illustrative multiplier, SCPHU surveillance evaluation 2025
16-49 2024-02-01 2024-02-29 1.6 1.3 2.0 Illustrative multiplier, SCPHU surveillance evaluation 2025
50-64 2024-01-01 2024-01-31 1.5 1.2 1.9 Illustrative multiplier, SCPHU surveillance evaluation 2025
50-64 2024-02-01 2024-02-29 1.4 1.1 1.7 Illustrative multiplier, SCPHU surveillance evaluation 2025
65+ 2024-01-01 2024-01-31 1.3 1.1 1.6 Illustrative multiplier, SCPHU surveillance evaluation 2025
65+ 2024-02-01 2024-02-29 1.2 1.0 1.5 Illustrative multiplier, SCPHU surveillance evaluation 2025

group_by and time_col tell corncrake() how to match rows in cases_monthly against rows in factors. Here time_col is left NULL and auto-detects as "month" from cases_monthly’s roost_meta:

cases_corrected <- corncrake(
  cases_monthly,
  factor_table = factors,
  group_by     = "age_group"
)

cases_corrected[, c("age_group", "month", "n", "ascertainment_factor",
                     "corrected_count", "corrected_count_lower", "corrected_count_upper")]
#> # A tibble: 10 × 7
#>    age_group month          n ascertainment_factor corrected_count
#>    <ord>     <date>     <int>                <dbl>           <dbl>
#>  1 0-4       2024-01-01     7                  3.2            22.4
#>  2 0-4       2024-02-01    13                  2.8            36.4
#>  3 5-15      2024-01-01    20                  2.6            52  
#>  4 5-15      2024-02-01    28                  2.3            64.4
#>  5 16-49     2024-01-01    87                  1.8           157. 
#>  6 16-49     2024-02-01    80                  1.6           128  
#>  7 50-64     2024-01-01    37                  1.5            55.5
#>  8 50-64     2024-02-01    28                  1.4            39.2
#>  9 65+       2024-01-01    48                  1.3            62.4
#> 10 65+       2024-02-01    52                  1.2            62.4
#> # ℹ 2 more variables: corrected_count_lower <dbl>, corrected_count_upper <dbl>
#> 
#> -- roost_meta --------------------------------------
#>  time_unit  : month 
#>  date_range : 2024-01-01 to 2024-02-29 
#>  group_cols : age_group 
#>  hemisphere : southern 
#>  n_rows_in  : 400

cases_corrected is still a roost_tbl — corncrake() preserves the input class, so it drops straight into any downstream code (or a future bowerbird::roost_plot()) written against roost() output.


Uncertainty: two different questions

ci_method controls what corrected_count_lower/corrected_count_upper actually represent, because “uncertainty” can mean two different things here:

  • "table" (default): the factor itself is uncertain (as in factors above); the observed count is treated as fixed. Bounds come from factor_lower/factor_upper.
  • "propagate": the factor is treated as fixed; the observed count is treated as a realisation of a Poisson process with its own sampling uncertainty (exact Poisson confidence interval), which is then propagated through the point factor.
  • "none": point estimate only.
corncrake(cases_monthly, factor_table = factors, group_by = "age_group",
          ci_method = "propagate")[
  , c("age_group", "month", "n", "corrected_count",
      "corrected_count_lower", "corrected_count_upper")
]
#> # A tibble: 10 × 6
#>    age_group month          n corrected_count corrected_count_lower
#>    <ord>     <date>     <int>           <dbl>                 <dbl>
#>  1 0-4       2024-01-01     7            22.4                  9.01
#>  2 0-4       2024-02-01    13            36.4                 19.4 
#>  3 5-15      2024-01-01    20            52                   31.8 
#>  4 5-15      2024-02-01    28            64.4                 42.8 
#>  5 16-49     2024-01-01    87           157.                 125.  
#>  6 16-49     2024-02-01    80           128                  101.  
#>  7 50-64     2024-01-01    37            55.5                 39.1 
#>  8 50-64     2024-02-01    28            39.2                 26.0 
#>  9 65+       2024-01-01    48            62.4                 46.0 
#> 10 65+       2024-02-01    52            62.4                 46.6 
#> # ℹ 1 more variable: corrected_count_upper <dbl>
#> 
#> -- roost_meta --------------------------------------
#>  time_unit  : month 
#>  date_range : 2024-01-01 to 2024-02-29 
#>  group_cols : age_group 
#>  hemisphere : southern 
#>  n_rows_in  : 400

Note the factor_table doesn’t need bounds at all for ci_method = "propagate" to work — the uncertainty here comes entirely from the count, not the factor.


Deriving a factor instead: the ratio (multiplier) method

If a more-complete secondary data stream is available — say, all positive laboratory results, independent of whether a notification was ever made — corncrake() can derive the factor directly rather than requiring you to supply one:

lab_positive_monthly <- cases_monthly
lab_positive_monthly$n <- round(cases_monthly$n * runif(nrow(cases_monthly), 1.3, 2.5))

cases_corrected2 <- corncrake(
  cases_monthly,
  method              = "ratio_estimate",
  group_by            = "age_group",
  secondary_data      = lab_positive_monthly,
  secondary_count_col = "n"
)

cases_corrected2[, c("age_group", "month", "n", "ascertainment_factor", "corrected_count")]
#> # A tibble: 10 × 5
#>    age_group month          n ascertainment_factor corrected_count
#>    <ord>     <date>     <int>                <dbl>           <dbl>
#>  1 0-4       2024-01-01     7                 1.86              13
#>  2 0-4       2024-02-01    13                 1.31              17
#>  3 5-15      2024-01-01    20                 1.95              39
#>  4 5-15      2024-02-01    28                 2.18              61
#>  5 16-49     2024-01-01    87                 1.59             138
#>  6 16-49     2024-02-01    80                 2.28             182
#>  7 50-64     2024-01-01    37                 1.81              67
#>  8 50-64     2024-02-01    28                 1.96              55
#>  9 65+       2024-01-01    48                 1.48              71
#> 10 65+       2024-02-01    52                 1.54              80
#> 
#> -- roost_meta --------------------------------------
#>  time_unit  : month 
#>  date_range : 2024-01-01 to 2024-02-29 
#>  group_cols : age_group 
#>  hemisphere : southern 
#>  n_rows_in  : 400

secondary_data must share the same time_col and group_by column names as the primary data. The derived factor is a point estimate only (ascertainment_factor_lower/_upper are NA); use ci_method = "propagate" if you want count-based bounds alongside a ratio-estimate factor.


Severity anchor: the CFR/IFR-anchor inversion

Both methods above need a second count data stream. Early in a novel outbreak — the scenario WHO pandemic-preparedness planning calls “Disease X”, where the pathogen is real but its identity, and therefore any tailored surveillance stream, doesn’t yet exist — that second count stream usually isn’t available yet. What often is available is an externally published severity estimate from a reference jurisdiction or a global body, together with your own locally observed severity ratio.

This is the case-fatality/infection-fatality-rate anchor inversion described in Smoll et al.’s Queensland COVID-19 under-ascertainment analysis. The logic: if surveillance ascertained every true infection, the observed case-fatality rate (deaths ÷ notified cases) would equal the true infection-fatality rate. Under-ascertainment inflates the observed rate above the true one — by exactly the ascertainment factor:

\[\text{UAF} = \frac{\text{CFR}_{\text{obs}}}{\text{IFR}_{\text{ref}}} = \frac{\text{deaths}_{\text{obs}} / \text{cases}_{\text{obs}}}{\text{IFR}_{\text{ref}}}\]

The same identity holds for any other severity outcome — a case-hospitalisation rate against a reference infection-hospitalisation rate works identically. corncrake() implements this generally as method = "severity_anchor", comparing severity_count_col / count_col in your data against a reference_rate.

A worked example, reproducing the paper’s own numbers

Suppose a jurisdiction early in a Disease X outbreak observes 100 registered deaths against 5,000 notified cases, and a reference IFR of 1.0% (with a plausible range of 0.5%–2.0%) is available from a high-ascertainment reference jurisdiction:

disease_x <- data.frame(
  month    = as.Date("2024-01-01"),
  n_cases  = 5000L,
  n_deaths = 100L
)

disease_x_corrected <- corncrake(
  disease_x,
  count_col             = "n_cases",
  method                 = "severity_anchor",
  time_col               = "month",
  severity_count_col     = "n_deaths",
  reference_rate         = 0.01,
  reference_rate_lower   = 0.005,
  reference_rate_upper   = 0.020,
  reference_source       = "WHO Disease X planning scenario, IFR 1.0% (0.5-2.0%)"
)

disease_x_corrected[, c("n_cases", "n_deaths", "ascertainment_factor",
                         "ascertainment_factor_lower", "ascertainment_factor_upper",
                         "corrected_count")]
#>   n_cases n_deaths ascertainment_factor ascertainment_factor_lower
#> 1    5000      100                    2                          1
#>   ascertainment_factor_upper corrected_count
#> 1                          4           10000

The observed CFR here is 100/5000 = 2.0%, double the 1.0% reference IFR, so UAF = 2.0 — implying the true infection burden was twice the notified case count, exactly matching the paper’s own worked example.

The bounds invert — and corncrake() handles that for you

This is the detail worth being deliberate about. Because UAF is divided by reference_rate, it’s a decreasing function of it: a higher reference rate implies less under-ascertainment, not more. That means the usual intuition — “lower bound in, lower bound out” — is backwards here:

  • ascertainment_factor_lower is computed from reference_rate_upper
  • ascertainment_factor_upper is computed from reference_rate_lower

In the example above, the 0.5%–2.0% reference range produces a UAF range of [1.0, 4.0] — and the lower UAF bound (1.0) comes from the higher reference rate (2.0%), not the lower one. Getting this backwards by hand is an easy mistake to make (the source paper calls it out as a dedicated remark), which is exactly why corncrake() encodes it once rather than leaving it as an instruction to re-derive on every use.

A stratified or time-varying reference rate

reference_rate doesn’t have to be a single scalar. Supply a factor_table-shaped data frame with a rate column (optionally rate_lower/rate_upper) for a reference rate that varies by group_by stratum or by time window, using exactly the same date_start/date_end validity-window mechanism as method = "user_supplied"’s factor_table:

disease_x_stratified <- data.frame(
  month     = as.Date(c("2024-01-01", "2024-01-01")),
  age_group = c("0-17", "18+"),
  n_cases   = c(1000L, 4000L),
  n_deaths  = c(1L, 99L)
)

reference_rates <- data.frame(
  age_group  = c("0-17", "18+"),
  date_start = as.Date(NA),  # open-ended: one reference rate per age group, all time
  date_end   = as.Date(NA),
  rate       = c(0.001, 0.02),
  source     = "Illustrative age-stratified reference IFR"
)

corncrake(
  disease_x_stratified,
  count_col           = "n_cases",
  method               = "severity_anchor",
  group_by             = "age_group",
  time_col             = "month",
  severity_count_col   = "n_deaths",
  reference_rate       = reference_rates
)[, c("age_group", "n_cases", "n_deaths", "ascertainment_factor")]
#>   age_group n_cases n_deaths ascertainment_factor
#> 1      0-17    1000        1               1.0000
#> 2       18+    4000       99               1.2375

severity_count_col and count_col need to sit in the same table at the same stratification/time — see vignette("flyway") for building exactly that shape from a linked cohort’s onset and fatality dates in one step.


Rates alongside corrected counts

If a population denominator is available, denominator_col adds corrected_rate (+ bounds) directly:

cases_with_pop <- cases_monthly
cases_with_pop$pop <- ifelse(cases_with_pop$age_group == "0-4", 8000,
                       ifelse(cases_with_pop$age_group == "5-15", 15000,
                       ifelse(cases_with_pop$age_group == "16-49", 45000,
                       ifelse(cases_with_pop$age_group == "50-64", 20000, 18000))))

corncrake(
  cases_with_pop, factor_table = factors, group_by = "age_group",
  denominator_col = "pop", rate_multiplier = 100000
)[, c("age_group", "month", "corrected_count", "corrected_rate")]
#> # A tibble: 10 × 4
#>    age_group month      corrected_count corrected_rate
#>    <ord>     <date>               <dbl>          <dbl>
#>  1 0-4       2024-01-01            22.4           280 
#>  2 0-4       2024-02-01            36.4           455 
#>  3 5-15      2024-01-01            52             347.
#>  4 5-15      2024-02-01            64.4           429.
#>  5 16-49     2024-01-01           157.            348 
#>  6 16-49     2024-02-01           128             284.
#>  7 50-64     2024-01-01            55.5           278.
#>  8 50-64     2024-02-01            39.2           196 
#>  9 65+       2024-01-01            62.4           347.
#> 10 65+       2024-02-01            62.4           347.
#> 
#> -- roost_meta --------------------------------------
#>  time_unit  : month 
#>  date_range : 2024-01-01 to 2024-02-29 
#>  group_cols : age_group 
#>  hemisphere : southern 
#>  n_rows_in  : 400

A note on missing factors

Not every stratum/time combination in your data needs to be covered by factor_table — but if one isn’t, corncrake() needs to know what to do about it. By default (on_missing = "warn_na") it leaves the row NA and issues a single warning naming how many rows were affected; on_missing = "error" stops outright, which is useful when you want to be certain your factor table has full coverage before proceeding.

sparse_factors <- factors[factors$age_group != "0-4", ]
corncrake(cases_monthly, factor_table = sparse_factors, group_by = "age_group",
          on_missing = "error")
#> Error:
#> ! (*)> mudnester::corncrake() — 2 row(s) had no matching ascertainment factor for their stratum/time and were left NA.

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.