Skip to content

Release 1.2: The Multipathogen Update

Choose a tag to compare

@Julian-Patzner Julian-Patzner released this 19 Aug 13:13
· 86 commits to main since this release
29133d2

Welcome to the 1.2 release! This release turns GEMS from a single-disease simulator into a multipathogen framework: an arbitrary number of pathogens can now circulate, interact, and co-infect the same population within a single run. On top of that core change it introduces a composable transmission system, a wider set of shipped transmission functions and modifiers (cross-immunity, viral interference, and seasonal forcing), customizable infectiousness and immunity profiles, per-pathogen seeding and recurring case importation, a full separation of host-level care and mortality from per-pathogen disease progression, and a contact-tracing index that makes TraceInfectiousContacts three to four orders of magnitude faster.

Backwards compatibility was a design goal throughout. Old-format configs and single-pathogen scripts are auto-detected and rerouted through a self-contained legacy layer, so existing work should keep running unchanged. The breaking surface that remains is narrow and listed immediately below.


Upgrading

Most code needs no changes. These are the exceptions:

  • Symptomatic was renamed to Mild. The old name still works, since Symptomatic is retained as a delegating compatibility type, but prefer Mild in new code.
  • Hospitalized is no longer a disease progression tier. Hospitalization is now a host-level decision made by a HealthProgression. The Hospitalized category survives in the legacy layer for old configs. If you are porting a config to the new format rather than leaning on that layer, care parameters have two destinations: the top-level [HealthProgression] section, or, for a single-pathogen simulation, embedded directly on the Severe and Critical categories, from where they are collected into the default policy automatically.
  • Per-individual care accessors were removed from the public API: hospital_admission, hospital_discharge, icu_admission, icu_discharge, ventilation_admission, and ventilation_discharge. Care is no longer a property of one infection; use health_episodes(rd) for reconstructed stays.
  • The infections() DataFrame changed columns. Seven columns were dropped: the six per-infection care ticks (hospital_admission, hospital_discharge, icu_admission, icu_discharge, ventilation_admission, ventilation_discharge) and death, since care and mortality are host-level now. Use health_episodes(rd) for stays. Five columns were added: pathogen_id, progression_category, critical_onset, critical_offset, and removed, the tick at which the infection ended, whether by recovery or by the host death that cut it short. Post-processing that read the dropped columns must be updated.
  • compartment_periods changed columns. The hospitalized, icu, and ventilated durations were replaced by a single critical duration, because those are host-level periods rather than per-infection ones. The same applies to aggregated_compartment_periods. Care durations now come from health_episodes(rd).
  • TestType and SeroprevalenceTestType take a pathogen_id::Int8 instead of a Pathogen object, and the accessor is pathogen_id(tt) rather than pathogen(tt). Note that pathogen is still exported as a Simulation accessor, so calling it on a test type fails at the call site rather than at import.
  • Start conditions must name their pathogen when a simulation has more than one. Pass pathogen = "Covid19", or pathogen = "all" to seed all of them.
  • Source files moved. structs/ and methods/ were dissolved into domain folders; see Repository Restructure.

Multipathogen Core

The simulation no longer tracks a single infection per agent. Each individual can carry multiple concurrent infections, one per distinct pathogen, and the entire disease timeline moved off the Individual struct into dedicated infection records.

  • Bits-type InfectionState: An immutable, isbits record holding one infection's complete timeline (exposure, infectiousness_onset, symptom_onset, …, recovery), its pathogen_id, a precomputed per-tick infectiousness::Int8, and an active flag. ImmunityState mirrors this for immunity records, tagged with IMMUNITY_SOURCE_NATURAL or IMMUNITY_SOURCE_VACCINE.
  • Pathogen bitmasks, bounded at 32: The scalar infected::Bool semantics were replaced with active_pathogens_mask::UInt32 and detected_mask::UInt32 on Individual, where bit id - 1 flags active infection and detection for a pathogen. This gives O(1) "is this agent infected with pathogen X?" queries without touching a registry, and the width of those masks is what sets MAX_PATHOGENS = 32 (shared with the test-registry key packing). Pathogen ids are validated and auto-assigned as a contiguous 1..N range by _finalize_pathogen_ids!, whether pathogens arrive from a config file or as a user-supplied tuple.
  • Single-active-infection invariant: All mask mutation is routed through infected!(ind, pid, val) and detected!(ind, pid, val), giving one place to enforce the rule that infect! upholds: at most one active infection per (individual, pathogen) pair, with a concurrent re-infection by the same pathogen warned and skipped rather than overwriting the existing one.

Record Storage & Lock-Free Concurrency

Those records live in two places: a small inline cache on the agent for the common case, and three new registry types holding what does not fit. InfectionRegistry and ImmunityRegistry share one design; TestRegistry is deliberately different.

  • Cache + overflow storage: Each Individual carries a small inline infection_cache::NTuple{INFECTIONS_CACHE_SIZE, InfectionState} (and a parallel immunity_cache) plus an infection_head::Int32 / immunity_head::Int32 pointer into a registry-backed linked list holding any additional concurrent infections. The common single-infection case stays allocation-free and cache-local, while co-infection is unbounded. Both cache sizes are tunable constants (default 1).
  • Flat storage, index-linked chains: An InfectionRegistry is a states::Vector{InfectionState} plus a free_slots::Vector{Int32} free list, and ImmunityRegistry mirrors it for ImmunityState. The overflow chain hanging off infection_head is a linked list of indices into that vector, not heap-allocated nodes, so walking a co-infected agent's records stays inside one contiguous array. Registries can be pre-sized at construction with overflow_fraction, which reserves capacity for the expected share of co-infected agents, divided across shards.
  • Recycled slots: When an infection ends, its index is returned to free_slots and reused by the next insertion instead of appending to states. Registry size therefore settles at the peak number of concurrent overflow records, not the cumulative number of infections over the run.
  • Cache promotion: When an infection living in an individual's inline cache ends, the head of that individual's overflow chain is moved into the freed cache slot and its node returned to the free list. The hot inline slot stays occupied while the agent has any active infection, so the fast path keeps working for agents who have been co-infected and partly recovered.
  • Thread-sharded ownership: Each registry type is stored on the Simulation as a Vector sharded by _owner_shard(id) = mod(id - 1, maxthreadid()) + 1, so every agent deterministically belongs to exactly one shard. Element-wise reads and writes during the threaded update and spread phases are therefore contention-free without locking.
  • Deferred structural mutation: Linked-list insertion and removal mutate registry structure, so they are staged as _PendingInfection / _EndedInfection records and applied in a per-shard threaded flush (flush_pending_infections!, flush_ended_infections!) between phases. The parallel hot loops stay purely element-wise. Overflow unlinks are ordered before cache promotions, so an ended overflow record can never be resurrected into a freed cache slot.
  • Lazy iterators: each_infection and each_immunity return type-constrained InfectionIterator / ImmunityIterator values that resolve shards lazily and walk an individual's inline cache and overflow chain, presenting all active records through one uniform interface regardless of where they live.
  • TestRegistry is a dictionary, not a sharded vector: Test state is sparse, since only tested agents have any, so it is held as a Dict{UInt64, TestState} keyed by a packed composite of individual id and pathogen id. This is the packing that MAX_PATHOGENS bounds alongside the bitmasks. get_test_state(ind, reg, pid) returns a default TestState() for a pair that was never tested rather than inserting one, so reads never grow the dictionary.
  • Immunity accessors: natural_immunity_recorded / vaccine_immunity_recorded, natural_immunity_active / vaccine_immunity_active, natural_immunity_pending / vaccine_immunity_pending, and immunity_active distinguish an immunity that exists, one that is currently protective, and one still building up.

Type-Stable Pathogen Dispatch

Supporting many pathogens naively would put a dynamic dispatch on every contact. Instead, the pathogen set is encoded in the type system.

  • Fully parameterized Pathogen{PRG, PA, TF, IP, IM}: A pathogen carries the concrete types of its progressions tuple, progression-assignment function, transmission function, and infectiousness and immunity profiles, so every dispatch on it resolves to a concrete method.
  • Simulation{P<:Tuple, HP<:HealthProgression}: Pathogens are stored as a statically-typed heterogeneous tuple and the health policy as a second type parameter, making the whole pathogen and care configuration part of the simulation's type.
  • @generated lookup: get_pathogen(sim, pid) unrolls the tuple search at compile time, recovering the concrete Pathogen{…} type even when indexed by a runtime pid::Int8. The spread and progression loops specialize on each pathogen's concrete transmission, infectiousness, and immunity types instead of paying dynamic dispatch per contact.
  • Net effect: Adding pathogens does not add per-tick dispatch overhead. Because the same rework also removed dynamic dispatch that existed in the single-pathogen code, a default one-pathogen run is slightly faster than on 1.1.

Pathogens: API & Configuration

  • Simulation(pathogens = ...) accepts a single Pathogen, a Vector, or a Tuple. A Vector requires all elements to share one concrete type; for heterogeneous pathogens pass a Tuple, and the error message says so explicitly.
  • Old calls still work: pathogen = ... (singular) is upgraded to pathogens automatically, so existing single-pathogen scripts are unchanged.
  • Accessors: pathogens(sim), get_pathogen(sim, id), first_pathogen(sim), and pathogen(sim) for the single-pathogen case.
  • TOML schema: The config defines named [Pathogens.<Name>] tables, each with its own transmission_function, infectiousness_profile, immunity_profile, progressions, and progression_assignment blocks.
  • Indexed calibration paths: assign_values_to_parameters! now understands name[i] segments in a parameter path, so calibration can target one element of a collection (pathogens[2].transmission_function.transmission_rate) rather than only named fields. pathogen also resolves against a multi-pathogen simulation for single-pathogen paths.

Composable Transmission

Transmission was decomposed into a base rate × stackable modifiers architecture, separating TransmissionFunction (produces a base probability) from TransmissionModifier (adjusts it through a transmission_factor contract). TransmissionFunction now returns only the base rate; infectiousness and immunity are applied exactly once, centrally, rather than being folded in by each function.

  • CompositeTransmissionRate{B, M}: Wraps a base TransmissionFunction together with a statically-typed tuple of modifiers, multiplying in each modifier's transmission_factor without re-multiplying base rates. The tuple is part of the type, so the composition stays allocation- and dispatch-free.
  • CrossImmunityModifier / CrossImmunityTransmissionRate: Prior exposure or immunity to one pathogen partially protects against another. Configured as a Dict{String, Float64} of cross_immunities plus a default_cross_factor. Pathogens are referenced by name and resolved to ids at construction, so the configuration cannot silently depend on id ordering.
  • ViralInterferenceModifier / ViralInterferenceTransmissionRate: An active infection with one pathogen suppresses transmission of another. Takes interferences and persistence dictionaries (with defaults) so interference can outlast the interfering infection itself. It reads interfering state from the registry rather than same-tick staged state, making it independent of intra-tick spread order.
  • SinusoidalSeasonalModifier / SinusoidalSeasonalTransmissionRate: Sinusoidal seasonal forcing with an amplitude and a peak_day, for flu- and RSV-style seasonality.
  • ConstantTransmissionRate and AgeDependentTransmissionRate are retained, adapted to the new interface, and concretely typed.

Each modifier also ships as a standalone …TransmissionRate for the common case of one base rate plus one modifier, so simple setups do not need to build a composite by hand.


Infectiousness & Immunity Profiles

Pathogen is now a composition root: progression, transmission, infectiousness, and immunity are each independently swappable per pathogen.

InfectiousnessProfile (calculate_infectiousness) computes an individual's shedding level over the course of an infection from its InfectionState. The result is cached per tick on the agent, so the spread phase reads it directly.

  • ConstantInfectiousness: One fixed level across the whole infectious window.
  • StagedInfectiousness: A separate level per disease stage (asymptomatic, presymptomatic, symptomatic, severe, critical), which is what makes presymptomatic transmission modelling possible.
  • BetaInfectiousness: A continuous shedding curve with a time_to_peak and a concentration shape parameter, scaled by level and optionally adjusted per stage via asymptomatic_factorcritical_factor. The default concentration is drawn from the SARS-CoV-2 viral load and shedding kinetics literature (Eckerle et al.).

ImmunityProfile (calculate_immunity) defines how immunity is acquired and decays.

  • FullImmunity (the previous behaviour) and NoImmunity.
  • ExponentialWaning and SigmoidalWaning: Configurable halflife (default 180), a protection floor, a vaccine_buildup_duration ramp from vaccination to peak immunity, and, for the sigmoidal form, a hill shape parameter.
  • Independent barriers: Natural and vaccine-derived immunity are combined as independent protective barriers rather than one overwriting the other, and susceptibility_factor exposes the resulting protection to the transmission path.
  • immunity_is_stable: A hook that lets the engine clear needs_immunity_update once an individual's immunity has settled, skipping per-tick recomputation while it is constant.
  • Immunity can shape progression, not only susceptibility: infect! calls two immunity-aware entry points, each falling through to its immunity-blind form so they are opt-in. assign(individual, pa_func, immunities, pathogen_id, rng) lets an assignment function choose the progression category from the infectee's pre-exposure immunity, and calculate_progression(individual, tick, dp, immunities, pathogen_id, rng) lets a category shape the timings within it. Both can read same-pathogen protection through immunity_level(individual, immunities, pathogen_id) or prior exposure to other pathogens through each_immunity(individual, immunities), which is what makes a milder course for previously immune hosts expressible.

Host Health Progression: Decoupling Disease from Care

This is the conceptually largest change in the release. Hospitalization, ICU, ventilation, and death used to be baked directly into the disease-progression categories. That is coherent for one pathogen and incoherent for two: if pathogen A schedules a hospital admission and pathogen B independently schedules another, neither stay is authoritative, and two independently scheduled deaths for one host is not an outcome at all.

Disease progression now describes severity only: AsymptomaticMildSevereCritical. Everything host-level is decided separately.

  • HealthProgression and HealthProfile: A HealthProgression is a host-level policy that looks at all of an individual's active infections at once. select_health_profile(hp, infection) picks the profile for one infection's peak tier, and calculate_health_profile(profile, individual, infection, rng) turns it into a (CareContribution, HealthOutcome) pair. calculate_health_progression! folds those contributions into a single forward plan for the host.
  • CareContribution and CareLevel: Each infection contributes care demand at one or more of CARE_HOSPITAL, CARE_ICU, CARE_VENTILATION, as an admission/discharge tick pair. HealthOutcome carries a proposed death and the pathogen_id to attribute it to; combine_outcome resolves competing proposals to the earliest death, so mortality has one unambiguous cause.
  • Per-level demand counters: Individual carries hospital_demands, icu_demands, and ventilation_demands counters rather than a single care episode. A host is at a care level while at least one infection demands it, which makes overlapping stays from different pathogens compose correctly instead of one overwriting the other.
  • HealthSchedule: A tick-keyed schedule of CareTransition records (host, level, admission-or-discharge) with a wake_ticks set for scheduled deaths. Transitions are emitted per level and applied when due. This replaces an earlier single-episode host timeline that merged two disjoint care episodes into one incorrectly widened stay; separate episodes now stay separate. Schedules are sharded per thread on the Simulation alongside the registries.
  • Forward plans are validated, not clamped: calculate_health_progression! produces a plan that _validate_health_plan checks against the current tick, instead of silently clamping inconsistent ticks into range. Corrupted care and mortality output under co-infection was a real bug this fixed.
  • DefaultHealthProgression{S, C}: The shipped policy, splitting care by peak tier into SevereHealthProfile (ward admission only) and CriticalHealthProfile (hospital → ICU → ventilation, plus a death probability that is ungated by hospital or ICU admission, so a host can die without ever being admitted).
  • HealthLogger and health_episodes(): Care is logged as transitions and reconstructed into clean, joinable hospital, ICU, and ventilator episodes, with pp_health_episodes.jl and pp_hospital_df.jl backing the post-processing.
  • Config: A top-level [HealthProgression] section configures the policy, with severe and critical parameter blocks (hospital_probability, hospital_to_icu_probability, icu_to_ventilation_probability, death_probability, and the corresponding delay distributions).
  • Embedded care for simple cases: Rather than writing a policy, care parameters can be attached directly to the Severe and Critical categories, either as flat keyword arguments or as a prebuilt HealthProfile, and are harvested into a DefaultHealthProgression when the simulation is assembled. Harvesting is global across the pathogen set and allows only one embedded profile per tier, so it raises an ArgumentError if two pathogens both embed care for Severe. It is a single-pathogen convenience; multi-pathogen setups configure care through [HealthProgression] or their own policy.

Breaking changes

  • SymptomaticMild, and Hospitalized is no longer a disease tier. Both remain available through the legacy layer.
  • The Critical progression's parameters changed shape: severeness_onset_to_hospital_admission, hospital_admission_to_icu_admission, icu_admission_to_icu_discharge, icu_discharge_to_hospital_discharge, icu_admission_to_death, and death_probability moved out of the category and into the health progression, replaced by severeness_onset_to_critical_onset, critical_onset_to_critical_offset, and critical_offset_to_severeness_offset.
  • The per-individual care accessors (hospital_admission, icu_admission, ventilation_discharge, …) are no longer exported. Use health_episodes(rd).
  • The HOSPITAL_STATUS_* constants were removed.

The legacy compatibility layer

src/pathogen/legacy_progressions/ reproduces the pre-decoupling behaviour for old configs and code, and is deliberately isolated rather than scattered through the codebase: the folder and its single include can be deleted outright to remove it.

  • Symptomatic is a delegating wrapper around Mild that accepts the same parameters and produces an identical progression.
  • Hospitalized resolves directly to its legacy type.
  • Critical is the hard case, because the name is reused rather than renamed. _is_legacy_critical detects an old-format Critical block by its parameters and _normalize_legacy_pathogen! reroutes just that block to LegacyCritical. A new-format Critical is left untouched even when the same config uses legacy names for other tiers.
  • LegacyHealthProgression / LegacyCriticalHealthProfile reproduce the old coupled care course: deterministic hospital → ICU admission with death anchored to the drawn ICU-admission tick, preserving the previous ICU/death coupling.

See Custom Health Progression for writing your own policy.


Seeding, Start Conditions & Imported Cases

  • Pathogen-aware start conditions: InfectedFraction, PatientZero, PatientZeros, and RegionalSeeds all take a pathogen name and seed that pathogen. Previously a pathogen argument was accepted and silently ignored with a warning that GEMS only supported single-pathogen simulations.
  • Explicit selection, resolved once: A condition names a pathogen, uses ALL_PATHOGENS ("all"), or leaves the name empty for a single-pathogen simulation. _expand_pathogens resolves the selection against the simulation's pathogens at construction: ALL_PATHOGENS expands into one condition per pathogen, a name is validated against the pathogen set, and an empty name is rejected with an actionable error when more than one pathogen exists. Because this happens at construction, the config-file and programmatic paths agree.
  • MultiStartCondition: A composite condition wrapping one sub-condition per pathogen, initialized in order, for example a different infected fraction per disease in one run. It is constructed automatically when multiple pathogens are supplied.
  • ImportedCases: A new start condition for recurring case importation, which previously had no equivalent. ticks accepts an explicit vector, a schedule, or a window with an interval and offset that may each be a distribution; count accepts a fixed number, a vector, or a distribution; and ags optionally restricts imports to a region or set of regions. An import only lands on a host that is alive, not hospitalized, not self-isolating, and not already carrying that pathogen. It makes multi-season and endemic scenarios expressible without driving importation from outside the simulation.

Contact Tracing: InfecterIndex

get_infections_between, which backs the contact-tracing measures, linearly scanned every thread's infection log shard on every call. It was the one hot path deliberately left out of the 1.1 intervention-optimization pass because it needed the most implementation work of that batch.

  • InfecterIndex: An inverted index mapping each infecter to a tick-ordered list of its infectees. get_infections_between walks the relevant slice of that list instead of scanning the log.
  • Built lazily, maintained incrementally: The index is constructed on first use, so runs without contact tracing pay nothing. New links are staged per thread as _PendingLink records and merged in with a k-way merge that preserves tick ordering, keeping the parallel logging path lock-free.

Measured impact

Average process_measure runtime for TraceInfectiousContacts, measured on a 16-thread desktop CPU using the default configuration and the shipped population files:

Population Speedup
Saxony 490×
North Rhine-Westphalia 1,750×
Germany ~4 orders of magnitude (extrapolated)

The scaling is expected: the old implementation's cost grew with the size of the infection log, while an indexed lookup grows with the number of infectees actually traced. The Germany figure is extrapolated from the Saxony and NRW measurements, not measured directly.


Pathogen-Aware Analysis, Reporting & Plotting

  • Pathogen-tagged loggers: Infection, test, vaccination, and health events carry a pathogen_id. The logging layer was reorganized into logging/{logger, event_loggers, infection_logger, tick_loggers, infecter_index}.jl.
  • Multipathogen post-processing: The full pp_* suite (incidence, attack rates, effective R, R0, serial and generation intervals, tick cases, deaths, tests, and detection) is pathogen-aware, splitting and grouping results per pathogen.
  • Per-pathogen aggregation: aggregate_by_pathogen was added, and ResultData and Batch aggregation were updated to render and store per-pathogen results across runs (also required because Simulation is now parameterized).
  • Plot filters: Every built-in plot gained a pathogen keyword accepting an id or a name, defaulting to showing all pathogens with an automatic per-pathogen subplot layout.

Repository Restructure

A large structural cleanup with no behavioural intent, reorganizing the source tree around its actual domains. Worth knowing about if you are looking for a file at its old path.

  • structs/ and methods/ were dissolved into pathogen/ (pathogens, progressions, health progression, transmission functions, profiles, vaccines), population/, settings/, registries/, simulation/ (simulation, infections, calibration, batch, initialization, termination), contacts/, analysis/ (formerly model_analysis/), and logging/ (formerly logger/).
  • agents.jl was split into population/individuals.jl (the struct) and population/individual_methods.jl (behaviour), with setting membership and lookup moved onto the individual.
  • Start conditions and stop criteria each moved into their own file under simulation/initialization/start_conditions/ and simulation/termination/stop_criteria/.
  • The public API surface was tightened further: internal helpers lost their exports and gained _ prefixes (e.g. _basefolder).

Reproducibility Notes

Seeded runs from 1.1.x will not reproduce bit-identically on 1.2.0. The infection model, the order in which random draws are made during progression and care assignment, and the separation of care from disease progression all change the draw sequence. Results remain equally valid; they are simply a different realization. Configs that relied on the coupled Critical behaviour and are rerouted through the legacy layer reproduce the old behaviour, not the old random stream.


Testing

  • test/registriestest.jl (new) covers cache and overflow promotion, sharded ownership, deferred flush ordering, and iterator traversal.
  • test/healthprogressiontest.jl (new) covers the contribution and schedule API, including a regression test that two disjoint care contributions produce two episodes rather than one widened stay.
  • test/pathogentest.jl was substantially expanded for the new transmission functions and modifiers, all three infectiousness profiles, the waning immunity profiles, and the single-active-infection guard.
  • test/testutils.jl is new, and agentstest, infectionstest, interventionstest, simulationtest, settingstest, postprocessortest, loggertest, and utilstest were updated for the per-(individual, pathogen) model, the bitmask accessors, and the new start conditions.