Skip to main content
🧠0% 🤖100%

odin: age-structured SEIR with a delay

A deterministic SEIR in 10 million people in five age groups (0–1, 1–5, 5–15, 15–30, and 30+ years), over one year in days, with births, deaths, aging, waning immunity, and a latent period of exactly 14 days written as a delay differential equation. It is odin’s own age-structured example (the README’s “more complicated example”, array_odin.R, with E and a delay added), and one of the six models odin’s tests check against hand-written deSolve code.

This is odin’s wave 3 example. model.py is the Starsim version, and compare.py checks it against reference/.

Comparison

Measure Original Starsim
Lines of model code 53 18 No line is over 100 characters
One run first 0.71 s, rerun 0.0014 s first 0.019 s, rerun 0.0032 s (Rust) The original’s first run is mostly compiling the C code. Its run of the delay and Starsim’s conveyor are both quick
Results, in total Agree. S, E, I, and R within 0.23%, 0.36%, 0.25%, and 0.15% of their largest values on every day (a quarter of that with a step of half a day), and I peaks at 6.19 million on day 125, as in the original. It needs the workaround of finding 1: without it the epidemic peaks 15 days late. compare.py allows 0.5%
Results, by age Agree, except in the youngest group, where E, I, and R differ by up to 15%, 5%, and 8%. The original’s delay does not let the latent age (finding 2); Starsim’s does. compare.py allows 1% in the other groups

Measured on 9 November 2026 with Starsim 4.0.0rc3 (the v4 branch, commit be38b9f6), on a 20-core Linux machine.

What was declared and what was not

Everything is declared: seven transitions and one stratification, whose aging between groups is generated, and the delay is a fixed dwell time. The stable age distribution is worked out by hand in two lines of NumPy, as the original works it out in two lines of odin (den[1] and den[2:N_age]). One transition, history, is a workaround (finding 1).

Original Starsim
dim(S) <- N_age and so on, for 9 arrays seir.stratify(age=[0, 1, 5, 15, 30], fractions=split)
age_width, age_rate, and the aging terms in all eight deriv() equations Generated by the age stratification, at one over each width
new_inf[] <- beta * S[i] * I_tot / N with I_tot <- sum(I) S >> E @ ss.perday(1) * I
lag_inf[] <- delay(new_inf[i] * surv, lat_hum) E >> I @ ss.constant(ss.days(14)), with deaths from E doing the work of surv
Births <- b * N into S[1], and - b * X[i] in every equation ss.birth(S, ...) and ss.death(...)
den[] split, worked out the same way
delay()’s value before day 14 history, an inflow to I for 14 days (finding 1)

Findings

  1. The history of a delay (blocks agreement; worked around). For the first 14 days, odin’s delay(new_inf * surv, 14) looks back to before day 0 and takes the value at day 0, so about one person a day reaches I as if infections had been going on before the start. A model whose latent state starts empty cannot do this, and with a reproduction number near 30, the 14 extra cases bring the epidemic forward by 15 days (without them, I differs by up to 66% of its peak). The workaround is a transition that adds them: history = ss.null >> I @ ss.timeseries({0: ss.perday(surv), 14: 0}). It adds 14 people that the original takes from E, a difference of one part in a million. Suggested change: let init_n fill a fixed delay as if people had been arriving at a constant rate (for example init_n=dict(E=ss.arrived(rate=1))), since the conveyor below has a natural place for them.

  2. The original’s latent do not age (disagreement by age, not a Starsim issue). In the original, E is aged like every state, but lag_inf[i] takes from E[i] everyone infected in group i 14 days earlier, including the 1 in 26 of the youngest who have aged out since. So E[1] falls short and I[1] gets more than it should; summed over ages nothing is lost. Starsim’s conveyor ages with everyone else. This was checked with a separate Python version of the original’s equations: as written, it matches the reference to within 0.2% in every group; with the latent aging, it differs from the reference in the youngest group by about what Starsim does. Both are reasonable, and the difference comes from how the original is written.

  3. Fixed: a fixed delay was many stages, which were slow and still approximate. The delay used to be 200 Erlang stages, which took 1.6 s a run and was 1.6% off in I. A fixed dwell time, ss.constant(ss.days(14)), is now a conveyor of whole steps under every solver: people who enter in a step wait in one slot for each step of the delay, thinned by deaths and moved along by aging. It is exact at the step, is not stiff, and takes milliseconds. The largest differences fall by about four times when the step is halved.

  4. Fixed: the first age group was labeled 0-0. An age group one year wide is now labeled by its single year, so the groups are 0, 1-4, 5-14, 15-29, 30+, and compare.py selects by label.

  5. describe() leaves out aging (half fixed). Inflows are now printed as inflows, with their rate. Still open: the aging transfers generated by the age stratification are not listed, though spec 05 says describe() lists them, so nothing in the output says that people age, or at what rates.

  6. Fixed: assigning a distribution to a transition after it was made broke the build. Assigning dur now sets the value of the transition’s parameter. One thing to know: ss.Sim(seir) keeps its own copy of the disease, so assigning seir.onset.dur after the sim is made changes nothing in that sim, silently. It would be clearer if the sim either used the disease it was given or said that it had copied it.

Why this model

odin’s other examples are a logistic curve, the Lorenz attractor, a Lotka–Volterra system, and the SIR models of the “discrete” vignette, one of which is the wave 1 example (stochastic_sir_arrays, whose arrays hold replicates). The suggested discrete_stochastic_seirds.R is a scalar SEIRDS with one multinomial branch, little beyond wave 1. This model is the one that uses what odin says it is for: arrays that are coupled (aging moves people from one element to the next, and everyone is infected by the sum over ages), equations written once for a range of indices with special cases for the first element, and delay(), which odin turns into a call to deSolve::dede.

What it exercises

  • Arrays over age, with equations for [1] and for [2:N_age], and sum() over an array
  • Aging between adjacent age groups at one over each group’s width, with nobody aging out of the last
  • Births into the youngest group at a constant total (b * N), and deaths at the same rate from every state, so the population stays at 10 million
  • An initial population in the stable age distribution, worked out in the model (den)
  • A latent period of exactly 14 days, as a delay: those leaving E are those who entered 14 days ago, less those who died (delay(new_inf[i] * surv, lat_hum))
  • Homogeneous mixing over all ages (I_tot / N), recovery, and waning immunity (R -> S)

Files

The model

For each age group i of 5:

  • Infection: new_inf[i] = beta * S[i] * I_tot / N, with beta = 1 per day, I_tot the sum of I over ages, and N = 1e7 a constant
  • Latency: lag_inf[i] = new_inf[i](t - 14) * exp(-b * 14): those infected 14 days ago who have not died leave E for I
  • Recovery at sigma = 1/30 per day, and waning from R to S at delta = 1/60 per day
  • Births of b * N a day into S[1], and deaths at b = 1/(365 * 50) per day from every state
  • Aging at 1 / (365 * width) per day, for widths of 1, 4, 10, and 15 years, and 0 from the last group
  • Initial conditions: S[i] = den[i] * (N - 1) and I[i] = den[i], where den is the stable age distribution given the aging and death rates; E and R start empty
  • Time: days 0 to 365 (mod$run(0:365)), solved by deSolve::dede

Two details of the delay matter for a translation:

  • Before day 14 the delay looks back to before the start. odin (like deSolve and Berkeley Madonna) then uses the value of the delayed expression at the start, so for the first 14 days lag_inf is the initial rate of infection, about one a day, as if infections had been happening at that rate for the 14 days before day 0. About 14 people reach I in the first two weeks who never passed through E (which runs 14 below what it would otherwise hold, for good). With a reproduction number near 30 this brings the epidemic forward by about 15 days.
  • The latent age but leave by the age they were infected at. E is aged like every other state, but lag_inf[i] takes from E[i] everyone infected in group i 14 days earlier, whether or not they have since aged into group i + 1. In total this makes no difference, but in the youngest group, where 1 in 26 of the latent age out within 14 days, E falls short of what it holds and I gets more than it would.

Results

With one initial case the epidemic grows slowly for two months and then fast: I peaks at 6.19 million on day 125, and on day 365 there are 0.33 million in S, 1.31 million in E, 2.79 million in I, and 5.57 million in R, heading for an endemic equilibrium.

Lines of model code

53: the 50 lines of seir_array_odin.R that are not blank or a comment, and the 3 lines of odin’s tests that compile and run it (odin::odin(), $new(), $run()). No line is over 100 characters.