Simplify#53
Merged
Merged
Conversation
Drop the jurisdiction dimension from the timeseries wrapper layer: as_matrix() now returns a plain date-indexed vector instead of a date x jurisdiction matrix, and requires target_infection_dates so every vector is reindexed onto a shared master date axis at construction time rather than having its date range inferred from the data (fill_date_gaps() is removed). The create_epiwave_*_timeseries() constructors drop their jurisdictions parameter accordingly. Part of the jurisdiction-dimension removal refactor (see plan at .claude/plans/logical-waddling-reef.md).
Drop the n_juris_ID indexing parameter now that obs_data/obs_prop are target_infection_dates-aligned vectors rather than jurisdiction-columns of a matrix. This also fixes a pre-existing bug: obs_prop[n_juris_ID] used linear indexing into a date x jurisdiction matrix and silently returned a single scalar (row n_juris_ID) instead of that jurisdiction's proportion series; obs_prop is now indexed consistently with obs_data by date. Also explicitly drops observable dates that fall outside target_infection_dates (e.g. a delay-shifted date walking off the front edge of the window), which previously could produce an NA in a subscripted assignment downstream.
…tion Drop target_jurisdictions and the associated column-matching/reordering logic entirely. Move DOW correction out of this function (it now needs n_jurisdictions, which this layer no longer knows about; DOW gets applied later by stack_jurisdictions()). Build the forward convolution matrix here, directly from this jurisdiction's own delay distribution, replacing the lapply-over-jurisdictions convolution-matrix construction that used to live in create_observation_model()/create_small_sero_model(). Collapse the per-jurisdiction inits loop and matrix-reassembly into a single direct inits_by_jurisdiction() call. Also fixes a field-name bug: seroprevalence data was stored as size_vec but create_small_sero_model() read size_mat (always NULL); size_vec is now properly coerced via as_matrix() and carried through to be renamed/stacked correctly at the jurisdiction-combining step. Drops the vestigial inform_inits field, which was never set by define_observation_data()/ define_sero_data() and never consumed downstream. Manually verified against a synthetic single-jurisdiction stream: all returned vectors are target_infection_dates-length and correctly aligned.
…ction
Drop target_jurisdictions and the dead x parameter (never had a default, no
call site ever supplied it, and the function body never referenced it).
Replace the pmax()-over-2D-arrays / abind-into-3D-array + apply(mean) logic
with plain 1-D Reduce('|', ...) and rowMeans(cbind(...)) now that each
stream's prepare_observation_data() output is a vector rather than a
jurisdiction-matrix. Drop the now-unused abind import.
define_observation_data()/define_sero_data() keep their existing argument
names/count unchanged; only their docs are updated to note that
total_pop/size_vec/proportion_infections are now scoped to a single
jurisdiction per call, since these are called once per jurisdiction going
forward.
Manually verified: a two-stream (cases + hospitalisations) bundle for one
jurisdiction produces correctly shaped, correctly aligned output.
This is the one place jurisdiction becomes a dimension again. It takes a named list of per-jurisdiction observation model bundles (each already 1-D, from define_observation_model()) and stacks their vectors into the date x jurisdiction matrices that create_infection_timeseries(), create_observation_model()/create_small_sero_model(), and create_dow_priors() already expect -- those functions' existing n_jurisdictions-parameterised, shared-hyperparameter behaviour (partial pooling) is unchanged, just fed from a different source. A length-1 list collapses to today's single-jurisdiction case. DOW correction moves here from prepare_observation_data(), since applying it requires knowing n_jurisdictions, which per-jurisdiction prep structurally doesn't. Requires every jurisdiction to supply the same set of streams and identical dow_model flags per stream, erroring otherwise -- both are intentional, documented limitations for now rather than permanent constraints. Also stacks total_pop/size_vec into total_pop/size_mat for seroprevalence streams. Because every per-jurisdiction vector coming out of layer 1/2 is already target_infection_dates-length and positionally aligned (the master-date-axis invariant enforced by as_matrix()), stacking here is a plain cbind() with no date-reconciliation logic needed. Manually verified with two jurisdictions on staggered, non-identical date ranges, a seroprevalence stream, and dow_model = TRUE on the cases stream: shapes are correct, the alignment check confirms row i means the same calendar date in every jurisdiction's column, the hierarchical DOW branch (create_dow_priors(2)) runs without error, and the mismatched-dow_model / unnamed-list error paths fire with clear messages.
create_small_sero_model() Both functions now read convolution_matrices, case_mat, prop_mat (and, for sero, size_mat/total_pop) directly from stack_jurisdictions()'s output instead of rebuilding a per-jurisdiction convolution matrix list themselves via lapply(unique(delays$jurisdiction), ...) -- that logic now lives once, in prepare_observation_data(), rather than being duplicated across these two functions. Also deletes the now-dead data_idx/expected_cases_idx trimming step (it inferred which rows of a full-length convolution output to keep by matching case_mat's rownames-derived date range; case_mat is now guaranteed target_infection_dates-length and positionally aligned by construction, so expected_cases is already the right shape). This incidentally fixes create_small_sero_model(), which was unreachable in practice: its convolution-matrix block called new_convolution_matrix(delays, x, n_dates), a 3-argument call that doesn't match the current 2-argument new_convolution_matrix(pmf, n) signature, and it read a size_mat field that prepare_observation_data() never actually set (it stored size_vec instead -- fixed in a prior commit). Both are resolved automatically now that convolution matrices and size_mat arrive pre-built and correctly named. Dropped the now-unused infection_days parameter from both functions (it was only used by the deleted convolution-matrix-building and data_idx logic). Manually verified: building the greta graph for both a cases stream and a sero stream from stacked two-jurisdiction data produces correctly-named, correctly-shaped (150 x 2) outputs with no errors.
…gthscale Rename the observations parameter to observations_by_jurisdiction and make the first line of the function stack_jurisdictions(observations_by_jurisdiction) -- this is what lets the common single-jurisdiction case stay a one-call ergonomic path (a length-1 named list) while keeping the stacking logic factored out into its own, independently-testable function. Fix the sero/count dispatch: previously create_small_sero_model()'s call was commented out, so every stream -- including seroprevalence -- was unconditionally run through create_observation_model()'s negative-binomial cases likelihood. Streams are now dispatched by the presence of total_pop. Fix greta::initials(gp_lengthscale = rep(0.5, n_jurisdictions)): gp_lengthscale is a scalar node in create_infection_timeseries() regardless of n_jurisdictions (it's a shared kernel hyperparameter across jurisdiction columns), so rep(..., n_jurisdictions) would produce a dimension mismatch on the first real n>1 fit; this was previously never exercised since the package has never been run with more than one jurisdiction. Manually verified end-to-end: a single-jurisdiction synthetic fit (flat_prior, small MCMC settings) runs to completion with correctly-shaped output.
Run devtools::document() after the R/ changes: exports stack_jurisdictions, drops importFrom(abind, abind)/importFrom(methods, is)/importFrom(rlang, .data)/importFrom(tidyr, pivot_wider) (all now unused after the jurisdiction- dimension refactor), and adds importFrom(greta, binomial) (now used directly in create_small_sero_model()'s roxygen tags). Drop abind/methods/rlang from DESCRIPTION Imports accordingly (confirmed via grep no remaining usage of any of the three anywhere in R/). Drop "jurisdiction" from R/epiwave-package.R's globalVariables() -- its only NSE use site (dplyr::filter(jurisdiction == ...) in the old inits_by_jurisdiction()) was removed in an earlier commit. Verified via devtools::check() that this doesn't reintroduce a NOTE (the remaining `~jurisdiction` reference in plot_infection_traj.R's facet_wrap() is a formula, not flagged by codetools the way NSE dplyr columns are). devtools::check() otherwise shows 0 errors; the remaining WARNINGs/NOTEs (undeclared cli/distributional imports, hidden .github, draw/value bindings in plot_infection_traj.R) all trace to files this refactor doesn't touch and are pre-existing on the simplify branch.
These aren't automated tests (no tests/testthat/ suite exists), but they're the primary way this pipeline gets manually exercised against real data, so keep them runnable: drop jurisdictions=/target_jurisdictions= args from create_epiwave_greta_timeseries()/create_epiwave_massfun_timeseries()/ define_observation_model() calls, and wrap each single-jurisdiction define_observation_model() bundle in a named list (setNames(list(...), jurisdictions)) before passing it to fit_waves(observations_by_jurisdiction = ...). testing.R's sero total_pop = c(8e6, 7e6) becomes a scalar (total_pop = 8e6), reflecting that define_sero_data() is now called once per jurisdiction. Also fixes a pre-existing typo in testing.R (fit_waves(..., infection_model = 'gp_growth_rate')) to the actual parameter name, infection_model_type -- unrelated to this refactor, but the line was already being touched.
… wrapper epiwave.params now has native discrete_pmf_series/discrete_weights_series objects (new_discrete_series(), with their own date-based subsetting, validation, print/summary methods). epiwave's own epiwave_massfun_timeseries wrapper predates these and is now purely redundant: prepare_observation_data() was coercing delay_from_infection into that wrapper tibble and then immediately reconstructing a discrete_pmf_series from it via new_discrete_series() -- a pointless round-trip. Remove create_epiwave_massfun_timeseries()/epiwave_massfun_timeseries entirely. prepare_observation_data() now accepts delay_from_infection as either a single discrete_pmf/discrete_weights object (replicated across target_infection_dates via new_discrete_series()) or an already time-varying discrete_pmf_series/discrete_weights_series (aligned via the series' own Date-based subsetting), and builds the convolution matrix directly from that series object. Widening to accept discrete_weights (not just discrete_pmf) matters for the seroprevalence pathway specifically: unlike case/hospitalisation notification (a one-time event, correctly modelled as a normalised discrete_pmf), seroconversion is typically persistent -- a person may test positive for many consecutive days -- so it's better represented as an unnormalised discrete_weights curve.
delays is now a discrete_pmf_series/discrete_weights_series object (see previous commit), not a data.frame/tibble, so dplyr::filter(delays, date %in% case_dates) no longer works -- replaced with the series' own Date-based subsetting, delays[case_dates]. Also fixes a live bug: expected_delay_vals was computed as sum(x$delay * x$mass), but discrete_pmf objects have columns step/prob, not delay/mass -- those fields don't exist, so this silently evaluated to 0 every time (the mean delay used to shift observed dates into inferred infection dates was always treated as zero, regardless of the actual delay distribution). Now uses epiwave.params's mean.discrete_pmf(), generalised to discrete_weights via epiwave.params::normalise() first (weights aren't a proper distribution, so they're normalised to a pmf before taking a mean). Verified directly: for a gamma(shape=3, rate=0.5) delay (true mean 6, vs 0 before this fix), expected_delay_vals now correctly computes ~6.
Mechanically identical to the discrete_pmf path (a day-difference matrix looked up against the object's step column), just using $weight instead of $prob. This is what makes the seroprevalence pathway able to use an unnormalised persistence curve (discrete_weights/discrete_weights_series) instead of being forced into a normalised discrete_pmf. Verified standalone: new_convolution_matrix() on a discrete_weights object produces row sums that don't sum to 1 (confirming weights aren't being force-normalised), and a discrete_weights_series built from a single replicated discrete_weights object produces an identical matrix to the single-object path.
…adoption
Document the discrete_pmf vs discrete_weights choice in
create_small_sero_model()/define_observation_data()/define_sero_data()'s
roxygen: sero streams should typically supply delay_from_infection as a
discrete_weights/discrete_weights_series persistence curve, not a normalised
discrete_pmf.
Fix tests/test_workflow/{testing.R,test2.R}: epiwave.params::add_distributions()
no longer exists (renamed to add_discrete()/the + operator) -- these calls
were already stale before this branch, unrelated to the jurisdiction
refactor, but directly relevant now that this pass touches the discrete
object plumbing throughout. Also drop testing.R's now-unnecessary
create_epiwave_massfun_timeseries() wrapping step, since
prepare_observation_data() accepts a raw discrete_pmf directly.
devtools::document() after the discrete-series adoption changes: drops the create_epiwave_massfun_timeseries export/man page, updates man pages for the touched functions. Also declares cli (Imports) and distributional (Suggests) in DESCRIPTION -- new_convolution_matrix.R calls cli::cli_abort()/cli::cli_warn() directly (including two new call sites added in this pass) and its @examples block uses distributional::dist_gamma(), neither of which were formally declared despite being used directly (previously only working because epiwave.params happens to depend on both transitively). devtools::check() now shows 0 errors, 0 warnings (down from 2 warnings); the remaining 3 NOTEs are pre-existing and trace to files this refactor doesn't touch.
Remove jurisdiction as a threaded dimension
…series The date-alignment check in prepare_observation_data() assumed any already- classed epiwave_timeseries object stores its dates in a flat $date column, but epiwave_greta_timeseries objects (e.g. the IHR-from-CHR pattern used for hospitalisation proportions in tests/test_workflow/testing.R and test2.R, built via create_epiwave_greta_timeseries()) are a list wrapping a greta array alongside the date tibble, with dates nested at $timeseries$date instead. The check was comparing as.Date(NULL) against target_infection_dates and always failing with "`proportion_infections` dates must match `target_infection_dates`". Found by actually running tests/test_workflow/testing.R's basic (cases + hospitalisations) fit against real data in ../data -- this is a core supported pattern (proportion driven by a greta array), not an edge case, and would have broken on the first real multi-stream fit using it. Verified: the same fit now runs to completion against real data.
Adds as_epiwave_timeseries(data): a new exported helper that takes a plain
data.frame/tibble with date/value columns (partial date coverage is fine)
and classes it as epiwave_fixed_timeseries -- replacing the fragile,
unvalidated class(x) <- c('epiwave_fixed_timeseries', 'epiwave_timeseries',
class(x)) pattern users previously had to hand-write for every observation
stream. Already-classed epiwave_timeseries objects and plain numeric
values/vectors pass through unchanged, so it's safe to call unconditionally.
prepare_observation_data() now calls this automatically on timeseries_data
and size_vec (mirroring the auto-coercion that delay_from_infection/
proportion_infections already had), so users can pass a plain data.frame
straight into define_observation_data()/define_sero_data() without ever
manually setting a class themselves. This directly targets the messy
"before define_observation_data()" data prep in tests/test_workflow/
testing.R, where every observation stream needed its own hand-rolled
class(x) <- ... line.
Verified: a raw unclassed data.frame produces identical output to the
equivalent pre-classed object; a data.frame missing the required columns
errors with a clear message instead of failing deep inside as_matrix().
define_observation_data() now accepts a plain date/value data.frame directly
(via as_epiwave_timeseries(), added earlier on this branch), so the
hand-rolled class(notif_dat) <- c('epiwave_fixed_timeseries',
'epiwave_timeseries', class(notif_dat)) lines are no longer needed.
Verified: re-ran testing.R's basic (cases + hospitalisations) fit against
real data in ../data with the classing removed -- notif_dat/hosp_dat stay
plain tbl_df objects throughout, and the fit runs to completion identically.
testing.R/test2.R depend on local, not-synced data (../data, simdata/) that isn't in the repo -- useful as personal real-data scripts, but not runnable by anyone else who clones it, and don't exercise the multi-jurisdiction path at all (the package has never had a multi-jurisdiction reprex before this branch). Adds two fully self-contained scripts (fabricated data, no external dependencies) that demonstrate the current API end to end: - single_jurisdiction_workflow.R: raw data.frames passed straight to define_observation_data() (no manual class(x) <- c(...)), delay distributions built via epiwave.params (discrete_pmf, combined with `+`), a greta-array-derived proportion (IHR-from-CHR), and an "advanced" section showing a time-varying discrete_pmf_series delay. - multi_jurisdiction_workflow.R: two jurisdictions with deliberately staggered/non-identical date coverage, dow_model = TRUE to exercise hierarchical DOW pooling, and an explicit alignment check confirming a date only one jurisdiction has data for stays correctly non-NA/NA in the right columns. Both run to completion (small MCMC settings for speed) with their stopifnot() checks passing.
Captures why the seroprevalence pathway (define_sero_data(), create_small_sero_model(), the discrete_weights widening in new_convolution_matrix()/evaluate()) is being pulled back out of this MVP, the recommended architecture for re-integrating it (merge into define_observation_data()/create_observation_model() rather than parallel functions with if-based dispatch), the real bugs found in the original sero code during this work, and what validation is still needed before it ships. dev/ is excluded from R CMD build via .Rbuildignore, matching the existing pattern for .claude/.positai.
Per discussion: don't ship if-statement dispatch for a not-yet-validated feature. Removes define_sero_data(), create_small_sero_model(), the total_pop/size_vec plumbing through prepare_observation_data() and stack_jurisdictions()'s stack_stream(), the model_fn dispatch in fit_waves(), and the discrete_weights/discrete_weights_series widening in new_convolution_matrix()/evaluate() (which had no other consumer once sero is removed). inits_by_jurisdiction()'s mean_step() helper simplifies back to a plain mean(x) call accordingly. The design for how to re-integrate this cleanly (merge into define_observation_data()/create_observation_model() rather than parallel functions with if-based dispatch, as agreed before removing it), the real bugs found in the original sero code, and what's still needed before it ships are captured in dev/sero-integration-notes.md (previous commit). fit_waves()/stack_jurisdictions()/prepare_observation_data() now only handle the validated cases/hospitalisations pathway.
testing.R's sero fit block called define_sero_data(), which no longer exists after the previous commit -- it was already non-functional (sero_dat/ sero_size_mat/sero_conversion were never assigned), so this just removes dead code rather than breaking anything that worked. Left short pointer comments to dev/sero-integration-notes.md at both spots (data prep, and where a third define_observation_data() block would go) for whoever picks this back up. devtools::document() drops create_small_sero_model/define_sero_data exports and man pages, and the now-unused importFrom(greta, binomial). Verified after all sero-removal changes: devtools::check() is 0 errors, 0 warnings (same 3 pre-existing NOTEs); testing.R's basic fit still runs against real data in ../data; both single_jurisdiction_workflow.R and multi_jurisdiction_workflow.R still run to completion with their stopifnot() checks passing.
Two changes, implemented together since they touch the same functions: 1. stack_jurisdictions() becomes the explicit, user-facing combining step, taking named ... args (mirroring define_observation_model()'s own cases=/hospitalisations= pattern) instead of a pre-built named list -- e.g. stack_jurisdictions(VIC = obs_vic, NSW = obs_nsw). A single jurisdiction needs no combining step at all: define_observation_model()'s output (now stamped with class epiwave_observation_model) goes straight to fit_waves(), which detects whether it already received a stacked object (class epiwave_stacked_observations) or a raw single-jurisdiction one and wraps the latter internally (auto-labelled, since there's no jurisdiction identity to preserve for n=1). This removes the setNames(list(...), jurisdictions) step that was previously required even for a single jurisdiction. 2. GAM-based initial values (inits_by_jurisdiction(), and the cross-stream union/mean that combines them) are only ever consumed by infection_model_type = 'flat_prior' -- every GP-based infection model ignores them entirely (confirmed in create_infection_timeseries(): observable_infection is referenced only inside the flat_prior branch). Previously this ran unconditionally for every stream x jurisdiction during data prep, running an mgcv::gam() fit that was silently discarded whenever a GP model was used instead. prepare_observation_data() no longer computes inits at all (keeps delays in its output instead, needed later); define_observation_model()/stack_jurisdictions() no longer eagerly combine them either. A new compute_flat_prior_inits() does this work on the fully-stacked object, called lazily by fit_waves() only inside its flat_prior branch, before create_infection_timeseries() (which needs the result) and reused for the MCMC initial values. Verified: a standalone check exercising both a single-jurisdiction fit (flat_prior AND gp_growth_rate, confirming the GP path never triggers inits computation) and a two-jurisdiction stack_jurisdictions() combination (plus the unnamed-args error path) all pass. devtools::check() is 0 errors, 0 warnings. testing.R's real-data fit (../data) verified against both flat_prior and gp_growth_rate. Both synthetic workflow scripts updated to the new API and re-verified end to end.
compute_flat_prior_inits() was reading prop_mat from the stacked object, but stack_stream() applies DOW correction to prop_mat before returning it -- so for any stream with dow_model = TRUE, inits_by_jurisdiction() was being handed a greta array (the DOW-corrected proportion) instead of a plain numeric one. It didn't error (there's already a greta_array branch in inits_by_jurisdiction()), but that branch draws from the *prior* predictive distribution of the DOW effect via greta::calculate(nsim = 100) -- no MCMC has run yet at this point -- making the "deterministic smoothed guess" inits computation silently stochastic, and adding an unnecessary 100-draw simulation per DOW-modelled stream x jurisdiction whenever flat_prior is used. Checked against the original, pre-refactor code to confirm this wasn't just a style difference: it explicitly computed inits from prop_mat *before* applying DOW correction, deliberately. My first refactor (PR #51) preserved this ordering (DOW correction lived in stack_jurisdictions(), which ran after prepare_observation_data()'s inits computation); deferring inits to run after stacking (previous commit) inverted it by accident. Fix: stack_stream() now also returns prop_mat_raw (captured before DOW correction is applied), and compute_flat_prior_inits() reads that instead of prop_mat. create_observation_model() is unaffected -- it still reads the DOW-corrected prop_mat for the actual likelihood, as it should. Verified: prop_mat_raw is a plain numeric column (not a greta_array) even when dow_model = TRUE, so inits_by_jurisdiction() no longer needs greta::calculate() at all for this case. Re-ran both synthetic workflow scripts and testing.R's real-data fit (flat_prior + dow_model = TRUE, the exact previously-affected combination) end to end. devtools::check() is still 0 errors, 0 warnings.
It's called from exactly one place (define_observation_model()) and never appears in any workflow script -- it's implementation detail, the same category as stack_stream()/compute_flat_prior_inits(), which are already @nord. Exporting it overstated it as public API and was contributing to a "why are there three equally-important functions here" feeling, when really only define_observation_data()/define_observation_model() (plus stack_jurisdictions()/fit_waves()) are meant to be called directly. No logic changed -- devtools::check() still 0 errors/0 warnings, both synthetic workflow scripts re-verified end to end.
Same situation as prepare_observation_data(): called from exactly one place (compute_flat_prior_inits()), never from a workflow script. Pure internal machinery. No logic changed -- devtools::check() still 0 errors/0 warnings, both synthetic workflow scripts re-verified end to end.
… simplify prop coercion create_epiwave_fixed_timeseries() had zero call sites anywhere in the package or workflow scripts, and did nothing as_epiwave_timeseries() doesn't already do better (it takes a real date/value table directly, rather than requiring separate dates + value arguments). create_epiwave_timeseries() had exactly one remaining caller, in prepare_observation_data()'s proportion_infections coercion -- but that coercion step was itself redundant: as_matrix() already has an as_matrix.numeric method that recycles a scalar or validates a same-length vector against target_infection_dates, identically to what going through create_epiwave_timeseries() first produced (verified: identical() output both ways, for both the scalar and vector cases). prepare_observation_data() now only branches on already-classed epiwave_timeseries objects (for date validation); a bare numeric prop passes straight to as_matrix(), which dispatches correctly on its own. create_epiwave_timeseries.R now holds only create_epiwave_greta_timeseries() (kept in epiwave rather than epiwave.params specifically because it depends on greta, which epiwave.params should not) and as_epiwave_timeseries(). Verified: both synthetic workflow scripts and testing.R's real-data fit produce identical output to before these changes (e.g. single-jurisdiction posterior median range unchanged at 0-7185). devtools::check() still 0 errors, 0 warnings.
…-export as_epiwave_timeseries() Renamed to pair naturally with as_epiwave_timeseries() now that the two create_*/dead functions removed from this file are gone -- as_* matches the coercion-style naming already used elsewhere (epiwave.params::as_discrete_pmf(), as_matrix()). The class it produces is renamed to match: epiwave_greta_timeseries -> greta_timeseries (still inherits from epiwave_timeseries, so existing inherits(x, 'epiwave_timeseries') checks that need to catch both plain and greta-backed timeseries objects continue to work unchanged). as_matrix.epiwave_greta_timeseries() is renamed to as_matrix.greta_timeseries() to match, and prepare_observation_data()'s inherits() check on proportion_infections updated accordingly. as_epiwave_timeseries() becomes internal (@nord): checked and confirmed no workflow script actually calls it by name -- users get its behaviour indirectly (passing a raw data.frame to define_observation_data(), which coerces internally). It remains available as an unexported helper rather than being deleted, since it's still the one general-purpose coercion path, just not meant as public API for now. Verified: class(ihr) is now `greta_timeseries epiwave_timeseries list`, dispatch to as_matrix.greta_timeseries() works correctly, both synthetic workflow scripts produce identical output to before the rename (e.g. single-jurisdiction posterior median range unchanged at 0-7185), and testing.R's real-data fit runs to completion. devtools::check() still 0 errors, 0 warnings.
Two minimal examples -- single jurisdiction, and multiple jurisdictions combined via stack_jurisdictions() -- mirroring tests/test_workflow/single_jurisdiction_workflow.R and multi_jurisdiction_workflow.R, trimmed to the essential shape for a README. Both verified to run to completion exactly as written before adding them; chunks are eval = FALSE in README.Rmd since running them needs a full greta/Python setup, which shouldn't be a precondition for a routine README re-knit. README.Rmd had been missing from the repo (only README.md existed) despite a pre-commit hook expecting both to move together; it's restored here from an old, partially-filled scaffold (its Installation section still suggested pak::pak() and its Example section was a never-filled-in placeholder). Updated Installation to match the remotes::install_github() already decided in README.md, and replaced the Example placeholder with the new Usage section. Left the empty Citation section and placeholder Support text untouched rather than inventing content, and left the informal draft paragraph after the srr-tags chunk untouched (out of scope for this change). Also gitignores README.html, a preview-render artifact that isn't source.
Clean up workflow ergonomics, defer flat_prior-only inits, drop sero for now
Consolidates the previous "one jurisdiction" / "multiple jurisdictions" examples into a single example showing the full parallel structure clearly: each jurisdiction bundles the same two streams (cases, hospitalisations) via define_observation_data() calls inside one define_observation_model() call, and the two jurisdictions (with deliberately non-identical date coverage) are combined explicitly via stack_jurisdictions(). No wrapper/helper function around define_observation_data() this time (the actual tests/test_workflow/multi_jurisdiction_workflow.R script uses one for brevity) -- every call is spelled out in full so a reader can see the real API shape directly, at the cost of some repetition between jurisdiction A and B. Verified end to end before adding: both jurisdictions' cases and hospitalisations streams fit correctly together.
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
1 participant
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.
Summary
This brings
mainup to date with everything built onsimplifyso far — the jurisdiction-dimension refactor and its follow-up work (PRs #51, #52, plus direct commits). Highlights:Jurisdiction is no longer a dimension threaded through every function. It used to be a
date x jurisdictionmatrix baked into nearly every function (convolution operators as lists of matrices, inits assembled viaabindinto 3D arrays,sweep/apply(MARGIN=2)broadcasting everywhere). Now:prepare_observation_data()/define_observation_model()are purely 1-D (one jurisdiction, one stream at a time);stack_jurisdictions()is the one explicit place jurisdiction becomes a dimension, taking named...args (stack_jurisdictions(VIC = ..., NSW = ...)) and combining them into the matrices the model-fitting functions need. A single jurisdiction needs no combining step at all —define_observation_model()'s output goes straight tofit_waves(). Partial pooling (shared GP kernel hyperparameters, shared hierarchical day-of-week prior) is preserved throughout.A master-date-axis invariant ensures every per-jurisdiction vector is aligned to a shared date sequence at construction time, so jurisdictions with different/staggered data coverage stack safely (verified with non-identical date ranges specifically to catch misalignment).
epiwave.params's native discrete series objects are adopted directly (discrete_pmf/discrete_pmf_series), replacing a now-redundant internal wrapper (epiwave_massfun_timeseries) that just round-tripped through them.Deferred, model-type-specific computation. GAM-based initial values (only ever used by
infection_model_type = 'flat_prior') no longer run unconditionally for every stream x jurisdiction during data prep — they're computed lazily, only when actually needed, viacompute_flat_prior_inits().Several real bugs found and fixed along the way: a linear-indexing bug in proportion lookups, an NA-subscript bug when a delay pushes an inferred date off the front of the date window, a date-alignment check that broke for greta-array-backed proportions, and — found via a colleague's review — inits being computed from a DOW-corrected proportion instead of the raw one (inits should be deterministic, not depend on a prior-predictive draw of the DOW effect).
Seroprevalence support was built, then deliberately pulled back out of this pass — it was only ever validated against fabricated data, and we didn't want unvalidated
if-dispatch logic in the MVP's core functions. The design for reintroducing it cleanly is documented indev/sero-integration-notes.md(excluded from the package build).API surface cleanup: several internal-only functions (
prepare_observation_data(),inits_by_jurisdiction()) are no longer exported; two fully dead functions (create_epiwave_timeseries(),create_epiwave_fixed_timeseries()) are deleted;create_epiwave_greta_timeseries()is renamed toas_greta_timeseries()to pair with the newas_epiwave_timeseries()auto-coercion helper (which replaces the need to hand-writeclass(x) <- c(...)on raw observation data).New self-contained example workflows (
tests/test_workflow/{single,multi}_jurisdiction_workflow.R, no external data needed) and a README Usage section with a full two-jurisdiction, two-stream walkthrough.Test plan
devtools::check()— 0 errors, 0 warnings at every stage (3 pre-existing NOTEs, unrelated to this work, remain)testing.R/test2.R, against local not-synced data) re-verified end to end after every changedow_model = TRUEcombined withflat_prior(the exact combination that exposed the DOW/inits bug), hierarchical DOW and GP pooling for n>1 jurisdictions, and GP-based models never triggering the GAM-based inits computation at all🤖 Generated with Claude Code