Function reference

daedalus(country, infections::Vector{InfectionData}; npi, log_rt, time_end, increment, n_threads)

Run the daedalus epidemic model with multiple infection parameters in parallel.

This is the vector-InfectionData dispatch of the daedalus function. It orchestrates multi-run ODE solving across multiple InfectionData objects (extracting R0 values) using shared computation and optional multi-threading.

Arguments

• country: Country to model (String or CountryData struct)
• infections::Vector{InfectionData}: Vector of infection parameter objects
• npi: Non-pharmaceutical intervention (optional)
• log_rt: Whether to log effective reproduction number (default: true)
• time_end: Final simulation time in days (default: 100.0)
• increment: Savepoint interval in days (default: 1.0)
• n_threads: Number of threads for parallel execution (default: 1)

Returns

A Vector{NamedTuple} where each element contains:

• sol: The ODE solution for that infection
• saves: Saved event values (if applicable)
• npi: The NPI specification used
• r0: The R0 value for this infection

daedalus_internal(n_runs, shared_data, param_sets, cb_set)::EnsembleSolution

Internal function that orchestrates multi-run ODE solving using EnsembleProblem.

Uses SciMLBase.EnsembleProblem with EnsembleThreads() solver for efficient parallel execution across multiple parameter sets. The base ODE problem is created once, then remade for each trajectory with its corresponding parameters.

Arguments

• n_runs::Int: Number of parameter sets / trajectories
• shared_data::NamedTuple: Pre-computed country data from prepare_shared_data()
• param_sets: Array of Params structs, one per trajectory
• cb_set::CallbackSet: Callback set (pre-constructed by caller)

Returns

An EnsembleSolution containing solutions for all trajectories at specified savepoints

Details

The function constructs an EnsembleProblem with:

• A single base ODE problem created from the first parameter set
• A prob_func that remakes the problem with the i-th parameter set
• EnsembleThreads() solver for multi-threaded execution (respects Julia's thread count)
• Savepoints from shared_data.savepoints

extract_infection_params(infection::DataLoader.InfectionData)::NamedTuple

Extract and process epidemiological parameters from an InfectionData object.

Converts base parameters to model-compatible form:

• Expands age-varying parameters (eta, omega, gamma_H) to work with both age groups (4) and economic groups (45), matching the 49-component state vector
• Computes working-age indices for proper parameter replication

Arguments

• infection::DataLoader.InfectionData: Infection parameters to extract

Returns

A NamedTuple with fields:

• r0::Float64: Basic reproduction number
• sigma::Float64: Incubation rate (1/latent period)
• p_sigma::Float64: Proportion symptomatic (presymptomatic transmission scaling)
• epsilon::Float64: Relative transmissibility of asymptomatic cases
• rho::Float64: Rate of infectiousness decay in hospitalized
• eta::Vector{Float64}: Hospitalization rate (expanded to 49 groups)
• omega::Vector{Float64}: Recovery rate for asymptomatic (expanded to 49 groups)
• gamma_Ia::Float64: Recovery rate from asymptomatic infection
• gamma_Is::Float64: Recovery rate from symptomatic infection
• gamma_H::Vector{Float64}: Recovery rate for hospitalized (expanded to 49 groups)
• nu::Float64: Fixed at 0.0 (not user-configurable)

prepare_shared_data(cd::DataLoader.CountryData, time_end::Float64,
    increment::Float64)::NamedTuple

Prepare country/structure data that is computed once and reused across all parameter sets.

Expensive operations (contact matrix processing, demographics aggregation) are performed once rather than repeated for each r0 value in ensemble runs.

Arguments

• cd::DataLoader.CountryData: Country demographic and economic data
• time_end::Float64: Final simulation time
• increment::Float64: Savepoint interval (typically 1.0 for daily output)

Returns

A NamedTuple with fields:

• init_state: Initial epidemic state (with appended Rt slot)
• contacts_unscaled: Raw contact matrix
• cw: Worker contact rates
• contacts_array: 3D contact array for ODE
• cm_temp: Temporary contact matrix workspace
• settings: Number of contact settings
• demog: Population vector
• timespan: Tuple (0.0, time_end)
• savepoints: Range of savepoint times

CountryData

Demographic and economic data for a single country, matching the structure of the country_data list in the R package {daedalus.data}.

This struct is mutable to allow reassigning values to examine different scenarios.

EconomicContacts

Per-sector social contact rates for economic sectors.

InfectionData

Epidemiological parameters for a single pathogen, matching the structure of infection_data in the R package {daedalus.data}.

This struct is mutable to allow reassigning values to examine different scenarios.

VaccinationScenario

Parameters for a single pre-pandemic vaccine investment scenario.

Compute the 4×4 population-weighted contact matrix from 16-bin CM values.

Replicates the weighted aggregation in the R data-raw/country_data.R script: contacts4x4[i, j] = sum{to in groupi, from in groupj}(CM[to,from] * pop[to]) / sum{to in groupi} pop[to]

Aggregate a 16-element vector of 5-year-bin populations into 4 DAEDALUS groups. Bins 1-16 correspond to ages 0-4, 5-9, ..., 75-79. Bins 14-16 map to 65+.

get_closure_strategy(name::String) -> Vector{Float64}

Return a 45-element vector of sector openness coefficients (0–1) for the named pandemic response strategy. Call list_closure_strategies for available names.

get_country(name::String) -> CountryData

Return demographic and economic data for the named country. Call list_countries to see available countries.

get_economic_contacts() -> EconomicContacts

Return within-sector and between-sector social contact rates for all 45 economic sectors. Data sourced from sectorcontacts.csv.

get_pathogen(name::String) -> InfectionData

Return infection parameters for the named pathogen. Call list_pathogens to see available names.

get_sector_names() -> Vector{String}

Return the names of all 45 economic sectors.

get_vaccination_scenario(name::String) -> VaccinationScenario

Return parameters for the named vaccination scenario. Call list_vaccination_scenarios to see available names.

list_closure_strategies() -> Vector{String}

Return names of available pandemic response strategies.

list_countries() -> Vector{String}

Return a sorted list of all available country names.

list_pathogens() -> Vector{String}

Return the names of all available pathogens.

list_vaccination_scenarios() -> Vector{String}

Return names of all vaccination scenarios: "none", "low", "medium", "high".

consumer_worker_contacts(demography; scaled=true)::Matrix{Float64}

Get the consumer-worker contact matrix (45 sectors × 4 age groups).

Represents contacts between consumer-sector workers and age groups.

Arguments

• demography::Vector{Float64}: Population vector for 4 age groups
• scaled::Bool: If true (default), contacts are scaled by age-group demographics

Returns

A Matrix{Float64} of size (45, 4) with consumer-worker contact rates

consumer_worker_contacts(cd::CountryData; scaled=true)

Get the 45×4 consumer-worker contact matrix scaled by the age-group demography in cd.

contacts3d(cd::CountryData) -> Array{Float64, 3}

Return the scaled contact matrices for a country as a single 3D array of size (N_TOTAL_GROUPS, N_TOTAL_GROUPS, K), where K is the number of closure settings (economic contact scenarios).

Calls prepare_contacts(cd) and normalises the result into 3D form:

• If the country has multiple contact matrices (one per closure setting), they are stacked along the third dimension to give shape (49, 49, K).
• If the country has a single contact matrix, it is reshaped to (49, 49, 1) so the return type is always a 3D array regardless of K.

The third dimension of the returned array is the one contracted by weighted_slice_sum! in the ODE force-of-infection calculation.

Arguments

• cd::CountryData: Country data struct from DataLoader

Returns

An Array{Float64, 3} of size (N_TOTAL_GROUPS, N_TOTAL_GROUPS, K).

expand_contacts(cm::Matrix{Float64})::Matrix{Float64}

Expand a 4×4 contact matrix to 49×49 covering all age groups and economic sectors.

Maps the 4×4 age-group contact matrix to the full 49×49 population, assuming economic sectors adopt the contact patterns of the working-age group.

Arguments

• cm::Matrix{Float64}: 4×4 contact matrix for age groups

Returns

A Matrix{Float64} of size (49, 49)

get_settings(cd::CountryData)::Int

Return the number of contact matrix settings (closure strategies) for a country.

Arguments

• cd::CountryData: Country data struct

Returns

The number of contact matrix settings: typically 1 (single matrix) or >1 (multiple scenarios)

initial_state(cd::CountryData) -> Array{Float64,3}

Construct the initial epidemic state (NTOTALGROUPS × NCOMPARTMENTS × NVACCINE_STRATA) for a country given as a DataLoader.CountryData struct. A fraction 1e-6 of each group is seeded as symptomatic infectious.

initial_state(country::String) -> Array{Float64,3}

Construct the initial epidemic state for a named country. Calls DataLoader.get_country and delegates to initial_state(::CountryData).

Example

state = initial_state("United Kingdom")

prepare_community_contacts(cm; scaled=true)::Matrix{Float64}

Get a 49×49 community-only contact matrix for all age-groups and economic sectors. Unlike prepare_contacts, this function does not add within-sector workplace contacts or consumer-worker contacts to the matrix. Those routes are kept separate for use in the ODE force-of-infection calculation (see plan_ode.md).

Arguments

• cm::Matrix{Float64}: 4×4 contact matrix for age groups
• scaled::Bool: If true (default), contact matrix is scaled by demographics

Returns

A Matrix{Float64} of size (49, 49) with community contacts only

prepare_community_contacts(cd::CountryData; scaled=true)

Get a 49×49 community-only contact matrix for all age-groups and economic sectors using cd. Workplace and consumer-worker contacts are excluded.

prepare_community_contacts(country::String; scaled=true)

Get the 49×49 community-only contact matrix for a named country.

prepare_contacts(cd::CountryData; scaled=true)

Get a 49×49 contact matrix for all age-groups and economic sectors using the demographic and contact data in cd.

prepare_contacts(country::String; scaled=true)

Get the 49×49 contact matrix for a named country.

prepare_demog(demog, workers)::Vector{Float64}

Get a 49-element population vector for all age-groups and economic sectors.

Concatenates 4 age groups with 45 economic sector worker counts to form the population vector used in force-of-infection calculations.

Arguments

• demog::Vector{Float64}: Demographics (4 age groups)
• workers::Vector{Int}: Worker counts (45 economic sectors)

Returns

A Vector{Float64} of length 49 (4 age groups + 45 workers)

prepare_demog(cd::CountryData) -> Vector

Get the 49-element population vector (4 age groups + 45 worker sectors) for any country from a DataLoader.CountryData struct.

Worker counts are clamped to a minimum of 1 to avoid division-by-zero when this vector is used as a denominator (e.g. in contact-matrix scaling and the Rt callback). This matches the +1 padding applied in initial_state.

prepare_demog(country::String) -> Vector

Get the 49-element population vector for a named country.

total_contacts(contacts) -> Matrix{Float64}

Reduce a contact matrix representation to a single aggregated 2D matrix.

When contacts is a Vector{Matrix{Float64}} (one matrix per closure setting), the matrices are summed element-wise to produce the total contact matrix. When contacts is already a single Matrix{Float64}, it is returned unchanged.

Arguments

• contacts: Either a Vector{Matrix{Float64}} or a Matrix{Float64}

Returns

A single Matrix{Float64} representing total contacts across all settings.

worker_contacts(workers; scaled=true)::Vector{Float64}

Get per-capita social contacts within each economic sector.

Data sourced from sectorcontacts.csv via DataLoader. When scaled=true (default), values are divided element-wise by sector workforce counts so that the result is contacts per worker (as used in the ODE force-of-infection).

Arguments

• workers: Worker counts for each of the 45 economic sectors
• scaled::Bool: If true (default), contacts are divided by worker counts

Returns

A Vector{Float64} of length 45 with per-capita within-sector contact rates

worker_contacts(cd::CountryData; scaled=true) -> SVector

Get per-capita social contacts within each economic sector using workforce counts from cd. Sectors with zero workers are treated as having 1 worker to avoid division by zero.

worker_contacts(country::String; scaled=true) -> SVector

Get per-capita within-sector contact rates for a named country.

daedalus_ode!(du, u, p, t)

Compute derivatives for the DAEDALUS epidemic ODE system.

The core compartmental model with 7 compartments (S, E, Is, Ia, H, R, D) across 49 demographic groups (4 age groups + 45 economic sectors) and 2 vaccine strata. Implements force-of-infection, disease progression, vaccination, and death.

Arguments

• du::Array: Pre-allocated output array (derivatives), same shape as u
• u::Array: State vector (length = NTOTALGROUPS × NCOMPARTMENTS × NVACCINE_STRATA + 1)
• p::Params: Parameter struct containing contact matrices, rates, and NGM
• t::Number: Current simulation time

State layout

The state array is reshaped internally as (49, 7, 2) where:

• Dimension 1: Demographic groups (4 age + 45 economic sectors)
• Dimension 2: Compartments (S, E, Is, Ia, H, R, D)
• Dimension 3: Vaccine strata (unvaccinated, vaccinated)
• Final element: Effective reproduction number Rt (updated by callbacks)

Details

This function is called internally by OrdinaryDiffEq's ODE solvers. It computes force-of-infection from the contact matrix and current infectious compartments, then updates all compartments according to the SEIR-like disease dynamics.

Npi

A container for parameter modifications (effects), each with trigger conditions.

Fields

• effects::Vector{Effect}: Vector of parameter modifications to apply

Constructor

Npi(effects::AbstractVector{<:Effect})

Example

trigger_on = ReactiveTrigger(5000.0, "H")
trigger_off = ReactiveTrigger(1.0, "Rt")
effect = ParamEffect(
    :beta, x -> x .* 0.4, x -> x ./ 0.4, trigger_on, trigger_off)
npi = Npi([effect])

ParamEffect

A parameter modification triggered by state or time conditions.

Fields

• target::Symbol: Parameter to modify (e.g., :beta, :omega)
• func::Function: Applied when activated
• reset_func::Function: Applied when deactivated
• trigger_on::Trigger: Activation condition (ReactiveTrigger or TimeTrigger)
• trigger_off::Trigger: Deactivation condition
• saved_values::SavedValues: Internal state tracking
• ison::Bool: Current activation state

Constructor

ParamEffect(target::Symbol, func::Function, reset_func::Function,
            trigger_on::Trigger, trigger_off::Trigger)

Example

trigger_on = ReactiveTrigger(5000.0, "H")
trigger_off = ReactiveTrigger(1.0, "Rt")
effect = ParamEffect(:beta, x -> x .* 0.4, x -> x ./ 0.4, trigger_on, trigger_off)

Params

A mutable struct that holds the parameters for the DAEDALUS model. It contains the contact matrices, contact weights, and various epidemiological parameters.

ReactiveTrigger

A trigger based on epidemic state (e.g., hospitalization count, Rt).

Fields

• name::String: Compartment name (e.g., "H", "Rt")
• value::Float64: Threshold value for triggering

Constructor

ReactiveTrigger(value::Float64, name::String)

TimeTrigger

A trigger based on simulation time.

Fields

• name::String: Always "time"
• value::Float64: Time point (days) for triggering

Constructor

TimeTrigger(value::Float64, name::String = "time")

make_events(npi::Npi, savepoints)::CallbackSet

Create a CallbackSet of callbacks for all effects in an NPI.

Arguments

• npi::Npi: Npi struct with effects
• savepoints: Time points for state checking (for reactive triggers)

Returns

A CallbackSet with all activation/deactivation callbacks

make_events(eff::ParamEffect, savepoints)::CallbackSet

Create a CallbackSet for an effect's activation and deactivation callbacks.

Arguments

• eff::ParamEffect: Effect with triggers
• savepoints: Time points for state checking (used if trigger is ReactiveTrigger)

Returns

A CallbackSet with on/off callbacks

make_param_changer(eff::Effect)::Function

Create a callback function that modifies a parameter when activated.

Handles both ReactiveTrigger (state-dependent) and TimeTrigger (time-based) activation.

Arguments

• eff::Effect: Effect with trigger conditions

Returns

A Function that modifies the target parameter in-place

make_param_reset(eff::Effect)::Function

Create a callback function that resets a parameter when deactivated.

Handles both ReactiveTrigger (state-dependent) and TimeTrigger (time-based) deactivation.

Arguments

• eff::Effect: Effect with trigger conditions

Returns

A Function that resets the target parameter in-place

make_rt_logger(savepoints)

Create a callback that computes and logs the effective reproduction number Rt.

Uses power iteration on the susceptibility-adjusted NGM to compute Rt at each savepoint. Maintains a warm-start eigenvector for improved convergence.

Arguments

• savepoints: Time points at which to compute and log Rt

Returns

A PresetTimeCallback that updates the Rt state variable at each savepoint

make_save_events(npi::Npi, savepoints)

Create SavingCallbacks for all effects that require state tracking.

Returns a vector of SavingCallbacks, one per effect with at least one ReactiveTrigger. Effects with only TimeTrigger conditions are skipped.

Arguments

• npi::Npi: Npi struct with effects
• savepoints: Time points at which to save state values

Returns

A Vector{SavingCallback} — one per effect that has reactive triggers

make_save_events(eff::ParamEffect, savepoints)::Union{SavingCallback, Nothing}

Create a SavingCallback for an effect that records trigger compartment values.

Skips effects with all TimeTrigger conditions (no state tracking needed). For effects with at least one ReactiveTrigger, captures trigger compartment values at each savepoint.

Arguments

• eff::ParamEffect: Effect with trigger conditions
• savepoints: Time points at which to save state values

Returns

A SavingCallback for reactive effects, or nothing if both triggers are time-based

dominant_eigenvalue(A; v_init=nothing, max_iter=100, tol=1e-6)

Compute the dominant (largest magnitude) eigenvalue using power iteration. This is significantly faster than computing all eigenvalues when only the largest eigenvalue is needed.

Arguments

• A::AbstractMatrix: The matrix whose dominant eigenvalue is to be computed
• v_init::Union{Nothing,AbstractVector}: Optional initial vector for iteration. If not provided, a random vector is used.
• max_iter::Int: Maximum number of iterations (default: 100)
• tol::Float64: Convergence tolerance (default: 1e-6)

Returns

• The dominant eigenvalue (largest in magnitude)

Algorithm

Uses the power iteration method: repeatedly multiplies a vector by the matrix and normalizes. The Rayleigh quotient converges to the dominant eigenvalue.

Performance

For n×n matrices, this method is O(kn²) where k is the number of iterations (typically k << n), compared to O(n³) for full eigendecomposition.

Example

A = rand(49, 49)
λ_dom = dominant_eigenvalue(A)

get_beta(cm, r0, sigma, p_sigma, epsilon, gamma_Ia, gamma_Is)::Union{Float64, Vector{Float64}}

Calculate transmission rate(s) beta given a contact matrix and epidemiological parameters.

Uses the next-generation matrix (NGM) method to compute beta values that achieve the target basic reproduction number(s) R₀.

Arguments

• cm::Matrix{Float64}: Contact matrix (ngroups × ngroups)
• r0::Union{Float64, Vector{Float64}}: Target basic reproduction number(s)
• sigma::Float64: Rate of progression from exposed to infectious
• p_sigma::Float64: Proportion of exposed who develop symptoms
• epsilon::Float64: Reduction factor for asymptomatic transmission
• gamma_Ia::Float64: Recovery rate for asymptomatic infectious
• gamma_Is::Float64: Recovery rate for symptomatic infectious

Returns

• If r0 is a scalar: returns a single Float64 value for beta
• If r0 is a vector: returns a Vector{Float64} of beta values

get_ngm(cm, beta, sigma, p_sigma, epsilon, gamma_Ia, gamma_Is)::Union{Matrix{Float64}, Vector{Matrix{Float64}}}

Calculate the next-generation matrix (NGM) given a contact matrix and transmission rate(s).

The NGM is used to compute the effective reproduction number Rt accounting for current susceptibility levels in the population.

Arguments

• cm::Matrix{Float64}: Contact matrix (ngroups × ngroups)
• beta::Union{Float64, Vector{Float64}}: Transmission rate(s)
• sigma::Float64: Rate of progression from exposed to infectious
• p_sigma::Float64: Proportion of exposed who develop symptoms
• epsilon::Float64: Reduction factor for asymptomatic transmission
• gamma_Ia::Float64: Recovery rate for asymptomatic infectious
• gamma_Is::Float64: Recovery rate for symptomatic infectious

Returns

• If beta is a scalar: returns a single Matrix{Float64}
• If beta is a vector: returns a Vector{Matrix{Float64}} of NGMs

get_ngm(cm::Matrix, beta::Vector, sigma, p_sigma, epsilon, gamma_Ia, gamma_Is)::Vector{Matrix{Float64}}

Vector-dispatch version: compute one NGM for each transmission rate in beta.

sum_by_age(state, index)::Vector{Float64}

Sum state array over economic groups into age groups for a given compartment.

Takes the full state array and sums all economic sector rows into the working-age group, producing a 4-element vector (one per age group).

Arguments

• state: 3D state array (NTOTALGROUPS × NCOMPARTMENTS × NVACCINE_STRATA)
• index: Compartment index (integer)

Returns

A Vector{Float64} of length 4 (one value per age group) with economic sectors summed into the working-age group.

weighted_slice_sum!(X, v, result)

Compute the weighted sum of 2D slices of a 3D array in-place, writing the result into result.

Mathematically: result[m, k] = Σᵢ X[m, k, i] * v[i]

This is a tensor contraction over the third dimension of X, equivalent to a matrix-vector product after reshaping. X is reinterpreted (without copying) as a (M*K) × N matrix whose columns are the flattened slices X[:,:,i]; BLAS gemv then computes the weighted sum in a single call.

Arguments

• X::Array{T,3}: Three-dimensional input array of size (M, K, N)
• v::Vector{T}: Weight vector of length N; element v[i] scales slice X[:,:,i]
• result::Array{T,2}: Pre-allocated output array of size (M, K); overwritten in-place

Notes

• reshape and vec are O(1) view operations — no data is copied.
• Efficient for large M, K. For small N (e.g. N_VACCINE_STRATA = 2), the BLAS dispatch overhead may dominate; a plain @inbounds @simd loop can be faster in that regime.

get_times(output)

Return the unique timepoints at which the model state was saved.

Arguments

• output: Named tuple returned by daedalus(), with field sol (the ODE solution).

Returns

A sorted Vector{Float64} of unique save times (e.g. [0.0, 1.0, 2.0, ..., tmax] for a daily-increment run).

get_values(output, comp::String, timebin::Int=90)

Extract values for a given compartment from a daedalus() model output, with optional aggregation into time bins.

Assumes the model was solved with daily savepoints (integer timesteps from 0).

Arguments

• output: Named tuple returned by daedalus(), with field sol (the ODE solution).
• comp: Compartment name. One of: "S", "E", "Is", "Ia", "H", "R", "D", or "Rt".
• timebin: Number of days per aggregation bin. Use 1 to return unaggregated daily values (default: 90).

Returns

A Vector{Float64} where each element is the sum of the compartment values across all age and vaccination groups:

• timebin == 1: one value per day, from day 0 to day tmax (length tmax + 1).
• timebin > 1: one value per non-overlapping window of timebin days. The final bin covers the remaining days and may be shorter than timebin.

Examples

output = daedalus(country="Australia", time_end=365.0)

# Daily hospitalizations
daily_H = get_values(output, "H", 1)

# Quarterly totals (90-day bins)
quarterly_H = get_values(output, "H", 90)

# Daily effective reproduction number
rt = get_values(output, "Rt", 1)

get_indices(compartment, groups=nothing)

Return the linear indices into the ODE state vector for compartment and groups.

The state vector is laid out as N_VACCINE_STRATA repetitions of N_COMPARTMENTS blocks, each block containing N_TOTAL_GROUPS consecutive entries (one per demographic/economic group). The Rt scalar parameters are appended after all epi-compartment entries.

Arguments

• compartment: compartment name — one of S, E, Is, Ia, H, R, D or "Rt".
• groups: groups to select within the compartment block:
  • nothing (default) — all groups (1:N_TOTAL_GROUPS)
  • Int — a single group
  • AbstractVector{Int} or UnitRange{Int} — an arbitrary subset of groups

Returns

A UnitRange{Int} when all groups are requested or a contiguous range is passed, otherwise a Vector{Int}.