Starsim v3 to v4 migration guide
Overview
Starsim 4 adds a way to write a disease as declared states and transitions, which can then be run as an ODE, as a stochastic compartmental model, or as agents. For a Starsim 3 model (custom ss.Infection subclasses with their own states and step methods, ss.MFNet and ss.RandomNet, ss.Births and ss.Deaths, custom interventions and connectors, and ss.Calibration), Starsim 4 is an upgrade rather than a switch: most of it still runs, and what is declared can then be described (sim.describe()), run as an ODE or a stochastic compartmental model when nothing in the sim needs agents, swept and compared with ss.parallel() and msim.diff(), and fitted with sim.fit(). This guide shows the same model in both, maps Starsim 3 code to Starsim 4 job by job, lists what changes meaning on the way, and ends with what works and what does not yet. It is written for both humans and LLMs.
Starsim 3 code is of two kinds, and they migrate differently:
- Modules written as code (
ss.Module,ss.Diseaseandss.Infectionsubclasses (ss.Infectionis nowss.AgentDisease, which the migration script renames) withdefine_states(),set_prognoses(), andstep_state(), custom networks,ss.Interventionandss.Analyzersubclasses,ss.Births,ss.Deaths,ss.RandomNet, andss.MFNet, nowss.SexualNet) are still part of Starsim 4. They run unchanged, as agents, and can sit in the same sim as declared diseases.ss.Diseaseis the base of every disease, as in Starsim 3, whether it is declared, written as code, or both in one class;ss.Infectionadds transmission written as code. Per-agent quantities that change continuously (a CD4 count, a level of immunity, time since an event) have no declared form, and keeping them as code is a supported choice. - Starsim 3 spellings that Starsim 4 replaces (the library
ss.SIR,ss.SIS, andss.SEIRas code;ss.routine_vx,ss.campaign_vx,ss.routine_screening,ss.treat_num, and the productsss.Vx,ss.Dx,ss.Tx;ss.seasonality) have been removed, as havess.Calibrationwith its components, andss.Samplesandss.Dataset. Use the table of jobs below to port them.
The migration script makes the mechanical changes for you and lists, by file and line, what needs a person: python -m starsim.converters.migrate mymodel.py --check lists what it would change, and without --check it rewrites the files (or folders) in place. Read its diff afterwards. A removed name raises an AttributeError in Starsim 4; the table of jobs says what replaced it.
Note: If you are an LLM reading this, migrate one module at a time, run the model after each step, compare it with the Starsim 3 results (allowing for the changes in the next section), and confirm with the user before going on.
The same model in both
An SIS on a random network, as Starsim 3 wrote it, with the library disease written as code and a beta for each network:
import starsim as ss # Starsim 3
sis = ss.SIS(beta=dict(randomnet=0.05), dur_inf=ss.years(1), init_prev=0.01)
sim = ss.Sim(diseases=sis, networks=ss.RandomNet(n_contacts=4), n_agents=5000, start=2000, stop=2010, dt=1/12)
sim.run()In Starsim 4, the same disease is declared as states and transitions. Two things change meaning on the way (see What changes the numbers). A bare number in a Starsim 3 beta dict was a probability per step, so here it is written ss.prob(0.05, 'month'); copied as a bare number, it would be a rate per year, 12 times too low at monthly steps (item 3). And a duration in a declared transition is exponential with that mean, unless a distribution is given, while Starsim 3’s ss.SIS drew a lognormal dur_inf; dur=ss.lognormal(mean=ss.years(1)) keeps the old shape (item 5). Here the textbook exponential is used:
import starsim as ss
S, I = ss.states('S I')
sis = ss.Disease(
infection = ss.Transition(S, I, rate=dict(randomnet=ss.prob(0.05, 'month')) * I),
recovery = ss.Transition(I, S, dur=ss.years(1)),
init_prev = 0.01,
)
sim = ss.Sim(sis, networks=ss.RandomNet(n_contacts=4), n_agents=5000, start=2000, stop=2010, dt=ss.months(1), verbose=0)
sim.run()
print('Infections:', round(sim.results.sis.infection.sum()))
sim.plot()This gives 22,860 infections with the default seed, against Starsim 3’s 21,557 (two seeds differ by about 6%). The Starsim 3 model itself, run unchanged as ssl.AgentSIS, is in What works and what does not yet.
The migration script
Most of the changes are mechanical, and a tool makes them: python -m starsim.converters.migrate rewrites Starsim 3 code for Starsim 4, in .py files and in the code cells of notebooks, Quarto, and Markdown files. It renames what was renamed (ss.MFNet is ss.SexualNet, ss.time_prob() is ss.probperyear(), the library’s modules are ssl.X), turns subclasses of Starsim 3’s ss.SIR, ss.SIS, and ss.SEIR written as code into subclasses of ssl.AgentSIR and the like, rewrites a beta pair on a sexual network as a probability per act in each direction, and renames the results that count people in a state. Everything else it lists by file and line, as something for a person to check. Give --check to see what it would change without changing anything. For this Starsim 3 script, mymodel.py:
import starsim as ss
class MySIS(ss.SIS):
def step_state(self):
super().step_state()
sis = MySIS(beta=dict(mf=[0.08, 0.04]), dur_inf=ss.years(1), init_prev=0.01)
net = ss.MFNet()
sim = ss.Sim(diseases=sis, networks=net, n_agents=5000, start=2000, stop=2010, dt=1/12, verbose=-1)
sim.run()
print(sim.results.mysis.n_infected[-1])it says:
python -m starsim.converters.migrate mymodel.py --checkmymodel.py: would change 6 line(s)
0: added import starsim.library as ssl
3: Starsim 3's ss.SIS, written as code, is ssl.AgentSIS (ss.SIS is now declared)
7: a beta pair on ss.MFNet (now ss.SexualNet; [p1 to p2, p2 to p1], with men in p1) is ss.peract(m2f=..., f2m=...), a probability per act in each direction
8: ss.MFNet is ss.SexualNet (between='mf' is its default)
9: verbose=-1 ("brief") is verbose=0
11: the count of a state is the result named after the state (e.g. results.alive, res.infected, results["sir.infected"])
6 change(s) in 1 file(s); 0 item(s) to do by handWithout --check, it makes the changes in place (so run it on code that is under version control, and read the diff afterwards); give it folders to upgrade a whole package. What it can’t do mechanically, such as a beta given as a bare number that was a probability per step, is in the next section, and an AI assistant with Starsim’s migrate-starsim3 skill can work through the list it leaves.
What changes the numbers
Read this first: each of these changes results without an error.
1. Stocks are recorded at the start of the step
In Starsim 4, each time point is a period: a step is taken from every time point, as in Starsim 3, and row ti of a result is the step from time point ti. A flow at ti (e.g. new_infections[ti], now infection[ti]) is the number during that step, as in Starsim 3. A stock (e.g. infected[ti]) is the number at its start, where Starsim 3 recorded it at the end of the step, in update_results(), but labeled it with the step’s start: Starsim 3’s infected[ti] was the number after step ti, which is Starsim 4’s infected[ti + 1]. So infected[0] is now the number at the start of the run, before any step, which Starsim 3 did not record, and the state after the last step is not recorded. The same script, with Starsim 3’s SIR (ssl.AgentSIR in Starsim 4; see 5), gives:
infected[:4] new_infections[:2]
Starsim 3.7.2 53, 99, 186, 356 26, 46
Starsim 4 31, 64, 128, 240 33, 64(ss.Sim(n_agents=3000, diseases=ss.SIR(beta=0.1), networks='random'), with ssl.AgentSIR(beta=0.1) in Starsim 4.) The counts of states are one index later, and the random draws differ after that, so results agree statistically (over several seeds), not number for number. For a model with yearly steps, the number infected in a row is the number on 1 January of that year, rather than on 31 December, while the infections in a row are still those during its year. To compare a stock with Starsim 3’s, compare Starsim 4’s x[1:] with Starsim 3’s x[:-1]. To reproduce published Starsim 3 numbers exactly, use pip install starsim==3.7.2.
With coarse steps, a model can also differ by more than one step. In Starsim 3, a person could pass through several states within one step (e.g. infected, then latent, then active, all in the same step), each taken for a whole step. In Starsim 4, people infected in a step are infected at its end, so they face the exits of the state they enter from the next step on (unless the exit is immediate, with no rate= or dur=), as in the compartmental solvers; someone moved by a rate during a step faces half a step of the rates out of the state they arrive in, so may move again in the same step (with agents by midsteps, the default; with solver=ss.solvers.abm(step='wholestep'), from the next step). So with long steps, people wait longer in their first state after infection. For tbsim with dt=0.1 years, latent TB is 36% higher and active TB 13% higher than in Starsim 3; at dt=0.005 they agree. If results differ at a coarse step, check them at a smaller dt (sim.check_dt() says whether the step is too long for the model’s rates).
So an agent can still move more than once in a step: through two states by rates, or out of a state by a rate and then infected again by the step’s transmission (e.g. asymptomatic, then recovered, then reinfected). Comparing each agent’s state at the end of one step with the next therefore misses some moves (about 1% of infections in tbsim at monthly steps). For flows, read the results of the transitions, which count every move (sim.results.tb.infection, by group with .by()), or, for who was infected, the infection log (disease.infection_log.to_df(ti=self.ti) in an analyzer; see Analyzers), rather than differences of states.
In a module written as code, split update_results() into record_states(), called at the start of each step, before any module steps, for the results that are the state then (a count of a state, a mean level), and record_events(), called at the end of the step, after the deaths, for those about what happened in it (doses, births, new cases). A class that still defines update_results() raises a TypeError when it is created, saying so. A result that needs both (births over the people alive) goes in record_events(), reading the stock from its result (self.sim.results.alive[self.ti]), or in finalize_results():
# Starsim 3
def update_results(self):
super().update_results()
self.results.n_on_art[self.ti] = self.on_art.sum()
self.results.new_doses[self.ti] = self.n_doses
# Starsim 4
def record_states(self): # At the start of the step
super().record_states()
self.results.n_on_art[self.ti] = self.on_art.sum()
def record_events(self): # At the end of the step
super().record_events()
self.results.new_doses[self.ti] = self.n_dosesCode that reads a stock during a step now reads the step’s own row: sim.results.sir.I[sim.ti] is the number infectious at the start of this step, which in Starsim 3 was [sim.ti - 1], the end of the last step. A flow of this step is not known until it ends, so [sim.ti - 1] of a flow is still the last step’s.
2. One time index in every phase of a step
During the step from time point k to k+1, sim.ti (and self.ti) is k in every module and every phase, from record_states() and step() to record_events() and the analyzers, and self.now is time point k, the time the step starts, at which it reads inputs that change in time, as in Starsim 3. The first step is ti == 0, as in Starsim 3, so code that does something if self.ti == 0 in step() still does it in the first step, and self.results.x[self.ti] = v is right wherever it is written. An analyzer runs at the end of each step, as in Starsim 3; the people it sees are as at the next time point.
A declared disease’s ti_<state> is the index of the step in which each person entered the state (a move is taken at the start of its step, and an infection happens by its end); its clocks (dis.P.time_in) are right, and for how long people stayed in a state, use dis.P.stays.
A disease written as code dates its initial infections by the step before the first (self.ti is -1 while ss.AgentDisease.init_post() calls set_prognoses() for them), as a declared disease dates its seeds, where Starsim 3 dated them by the first step. So in a Starsim 3 disease that schedules recovery at self.ti + dur and takes it when ti_recovered <= self.ti, the seeds are infectious for one step less than in Starsim 3, as long as later cases are. Both are about half a step short of dur on average; to give every case dur steps on average, as a declared disease does, schedule from the step after the infection, since an infection happens at the end of its step, and round as the engine does: ss.TimeArr('ti_recovered'), self.ti_recovered[uids] = self.ti + 1 + dur, and self.ti_recovered.due() (see ss.Infection; ssl.AgentSIR and the others do this).
3. A beta dict of bare numbers changes meaning
In Starsim 3, beta=0.05 is a rate per year, but beta={'randomnet': 0.05} is a probability per step, used as it is. In Starsim 4 both are rates per unit time, so a dict copied verbatim is too low by a factor of about 1/dt: 12 times at monthly steps, 52 times at weekly ones (5 infections instead of 22,000 in an SIS on a random network at monthly steps). Give each value its Starsim 3 meaning explicitly:
| Starsim 3 (monthly steps) | Starsim 4 |
|---|---|
beta={'randomnet': 0.05} |
beta={'randomnet': ss.prob(0.05, 'month')} |
beta={'mf': [0.008, 0.004]} |
beta={'sexual': ss.peract(m2f=0.008, f2m=0.004)} |
Use ss.prob(..., 'week') or ss.prob(..., 'day') for weekly or daily steps. Starsim 4 does not warn about a dict of bare numbers; the migration script flags each one, but cannot know your step, so change them by hand. On ss.SexualNet (Starsim 3’s ss.MFNet), which has acts, a beta is a probability per act, the same at any step: Starsim 3’s pair, [men to women, women to men], is ss.peract(m2f=..., f2m=...), and the network is sexual in a beta dict, not mf. Over any network, rel_trans, rel_sus, an edge’s weight and the modifiers now scale the hazard of a contact rather than its chance: Starsim 3’s p·rel_sus is 1 − (1 − p)^rel_sus in Starsim 4. The two are close when p is small; with a large chance per contact, a 50% reduction in rel_sus reduces infections over that contact by less than half.
4. “Routine” delivery means something else
Starsim 3’s ss.routine_vx(prob=0.05, start_year=2005) gives each eligible person a 5% chance of vaccination every year. Starsim 4’s routine delivery, ss.vaccinate(..., coverage=0.05, start='2005'), offers the vaccine once, when a person first becomes eligible. In an SIS model with 5,000 agents from 2000 to 2020, the first reaches about half the population and the second about 5%. The Starsim 3 behavior is delivery in yearly rounds:
ss.vaccinate(efficacy=0.8, when=ss.schedule(every=ss.years(1), start='2005'), coverage=0.05)5. A fixed duration becomes exponential
In a disease you declare with ss.Disease(...) (or a class with states and transitions), a duration given as a number (dur=ss.years(1)) is the mean of an exponential distribution, as in an ODE; describe() marks this as a guess. Starsim 3’s library diseases used a lognormal (ss.lognorm_ex(), now ss.lognormal()). To keep a distribution, give it: I >> R @ ss.lognormal(mean=ss.years(6)). (A distribution that is not exponential runs as Erlang stages in a compartmental run, with a note.)
This includes the library ss.SIR, ss.SIS, and ss.SEIR, which are now declared, with Starsim 3’s keywords but the textbook defaults, of a respiratory infection, in days: ss.SIR() has beta=ss.perday(0.3) and dur_inf=ss.days(10) (R0 = 3 without a network), where Starsim 3 had 0.1 a year and 6 years, and a sim of it runs for 100 days by default, not 50 years. A bare number is in the unit of the model’s typed values, so ss.SIR(beta=0.1, dur_inf=10) is now per day and 10 days; give the units (beta=ss.peryear(0.1), dur_inf=ss.years(10)), which the migration script does for dur_inf and asks for beta. dur_inf (and dur_exp) given as a number or an ss.Duration is the mean of an exponential, where Starsim 3 had a lognormal with a standard deviation of 1 in the same unit (nearly fixed); and nobody dies of the infection unless p_death is given, where Starsim 3 had p_death=0.01. The shape matters: an SIR with R0 = 2.5 and a 7-day infectious period, with beta cut by 40% from day 30 to 90, infects 715,000 of a million with the exponential and 786,000 with Starsim 3’s nearly fixed period. Starsim 3’s models, with their defaults, are ssl.AgentSIR, ssl.AgentSIS, and ssl.AgentSEIR (import starsim.library as ssl): use them to keep Starsim 3’s model, or give ss.SIR(dur_inf=ss.lognormal(mean=ss.years(6), std=ss.years(1)), p_death=0.01).
6. Result names follow the states
A declared disease’s results are named after its states and transitions: I, not n_infected (a module written as code also counts each state under its own name: infected, and sim.results.alive, not n_alive); the flow of a transition is named after it (infection, recovery, with the cumulative total res.infection.cum); and new_infections and cum_infections are infection and infection.cum, for a disease written as code too. Map the names before comparing with Starsim 3 results.
The results of demographics follow the same rule: ss.Births and ss.Deaths have the flows birth and death (were new), with totals birth.cum and death.cum (were cumulative); ss.Pregnancy has preterm, very_preterm, and unborn (were n_preterm, n_very_preterm, and n_unborn), and maternal_death and neonatal_death (were plural).
7. Write dates as strings
Write calendar dates as strings (start='2000', stop='2030-06-30'), and step sizes with a unit (dt=ss.months(1), not dt=1/12). Whether a sim has a calendar is decided once, from its start: a date or a year (a bare number from 1800 on, start=2000) gives it one, and then every bare number given as a time (stop=, when=, the times of ss.schedule and ss.timeseries) is a year; otherwise a bare time is a time in the sim’s unit from 0. A duration given as a time is that long after the start, as in Starsim 3 (when=ss.days(30), start=ss.days(30), sim.run(until=ss.years(5))). stop= is a time, read as the start is, and the length of the run is dur= (ss.Sim(start=2000, dur=10), or dur=ss.years(10), or ss.Sim(stop=2010, dur=10), which starts in 2000); with a calendar, stop=ss.years(10) is an error that suggests dur=ss.years(10).
The keywords of the distributions are their moments, where they have them: ss.normal(mean=, std=) (was loc=, scale=), ss.truncnorm(mean=, std=, low=, high=), ss.expon(mean=) (was scale=), ss.poisson(mean=) (was lam=), ss.weibull(shape=, scale=) (was c=), and ss.gamma(shape=, scale=) (was a=) or ss.gamma(mean=, std=); the order by position is unchanged. ss.datedur(months=1) is ss.months(1). The migration script renames these.
Jobs and spellings
Each row is one job, with the Starsim 3 spelling and the Starsim 4 one. In the Starsim 4 column, sir is a declared disease.
Diseases
| Job | Starsim 3 | Starsim 4 |
|---|---|---|
| Starsim 3’s SIR, with its numbers | ss.SIR(beta=0.1, dur_inf=10) (a lognormal dur_inf, p_death=0.01) |
ssl.AgentSIR(beta=0.1, dur_inf=10), and ssl.AgentSIS, ssl.AgentSEIR (agents only); ss.SIR(...) is now the textbook SIR, with an exponential dur_inf and no deaths |
| An SIR | ss.SIR(beta=0.1, dur_inf=10, init_prev=0.01) |
S, I, R = ss.states(3), then infection = S >> I @ ss.peryear(0.1) * I, recovery = I >> R @ ss.years(10), and sir = ss.Disease(infection, recovery, init_prev=0.01) |
| An SIS | ss.SIS(beta=0.05, dur_inf=10) |
S, I = ss.states(2), then ss.Disease(infection=S >> I @ ss.peryear(0.05) * I, recovery=I >> S @ ss.years(10), init_prev=0.01, name='sis') (no waning immunity: add an R state and R >> S for that) |
| An SEIR | ss.SEIR(dur_exp=ss.days(5), dur_inf=ss.days(10)) |
S >> E @ beta * I, E >> I @ ss.days(5), I >> R @ ss.days(10) (the state in the rate, I, says who infects) |
| A disease given as a dict | diseases=dict(type='sir', beta=0.1) |
the module itself: ss.Sim(sir) |
| Parameters | define_pars(beta=ss.peryear(0.1)), pars=dict(beta=...) |
beta = ss.par(ss.peryear(0.1), bounds=[0.01, 0.5]), used in the rate; bounds make it fittable |
| The sim’s settings | ss.Sim(pars=dict(n_agents=5000, pop_scale=100, networks='random')) |
ss.Sim(n_agents=5000, total_pop=500_000, networks='random'): settings are keywords (or a dict in place of the model, ss.Sim(dict(...))), or in pars= as in Starsim 3, beside the model’s parameters by name, e.g. pars=dict(n_agents=5000, beta=0.1) (a name is a setting only if no parameter of the model has it); sim.pop_scale is total_pop/n_agents |
| Initial infections | init_prev=ss.bernoulli(0.01) |
init_prev=0.01 (always a fraction), or init_n=100 for a number of people seeded the same way (ss.choose_n(100)), or init_n=dict(I=100) by state (init=dict(I=0.01) is a fraction); init_prev= seeds the state a new infection enters (E in an SEIR), and init= and init_n= name the state |
| Where people start | a state whose default is True (ss.BoolState('susceptible', default=True)) |
one named with default_state='N' (in a class, default_state = 'N' in its body), else the state named S or susceptible, else the one state transmission infects people from, else the first state made |
| Population size | n_agents=10_000 (with total_pop= to scale) |
total_pop=1e6 for any solver; n_agents= makes solver='auto' run as agents (and, without total_pop=, is an error in a compartmental run, which has no agents) |
| A new state, by subclassing | an ss.Infection subclass with define_states(), set_prognoses(), and step_state() |
a class with states and transitions in its body (below), or keep it as code |
Moves scheduled with ti_<event> |
set_prognoses() draws dur_inf; step_state() moves people whose ti_recovered <= ti |
I >> R @ dur_inf (a duration after @ is a time in the state) |
| Competing exits | self.transition(uids, to={A: rate_a, B: rate_b}) |
one transition per exit: X >> A @ rate_a, X >> B @ rate_b |
| Death from infection | p_death=0.01, drawn at infection |
two transitions with the same dur and p=: ss.Transition(I, R, dur=dur, p=0.99) and ss.Transition(I, ss.dead, dur=dur, p=0.01) |
| A rate scaled per agent | rel_sus[uids] *= 2 in step() |
a stratification with a value by level (ss.modify(sir.rel_sus, multiply=ss.by(risk=[1, 2]))), or ss.modify(sir.rel_sus, multiply=2, who=...) |
| A lognormal | ss.lognorm_ex(mean=10, std=3); ss.lognorm_im(mean=2, sigma=0.5) (of the logarithm) |
ss.lognormal(mean=10, std=3); ss.lognormal(log_mean=2, log_std=0.5) |
| A beta distribution by its mean | ss.beta_mean(mean=0.3, var=0.01) |
ss.beta_dist(mean=0.3, std=0.1) (a standard deviation, not a variance) |
| A constant | ss.constant(v=5) |
ss.constant(5), or ss.constant(value=5) |
| The type of a duration | isinstance(x, ss.dur), ss.dur(3) |
isinstance(x, ss.Duration); to make one, ss.years(3) |
| Internal classes | ss.Base, ss.Loop, ss.IndexArr |
in their submodules: ss.modules.Base, ss.loop.Loop, ss.arrays.IndexArr |
A disease written as a class, when it needs a method as well as transitions:
class SEIR(ss.Disease):
S, E, I, R = ss.states(4)
beta = ss.par(ss.perday(0.3), bounds=[0.1, 1])
exposure = S >> E @ beta * I
onset = E >> I @ ss.days(5)
recovery = I >> R @ ss.days(10)
seir = SEIR(init_prev=0.01) # Named "seir", after the classMixing and demographics
| Job | Starsim 3 | Starsim 4 |
|---|---|---|
| Homogeneous mixing | a ss.RandomNet() with many contacts |
no network (ss.Sim(sir)); this runs as an ODE |
| A contact matrix by age | a custom network or mixing pools | sir.stratify(age=[0, 15, 65]) and networks=ss.contacts(ages=[0, 15, 65], matrix=matrix) (ages= gives the lower edge of each of the matrix’s age groups), or ss.contacts('Kenya') |
| A network | ss.RandomNet(), ss.MFNet(), networks='random' |
the same, with ss.MFNet() now ss.SexualNet() (between men and women by default; ss.MSMNet() is ss.SexualNet(between='mm')); a network makes the sim run as agents |
| Transmission on a route, and the weight of an edge | ss.RandomNet(beta=0.5), ss.HybridNet(beta=dict(h=3, ...)), ss.MixingPool(beta=..., diseases=...); self.append(p1=, p2=, beta=...), net.edges.beta |
the disease’s beta by route, ss.SIS(beta=dict(randomnet=0.05)) (0 for a route it does not spread on), and for Covasim’s layers ss.SIS(beta=ssl.HybridNet.betas()); weight=, net.edges.weight. A route takes no beta= or diseases= |
| Contacts that last | ss.RandomNet(dur=...) (dur=0, the default, redrew them each step) |
ss.RandomNet(dur=ss.days(5)), with the unit, since a bare number is in the sim’s unit; ss.DegreeNet(n_contacts=..., dur=...) for partnerships formed toward a target number of partners |
| A matrix of mixing pools | ss.MixingPools(src=, dst=, contacts=M), with M[source, target] |
ss.contacts(M.T, ages=[...]), or ss.networks.MixingPools(source=, target=, n_contacts=M.T): row i is the contacts of a person in target group i, as a contact survey has it. A number in the matrix is now the same for everyone in the group (Starsim 3’s was the mean of a Poisson draw for each person): for Starsim 3’s contacts, give n_contacts=ss.poisson(M.T) |
| Births and deaths | ss.Births(birth_rate=20), ss.Deaths(death_rate=15), or ss.Sim(birth_rate=20, death_rate=15) (per 1,000 per year) |
for the population: ss.Sim(demographics=[ss.birth(rate=ss.peryear(0.02)), ss.death(rate=ss.peryear(0.015))]) (under every solver; a rate above 1 a year, bare or typed, is per 1,000, so ss.birth(rate=20) and ss.peryear(20) are the same); or keep the modules (agents only), which read their rates the same way: per 1,000 people if above 1 (20 or ss.peryear(20), as in Starsim 3), else per person (ss.peryear(0.02); a rate of 1 or below was per 1,000 in Starsim 3), with rate_per= to say which. birth_rate= and death_rate= on the sim are removed |
| Pregnancy | ss.Pregnancy(), written as code, with methods to override (make_embryos(), …) and arrays (pregnant, ti_pregnant, ti_delivery) |
ss.Pregnancy() is now declared, with Starsim 3’s arguments: the states NP and P and transitions (preg.conception, preg.delivery, preg.birth, …), under every solver. preg.pregnant is preg.P, ti_pregnant is ti_P, the scheduled delivery is ti_exit_P, and the results pregnancies, births, miscarriages, and stillbirths are its flows conception, birth, miscarriage, and stillbirth. A subclass’s overrides of Starsim 3’s methods are never called (the migration script lists them): use its seams instead: its three rates, fertility_rate=, loss_rate=, and infant_death_rate= (each a rate, or a function of (sim, uids); fertility_rate= also takes a table of rates or a place, which conception_rate() turns into the hazard of conception), ss.modify() of its transitions, since_birth and preg.lactating (breastfeeding) after a birth, and ss.on(preg.birth) or, from code, preg.add_callback(preg.birth, fn) (Starsim 3’s add_conception_callback(fn) and add_delivery_callback(fn) are removed). ss.PregnancyPars is removed too: ss.Pregnancy().pars has the parameters |
| The number of people alive | sim.people.n_alive |
sim.people.n: as agents, the people the agents stand for, so the number of agents only with no population scale (sim.people.alive.count() counts agents); the same under every solver. len(sim.people) is the number of rows: the agents, or in a compartmental run, the groups |
| The children of some parents | sim.people.find_children(uids) |
ss.uids(ppl.auids[np.isin(ppl.parent, uids)]) (People.find_children() is removed) |
| An age range | ss.parse_age_range('[15, 25)'), ss.apply_age_range() |
a condition, (ppl.age >= 15) & (ppl.age < 25), or a stratification, sim.stratify(age=[15, 25, 50]) (both functions moved to STIsim) |
Interventions and connectors
| Job | Starsim 3 | Starsim 4 |
|---|---|---|
| A vaccine | ss.simple_vx(efficacy=0.8), ss.Vx(...) |
the product’s keywords on its delivery, ss.vaccinate(efficacy=0.8, ...), or the product itself, ss.Vaccine(efficacy=0.8), when one product is given by several deliveries. It is leaky, as ss.simple_vx was by default (everyone vaccinated has their chance of infection cut by 80%); for all-or-nothing (Starsim 3’s leaky=False: 80% of those vaccinated are fully protected, the rest not at all), leaky=False |
Routine vaccination, p a year |
ss.routine_vx(product=vx, prob=0.05, start_year=2005) |
ss.vaccinate(efficacy=0.8, coverage=0.05, when=ss.schedule(every=ss.years(1), start='2005')) (see item 4 above) |
| A campaign | ss.campaign_vx(product=vx, years=[2020], prob=0.6) |
ss.vaccinate(efficacy=0.8, when=['2020'], coverage=0.6) |
| Who is eligible | eligibility=lambda sim: sim.people.age < 15 |
who=lambda sim: sim.people.age < 15, or a condition: who=sir.age < 15 when sir is stratified by age |
| Screening and treatment | ss.routine_screening(product=ss.Dx(...), prob=0.1) with a treatment intervention |
screening = ss.screen(sensitivity=0.9, specificity=0.98, coverage=0.1, when=ss.schedule(every=ss.years(1))) and ss.treat(cure=0.95, when=screening.positive, coverage=0.8) |
| Treatment with limited capacity | ss.treat_num(product=ss.Tx(...), max_capacity=10) |
ss.treat(cure=0.95, when=sir.infection, capacity=10) |
| A parameter change at a time | an ss.Intervention that edits sim.diseases.sir.pars.beta in step() |
ss.modify(sir.pars.beta, multiply=0.5, start='2010') |
| A change while something holds | the same, with a test in step() |
ss.modify(sir.pars.beta, multiply=0.6, start='2020-03-15', stop='2020-06-01') |
| Seasonality | ss.seasonality(...) |
beta=ss.seasonal(mean=ss.peryear(0.5), amplitude=0.4, peak='January') |
| One disease changing another | an ss.Connector that writes sis.rel_sus from hiv’s states |
connectors=ss.modify(sis.rel_sus, multiply=3, who=hiv.I) |
| A custom intervention | an ss.Intervention subclass with step(), with eligibility= and check_eligibility() |
keep it (agents), with who= for eligibility= (check_eligibility() is unchanged), or ss.vaccinate(..., who=fn) with a function returning who receives it |
A connector that sets an array (sis.rel_sus[:] = 1.0, then raises it for some people) also erases what other modules set, such as a disease’s own immunity. ss.modify() multiplies a base value that it never overwrites, so modifiers from several modules combine.
Runs, results, and calibration
| Job | Starsim 3 | Starsim 4 |
|---|---|---|
| Many runs of one sim | ss.MultiSim(sim, n_runs=10).run(), then msim.mean() |
msim = ss.parallel(sim, n_runs=10), then msim.summarize() (a table of the median, a 95% interval, the mean, and the standard deviation of each result over the runs), msim.plot() (the median and 95% interval), or msim.results.sir.I.median() and .quantiles() at each time, for any sim; msim.reduce(), mean(), and median() are gone. ss.MultiSim(sim) alone is now one run (Starsim 3: four), and msim.show(output=True) is text = msim.show() |
| Scenarios | ss.parallel(sim1, sim2) |
ss.parallel({'Baseline': sim1, 'Masks': sim2}), labeled as you like, and msim.diff('Baseline') |
| A parameter sweep | a loop over ss.Sim(pars=...) |
ss.parallel(sim, pars={'sir.beta': [0.2, 0.3, 0.4]}) |
| A result | sim.results.sir.n_infected |
sim.results.sir.I |
| A result by group | an analyzer | sim.results.sir.I.by('age'), .select(age='0-14') |
| The modules given to a sim | ss.Sim(..., copy_inputs=False), so that sir.results is the sim’s |
a sim always copies the modules it is given (copy_inputs= is removed; the migration script removes it and lists the line): read them, and their results, from the sim (sim.diseases.sir, sim.results.sir) |
| A copy with changes | sc.dcp(sim), then edit |
sim.new(n_agents=5000), a new sim from the same inputs with what is given changed; to run the same sim again with another solver, sim.new(solver='abm').run() (a sim that has run is run again from the start) |
| Calibration | ss.Calibration(sim, calib_pars=dict(beta=dict(low=0.01, high=0.3, guess=0.1)), build_fn=..., components=[ss.Normal(...)], total_trials=100), then calib.fit() |
beta=ss.par(0.1, bounds=[0.01, 0.3]) in the model, ss.Sim(sir, data=df) with a column named after a result (prevalence), then fit = sim.fit() (method=ss.fitting.method('optuna', n_trials=100) for Optuna, used for agents) |
| A likelihood | ss.Normal(...), ss.GammaPoisson(...), ss.BetaBinomial(...) components |
sim.fit(likelihood=dict(prevalence=ss.normal(std=0.02))), or likelihood= on an ss.outcome(); ss.GammaPoisson is ss.nbinom(n=...), and a beta-binomial is a function of your own (Fitting) |
| Data shared by two diseases | an extract_fn per component |
name the column after the disease too: sis.prevalence |
| Data on a flow by period (weekly cases) | lined up by a component’s extract_fn |
each row is labeled by the start of its period, and summed over the steps to the next row, so no extract_fn is needed; move data labeled by the end of each period (“week ending”) back one period first, e.g. df['date'] -= pd.Timedelta(days=7) |
Migrating one disease at a time
A published model does not have to move all at once. Keep the hard disease as code and declare the easy one: the mixed sim runs as agents, and sim.new(), ss.parallel() scenarios and sweeps, msim.diff(), and sim.fit() work on the whole sim.
import starsim as ss
import starsim.library as ssl
hiv = ssl.HIV(beta={'sexual': ss.peract(m2f=0.008, f2m=0.004), 'random': 0}, init_prev=0.05) # Starsim 3's module, written as code (see 3 for its beta)
S, I = ss.states('S I')
sis = ss.Disease( # Newly declared
infection = ss.Transition(S, I, rate=dict(sexual=0, random=ss.par(ss.prob(0.005, 'month'), bounds=[0.001, 0.05])) * I),
recovery = ss.Transition(I, S, dur=ss.years(2)),
init_prev = 0.05,
name = 'sis',
)
susceptible = ss.modify(sis.rel_sus, multiply=3, who=hiv.infected)
sim = ss.Sim([hiv, sis], networks=[ss.SexualNet(), ss.RandomNet()], connectors=susceptible,
n_agents=2000, start='2000', stop='2020', dt=ss.months(1))
vaccination = ss.vaccinate(efficacy=0.8, disease=sis, when=ss.schedule(every=ss.years(1), start='2010'), coverage=0.05)
msim = ss.parallel({'Baseline': sim, 'Vaccine': sim.new(interventions=vaccination)}, n_runs=3)
msim.diff('Baseline')Two things to know in a mixed sim:
- A condition can name the arrays of a module written as code (
who=hiv.infected,who=hiv.cd4 < 200), or be a function of the sim (who=lambda sim: sim.diseases.hiv.infected). ss.modify()andss.vaccinate()act on a module written as code too, as agents, if it has what they change (such asrel_sus); givedisease=when there are several diseases. Its Starsim 3 connectors and interventions also still run.
To validate each step, run the Starsim 3 version and the ported one over several seeds (ss.parallel(sim, n_runs=10)) and compare means and spreads, after mapping the result names and allowing for the one-step shift.
What needs a declared disease
| Works on any module | Needs a declared disease |
|---|---|
Running as agents, networks, ss.Births/ss.Deaths, analyzers, ss.parallel(sim, n_runs=...) with msim.results across the runs, sim.new(), scenarios and sweeps with ss.parallel(), sim.fit() (with bounds on a parameter), ss.modify() and ss.vaccinate() targets (as agents), sim.describe() (in brief) |
The ODE and stochastic compartmental solvers, describe() in full, equations() |
The compartmental solvers need every module of the sim to have a compartmental form, so a sim with a disease written as code runs as agents.
What works and what does not yet
These are from Starsim’s switching demo for Starsim 3 users (tests/switching/test_starsim3.py in the lingua_franca repo), and from the review it was built from. All four of the demo’s checks pass today.
| Task | Today |
|---|---|
| Run a Starsim 3 model with the disease as code | Works: the SIS above as ssl.AgentSIS (one changed name) gives the Starsim 3 answer within seed noise (22,484 infections against 21,557, where two seeds differ by about 6%; a little higher, since each case is infectious for half a step longer than in Starsim 3, as a declared disease’s is) |
Use Starsim 4 verbs on a disease written as code: ss.modify(hiv.rel_trans, multiply=0.5, start=2005), conditions on its arrays (who=hiv.cd4 < 200), scenarios with ss.parallel() and msim.diff('baseline') |
Works, on the library HIV written as code alongside a declared SIS |
Fit with sim.fit(), with bounds given on a parameter (ss.SIS(beta=ss.par(ss.peryear(0.1), bounds=[0.05, 1]))) |
Works, for the declared ss.SIS and for ssl.AgentSIS |
Port a disease with its beta dict copied verbatim |
Silently wrong: the result is about 12 times too low at monthly steps, with no warning from Starsim, though the migration script flags it (see item 3) |
| Port one disease at a time | Works: a sim can mix diseases written as code and declared ones, and ss.parallel() sweeps and diff() work on the whole sim (see Migrating one disease at a time) |
| Speed of a model written as code | Faster: an unchanged Starsim 3 model with HIV, an SIS, an SIRS, two networks, demographics, and interventions ran about twice as fast under Starsim 4 (2.2 s against 4.9 s at 3,000 agents) |
| Speed of a declared disease run as agents | Slower: 1.4 to 1.6 times slower than the same model written as code, at 3,000 to 30,000 agents; making them equal is planned |
| Per-agent quantities that change continuously (a CD4 count, a level of immunity, time since an event) | No declared form: keep them as code, which is a supported choice |
Further reading
- A tour of Starsim: a declared model built step by step, with age groups, interventions, calibration, and uncertainty
- Getting started: the tutorials, which build up one model a step at a time
- Library: diseases: the Starsim 3 models written as code,
ssl.AgentSIR,ssl.AgentSIS, andssl.AgentSEIR, and the library’s other diseases - Starsim: releases: the changes in each version, including the last of Starsim 3 (3.7.2)
- Starsim 3 on PyPI: the version to install to reproduce Starsim 3 results exactly