feat(patch): meta-population model over the three affected provinces - #412
Merged
Merged
Conversation
seabbs-bot
commented
Jul 11, 2026
seabbs-bot
commented
Jul 11, 2026
seabbs-bot
commented
Jul 11, 2026
seabbs-bot
commented
Jul 11, 2026
Collaborator
Author
|
Review pass done (4 inline comments above). CI status: One further observation from watching the fit matrix fan out: So the patch model has unit tests ( This was opened by a bot. Please ping @seabbs for any questions. |
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## main #412 +/- ##
==========================================
+ Coverage 95.12% 96.73% +1.61%
==========================================
Files 15 16 +1
Lines 3731 4536 +805
==========================================
+ Hits 3549 4388 +839
+ Misses 182 148 -34 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
Contributor
|
📖 Documentation preview is ready! View the docs for this PR at: http://epiforecasts.io/BVDOutbreakSize/previews/PR412/ This preview will be updated automatically when you push new commits. |
seabbs
marked this pull request as ready for review
July 12, 2026 08:10
Collaborator
|
This needs to be the new headline model ie joint and be wired into the analysis and sensitivity ie for one week ahead forecasts. I don't see evidence it works here |
This comment was marked as outdated.
This comment was marked as outdated.
This comment was marked as outdated.
This comment was marked as outdated.
Collaborator
|
There are some conflicts. The much larger outbreak size is worrying and I am not sure plausible. You don't show estimates or forecats for any of the spatial patches or Rt by patch |
A meta-population / patch model formulation for the 2026 DRC BVD outbreak that splits the single-patch national model into three spatial patches (Ituri, Nord-Kivu, Sud-Kivu). Each patch has a province-specific reproduction number (hierarchical modifier on the national Rt), within-patch renewal transmission, and a between-patch importation kernel. National data streams are fitted to the sum of patch trajectories; per-province spatial-table data (confirmed cases, lab analysed, isolation) fit individual patches. See notes/patch-model.md for the full formulation. Relates-to: patch-model exploration
- patch_rt_model: hierarchical Rt (national walk + per-patch modifiers) - patch_infection_model: multi-patch renewal with importation - patch_infections() and importation_from_kernel() in renewal.jl - Verified: module compiles, prior predictive sampling runs Relates-to: patch-model exploration
- _patch_latent: wraps patch_infection_model and sums patch trajectories - bvd_patch_joint: full joint model using patch latent process with national-level observation submodels - All models compile; patch_infection_model prior predictive verified Relates-to: patch-model exploration #5
Use the national-level cases_state.bvd_reports_daily and cases_state.bg_daily for the treatment flow model, matching the signature of treatment_flow_model. Fix field name from deaths_state.onsets_to_deaths (does not exist) to the correct cases_state.bvd_reports_daily + cases_state.bg_daily pair. Remove duplicate := deterministics for parameters already traced inside the patch submodel (importation_epsilon, sigma_region, C_T_total).
- Expose C_T_patch_1/2/3, R_T_patch_1/2/3, infections_T_patch_1/2/3 as := deterministics in bvd_patch_joint - Add patch_summary_table() to summaries.jl — per-patch posterior summary of C_T, R_T, daily infections, and δ modifier - Export patch_summary_table from BVDOutbreakSize.jl
…onal Rt Replace the constant-modifier patch_rt_model with patch_rt_mvwalk_model: - Each patch gets its own log-Rt trajectory via an MVN random walk at weekly knots - Innovations are correlated across patches via LKJ(2) sampled correlation matrix - Each patch has its own step SD sigma_rw_patch and initial R0 variation - Allows Ituri (high Rt, sustained transmission) and Sud-Kivu (near-zero local Rt) to diverge arbitrarily Add implied_national_Rt() helper to renewal.jl: - Derives the national Rt from summed patch infections via the renewal equation - No separate national Rt parameter — it emerges from the aggregated dynamics Update patch_infection_model to use the new MV walk by default: - Removes rt_state.delta_patch and rt_state.sigma_region dependencies - Adds implied_Rt_national to the return tuple - Exposes sigma_rw_patch and Omega from the MV walk - Fix docstring to document new return values Add LinearAlgebra as a project dependency for cholesky factorization.
Replace delta_patch references with sigma_rw patch step SDs from the multivariate normal random walk.
- Fix intervention ramp in patch_rt_mvwalk_model: add intervention_effect parameter and apply ramp properly (was multiplied by 0.0) - Add @inbounds to MV walk core loops (Sigma_half, knot innovations, daily interpolation) - Add effect_prior keyword arg to patch_rt_mvwalk_model - Fix background_re code path: guard against empty histories with bg_lead pattern matching bvd_joint - Update design doc to document MV walk model as Option C
- Add confirmed_cases_patch_model to observations.jl: shares receipt delay and test sensitivity with the national confirmed stream, fits against per-province data from spatial tables (Tableau 1) - Wire into bvd_patch_joint: per-province fitting loop routes each patch's onsets through the shared lab pipeline - Add province_confirmed_history to load_observations() return (empty Dict by default, ready for data population) - Export s_test and spec from confirmed_cases_model so per-province models can access them - Export confirmed_cases_patch_model from BVDOutbreakSize.jl
- Add province_confirmed_history TOML block with per-province spatial
table data (Ituri, Nord-Kivu, Sud-Kivu confirmed cases from Tableau 1)
- Add province_history() parser in data.jl that reads per-province arrays
from TOML blocks into Dict{String, NamedTuple}
- Wire province_history() output into load_observations return tuple
- Fix MustNotOverwriteError: use = instead of := for
expected_patch_confirmed (loop runs multiple times)
The per-province confirmed counts are an exact partition of the national confirmed counts (verified: Ituri + Nord-Kivu + Sud-Kivu equals the national total at all 17 shared vintages). Fitting them with their own count likelihood alongside the national confirmed stream put the same observations into the joint density twice, double-weighting the confirmed stream against every other stream. Replace confirmed_cases_patch_model with province_composition_model, which factorises P(y) = P(N) x P(shares | N): the national stream keeps the total term, and the composition term scores only the spatial split, by stick-breaking over patches with an overdispersed Binomial. The vintage totals are conditioned on, never scored. Rt across space: the provincial shares are flat over the whole window (Ituri 91.4% -> 91.1%; Ituri grew 1.82x, Nord-Kivu 1.94x), so the data carry no signal for a time-varying divergence between provinces. Drop the multivariate-normal random walk on joint per-patch log-Rt, which fitted ~30 parameters to that absent signal, discarded the national Rt walk the headline model depends on, and had an inverted intervention sign (half normal on lower=0, raising Rt post-intervention, where rt_walk_model uses upper=0). patch_rt_model now layers a constant per-patch modifier on the unchanged national walk, reference-coded on the primary patch (delta_1 = 0) so the walk level and the mean of delta are not confounded. Importation: epsilon is now sampled only when a non-zero kernel is supplied. There is no mobility data, and under flat shares importation is confounded with the secondary-patch seeds, so the default kernel is zero and no epsilon enters the parameter space. Previously epsilon was sampled against an all-zero kernel, giving a parameter the likelihood never touched. Also: surface C_T, R_T, r, r0, T, CFR, R0 and doubling_time under the same names as bvd_joint so a patch chain drops into summary_table and the existing reporting unchanged; rewrite patch_summary_table, which reported the midpoint of the 20-80 percentile range as a "median"; export the patch API, which was not actually exported. Add test/test_patch_model.jl (198 assertions), covering the renewal equivalence, the importation kernel, the implied-national-Rt identity, the composition's invariance to the modelled level, the epsilon gating, and the headline quantities on a real chain. Co-authored-by: Sam Abbott <contact@samabbott.co.uk>
Rewrite the model methods so the meta-population renewal reads as the model rather than as an addendum. The bolt-on spatial section is gone and its content sits in the sections it belongs to: provincial reproduction numbers under Reproduction number, the per-patch cryptic seed under Seeding and growth, a new Mixing and importation subsection for the gravity kernel and the importation intensity, the patch renewal and the national read-back under Infection process, the export weights under Exported cases, and a new Province compositions subsection under Observation models. Every section now states its maths and its priors, the overview carries prose only, and each notes how the model collapses at one patch. The pooling target is renamed R^trend, the response effect becomes beta_R to free delta for the patch deviations, and the equation tags are renumbered. Also corrects prose that had gone stale against the source: the seed is C_T rather than 2^m, the exports are the weighted patch sum rather than Ituri alone, the patch count is four, and the quoted capital distances are the current patches'. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
…ries The per-province summary figure and its folded table carry the log-Rt deviation, its walk scale and the contrast against Ituri, so they sit with the reproduction-number trajectories rather than with the size tables. The cross-province overview table stays where it was. Cut the prose around both to the fact a reader needs, that the reproduction number and the relative case ascertainment are identified only as a product with the deaths breaking the tie, and drop the descriptions of the layout, the axes and the folded table. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
Drop the issue reference and the projection caveat, leaving what the forecast split is and how each province's count is built. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
`exponential_growth_model` takes the generation interval so it can count the cryptic phase in generations. The merge that brought that in updated `infection_model`'s call site but not the patch composer's, so every `bvd_joint` build on this branch threw a MethodError and no test that builds the model could run.
…ches The executive summary now describes the meta-population renewal and the trend the provinces pool toward, rather than a national single-population process. The reproduction-number results name the pooled patch instead of a Sud-Kivu panel that no longer exists, and the correlation note counts four patches. Drop the y-axis explanation on the modelled-infections figure. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
…urns The helper recomputed the national cryptic curve from `2^m`, which stopped matching when `m` became a count of generations. It now reads the anchor `patch_infection_model` surfaces, so it cannot drift from the growth prior's parameterisation again.
`plot_province_forecast` draws the per-province one-week-ahead forecast as one panel per stream, provinces side by side on a shared axis, each a median dot over nested 30/60/90% credible bars, in the style of `plot_patch_summary`. `province_forecast_archive` writes the same split in the `forecast_archive` long schema plus the province each row is a share of, so a release records the provincial forecast it made alongside the national one. Both read the split through `_province_forecast_draws`, which `province_forecast_table` now uses too, so the table, the figure and the archive multiply the national draw by the province share draw by draw in one place. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
The one-week-ahead split by province is a figure alongside the other forecast figures, with the table kept behind the drop-down for the numbers. The same draws are written to `output/province_forecast.csv`, in the national forecast's schema plus the province each row is a share of, so every release records the provincial forecast it made.
…n the step The dashboard carried the national and per-province reproduction numbers but nothing showing where the epidemic sits. It now also shows modelled infections by province and the per-province summary of size, reproduction number and relative ascertainment. The headline fits drop to a target acceptance of 0.80. The patch fit's cost is trajectory length rather than gradient cost, so a longer step is what buys the adaptation the effective sample size needs.
`interpolate_knots` only reads its knots, so slicing a row out of the knot matrix copied it for nothing and put one `getindex` per knot on the gradient tape. A view removes the copy. The log density and its gradient are bit-identical at both one and three patches. Measured on the isolated kernel (Mooncake gradient of the three-patch knot interpolation at the fitted sizes), 0.135 ms before and 0.103 ms after. The saving is too small to resolve in the full joint gradient. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
Tableau 1 gained a seventh province at SitRep 119 (10 September): Sud
Ubangi, one confirmed case that is already a death, health zones 1/16.
The scanner knew six provinces, so its rows fell one case and one death
short of the national totals at 119, 120 and 121 and it refused all
three vintages. The province blocks have been three vintages behind the
national series since.
Sud Ubangi is a source province, not a patch. It joins the pooled
`other` patch with Tshopo, Sud-Kivu and Bas-Uele, whose population and
capital are already derived from their members, so the source entry is
all that is needed. One confirmed case would inform nothing given its
own reproduction number and ascertainment.
Its population is 2 755 000, the 2019 figure from the Democratic
Republic of the Congo's Institut National de la Statistique, Annuaire
statistique RDC 2020 (March 2021), read from the same Wikipedia
tabulation the other six provinces come from. The six existing figures
were checked against that table and match to the unit, so all seven are
one source. Its capital is Gemena at 3.25651 N, 19.77234 E, the
GeoNames record for the seat of the first-order administrative division
of Sud-Ubangi (geonameid 2315728), the same source as the other six
capitals.
Gemena sits well west of the pooled patch's other members, so the
population-weighted mean capital moves about 1.6 degrees west, from
27.46 E to 25.85 E. That is the kernel reading where the pooled
population is, which is what it is for.
Every vintage from 15 June to 12 September now reconciles exactly with
the national confirmed case and death totals, 82 dates against the
previous 75. Four of the seven new dates (29 and 30 June, 6 and 7
August) are not new vintages but dates the national series has gained
since the province blocks were last written.
The laboratory series does not advance past 9 September. SitReps 119 to
121 word their negative and positive results as "tous se sont reveles
negatifs/positifs", which the bullet parser does not read, and Sud
Ubangi's own 119 bullet gives no sample count at all ("sur l'echantillon
analyse"). Those vintages are left out as before. The block still gains
Sud Ubangi's columns and two earlier dates.
The header quoted a steady 8.5-9% Nord-Kivu case share from an early window. Over the full scanned range it is 8.0-17.6%, against a death share of 13.6-22.6%, so the gap the deaths column identifies narrows rather than holding.
The patch renewal built the whole implied national reproduction number trajectory and read one entry from it. Nothing read the rest: the `Rt_national_implied` deterministic had no consumer in the analysis page, the sensitivity page, the summaries, the plots, the forecast, the tests or the docs prose. `implied_national_Rt_at` returns the single day, and the trajectory form is now written in terms of it. The point is the released chain rather than the gradient. Every draw carried a 212-day vector that nothing ever read, and that vector is now gone from the chain. The quantity never shipped, so no release loses it. The log density and its gradient are bit-identical at one and three patches, since the cut-off day is computed by the same arithmetic in the same order, and the trajectory only ever fed deterministics. Gradient 17.233 -> 16.819 ms at one patch and 20.159 -> 19.933 ms at three, against a run-to-run band of about 0.15 ms. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01DfoEJUruFPcNrWxtLUBBkz
The half-life is documented as the persistence of a provincial divergence but the spatial hyperparameter table did not carry it.
Takes v2.0.0 from main. `docs/examples/analysis.jl` keeps this branch's restructured model section and applies main's form on top: the remaining `@cite` site drops its author prefix for `@citet`, main's abscond competing-risk maths comes in as equation 37, and the equation tags renumber to a contiguous 1-59. `docs/src/news.md` keeps main's released v2.0.0 section and puts this branch's entries under a fresh `Unreleased` heading above it.
The spatial model ships in v2.0.0, so its notes belong in that section rather than under a heading of their own. The generations reparameter- isation moves to Model and the per-fit convergence summary to Report. Corrects three entries: the parameter set is compared against v1, the target acceptance is 0.80, and the seeding is no longer a partition. Completes the analysed-volume entry, which stated the old behaviour and stopped, and adds the provincial forecast figure, its release asset and the two new dashboard figures.
Comments had accumulated change history, arguments against alternatives and restatement of the docstring above them. They now document the current state once. Docstrings keep the summary, the contract and the maths the code does not make clear. The package loses 1,035 lines of prose and no code: the parsed syntax is unchanged in every file but plots.jl, where the estimate-evolution block parsed as a bare string between two helpers and so documented nothing. It now sits on plot_estimate_evolution, which is exported and had none. Corrects two stale claims. bvd_joint said Uganda exports are driven by Ituri alone, where export_pressure weights every patch, and plot_posterior_predictive said four streams lay out as a 2x2 grid, where four or more panels take three columns. Provenance, numerical guards, AD shape constraints and indexing conventions are kept.
The break-day gross check and the centred dispersion default are both chosen on measured sampler behaviour. Stating the choice without the numbers leaves no way to tell whether it still holds, and the runs they came from are not reproducible from the code.
The methods and the sensitivity page had grown a sentence after most equations re-explaining what the equation shows, sentences announcing the subsections that follow them, and arguments against alternatives that were never implemented. About 200 lines of prose come out. Corrections carried with it: the pooled patch has four members since Sud Ubangi joined, the onset curve does carry symptom-onset dates so the report-dating limitation was wrong, province populations and capital coordinates now name their sources, and the seeding says the outbreak is assumed to have begun in Ituri. The per-province summary bullets list one province per line rather than joining them with semicolons.
`province_forecast.csv` has been archived with every release and read by nothing, so the provincial forecast was never scored. Each row is now labelled `<stream> [<patch>]` and fed through the existing pipeline, so the persistence baseline, the CRPS decomposition, coverage and bias all apply to a province as they do to a national stream. A patch's truth is its member provinces' cumulative counts summed vintage by vintage. There is no per-province harmonisation data, so a window holding a break day is skipped rather than corrected, both for the truth and for the baseline's step pool. Releases predating the asset carry no province rows and are counted, not failed.
seabbs
approved these changes
Sep 15, 2026
seabbs-bot
added a commit
that referenced
this pull request
Sep 17, 2026
#717 and #713 each added their own news entry under the v2.0.0 heading when they merged. Re-filing them under v2.1.0 described the same two changes twice on the same page. The gradient entry already sits under v2.0.0 as "The joint gradient costs about 30% less", carrying the same 87,000 to 8,500 hazard-evaluation figures. The suspect-series entry already sits there as the 7 August resumption running to 13 September, cited against issue #708 rather than the pull request. #412 was left out of v2.1.0 on exactly these grounds. This applies the same test to the other two.
seabbs
pushed a commit
that referenced
this pull request
Sep 17, 2026
#717 and #713 each added their own news entry under the v2.0.0 heading when they merged. Re-filing them under v2.1.0 described the same two changes twice on the same page. The gradient entry already sits under v2.0.0 as "The joint gradient costs about 30% less", carrying the same 87,000 to 8,500 hazard-evaluation figures. The suspect-series entry already sits there as the 7 August resumption running to 13 September, cited against issue #708 rather than the pull request. #412 was left out of v2.1.0 on exactly these grounds. This applies the same test to the other two.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
3 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Makes the headline model a meta-population: one renewal equation per province, coupled by importation, with every national stream fitted against the summed provinces.
n_patches = 1collapses it exactly ontomain's single-population model, so there is one model rather than two.What the branch adds over
mainFour patches. Ituri, Nord-Kivu, Haut-Uélé, and
otherpooling Tshopo, Sud-Kivu and Bas-Uélé, which is the province set the situation reports now carry.Reproduction number.
R_{p,t} = R_t exp(δ_{p,t}), withR_tmain's national weekly-knot walk. The deviations share those knots, are correlated across provinces through an LKJ factor, revert toward zero with a shared half-life, and sum to zero at every knot, so no province is the reference and the national level is left toR_t.Seeding. One cryptic epidemic partitioned across provinces:
s_1 = 1/(1+Σf),s_p = f_{p-1}/(1+Σf), withf_p ~ LogNormal(log 0.05, 1). The shares sum to one, so2^mkeeps its elicited meaning as the country's cryptic size andC_Tstays comparable withmain.Importation. A gravity kernel weights each origin's exported transmission by destination population and inverse distance between provincial capitals, holding each origin's total outflow at
1 - N_q/N. The intensity is one level per origin, partially pooled, and changes at detection on the ramp the reproduction number already uses. Coupling is a same-day transfer: the origin is debited exactly what the destinations are credited.National quantities read back, not imposed.
I_t = Σ_p I_{p,t},C_T = Σ_t I_t, and the national reproduction number is recovered by inverting the renewal on the summed infections. The provinces are not rescaled to reproduce the trend, so the country runs at the force-weighted mean of the provincial reproduction numbers.Composition likelihoods. Per-vintage provincial shares of confirmed cases and of confirmed deaths, scored as stick-breaking BetaBinomials on the totals the national streams already carry. Cases identify only the product of relative ascertainment and incidence; the deaths column, under a shared case-fatality ratio and a tight death-ascertainment prior, is what separates them.
Data. Per-province confirmed cases and deaths from Tableau 1, and per-province laboratory throughput, both scanned from the situation-report PDFs and gated on summing exactly to the national totals on every date.
Report. Per-province summary and overview tables, per-province reproduction numbers, the provincial case-fatality comparison, and one-week-ahead forecast scoring by province.
Consequences worth knowing
The provinces are free to disagree with the single-population model, and the national streams are what hold the summed trajectory down. Because the faster province keeps gaining share of the force, the force-weighted mean converges on the fastest province, so the prior-implied national size depends on the deviation scale and the molecular-clock prior sets the trend rather than the country.
With no deviations the provinces reproduce the single-population trajectory to machine precision, which is what the partitioned seed and the transferred importation buy.
Verified
test/test_patch_model.jlcovers the renewal against the single-patch reference, importation transfer, the sum-to-zero deviations, the composition likelihoods on both columns, the per-province tables, and the headline quantities on a patch chain.This was opened by a bot. Please ping @seabbs for any questions.