Skip to content

seeds: derive each row's seed from the base, not base + i (#1365) - #1382

Open
Montana wants to merge 1 commit into
lanl:mainfrom
Montana:fix/1365-independent-replicate-seeds
Open

Montana wants to merge 1 commit into
lanl:mainfrom
Montana:fix/1365-independent-replicate-seeds

Conversation

@Montana

@Montana Montana commented Oct 29, 2026

Copy link
Copy Markdown
Contributor

Closes #1365.

What was wrong

Row i of run_batch, run_replicates and parameter_scan was seeded base_seed + i. So run_replicates(50, seed=1) was run_replicates(50, seed=0) shifted by one replicate. Each ensemble on its own was fine, but ensembles from nearby bases were nearly the same samples, with no warning. On main c0877fa, with the issue's birth-death model:

main
seeds of seed=0 / seed=1 [0, 1, 2, 3, …] / [1, 2, 3, 4, …]
trajectories shared by the two 50-replicate ensembles 49
shared by two 20-row run_batch calls, seeds 0 and 1 19

The fix

python/bngsim/_seed.py: derived_seeds(base_seed, n). It returns the seeds of rows 0..n−1, drawn from numpy's SeedSequence(base_seed). That stream is a hash of the base and the position, so nearby bases give unrelated seeds. Each seed is masked to 31 bits, the width _resolve_seed already draws for the NFsim, RuleMonkey and SSA backends.

The seeds of one call are distinct. The first n distinct values of the stream are taken, in order. A plain 31-bit draw repeats inside an ensemble of 10,000 about once in 40, and a repeat is a duplicated trajectory, which base + i never produced. The stream is prefix-stable, so row i's seed does not depend on n.

Where it is used. Each call site computes its row seeds once, before any row runs, so the threaded paths draw what the sequential paths do:

  • run_batch: row_seeds = derived_seeds(base_seed, n_sims). _run_single_batch now takes the row's seed (row_seed), not the base plus its index.
  • run_replicates: rep_seeds, on the sequential path and per thread.
  • parameter_scan: scan_seeds.

Result.seed is the derived seed, so run(seed=r.seed) still repeats that row alone. The docstrings of all three methods and docs/quickstart.md now describe this schedule.

Behaviour change

A seeded batch, ensemble or scan gives different trajectories than before. That is the point of the change, but numbers recorded from a seeded call under base + i will not be reproduced. seed=None (a fresh base from entropy), single run(seed=...) calls, and ODE (whose batch seed only breaks event ties) are unchanged. The changelog fragment is changed, not fixed, for this reason.

Tests

test_derived_seeds.py (8; on main the file fails at import, since derived_seeds does not exist):

  • ensembles, batches and scans from bases 0 and 1 share no trajectory (on main, ensembles and batches share 49 and 19; scans were not measured on main, but parameter_scan used the same base + i);
  • a replicate reruns alone from its Result.seed;
  • threads draw the seeds, and the trajectories, of the sequential path;
  • 200,000 seeds are distinct (about nine repeats would be expected without the dedupe) and fit in 31 bits;
  • row i's seed does not depend on n;
  • nearby bases share no seed.

Existing tests that encoded base + i, updated to derived_seeds while keeping what each checks:

  • test_seed_semantics.py: test_per_sim_seeds_are_consecutive becomes …_are_derived_from_the_base; and "seeds within a batch are consecutive" becomes "are distinct".
  • test_ssa_run_replicates.py: the clone-per-replicate reference, against which the sequential and threaded paths are compared bit for bit, seeds each replicate as run_replicates now does; and the exposed seeds are derived_seeds(100, N).
  • test_ssa_reaction_stats.py::test_batch_rows_carry_the_blocks: compares row 0 with a run at its seed (rows[0].seed), not at the base.

Three PSA mean tests had a standard error 4.7× too small. These are test_sbml_psa.py::test_psa_sbml_matches_analytical, test_ssa_psa_volume.py::test_psa_v_neq_1_mean_matches_analytical and ::test_psa_v_eq_1_unchanged. They bounded the PSA ensemble mean by 5 SEs of the exact SSA's binomial spread, √(N·p·q) = 46.5. PSA's leaps spread the population far wider: I measured a standard deviation of 218 over 4 bases × 500 replicates. So "5·SE" was about one of PSA's SEs, and each test passed or failed with its seeds. Reseeded, all three missed by 17.2 against a bound of 16.4. PSA is not biased: over those 2,000 replicates its mean is 6843.3 against 6839.4, z = 0.8. The tests now take the SE from the sample (std(ddof=1)/√n), with the measurement in a comment.

Mutations

mutant fails
D1: derived_seeds returns base + i (main's schedule) 4: the ensemble, batch and scan tests, and nearby bases
D2: no dedupe the 200,000-distinct test
D3: parameter_scan left on base + i the scan test
D4: run_batch left on base + i the batch test

Measured

Full suite (macOS arm64, main c0877fa + this branch): 16,420 passed and 4 failed. The 4 were the reaction-stats row test and the three PSA tests above, which now pass (those three files: 26 passed). pre-commit run --all-files (new files staged) and sphinx -W with docs/requirements-ci.txt are clean.

Row i of run_batch, run_replicates and parameter_scan was seeded base + i,
so an ensemble from seed 1 was the one from seed 0 shifted by a row: 49
of 50 trajectories bit-identical. Each row's seed is now drawn from
numpy's SeedSequence(base) (derived_seeds), distinct within a call and
31 bits wide, and stamped on Result.seed as before. Seeded ensembles
change; seed=None, single runs and ODE do not.

Three PSA mean tests bounded the mean by 5 binomial SEs, a fifth of
PSA's real spread (sd 218 for 46.5); they now use the sample's SE.

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

run_replicates / run_batch: ensembles from consecutive base seeds share almost every trajectory

1 participant