🧠0% 🤖100%

Implementation plan

Status of the Starsim v4 implementation of the design spec. The code is on the v4 branch of the Starsim repository; section numbers refer to the files in design/. Checked items run and have a test in tests/test_transitions.py, tests/test_deliveries.py, or tests/test_dimensions.py.

Build set

The nine build steps from models/README.md.

1. SIR, in all spellings
2. SEIR with an Erlang dwell time
3. A branching model with a competing death hazard
4. Cholera, with a water reservoir
5. Age-stratified SIR with a contact matrix
6. Vaccination as a stratification plus a transfer
7. A model with outcomes, data, and bounds
8. A model with custom code
9. A sexual-network model

Departures from the spec

Sims with no declared states keep Starsim 3’s time defaults.
The unit of a model’s bare numbers comes from the model (unit=, else its typed values, else dates, else a typed dt, stop, or dur), as the review’s Decision lang 1 recommends, with two departures from that recommendation: a year-like number for start or stop is still a calendar year (stop=2020), since Starsim 3 sims are written that way; and the typed defaults of modules written as code do not set the unit, since they convert themselves and would otherwise disagree with most declared models.
Rates are written as strings in the model spec, and the spec is sim.to_spec(); sim.to_json() is unchanged from Starsim 3. The spec is for reading by people and other tools; it is not loaded back (ss.load() of a spec and ss.from_spec() were removed).
The default ODE integrator is fixed-step RK4 with automatic substeps, in the Numba kernels, with or without gradients being asked for; JAX is used only for integrator='jax', sim.gradient(), and gradient-based fitting. integrator='jax' takes the same RK4 steps (not an adaptive JAX integrator), so its results are the kernels’.
What JAX differentiates is the run, from the numbers in the program and the compartments at the start to the compartments and flows at every step. The steps either side are Python (the parameters to those numbers, and the run to the results, outcomes, and goodness of fit), and are differentiated by central differences along the direction of each parameter. So a gradient is exact to about 1e-10, not to the last bit, and costs two restarts and two runs’ worth of results for each parameter.
Superseded by v4: an age stratification with no split= takes the named world population’s ages (or those of the people given), not ss.People’s defaults
Every agent’s dwell time is sampled, including exponential ones.
A clock of a set of states ((A | B).time) has to be made before the sim is built (the time each person entered the set is kept only for the sets known then); one made later raises an error saying so. X.enter.uids in custom code is read from the states declared on the module in the sim (sim.diseases.d.graph.states.X), since sim.diseases.d.X is the state’s array as agents.
The proposal’s initial placement firing events (seeds counted as having gone through the transitions to their state) is not done: people placed at the start have entered their state, and their clocks start then. Moves at the start from there (the cascade) are events and are counted at the first time point (agents).
A clock in a compartmental run is a conveyor of one-step levels by default (exact on the step grid), not 10 levels as products have (proposal 01 §2, Q3), and not Erlang aging (Q1); it does not reuse a dwell time’s Erlang stages (Q2), except that a dwell time kept as a conveyor is its own clock. Only ss.constant dwell times are conveyors by default, not every dwell time with a CV below 0.1 (03 Q4).
The property of people that says they hold a product is named has_<product> (e.g. has_bcg), and the record that a routine delivery has been offered is offered_<delivery>. The product’s own name is taken by the product module.
A product whose waning= is a duration is lost at an exponentially distributed time, under every solver, after which the person can receive it again. Only a function of time since delivery fades the efficacy while the person still holds the product.
A delivery never gives a product to someone who already holds it.
A bare number for a time (ss.change(at=...), a delivery’s start=) is a time in the sim’s unit, counted from the start, even in a sim with dates; it is a calendar year only in a sim whose own start/stop are year-like numbers. A time outside the sim raises E1101.
A modifier acts where the engine uses the quantity; reading sir.rel_sus from custom code gives the base value, and sir.rel_sus.effective the value with the modifiers applied (agents only). So .base is not needed.
In a compartmental run, a product whose protection fades by a function of time is held at levels of time since delivery (stage 1, …), so has_<product> has more than two levels.
A network’s parameter is the target of a modifier as a path ('school.pars.n_contacts'), since school.pars.n_contacts is a plain value.
W1203 and W1204 are not used: tau leaps adapt to the rates, and W102 covers large steps.
ss.modify of an undriven transition multiplies the efficacy of the treatments that move people along it, for the people the modifier applies to (proposal 02 §2.1), rather than raising E504 (proposal 03 §2.4): 02 owns treatment. A modifier of a transition with a rate does not change a treatment’s efficacy.
A transition with no rate, duration, or p= out of a state that has p= branches is undriven, not their remainder (03 §2.1). An immediate remainder is written as a split ('E -> P | A', p=[0.3, None]) or with p=None and a driver.
To the engines, an undriven transition is a rate of 0, as @ 0 was, so as agents a state whose only exits are undriven still draws for them each step.
As agents, each disease of a treatment draws its own number per person, so cures of two diseases by one product are independent.
Under gillespie (and tau), a value that varies smoothly in time is held for a quarter of a step at a time, at its middle, with a note, rather than by thinning, which needs a bound on every rate over the step; a rate can be any expression of the parameters.
A module written as code that works on groups as well as on agents says so with the class attribute compartmental = True; otherwise it makes the model agent-based.
The peak of ss.seasonal given as a month is the middle of that month. In a sim with no dates, time starts at the beginning of a year.
ss.gravity gives k × pop[i]^(a−1) × pop[j]^b / distance^c: for commuting, the share of time a person from i spends in j (capped so that it sums to at most 1, or rescaled to away=), and for migration the rate of moving.
An inflow whose rate uses a count (xi * I.n, or ss.birth) arrives in each place according to the count there; an inflow at a fixed rate is shared between places by their split.
A new infection is recorded under the strain acquired, and a state that a dimension does not apply to is counted at that dimension’s first level.
As agents, the number who start in each state is exact (rounded to whole agents), not a chance for each agent.
A level of a stratification named like an attribute of a result (low, high) is a column called <dimension>_<level> in result.by().
A test with several results draws once per person from the distribution of the worst result (or of the highest, combine='max'), not once per disease; the same distribution, and a binary test keeps the draws it had.
outcomes= cannot be given with sensitivity= or specificity=; false positives are its None key.
An imported rate with another state and no division by the total (beta * S * I, rule 3 of 16) is kept as counts (beta * I.n), which is exact, rather than given dependence='density' and an n_ref: there is no density-dependent mixing.
Results follow the convention of other frameworks under every solver: n_<state>[ti] is the number at time point ti (after ti steps, so n_<state>[0] is the number at the start), and a flow at ti (recovery[ti]) is the number during the step that ends there (0 at the start), so n_R[ti] - n_R[ti-1] == recovery[ti]; a sim with stop=T ends at T.
E215 is already PeopleAndCounts, so the pregnancy-structure code of proposal 04 is E216 (reserved, not raised yet).
Maternal death at delivery (p_maternal_death) is code in ss.Pregnancy (agents), not a declared branch: a split to ss.null is not a destination.

Known issues

The kernels are compiled the first time they are used after Starsim is installed or changed, which takes one or two minutes (a line says so); after that they are loaded from a cache.
With sim.run(n_runs=) in one call to the kernel, only the states and transitions are kept by group (not outcomes or results given as expressions), and reported outcomes are drawn for all the runs at once, so differ (in the draw, not the distribution) from those of a single run with the same seed.
A typed rate cannot be used as a plain number in an expression (1 - ss.peryear(0.008)); the error says to convert it first (rate.to_prob(dt)).
A bare number given for a parameter of a module written as code (ss.MFNet(duration=5), ss.SIS(dur_inf=10)) keeps the unit of its default, as in Starsim 3, while one in a declared model is in the sim’s unit. It warns when the two units differ (duration read in years in a sim in days).
As agents, a dwell time is sampled when someone enters a state, so an adjust on a duration does not follow a person who later moves to another level (through an ss.transfer, or by aging). A rate does.
A strain’s initial share is the split of the dimension (equal by default), which applies to everyone who starts in an infected state.
As agents, someone tested on the step they were infected may already have recovered (about dt/2 of a dwell time), so event-triggered tests find slightly fewer cases than in a compartmental run.
An ODE and an agent-based run of the same model differ by an amount of order dt (about 6% in the prevalence of an SIS model with a daily step).
As agents, people added by a movement register in a step are made at its end, so the other events of the same step cannot choose them; and people who moved into a place in a step are chosen after those who were already there.
A movement register cannot add people to a state that splits on entry into branches.
A modifier’s cached value is worked out once per row; people who enter a state later in the row (by the cascade) get the value from before they entered for the rest of it.
[B] In a compartmental run, a trigger on states and strata is exact for arrivals by transitions, but people who come to meet it by aging, a transfer, or another disease’s move wait until the next loop time (the proposal routes those flows too); a trigger on an event is agents only.
An ss.by of distributions of several families, or over several dimensions, is not supported; nor is arithmetic between ss.by of different dimensions, or a branch modifier through tr.to(B).
Arithmetic on stratified results keeps the values by group only for sums, differences, and scaling by a number; for anything else (a ratio), select() says the result is not stratified. (A ratio declared on the disease, prev=I/N, is stratified: it keeps its numerator and denominator.)
With immunity= between strains, a compartmental run keeps only each person’s last strain (noted), so with partial protection people can go from A to B and back to A; as agents, every strain is kept, but not yet over networks.
Code that changes a property of people that a dimension is made from (people.patch) has to call sim.groups_changed(uids), or the agents’ groups are not worked out again (dimensions binned from a property that changes on its own, such as age, are checked every step).
A state that is not people is deterministic under the stochastic solvers, so a vector population whose own randomness matters cannot be drawn.
Under the exact solver, places are run one at a time between steps only when people mix within their place and nothing moves them between groups within a step; otherwise every rate is worked out again after each event.
A plain duration= on an outcome is an exponential stay, while delay= is an exact shift.
As agents, TB (LSHTM) comes out 5–10% higher than the ODE with weekly steps, and 20–30% higher with monthly ones; this is the difference of order dt above, made larger by an epidemic that is near R = 1.
rel_sus and rel_trans on a state are numbers, not ss.pars, so they cannot be fitted; the TB models write them into the rates instead.
Once JAX has run in a process (sim.gradient(), a gradient fit, or integrator='jax'), running sims in parallel by forking processes (ss.parallel() of agent-based sims, Optuna with several workers) can deadlock, as JAX warns. Spawning instead would be safe but takes about 1.5 s longer to start.
Gradients need a sim that the kernels run from start to finish (no interventions other than deliveries the kernels take, no modifiers with when=, transfers whose rate varies, or movement register), and a rate that is a Python function of the sim has none (one written with NumPy functions of counts does).
As agents, a stock target is met to within an agent or two: a dwell in the end state drawn shorter than a step ends at once (the same-step cascade).
A treatment’s transitions given by name across diseases (moves='ct.cure') do not make the diseases tracked jointly in a compartmental run; given as objects (ct.cure), they do.
In a compartmental run, on_fail= given as a delivery receives the failures as a share of each state and group of the delivery’s disease, so it can act only on what that disease’s engine can select.
Under solver='binomial', people who arrive in a state during a step face none of its hazards until the next step (the chain-binomial’s whole steps); with a hazard that reads the clock it is the largest they miss, so BV recurrence (decaying over 3 months) is 6% below the ODE at weekly steps and 24% at monthly ones. Agents had the same miss for people entering in the state-change row, now fixed (see 03); substeps (substeps=) shrink it
As agents, back-dated seeds who die before the start are removed, so the population is that much smaller; with no over=, how far back they go is ten times the mean time to leave the state (from its dwell times and constant rates), or 20 years
A product keyed by diseases does not make them tracked jointly as compartments, so a delivery of it whose selection depends on one disease’s state is E804 for the others (as for any delivery of a multi-disease sim).
ss.probperyear(ss.timeseries(...)) (proposal 02 §4.4’s scale-up of a routine rate) is not supported; a time series of bare coverage is.
take= (an all-or-nothing vaccine) is not built.
[B] A recorded attribute’s levels are read at the middle of each bin (the open last bin at the unweighted mean of what can be recorded into it), so coarse bins bias what is read from them: the HIV test gives cd4_stop finer bins at the top than the CD4 stratification. Levels weighted by what is actually recorded into them would need the flows of each step
[B] With ss.spread in a compartmental run, everyone in a clock or conveyor moves up a level at the start of the first step (the back-dated people have been there a step), including anyone placed there at the start by init=
[B] A compartmental model with attributes works out each compartment’s rates in Python: the HIV CD4 model (6,720 compartments) takes about 7 s for 10 years of monthly steps
The JAX run is compiled the first time it is used for a model (about 1 s), and again if a step needs more substeps than it has allowed for.

Implementation plan

Build set

Held-out set

By section of the spec

01 Objects

02 Types

03 States and transitions

04 Parameters

05 Dimensions

06 Routes

07 Effects

08 Interventions

Pregnancy and births (proposal 04)

09 Results

10 to 12 Execution, random numbers, and solvers

13 Inference

14 to 16 Model spec, errors, and interop

Departures from the spec

Known issues

Review of 2026-10-17

Phase 1: unblock

Phase 2: contract tests

Phase 3: the silent wrong answers

Phase 4: the guessing rules

Phase 7, started early: inference

Phase 5: one system (first round)

Phase 6: docs

Phase 7: features for the seats

Phase 8: performance

Phase 9: surface and structure

Extending a model with code (2026-10-25)