test(imputation): covariate-path R parity anchor vs didimputation

CHANGELOG.md

@@ -62,6 +62,16 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
   deviation (~0.3% unbalanced); TWFE `hetero` has no public unclustered surface
   (auto-cluster-at-unit convention), so its scenario locks iid — which also pins the D4
   full-K rescale on an UNBALANCED panel for the first time.
+- **`ImputationDiD` covariate-path R parity anchor.** The no-covariate staggered panel was
+  the only `didimputation` R anchor; the covariate branch (first-stage imputation model
+  `y ~ x | unit + time` on the untreated sample, R `first_stage = ~ 0 + x | unit + time`
+  == diff-diff `covariates=["x"]`) now has its own golden: a time-varying, unit-correlated
+  covariate panel appended to `generate_didimputation_golden.R` (the base scenario's RNG
+  draws precede the new block, so the committed base panel and golden values reproduce
+  byte-identically) and a `TestImputationDiDCovariateParityR` class pinning overall +
+  per-horizon event-study ATT and SE. Observed agreement on the reference platform:
+  SE ~2e-10 (the covariate-augmented untreated `v_it` projection + clustering machinery),
+  ATT ~2e-7; asserted at abs=1e-6/1e-7 for cross-platform robustness.
 - **`CallawaySantAnna` ipw R-parity yardsticks folded into the golden fixture + no-covariate
   ipw structural-parity decision recorded.** `csdid_golden_values.json` regenerated (R 4.5.2,
   did 2.5.1, DRDID 1.3.0): all pre-existing data and result blocks reproduced byte-identically;

TODO.md

@@ -62,7 +62,6 @@ generic sparse-FE, QR+SVD rank-detection redundancy, `check_finite` bypass — m
 |-------|----------|--------|--------|----------|
 | SE-audit CI-lock — remaining "fiddly bits" after the second coverage batch landed (that batch pinned C2 `dof_hc2_bm`/`dof_per_coef` via CI-inversion, C3 LOO `df`, C4 estimatr HC1/CR1 intercept SE, C5 Yatchew `p`/`sigma2_lin`/`sigma2_diff`, and the G2 fixest cluster-SE band). Still deferred, each needing a golden regeneration, new computation, or a documented-deviation call — **(a) G2 hetero/cluster is DONE (2026-07-07)**: the unbalanced/heteroskedastic-DGP regen landed, the DiD-path hetero AND cluster CR1 SEs are machine-precision-locked (the plain-OLS CR1 matches fixest exactly — the DOF-convention deviation is absorbed-FE-only), and the TWFE cluster band-pin is retained for the documented non-nested-FE ssc deviation (tracked under "Needs external reference"); TWFE has no public unclustered-hetero surface (auto-cluster convention). Remaining: **(b) PlaceboTests `boundary_gap`** — a permutation randomization-inference margin NOT computed anywhere in code (a new feature + result field, not a coverage lock); **(c) StackedDiD intercept SEs** (`se_cr1/cr2_intercept`, C1) — MEASURED to diverge ~0.3% from R: a nuisance-parameter reference-cell/parameterization gap, NOT machine-precision lockable (the event-study interaction SEs already match ~2e-13; surfacing it would add an unasserted, R-divergent public field); **(d) estimatr `classical` intercept SE** — same documented `O(1/n)` projection/DOF deviation as the slope (reference-only, excluded from parity). Plus the tolerance-tightenings C6-C8. | `benchmarks/R/generate_fixest_did_twfe_golden.R`, `tests/test_fixest_did_twfe_parity.py`, `tests/test_methodology_stacked_did.py`, `tests/test_methodology_placebo.py` | SE-audit | Mid | Low |
 | Render `docs/methodology/REPORTING.md` and `REGISTRY.md` as in-site Sphinx pages so cross-refs can use `:doc:` instead of off-site `blob/main` URLs (stable-docs readers can otherwise land on a different revision than their package version). Two paths: (a) add `myst-parser` to `conf.py` + docs extras and link with `:doc:`, or (b) convert both to `.rst`. **Note:** REGISTRY.md is ~4.5k lines of LaTeX-heavy markdown — high risk under the `-W` (warnings-as-errors) Sphinx build; budget multiple rounds. | `docs/conf.py`, `docs/api/business_report.rst`, `docs/api/diagnostic_report.rst`, tutorials 18 & 19 | follow-up | Mid | Low |
-| `ImputationDiD` covariate-path variance lacks a dedicated parity anchor — only the no-covariate staggered panel is R-parity'd, though the covariate path shares the same validated projection code. Add a small dense-design **hand-calc** for the covariate projection (no external tooling), or a covariate (time-varying X) R `didimputation` golden asserting overall/ES SE parity (the golden variant needs local R). | `tests/test_methodology_imputation.py`, `benchmarks/R/generate_didimputation_golden.R` | imputation-validation | Mid | Low |
 | Add true half-sample BRR replicate-weight regressions per estimator family (current tests use Fay-like 0.5/1.5 perturbations; `test_survey_phase6.py` covers true BRR at the helper level). | `tests/test_replicate_weight_expansion.py` | #253 | Mid | Low |
 | Port the CI `<notebook-prose>` extraction into the reviewer-eval harness so `docs/tutorials/*.ipynb` cases (currently guarded out of `verify-corpus`/`run`) can be reviewed with CI-equivalent context. | `tools/reviewer-eval/adapters/ci_prompt.py` | local-review | Mid | Low |
 

benchmarks/R/generate_didimputation_golden.R

@@ -111,3 +111,68 @@ golden <- list(
 golden_path <- file.path("benchmarks", "data", "didimputation_golden.json")
 write_json(golden, golden_path, auto_unbox = TRUE, pretty = TRUE, digits = 12)
 message(sprintf("Wrote golden: %s", golden_path))
+
+# ---- Covariate scenario (time-varying X in the first stage) ----
+# Appended AFTER the base scenario so the base panel's RNG draws (and the
+# committed base CSV/JSON values) are unchanged on regeneration. Validates
+# the covariate branch of the untreated v_it projection: the first-stage
+# imputation model becomes y ~ x | unit + time on the untreated sample
+# (fixest-style first_stage formula), matching diff-diff's
+# ImputationDiD(covariates=["x"]).
+rows_x <- list()
+uid <- 0L
+for (g in cohorts) {
+  for (j in seq_len(n_per_cohort)) {
+    c_i <- rnorm(1)
+    for (t in seq_len(n_periods)) {
+      beta_t <- 0.5 * t
+      x_it <- 0.4 * c_i + 0.2 * t + rnorm(1)   # time-varying, unit-correlated
+      u <- 0.2 * rnorm(1)
+      treated <- (g > 0L) && (t >= g)
+      eff <- if (treated) tau_h(t - g) else 0.0
+      y <- c_i + beta_t + 0.7 * x_it + eff + u
+      rows_x[[length(rows_x) + 1L]] <- list(
+        unit = uid, time = t, first_treat = g, y = y, x = x_it
+      )
+    }
+    uid <- uid + 1L
+  }
+}
+panel_x <- rbindlist(rows_x)
+panel_x[, unit := as.integer(unit)]
+panel_x[, time := as.integer(time)]
+panel_x[, first_treat := as.integer(first_treat)]
+
+panel_x_path <- file.path("benchmarks", "data", "didimputation_covariate_panel.csv")
+fwrite(panel_x, panel_x_path)
+message(sprintf("Wrote covariate panel: %s (%d rows)", panel_x_path, nrow(panel_x)))
+
+overall_x <- did_imputation(
+  data = panel_x, yname = "y", gname = "first_treat",
+  tname = "time", idname = "unit", cluster_var = "unit",
+  first_stage = ~ 0 + x | unit + time
+)
+overall_x_att <- as.numeric(overall_x$estimate[1])
+overall_x_se <- as.numeric(overall_x$std.error[1])
+message(sprintf("Covariate overall ATT = %.8f (SE %.8f)", overall_x_att, overall_x_se))
+
+es_x <- did_imputation(
+  data = panel_x, yname = "y", gname = "first_treat",
+  tname = "time", idname = "unit", horizon = TRUE, cluster_var = "unit",
+  first_stage = ~ 0 + x | unit + time
+)
+es_x_h <- as.integer(gsub("tau", "", es_x$term))
+ord_x <- order(es_x_h)
+es_x_h <- es_x_h[ord_x]
+es_x_att <- as.numeric(es_x$estimate)[ord_x]
+es_x_se <- as.numeric(es_x$std.error)[ord_x]
+
+golden$covariate <- list(
+  panel = "didimputation_covariate_panel.csv",
+  first_stage = "~ 0 + x | unit + time",
+  overall = list(att = overall_x_att, se = overall_x_se),
+  event_study = list(horizons = es_x_h, att = es_x_att, se = es_x_se)
+)
+
+write_json(golden, golden_path, auto_unbox = TRUE, pretty = TRUE, digits = 12)
+message(sprintf("Re-wrote golden with covariate block: %s", golden_path))