Skip to content

API

Public functions and constants exported by BVDOutbreakSize.

BVDOutbreakSize.CUMULATIVE_INTEGRAL_ALG Constant
julia
CUMULATIVE_INTEGRAL_ALG

Gauss-Legendre quadrature scheme (n = 32) used for the at-risk person-time export integral and the deaths-among-exports convolution (outer and inner integrals).

source
BVDOutbreakSize.DEATH_INTEGRAL_ALG Constant
julia
DEATH_INTEGRAL_ALG

Gauss-Legendre quadrature scheme (n = 64) used for the deaths onset-to-death convolution, the no-onward-transmission counterfactual, and the forecast deaths integral.

source
BVDOutbreakSize.DELAY_SUPPORT_K Constant
julia
DELAY_SUPPORT_K

Number of standard deviations beyond the mean used as the clustering scale for a delay distribution in the onset-to-death convolution integrals. mean + DELAY_SUPPORT_K · std is the width near the cut-off over which the clustered integrate packs roughly half its nodes, so the quadrature tracks the delay's scale as it is sampled.

source
BVDOutbreakSize.EXPORT_DELAY_GRID_POINTS Constant
julia
EXPORT_DELAY_GRID_POINTS

Number of evenly spaced grid points used to precompute the onset-to-death CDF in ExportDeathDelay.

source
BVDOutbreakSize.ITURI_DAILY_TRAVEL Constant
julia
ITURI_DAILY_TRAVEL

Default prior mean for the daily outbound traveller volume from Ituri Province across seven points of entry.

source
BVDOutbreakSize.ITURI_DAILY_TRAVEL_SD Constant
julia
ITURI_DAILY_TRAVEL_SD

Default prior SD for the daily outbound traveller volume, covering point-of-entry-to-point-of-entry variation and reporting uncertainty in the underlying mobility survey.

source
BVDOutbreakSize.ITURI_POPULATION Constant
julia
ITURI_POPULATION

Source population for the Ituri Province (McCabe et al., Table 1).

source
BVDOutbreakSize.REPORT_SCENARIOS Constant
julia
REPORT_SCENARIOS

Published point estimates of cumulative cases C_T from McCabe et al. (Imperial College London, 20 May 2026 update), as (label, value) tuples in the order they appear in Tables 1 and 2.

source
BVDOutbreakSize.ExportDeathDelay Type

Onset-to-death delay carrying its CDF F_d precomputed on an evenly spaced grid over [0, gmax], built once and reused across every outer node and bin edge of the deaths-among-exports convolution rather than re-integrating the density at each (see integrate_exports_deaths). gmax must cover the largest delay argument reached, the detection window w. Pass one in place of the raw delay distribution to select this path by dispatch; the distribution methods are unchanged and remain the reference.

F_d is a cumulative trapezoid of the density, so the convolution still differentiates through the density alone (the AD backend lacks the Gamma CDF shape derivative). The density is never evaluated at 0, whose Gamma shape derivative is 0·log 0 = NaN under AD; for f_d(0) = 0 (Gamma shape > 1) treating it as zero is exact.

source
BVDOutbreakSize._gamma_cdf Method
julia
_gamma_cdf(α, θ, x) -> Any

Wrapper around cdf(Gamma(α, θ), x) as a 3-argument scalar scalar function to attach reverse-mode rule.

Instead we attach our own analytic rrule:

xF=f(x;α),θF=xθf(x;α),αF=ψ(α)P(α,y)+1Γ(α)0ytα1etlogtdt,

with y = x/θ, P the regularized lower incomplete gamma and ψ is the digamma function.

source
BVDOutbreakSize._gamma_cdf_partials Method
julia
_gamma_cdf_partials(α, θ, x) -> Tuple{Any, Any, Any}

Compute the partial derivatives of the gamma CDF with respect to α, θ, and x.

source
BVDOutbreakSize._grad_p_a_series Method
julia
_grad_p_a_series(a, z; rtol, maxiter) -> Any

Series sum of term derivatives for ∂_α P(α, z), using the absolutely-convergent Kummer expansion

P(α,z)=zαezn=0znΓ(α+n+1),αP(α,z)=log(z)P(α,z)zαezn=0ψ(α+n+1)znΓ(α+n+1).

The digamma factor advances by the recurrence ψ(α + n + 1) = ψ(α + n) + 1/(α + n), so each iteration costs only a division, an add and two multiplies — no per-iteration gamma or digamma calls. Stan's grad_reg_inc_gamma uses this same series as its small-x branch; the Julia formulation here is the direct port from EpiAware/CensoredDistributions PR #250.

source
BVDOutbreakSize.comparison_table Method
julia
comparison_table(
    C_draws::AbstractVector;
    scenarios
) -> DataFrames.DataFrame

For each published C_T scenario, the narrowest joint posterior credible interval (30, 60 or 90%) that contains it, or "outside 90%".

source
BVDOutbreakSize.default_adtype Method
julia
default_adtype() -> ADTypes.AutoMooncake{Mooncake.Config}

Mooncake reverse-mode AD with default Mooncake.Config(). Used as the NUTS adtype keyword.

source
BVDOutbreakSize.diagnostics_table Method
julia
diagnostics_table(fits::Pair{String}...) -> Any

DataFrame of fit-quality diagnostics with one row per fit. Pass each fit as "label" => chain. Columns :fit, :max_rhat, :min_ess_bulk, :divergences.

source
BVDOutbreakSize.expected_deaths Method
julia
expected_deaths(CFR, r, T, delay_dist; alg) -> Any

Expected cumulative deaths by time T from a single seeding case under exponential growth at rate r:

E[DT]=CFR0Tersf(Ts)ds,

with f the delay_dist onset-to-death density. The integrand returns zero past T so the convolution support is respected. Integrated with the clustered integrate so the quadrature nodes pack near T, where f(T − s) has mass, over a window set by the delay scale (see DELAY_SUPPORT_K). Uses DEATH_INTEGRAL_ALG.

source
BVDOutbreakSize.expected_deaths Method
julia
expected_deaths(
    CFR,
    r,
    T,
    delay_dist::Distributions.Gamma
) -> Any

Expected cumulative deaths by time T from a single seeding case under exponential growth at rate r, using analytic result specific to the Gamma distribution:

E[DT]=CFR0Tersf(Ts)ds=CFRerTM(r)F(T(1+θr)),

where f is the Gamma(α, θ) density, M is its moment-generating function, and F is its CDF. The CDF is routed through _gamma_cdf so Mooncake reverse-mode AD picks up a shape-parameter gradient (computed internally with ForwardDiff).

source
BVDOutbreakSize.expected_exports Method
julia
expected_exports(cumulative, p, q, t, window; alg) -> Any

Expected cumulative detected exports by elapsed time t, clamped to be strictly positive and finite:

E[exports(t)]=pqtwtC(s)ds,

with cumulative the trajectory C(s), p the detection probability, q the per-capita travel rate, and window the detection window w. Backs both the exports count likelihood (evaluated at t = T) and the first-export-detection survival term (evaluated at an earlier t). Uses CUMULATIVE_INTEGRAL_ALG.

source
BVDOutbreakSize.expected_exports_deaths Method
julia
expected_exports_deaths(
    cumulative,
    delay_dist,
    CFR,
    p,
    q,
    t,
    window;
    alg
) -> Any

Expected cumulative deaths among detected exports by elapsed time t, clamped to be strictly positive and finite:

E[Duganda(t)]=CFRpqtwtC(s)Fd(ts)ds,

with cumulative the trajectory C(s), delay_dist the onset-to-death distribution (CDF Fd), p the detection probability, q the per-capita travel rate, and window the detection window w. Backs the binned-Poisson export-death likelihood, evaluated at the daily bin edges. Uses CUMULATIVE_INTEGRAL_ALG.

source
BVDOutbreakSize.fit_diagnostics Method
julia
fit_diagnostics(
    chn
) -> NamedTuple{(:max_rhat, :min_ess_bulk, :n_divergent), <:Tuple{Float64, Float64, Any}}

NUTS fit-quality summary for one chain: the worst (maximum) R-hat and the smallest bulk effective sample size across parameters, and the number of divergent transitions.

source
BVDOutbreakSize.forecast_reported Method
julia
forecast_reported(chn; horizon = 7, daily_travellers, source_population,
                  obs_cases, obs_deaths, obs_exports, seed = 20260520)

One-week-ahead (default horizon = 7 days) posterior-predictive forecast. For each draw, continue exponential growth to T + horizon and apply the observation models, returning a DataFrame with one row per draw and columns:

  • :cases_cum, :deaths_cum, :exports_cum — replicated cumulative counts reported by T + horizon.

  • :cases_new, :deaths_new, :exports_new — new counts over the coming week (*_cum minus the corresponding observed count at T, floored at zero).

Reads :r, :T, :CFR, :α, :θ, :w, :p_drc, :p_uganda, :k from chn. DRC reported cases use the DRC ascertainment fraction p_drc; exports use p_uganda · q with q = daily_travellers / source_population. Assumes growth continues unchanged over the horizon (no interventions, no saturation).

source
BVDOutbreakSize.forecast_table Method
julia
forecast_table(fc::DataFrame; digits = 0)

Summarise a forecast_reported result into a DataFrame with one row per stream (cases, deaths, exports) and quantity (cumulative total by T + 7, or new this week), reporting the same equal-tailed 30/60/90% credible interval endpoints (lower_90 … upper_90) as the other summary tables.

source
BVDOutbreakSize.forecast_vs_truth Method
julia
forecast_vs_truth(fc::DataFrame; cases, deaths, exports, digits = 0)

Validate a forecast_reported projection against the counts that were later observed. cases, deaths and exports are the observed cumulative DRC reported cases, DRC deaths and Uganda exports at the forecast target date. Returns a DataFrame with one row per stream giving the observed count, the equal-tailed 30/60/90% predictive intervals (the same endpoints as the other summary tables), and whether the observed count falls inside the 90% interval. Use it to forecast from an earlier data snapshot and check the now-known truth against the projection.

source
BVDOutbreakSize.integrate Method
julia
integrate(f, lo, hi, scale; alg) -> Any

Integrate a scalar function f over [lo, hi] by Gauss-Legendre quadrature with the nodes clustered towards hi, resolving features of size scale near that limit. Use for onset-to-death convolutions, whose integrand mass piles up against the cut-off hi over a window set by the sampled delay scale: the uniform integrate spread across a wide [lo, hi] would resolve that peak too coarsely. The whole domain is still covered (no truncation), so a delay wider than [lo, hi] is integrated exactly; when scale ≥ hi - lo the map reduces to the uniform method. Returns zero when hi <= lo. The default scheme is DEATH_INTEGRAL_ALG.

source
BVDOutbreakSize.integrate Method
julia
integrate(f, lo, hi; alg) -> Any

Integrate a scalar function f over [lo, hi] by Gauss-Legendre quadrature. The reference domain [-1, 1] is mapped onto [lo, hi], so any forward integral in the package can share one integrator. Returns zero when hi <= lo. The default scheme is CUMULATIVE_INTEGRAL_ALG.

source
BVDOutbreakSize.integrate_cumulative Method
julia
integrate_cumulative(cumulative, lo, hi; alg) -> Any

Integrate a cumulative-incidence trajectory cumulative (a callable C(s)) over [lo, hi]. Backs the at-risk person-time export integral TwTC(s)ds. Uses CUMULATIVE_INTEGRAL_ALG.

source
BVDOutbreakSize.integrate_exports_deaths Method
julia
integrate_exports_deaths(
    cumulative,
    delay_dist,
    lo,
    hi,
    T;
    alg
) -> Any

Deaths-among-exports convolution lohiC(s)Fd(Ts)ds with F_d the delay_dist onset-to-death CDF. The CDF is itself written as the inner integral of the density, Fd(x)=0xfd(u)du, so the whole expression differentiates through the density alone. Previously, the reverse-mode AD backend did not support the gamma CDF shape-parameter derivative). Outer and inner integrals both use CUMULATIVE_INTEGRAL_ALG.

source
BVDOutbreakSize.integrate_exports_deaths Method
julia
integrate_exports_deaths(
    cumulative,
    delay_dist::Distributions.Gamma,
    lo,
    hi,
    T;
    alg
) -> Any

Deaths-among-exports convolution specialised for Gamma delay distributions. Same expression as the generic method — lohiC(s)Fd(Ts)ds — but the onset-to-death CDF is evaluated in closed form via _gamma_cdf rather than re-integrated from the density at each outer quadrature node. Drops one entire quadrature level (just the outer integral remains) and sidesteps the singularity that fixed-node Gauss-Legendre hits at the head-form 0ytα1etdu for α < 1. The shape-parameter derivative the AD backend needed for cdf(::Gamma, ·) is supplied by _gamma_cdf's rrule (see src/gamma_cdf.jl).

source
BVDOutbreakSize.integrate_exports_deaths Method
julia
integrate_exports_deaths(
    cumulative,
    ed::ExportDeathDelay,
    lo,
    hi,
    T;
    alg
) -> Any

Deaths-among-exports convolution using a precomputed ExportDeathDelay: lohiC(s)Fd(Ts)ds with F_d interpolated off the grid rather than re-integrated at each node. Mathematically identical to the distribution method up to the grid resolution.

source
BVDOutbreakSize.load_observations Function
julia
load_observations(

) -> NamedTuple{(:as_of_date, :exported_cases, :exports_deaths, :total_deaths, :reported_cases, :daily_outbound_travellers, :daily_outbound_travellers_sd, :source_population, :export_deaths_daily, :first_export_detection_delta, :genetic_tmrca_days, :genetic_tmrca_days_sd, :genetic_tmrca_alt_days, :genetic_tmrca_alt_days_sd, :sources), <:Tuple{Any, Any, Any, Any, Any, Any, Any, Any, Any, Union{Missing, Int64}, Union{Missing, Int64}, Any, Union{Missing, Int64}, Any, NamedTuple{(:exported_cases, :exports_deaths, :total_deaths, :reported_cases, :daily_outbound_travellers, :daily_outbound_travellers_sd, :source_population, :genetic_tmrca), <:NTuple{8, Any}}}}
load_observations(
    path::AbstractString
) -> NamedTuple{(:as_of_date, :exported_cases, :exports_deaths, :total_deaths, :reported_cases, :daily_outbound_travellers, :daily_outbound_travellers_sd, :source_population, :export_deaths_daily, :first_export_detection_delta, :genetic_tmrca_days, :genetic_tmrca_days_sd, :genetic_tmrca_alt_days, :genetic_tmrca_alt_days_sd, :sources), <:Tuple{Any, Any, Any, Any, Any, Any, Any, Any, Any, Union{Missing, Int64}, Union{Missing, Int64}, Any, Union{Missing, Int64}, Any, NamedTuple{(:exported_cases, :exports_deaths, :total_deaths, :reported_cases, :daily_outbound_travellers, :daily_outbound_travellers_sd, :source_population, :genetic_tmrca), <:NTuple{8, Any}}}}

Load the observation block from data/observations.toml and return it as a NamedTuple. Each observation in the TOML is a subtable with value = … and source = "…"; this function returns both the parsed numeric values and a parallel sources::NamedTuple of citation strings so they can be printed alongside the data.

Fields returned:

  • exported_cases::Int

  • exports_deaths::Int

  • total_deaths::Int

  • reported_cases::Int

  • daily_outbound_travellers::Real

  • daily_outbound_travellers_sd::Real

  • source_population::Int

  • genetic_tmrca_days::Union{Real, Missing} — estimated time to the most recent common ancestor (TMRCA) in days before as_of_date, a soft lower bound on the seeding time T; missing when no genetic_tmrca block is present.

  • genetic_tmrca_days_sd::Union{Real, Missing} — SD (days) on the location of that floor; missing when absent.

  • genetic_tmrca_alt_days::Union{Real, Missing} — TMRCA (days before as_of_date) under the alternative clock rate, for the clock-rate sensitivity; missing when no alt_date is present.

  • genetic_tmrca_alt_days_sd::Union{Real, Missing} — SD (days) on the alternative-clock floor; missing when absent.

  • sources::NamedTuple{(:exported_cases, :exports_deaths, :total_deaths, :reported_cases, :daily_outbound_travellers, :daily_outbound_travellers_sd, :source_population, :genetic_tmrca), NTuple{8, String}} — citation per field.

source
BVDOutbreakSize.nuts_sample Method
julia
nuts_sample(
    model;
    samples,
    chains,
    target_accept,
    seed,
    progress,
    adtype
) -> Any

NUTS on model, parallel chains via MCMCThreads. Chains initialise from the prior to keep the sampler away from the boundary of constrained variables.

source
BVDOutbreakSize.plot_cfr_prior Method
julia
plot_cfr_prior(
    prior::Distributions.Distribution
) -> Makie.Figure

Density of a prior over the case-fatality ratio (CFR) on [0, 1], plotted on the sub-range [0, 0.7]. The CDC central estimate of 55/169 ≈ 0.33 is drawn as a solid vertical rule, and the report's 26% and 40% scenario bounds as dashed rules, so the prior can be read against the published CFR scenarios.

source
BVDOutbreakSize.plot_cumulative_cases Method
julia
plot_cumulative_cases(
    streams::Pair{String, <:AbstractVector}...;
    scenarios,
    xmax
) -> AlgebraOfGraphics.FigureGrid

Overlaid posterior densities of C_T from one or more fits, built through AlgebraOfGraphics. The 15 published scenario point estimates are drawn as faint dashed Makie vlines on top of the AoG figure.

source
BVDOutbreakSize.plot_density_overlay Method
julia
plot_density_overlay(
    streams::Pair{String, <:AbstractVector}...;
    xlabel,
    title
) -> AlgebraOfGraphics.FigureGrid

Overlaid posterior densities of an arbitrary scalar quantity from one or more fits, built through AlgebraOfGraphics. Pass each fit as "label" => draws; xlabel and title set the axis text.

source
BVDOutbreakSize.plot_estimate_comparison Method
julia
plot_estimate_comparison(
    rows::AbstractVector;
    xlabel,
    xmax
) -> Makie.Figure

Horizontal point-and-interval comparison of cumulative-case estimates from several sources. rows is a vector of (label, central, lower, upper) tuples, drawn top to bottom with the central estimate as a point and [lower, upper] as a bar. Use it to place model posteriors next to published point estimates and their intervals.

source
BVDOutbreakSize.plot_forecast Method
julia
plot_forecast(fc::DataFrame)

Three-panel histogram of the new-this-week forecast counts (cases, deaths, exports) from forecast_reported.

source
BVDOutbreakSize.plot_forecast_vs_truth Method
julia
plot_forecast_vs_truth(fc::DataFrame; cases, deaths, exports,
                       baseline_cases = 0, baseline_deaths = 0,
                       baseline_exports = 0)

Validation figure for a forecast_reported projection, laid out as a 2×3 grid. The top row shows the cumulative forecast distribution per stream (DRC reported cases, DRC deaths, Uganda exports); the bottom row shows the new counts forecast over the horizon, mirroring the one-week-ahead forecast. Each panel is a histogram with the 90% predictive interval shaded and the later-observed count drawn as a dashed black rule. cases, deaths and exports are the observed cumulative counts; baseline_* are the counts at the forecast origin, so the observed new count is the cumulative truth minus the baseline.

source
BVDOutbreakSize.plot_no_onward_deaths Method
julia
plot_no_onward_deaths(df; obs_deaths)

Two-panel density of the no-onward-transmission counterfactual from predict_no_onward_deaths. The left panel shows the still expected deaths (:delta_deaths, the future deaths in cases already infected by T, net of the obs_deaths already observed). The right panel shows the projected total (:total_projected = obs_deaths + delta_deaths) with a dashed black rule at obs_deaths. Both are lower bounds: they assume every onward transmission stops at time T.

source
BVDOutbreakSize.plot_pair Method
julia
plot_pair(
    chn,
    params::AbstractVector{Symbol};
    thin,
    prior
) -> Makie.Figure

PairPlots.jl corner plot over the named posterior parameters, thinned by thin. Pass prior (another chain holding the same parameters) to overlay the prior as a second series with a legend, so the data's contribution to each marginal is visible.

source
BVDOutbreakSize.plot_posterior_predictive Method
julia
plot_posterior_predictive(
    pp_exports::Union{Nothing, AbstractVector},
    pp_deaths::Union{Nothing, AbstractVector},
    obs_exports::Union{Nothing, Real},
    obs_deaths::Union{Nothing, Real};
    pp_cases,
    obs_cases,
    pp_exports_deaths,
    obs_exports_deaths,
    predictive_label
) -> Makie.Figure

Posterior predictive histogram with one panel per supplied data stream. Pass pp_exports/pp_deaths as nothing to suppress either of the first two panels, and supply pp_cases and/or pp_exports_deaths to add the reported-cases and deaths-among-exports panels. Observed values are drawn as red vlines. With four streams the panels are laid out as a 2×2 grid (exports cases, exports deaths, DRC deaths, DRC reported cases); fewer streams are placed in a single row.

source
BVDOutbreakSize.plot_posterior_predictive_grid Method
julia
plot_posterior_predictive_grid(
;
    individual,
    joint,
    observed
)

Two-row × four-column comparison of posterior-predictive distributions. Top row: replicates from the per-stream fits. Bottom row: replicates from the joint fit, conditioning on all observed streams. Observed values shown as red vertical lines.

Each NamedTuple carries (; exports, exports_deaths, deaths, cases). Each panel is a histogram of replicated counts; rows share the same x-axis (the stream's count) so the per-stream and joint predictives are directly comparable.

source
BVDOutbreakSize.plot_prior_predictive Method
julia
plot_prior_predictive(
    pp_exports::Union{Nothing, AbstractVector},
    pp_deaths::Union{Nothing, AbstractVector},
    obs_exports::Union{Nothing, Real},
    obs_deaths::Union{Nothing, Real};
    pp_cases,
    obs_cases
) -> Makie.Figure

Prior predictive variant of plot_posterior_predictive, with the panel labels switched to "Prior".

source
BVDOutbreakSize.plot_start_date_pair Method
julia
plot_start_date_pair(chn; as_of_date, thin)

One-row, two-panel figure summarising when the outbreak began. The left panel is the posterior density of the outbreak start date, obtained by rescaling the days-since-seeding T to a calendar date (as_of_date minus T). The right panel is the joint (τ, T) posterior pair plot, which is positively correlated: slower growth (larger τ) needs a longer elapsed T to reach the same counts.

source
BVDOutbreakSize.posterior_summary Method
julia
posterior_summary(
    xs
) -> NamedTuple{(:lo90, :lo60, :lo30, :hi30, :hi60, :hi90), <:NTuple{6, Any}}

Return (lo90, lo60, lo30, hi30, hi60, hi90) equal-tailed credible interval endpoints from a vector of draws.

source
BVDOutbreakSize.predict_no_onward_deaths Method
julia
predict_no_onward_deaths(chn; obs_deaths,
                         alg = GaussLegendre(n = 64))

Per-draw projection of cumulative deaths under the counterfactual that every onward transmission stops at time T. Reads :r, :T, :α, :θ, :CFR from the posterior chn and integrates

ΔD=CFR0Trexp(rs)(1Fd(Ts))ds,

with F_d the Gamma(α, θ) onset-to-death CDF, returning a DataFrame with one row per draw:

  • :delta_deaths additional future expected deaths beyond obs_deaths

  • :total_projected obs_deaths + delta_deaths

obs_deaths is the number of deaths already observed at time T (e.g. obs.total_deaths from the bundled observations). alg is the quadrature scheme used for ΔD; defaults to GaussLegendre(n = 64), matching the rest of the package.

source
BVDOutbreakSize.streams_table Method
julia
streams_table(
    streams::Pair{String, <:AbstractVector}...;
    digits
) -> DataFrames.DataFrame

Side-by-side credible intervals for C_T from several fits. Pass each fit as "label" => draws_vector.

source
BVDOutbreakSize.summary_table Method
julia
summary_table(
    chn,
    params::AbstractVector{Symbol};
    digits
) -> DataFrames.DataFrame

DataFrame with one row per posterior parameter and the columns Quantity, Lower 90%, Lower 60%, Lower 30%, Upper 30%, Upper 60%, Upper 90% giving the lower and upper endpoints of the equal-tailed 30%, 60% and 90% credible intervals.

source
ChainRulesCore.rrule Method
julia
rrule(
    _::typeof(BVDOutbreakSize._gamma_cdf),
    α::Real,
    θ::Real,
    x::Real
) -> Tuple{Any, BVDOutbreakSize.var"#_gamma_cdf_pullback#_gamma_cdf_pullback##0"}

ChainRulesCore rrule for the gamma CDF, using the above partials.

NB: the NoTangent() is for the function argument itself, which is not a callable/functor. Using standard pullback convention where seed combines with the transposed Jacobian, although in this case all the partials are real scalars.

source