Stabilize multi-way cluster-robust standard errors

[Copilot]

Note

This branch was previously stacked on #654 (now merged --- see
c67c4d6de, which synced main back in). The diff below now reflects
only this PR's own scope again.

Description

• Cluster-robust standard errors could shrink below the model-based SE and cluster_var = c("commune", "household_id") collapsed to the most granular interaction instead of honoring multi-way clustering.
• This change makes multi-column cluster_var use multi-way cluster combinations, applies CR1 small-sample correction to each one-way subset term by default, and makes the floor to the model-based variance optional and off by default.
• It also adds debug/decomposition output so the multi-way variance terms can be inspected directly.

Implementation

• Refactor cluster-robust variance computation into:
  • one-way cluster score aggregation
  • multi-way combination logic over all cluster-variable subsets
  • term-level decomposition storage for subset variances and signed inclusion-exclusion terms
• Replace the previous single interaction()-only behavior for multi-column clustering.
• Add small_sample = c("none", "CR1") to the one-way subset calculation and thread small_sample, floor_to_standard, and debug_cluster through clustered summary paths, including summary() on stratified fits.
• Attach decomposition details as a cluster_decomp attribute and print concise V_commune / V_household / V_intersection / V_raw / V_final diagnostics when debug output is enabled.
• Preserve existing point estimates while making clustered variance behavior more transparent and configurable.

Regression coverage

• Add focused tests for:
  • CR1 increasing variance when the number of clusters is small
  • decomposition metadata existing and summing back to the raw robust variance
  • optional floor behavior with and without floor_to_standard = TRUE
  • debug output reporting the two-way variance terms
  • nested commune + household_id clustering using tolerance-based diagnostic checks instead of over-strict exact equality
• Add a NEWS entry describing the bug fix and follow-up stabilization.

Example

summary(
  est_seroincidence(
    pop_data = data,
    sr_param = curves,
    noise_param = noise,
    antigen_isos = c("HlyE_IgG", "HlyE_IgA"),
    cluster_var = c("commune", "household_id")
  ),
  debug_cluster = TRUE,
  floor_to_standard = FALSE
)

[codecov]

Codecov Report

❌ Patch coverage is 98.92473% with 2 lines in your changes missing coverage. Please review.
✅ All tests successful. No failed tests found.

Files with missing lines	Patch %	Lines
R/compute_cluster_var_oneway.R	98.11%	1 Missing ⚠️
R/summary.seroincidence.by.R	83.33%	1 Missing ⚠️

Files with missing lines	Coverage Δ
R/combine_cluster_decomp.R	100.00% <100.00%> (ø)
R/compute_cluster_robust_var.R	100.00% <100.00%> (+2.04%)	⬆️
R/est_seroincidence.R	81.48% <ø> (ø)
R/print_cluster_decomp.R	100.00% <100.00%> (ø)
R/summary.seroincidence.R	79.66% <100.00%> (+2.73%)	⬆️
R/compute_cluster_var_oneway.R	98.11% <98.11%> (ø)
R/summary.seroincidence.by.R	88.37% <83.33%> (+1.19%)	⬆️

[kaiemjoy]

@copilot I reviewed the multi-way cluster variance implementation and found a likely issue in R/compute_cluster_robust_var.R inside the subset-combination loop.

Current code uses:

R
(-1)^(n_vars_per_subset[[i]] + 1) * subset_var_log_lambda
But n_vars_per_subset is created via vapply(..., integer(1)), i.e., an atomic integer vector. Using [[i]] here is fragile; this should use single-bracket vector indexing:

R
(-1)^(n_vars_per_subset[i] + 1) * subset_var_log_lambda
Because this term controls inclusion–exclusion signs, wrong indexing here can misassemble the multi-way robust variance.

Suggested patch:

R
robust_var_log_lambda <- robust_var_log_lambda +
(-1)^(n_vars_per_subset[i] + 1) * subset_var_log_lambda
Optional hardening in R/compute_cluster_var_oneway.R:

coerce Hessian to scalar numeric (hessian <- as.numeric(fit$hessian))
guard against non-finite or non-positive Hessian before dividing.

[Copilot]

@copilot I reviewed the multi-way cluster variance implementation and found a likely issue in R/compute_cluster_robust_var.R inside the subset-combination loop.

Current code uses:
...

Addressed in 1a361a9.

[kaiemjoy]

@copilot
Please update this PR to stabilize and debug the cluster-robust SE calculation for multi-way clustering.

Goals
Add term-level decomposition output for multi-way clustering.
Add CR1 small-sample correction to each one-way subset term.
Make the current floor-to-standard behavior optional (default OFF).
Update tests to avoid over-strict equality assumptions.
Required code changes

R/compute_cluster_var_oneway.R
Add argument small_sample = c("none","CR1") (default "CR1").
After computing one-way variance v = score_variance / (hessian^2), apply:
if small_sample == "CR1" and G > 1, v <- v * (G / (G - 1)), where G = length(unique(cluster_ids)).
Return corrected v.
R/compute_cluster_robust_var.R
Add args:
small_sample = c("none","CR1") (default "CR1")
floor_to_standard = FALSE
debug_cluster = FALSE
For each subset term, call .compute_cluster_var_oneway(..., small_sample = small_sample).
Store decomposition rows with:
subset name (e.g. "commune", "household_id", "commune + household_id")
order
sign
subset variance
signed term
Compute:
robust_raw = sum(signed terms)
robust_final = if (floor_to_standard) max(standard_var, robust_raw) else robust_raw
If debug_cluster = TRUE, print a concise decomposition summary and indicate if floor was triggered.
Attach decomposition as attribute on returned value:
attr(ret, "cluster_decomp") <- list(standard_var=..., robust_raw=..., robust_final=..., terms=decomp_df, floor_applied=...).
Thread arguments through callers
In functions that trigger cluster-robust SE (including summary() path and est.incidence.by path), thread through:
small_sample, floor_to_standard, debug_cluster.
Keep backward compatibility (defaults should preserve typical behavior except floor now off by default).
Tests to add/update (tests/testthat/test-cluster_robust_se.R)
CR1 increases variance when clusters are few (synthetic small G case).
Decomposition attribute exists and robust_raw == sum(signed terms) within tolerance.
No exact-equality assumptions like “nested always equals broader exactly”; replace with tolerance/diagnostic checks.
Optional floor behavior:
with floor_to_standard=TRUE, final variance >= standard variance;
with FALSE, return raw robust variance unchanged.
Acceptance criteria
Two-way runs can print V_commune, V_household, V_intersection, V_raw, V_final when debug is enabled.
CR1 is applied per subset term.
Floor is optional and explicitly reported when used.
Tests pass.
Please implement as a clean commit series and summarize what changed and why.

[Copilot]

@copilot
Please update this PR to stabilize and debug the cluster-robust SE calculation for multi-way clustering.

...

Addressed in 76a755d. This adds per-subset CR1 corrections, optional floor_to_standard/debug_cluster handling with decomposition output, threads the new summary arguments through the .by path, and updates the cluster-robust tests to use tolerance/diagnostic checks instead of exact equality assumptions.

[dem-extra1]

@claude review

[d-morrison]

CI fully green (all 17 checks pass) and the review at head c67c4d6de is clean: Ready for merge, confirming review, zero new findings. mergeStateStatus now reads CLEAN.

This is genuinely done now: both the code (findings #1-#4 all fixed and independently re-verified across every review round) and the #654 stacking dependency (merged, synced, diff confirmed back down to this PR's own 19-file scope) are resolved. Nothing else outstanding.

Done for now -- unclaiming. Leaving the merge decision to the maintainer.

[d-morrison]

Reopening work on this -- paws off until I'm done.

[d-morrison]

Per the maintainer's cai correction (and consistent with #654's formal div structure for the rest of the vignette), reorganized this PR's own new content -- the Multi-way clustering, Small-sample correction, and Safeguards sections -- into formal, independently-citable definition/example divs:

• #def-multiway-clustering-variance (the CGM inclusion-exclusion sum)
• #exm-two-way-clustering (the $p=2$ worked example -- matches the existing definition-then-example pattern already used for #exm-hessian-exponential)
• #def-finite-difference-score (the $U_c$ numerical approximation)
• #def-cr1-correction (the CR1 small-sample adjustment)
• #def-variance-safeguards (the floor-at-0 / floor-to-standard / missing-value fallback procedure)

None of this content had a stable id before. As a side effect, replaced a vague "defined above" backward reference with a real (@def-cluster-robust-variance) crossref.

Verified: every crossref in the vignette still resolves (mechanical check across all files), rendered vignettes/methodology.qmd to HTML locally -- clean, zero broken references, each new div renders with its own sequential number (Definition 8-11, Example 3). Spellcheck clean, no new WORDLIST entries needed. Semantic-line-break and punctuation checks clean, run post-commit with the correct three-dot range.

Pushed as ff79091e2.

[d-morrison]

/review

[d-morrison]

Per a follow-up cai correction ("keep divs focused"), moved commentary out of two of the new definition divs and into trailing ::: notes blocks, matching this document's own established convention:

• #def-multiway-clustering-variance's closing sentence about computational cost ($2^p - 1$ evaluations, small $p$) -- a practical scope note, not part of what $V_{\text{multiway}}$ is.
• #def-variance-safeguards's "this is a documented property... not a bug" -- reassurance rather than specification.

Re-rendered, clean, all crossrefs still resolve.

Pushed as b87938747.

[d-morrison]

/review

[d-morrison]

CI fully green (all 17 checks pass) and the review at head b87938747 is clean: Ready for merge, confirming review, no new findings.

This PR has now been through three closing rounds as follow-up cai corrections landed (multi-way clustering math, formal div restructuring, then focused-divs cleanup) -- all resolved. Done for now -- unclaiming. Leaving the merge decision to the maintainer.

[d-morrison]

Stepping back for now -- this PR is code-quality clean (CI green, review "Ready for merge" at head b87938747), but main has since advanced (#659 merged after this PR's last sync). Whoever picks this back up should run a routine git merge origin/main before requesting a final review -- no known conflicts expected, since #659 only touched vignettes/methodology.qmd and this branch's own content lives in vignettes/methodology/_cluster-robust-se.qmd and R/.

[d-morrison]

Merged main in again to resolve a fresh conflict from PR #661 ("rewrite dev-targeted R comments in rendered methodology chunks"). Checked file overlap before merging: #661 only touches DESCRIPTION/NEWS.md/inst/WORDLIST plus three unrelated vignette fragments (_antibody-response-model.qmd, _estimating-incidence-in-each-cluster.qmd, _simulating-a-cross-sectional-survey.qmd) -- no overlap with this PR's own content.

One real conflict, in DESCRIPTION (both branches had independently bumped to different versions past main's prior tip) -- resolved to 1.4.1.9042, exceeding both. NEWS.md and inst/WORDLIST auto-merged cleanly (splice-delta and sort both verified). Re-rendered vignettes/methodology.qmd to HTML -- clean, zero broken crossrefs, since #661 touched other included fragments. No R files touched by this merge.

Pushed as 9b4adc999.

[d-morrison]

/review

[github-actions]

Claude finished review — View run

---

Review of PR #545 — Stabilize multi-way cluster-robust standard errors

Reviewed at head commit 9b4adc9 (9b4adc9998cd5b7b8c8e151b1f425d1822370370).

This is a confirming review. The most recent verdict was Ready for merge, posted at 2026-08-09T22:29:08Z at commit b87938747, itself confirming the chain of prior rounds where all originally-tracked findings (CR1 default, negative-multiway-variance flooring, degenerate-Hessian warning, floor_applied NA-safety, and the stale cluster_var doc paragraph) were verified fixed with regression coverage.

The only change since b87938747 is one merge commit, 9b4adc999 ("Merge remote-tracking branch 'origin/main' into copilot/fix-cluster-robust-standard-errors"), made to resolve a fresh conflict after main advanced (PR #661, "rewrite dev-targeted R comments in rendered methodology chunks"). I independently verified the author's description of this merge rather than trusting the commit message, by pulling the PR's full diff against main (gh pr diff 545) and confirming its contents match exactly what was already reviewed and approved in prior rounds:

• File scope unchanged: still the same 19 files (DESCRIPTION, NEWS.md, five R/ files, five man/*.Rd, inst/WORDLIST, the test file, one vignette fragment, vignettes/references.bib) — no new files, no files dropped.
• DESCRIPTION: the only content change anywhere in this delta is the version bump 1.4.1.9041 → 1.4.1.9042 (resolving the conflict from two branches independently bumping past main's prior tip). The testthat (>= 3.2.0) bump was already reviewed and justified in an earlier round (needed for testthat::local_mocked_bindings()) and is unchanged here.
• NEWS.md: identical single bullet describing the floor-at-0 / missing-SE-on-degenerate-Hessian behavior, already verified accurate against the shipped code in prior rounds. No duplicate or garbled entries from the additional main merge.
• All five R/ files (combine_cluster_decomp.R, compute_cluster_robust_var.R, compute_cluster_var_oneway.R, est_seroincidence.R, print_cluster_decomp.R) and summary.seroincidence.R/summary.seroincidence.by.R: byte-for-byte the same content already independently re-verified in the immediately preceding rounds — negative-variance floor with cli::cli_warn() independent of floor_to_standard, degenerate-Hessian warning before returning NA_real_, structurally-NA-safe floor_applied, and the corrected (non-stale) multi-way cluster_var roxygen doc.
• inst/WORDLIST and all man/*.Rd: unchanged from the previously-reviewed, already-clean state.
• tests/testthat/test-cluster_robust_se.R: the same 84-test suite reviewed and confirmed in prior rounds, unchanged.
• vignettes/methodology/_cluster-robust-se.qmd and vignettes/references.bib: unchanged from the div-restructured, "keep divs focused" version confirmed clean in the last two rounds.

No <<<<<<< conflict markers or other merge artifacts anywhere in the diff. No R, test, or dependency-relevant content changed in this round — only a routine version-number conflict resolution.

I found no new issues in this delta.

Verdict

Ready for merge. The prior Ready for merge verdict (commit b87938747) still stands. The only change since then is a routine git merge origin/main (commit 9b4adc999) resolving a DESCRIPTION version-number conflict, which I independently verified is clean and content-neutral — no R, test, vignette, or dependency changes beyond the version bump itself. All previously-tracked findings remain fixed and tested in the current code.

[github-actions]

💰 Cost: $1.5020 (review) — run

[github-actions]

PR Preview Action v1.8.1
🚀 View preview at https://UCD-SERG.github.io/serocalculator/pr-preview/pr-545/
Built to branch gh-pages at 2026-08-10 07:04 UTC. Preview will be ready when the GitHub Pages deployment is complete.