From ad2454241aded653b23e3235b8a4f3460a2493e6 Mon Sep 17 00:00:00 2001 From: igerber Date: Sun, 17 May 2026 06:44:44 -0400 Subject: [PATCH 1/4] estimators: lift MultiPeriodDiD-absorb HC2/HC2-BM via auto-route to fixed_effects= Mirrors PR #458 (DiD-absorb auto-route) on MultiPeriodDiD: when absorb= is paired with vcov_type in {hc2, hc2_bm}, the fit promotes the absorb columns to fixed_effects= internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (treated + period_X dummies + treated:period_X interactions + factor(unit)). Verified at ~1e-15 vs lm() + sandwich::vcovHC(type="HC2") and lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2") on a new 5-cohort x 5-period mpd_absorbed_fe_did fixture. Includes three-guard reorder so the auto-route sits BETWEEN the absorb + fixed_effects mutual-exclusion check (above) and the multi-absorb + survey-weights reject (below), matching the DiD ordering. The survey-replicate absorb-refit branch at estimators.py:1689 is short- circuited under the auto-route (the standard compute_replicate_vcov path applies on the fixed full-dummy design; no per-replicate refit needed). Tests: new TestMPDAbsorbedFERParity class (7 tests) mirrors PR #458's TestDiDAbsorbedFERParity, pinning parity targets on per-period interaction coefficients (treated:period_4) to avoid the treated x unit collinearity baked into MPD's time-invariant ever-treated indicator. Existing test_multi_period_absorb_rejects_hc2_and_hc2_bm deleted. REGISTRY.md per-estimator status block updated (MPD moves REJECT -> SUPPORTED; TWFE remains the only REJECT case). TODO row 99 narrowed to TWFE-only. Co-Authored-By: Claude Opus 4.7 (1M context) --- CHANGELOG.md | 1 + TODO.md | 2 +- benchmarks/R/generate_clubsandwich_golden.R | 50 ++++ benchmarks/data/clubsandwich_cr2_golden.json | 17 +- diff_diff/estimators.py | 68 ++++-- docs/methodology/REGISTRY.md | 4 +- tests/test_estimators_vcov_type.py | 230 ++++++++++++++++--- 7 files changed, 319 insertions(+), 53 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 39888e0d..2ab0c506 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,6 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] ### Added +- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that is collinear with unit dummies, so the `treated` main-effect coefficient becomes NaN under solve_ols's rank-deficiency handling — this is expected; per-period interaction parity targets (`treated:period_X`) are unaffected. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. - **`DifferenceInDifferences(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:382`). Previously raised `NotImplementedError` because the HC2 leverage correction and CR2 Bell-McCaffrey DOF depend on the FULL FE hat matrix, while within-transformation (FWL) preserves coefficients and residuals but not the hat. Lift via internal auto-route: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, the fit promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov. Empirically matches `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=..., type="CR2")` at ~1e-10 (verified via new `tests/test_estimators_vcov_type.py::TestDiDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `absorbed_fe_did`, with the R generator using the singleton-cluster CR2 trick for one-way HC2-BM Satterthwaite DOF). HC1/CR1 paths unchanged. `MultiPeriodDiD(absorb=...)` and `TwoWayFixedEffects` rejections remain as follow-ups (different fit-path structure). **Behavioral note (full `DiDResults` surface change under auto-route):** under the auto-route, the entire returned `DiDResults` reflects the full-dummy fit rather than the within-transformed fit. Specifically, `result.coefficients` and `result.vcov` include the FE-dummy entries (matching the `fixed_effects=` path), `result.residuals` and `result.fitted_values` are on the un-demeaned outcome scale, and `result.r_squared` is computed on the un-demeaned outcome (so it absorbs the FE variance and will typically be higher than the within-R²). `result.att` is invariant to this routing (FWL guarantee). Downstream consumers reading `result.att` are unaffected; consumers reading the broader result surface should expect the full-dummy values. **Survey-design scope:** the auto-route changes the FE handling (and removes the prior absorbed-FE rejection), but `survey_design=` continues to drive its own variance path (Taylor-series linearization or replicate-weight variance, per the existing survey contract) rather than the analytical HC2/HC2-BM sandwich. The auto-route is therefore methodologically meaningful for non-survey fits and for the FE-handling side of survey fits; analytical small-sample inference under `vcov_type in {"hc2","hc2_bm"}` is bypassed when a survey design is supplied. - **BaconDecomposition R parity goldens.** Closes the PR-B deferral row in `TODO.md`. JSON goldens at `benchmarks/data/r_bacondecomp_golden.json` generated from the committed `benchmarks/R/generate_bacon_golden.R` script (3 fixtures: `uniform_3groups_with_never_treated`, `two_groups_no_never_treated`, `always_treated_remapped`) against `bacondecomp 0.1.1` on R 4.5.2. `tests/test_methodology_bacon.py::TestBaconParityR` now active (4 tests, no skips): TWFE coefficient parity at `atol=1e-6` across all 3 fixtures; weights-sum parity at `atol=1e-6` across all 3 fixtures; per-component estimate + weight parity at `atol=1e-6` on the 2 non-remap fixtures **and on the 6 timing-vs-timing rows of `always_treated_remapped`** (carve-out narrowed to U-bucket rows only); plus a dedicated fold-back test (`test_always_treated_remapped_fold_back_matches_r`) that pins the **documented convention divergence** on `always_treated_remapped` (R keeps `first_treat=1` as a distinct timing cohort and emits `Later vs Always Treated` comparisons; Python's paper-footnote-11 convention remaps those units to `U` and folds them into a single `treated_vs_never` cell per treated cohort) by aggregating R's split rows per cohort and asserting they match Python's single fold at `atol=1e-6`. The aggregate is invariant per Theorem 1; the per-component breakdown differs structurally between conventions but the fold-back is now directly asserted. New `**Note (R parity convention divergence on always-treated)**` and `**Deviation (first-period boundary extension on always-treated remap)**` in `docs/methodology/REGISTRY.md`. **First-period boundary deviation:** the paper uses strict `t_i < 1` for the always-treated bucket; the library uses the inclusive `first_treat <= min(time)` rule and folds `first_treat == min(time)` cohorts into `U`. R does NOT apply this fold (it keeps such cohorts as their own bucket). When `min(time) > 1` the rules coincide. Explicitly labeled in REGISTRY's Deviations block and mirrored in `METHODOLOGY_REVIEW.md` and `bacon.py`. METHODOLOGY_REVIEW.md tracker row promoted `**Complete** (R parity goldens pending)` → `**Complete**`. - **`generate_ddd_panel_data` — panel-structured DGP for Triple-Difference power analysis** (`diff_diff/prep_dgp.py`). New public function exported from `diff_diff` and `diff_diff.prep` for panel DDD simulations. Cross-sectional `generate_ddd_data` remains available unchanged. Produces a balanced panel of `n_units × n_periods` with two unit-level binary dimensions (`group`, `partition`) and a derived `post = 1[period >= treatment_period]` indicator; columns: `unit, period, outcome, group, partition, post, treated, true_effect` (+ `x1, x2` when `add_covariates=True`). DDD-CPT identification holds because the `group * partition` interaction enters as a unit-level (time-invariant) term, leaving the triple-interaction `treatment_effect * group * partition * post` as the sole source of differential group × partition trend. Compatible with `TripleDifference(cluster="unit").fit(..., time="post")` (the cluster kwarg is required because `TripleDifference` is the repeated-cross-section `panel=FALSE` estimator and unclustered SE on panel-generated rows understates variance under within-unit serial correlation; the point estimate `att` is invariant to clustering — see the new `TripleDifference` REGISTRY note on panel-shaped input). Users get panel-realistic unit fixed effects and within-unit serial correlation while the binary 2×2×2 estimator surface is unchanged. **Stratified allocation:** the partition split is drawn stratified-by-group at the requested `partition_frac` so every `(group, partition)` cell receives at least one unit; a targeted `ValueError` is raised at fit-time when the rounded cell counts (`n_units`, `group_frac`, `partition_frac`) would leave any cell empty. This guarantees the 2x2x2 DDD surface is populated for any valid input — independent marginal sampling (the cross-sectional `generate_ddd_data` convention) could collapse cells when marginals are small (e.g., `n_units=4, group_frac=partition_frac=0.25`). Validates `1 <= treatment_period < n_periods`, `group_frac` and `partition_frac` strictly in `(0, 1)`, and `n_units >= 4`. Deterministic recovery (`noise_sd=0`) matches `treatment_effect` to ~1e-15 (covered by `tests/test_prep.py::TestGenerateDddPanelData`, 16 tests including infeasible-config rejection and smallest-feasible-config round-trip through `TripleDifference.fit`). `power.simulate_power` is NOT yet auto-routed to the panel DGP for `TripleDifference` (the existing `_ddd_dgp_kwargs` registry entry still ignores `n_periods` and the existing `_check_ddd_dgp_compat` warning still fires on non-default kwargs) — that wiring is tracked as a follow-up in TODO.md. diff --git a/TODO.md b/TODO.md index b1079ad3..667ff7bd 100644 --- a/TODO.md +++ b/TODO.md @@ -96,7 +96,7 @@ Deferred items from PR reviews that were not addressed before merge. | WooldridgeDiD: Stata `jwdid` golden value tests — add R/Stata reference script and `TestReferenceValues` class. | `tests/test_wooldridge.py` | #216 | Medium | | Thread `vcov_type` (classical / hc1 / hc2 / hc2_bm) through the 8 standalone estimators that expose `cluster=`: `CallawaySantAnna`, `SunAbraham`, `ImputationDiD`, `TwoStageDiD`, `TripleDifference`, `StackedDiD`, `WooldridgeDiD`, `EfficientDiD`. Phase 1a added `vcov_type` to the `DifferenceInDifferences` inheritance chain only. | multiple | Phase 1a | Medium | | Weighted one-way Bell-McCaffrey (`vcov_type="hc2_bm"` + `weights`, no cluster) currently raises `NotImplementedError`. `_compute_bm_dof_from_contrasts` builds its hat matrix from the unscaled design via `X (X'WX)^{-1} X' W`, but `solve_ols` solves the WLS problem by transforming to `X* = sqrt(w) X`, so the correct symmetric idempotent residual-maker is `M* = I - sqrt(W) X (X'WX)^{-1} X' sqrt(W)`. Rederive the Satterthwaite `(tr G)^2 / tr(G^2)` ratio on the transformed design and add weighted parity tests before lifting the guard. | `linalg.py::_compute_bm_dof_from_contrasts`, `linalg.py::_validate_vcov_args` | Phase 1a | Medium | -| HC2 / HC2 + Bell-McCaffrey on absorbed-FE fits — REMAINING sub-gates: `TwoWayFixedEffects` (`twfe.py:154` rejects unconditionally); `MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` (`estimators.py:1458` rejects). The DiD sub-gate (`DifferenceInDifferences(absorb=..., vcov_type in {"hc2","hc2_bm"})`) was lifted via auto-route to `fixed_effects=` internally; clubSandwich-parity at 1e-10 verified. The same auto-route pattern can apply to MPD-absorb; TWFE is its own class and may need different surgery (TWFE always within-transforms with no equivalent `fixed_effects=` path). Within-transformation preserves coefficients and residuals under FWL but not the hat matrix; HC1/CR1 are unaffected (no leverage term). | `twfe.py::fit`, `estimators.py::MultiPeriodDiD.fit` | follow-up | Medium | +| HC2 / HC2 + Bell-McCaffrey on absorbed-FE fits — REMAINING sub-gate: `TwoWayFixedEffects` (`twfe.py:154` rejects unconditionally). The DiD sub-gate and the MultiPeriodDiD sub-gate were both lifted via auto-route to `fixed_effects=` internally (DiD: PR #458, ~1e-10 vs clubSandwich; MPD: this release, ~1e-10 vs sandwich::vcovHC and clubSandwich::vcovCR). TWFE has no equivalent `fixed_effects=` code path (always within-transforms), so the same auto-route surgery is not directly applicable — lifting requires either building the full-dummy design inline or refactoring TWFE to delegate to DiD. Within-transformation preserves coefficients and residuals under FWL but not the hat matrix; HC1/CR1 are unaffected (no leverage term). | `twfe.py::fit` | follow-up | Medium | | Weighted CR2 Bell-McCaffrey cluster-robust (`vcov_type="hc2_bm"` + `cluster_ids` + `weights`) currently raises `NotImplementedError`. Weighted hat matrix and residual rebalancing need threading per clubSandwich WLS handling. | `linalg.py::_compute_cr2_bm` | Phase 1a | Medium | | Unify Rust local-method `estimate_model` solver path to `solve_wls_svd` (the same SVD helper used by the global-method since PR #348) for sub-1e-14 bootstrap SE parity. Current local-method bootstrap parity test (`tests/test_rust_backend.py::TestTROPRustEdgeCaseParity::test_bootstrap_seed_reproducibility_local`) passes at `atol=1e-5` — the residual ~1e-7 gap is roundoff between Rust's `estimate_model` matrix factorization and numpy's `lstsq`, which accumulates differently across per-replicate bootstrap fits. Main-fit ATT parity is regime-dependent (`atol=1e-14` for `lambda_nn=inf`, `atol=1e-10` for finite `lambda_nn` — see `test_local_method_main_fit_parity`); the bootstrap gap is a same-solver-path roundoff concern and not a user-visible correctness bug. | `rust/src/trop.rs::estimate_model`, `rust/src/linalg.rs::solve_wls_svd` | follow-up | Low | | Rust multiplier-bootstrap weight RNG (`generate_bootstrap_weights_batch` in `rust/src/bootstrap.rs:9-10, 57-75`) uses `Xoshiro256PlusPlus::seed_from_u64(seed + i)` per row for Rademacher/Mammen/Webb generation. If any Python caller (SDID / efficient-DiD multiplier bootstrap) has a numpy-canonical equivalent, the two backends likely diverge under the same seed. Audit Python callers (`diff_diff/sdid.py`, `diff_diff/efficient_did_bootstrap.py`, `diff_diff/bootstrap_utils.py::generate_bootstrap_weights_batch_numpy`) for parity-test gaps. Same fix shape as TROP RNG parity (PR #354): pre-generate weights in Python via numpy and pass them to Rust through PyO3. | `rust/src/bootstrap.rs`, `diff_diff/bootstrap_utils.py` | follow-up | Medium | diff --git a/benchmarks/R/generate_clubsandwich_golden.R b/benchmarks/R/generate_clubsandwich_golden.R index a4e748cb..cf11f957 100644 --- a/benchmarks/R/generate_clubsandwich_golden.R +++ b/benchmarks/R/generate_clubsandwich_golden.R @@ -122,6 +122,56 @@ output$absorbed_fe_did <- list( dof_cr2 = as.numeric(ct_did_cr2$df_Satt) ) +# --- Absorbed-FE MultiPeriodDiD event-study scenario (gate lift PR) ---------- +# Mirrors MPD(fixed_effects=["unit"]) destination of the absorb auto-route on +# MultiPeriodDiD. MPD parameterization: const + treated + period_f (non-ref) +# + treated:period_X (non-ref) + factor(unit). Build the interaction columns +# explicitly so the R fit's coefficient names match MPD's `treated:period_X`. + +make_mpd_panel <- function(n_total, units_per_cohort, n_periods, seed) { + set.seed(seed) + d <- expand.grid(unit = seq_len(n_total), period = seq_len(n_periods)) + d$cohort <- ((d$unit - 1L) %/% units_per_cohort) + 1L + n_cohorts <- n_total %/% units_per_cohort + # Last cohort is never-treated control; preceding cohorts ever-treated. + d$treated <- as.integer(d$cohort < n_cohorts) + d$y <- 1 + 0.5 * d$treated * (d$period >= 3) + + rnorm(nrow(d), sd = 0.5) + + 0.1 * d$unit + 0.2 * d$period + d +} + +d_mpd <- make_mpd_panel(n_total = 25, units_per_cohort = 5, n_periods = 5, + seed = 12345) +d_mpd$period_f <- relevel(factor(d_mpd$period), ref = "1") +# Explicit interaction columns to match MPD's parameterization exactly. +for (p in 2:5) { + d_mpd[[paste0("treated_period_", p)]] <- d_mpd$treated * (d_mpd$period == p) +} +fit_mpd <- lm(y ~ treated + period_f + + treated_period_2 + treated_period_3 + + treated_period_4 + treated_period_5 + + factor(unit), + data = d_mpd) +vcov_mpd_hc2 <- sandwich::vcovHC(fit_mpd, type = "HC2") +vcov_mpd_hc2_bm <- vcovCR(fit_mpd, cluster = seq_len(nrow(d_mpd)), type = "CR2") +ct_mpd_hc2_bm <- coef_test(fit_mpd, vcov = vcov_mpd_hc2_bm) +output$mpd_absorbed_fe_did <- list( + unit = d_mpd$unit, + period = d_mpd$period, + treated = d_mpd$treated, + y = d_mpd$y, + coef = as.numeric(coef(fit_mpd)), + coef_names = names(coef(fit_mpd)), + vcov_hc2 = as.numeric(vcov_mpd_hc2), + vcov_hc2_shape = dim(vcov_mpd_hc2), + vcov_hc2_bm = as.numeric(vcov_mpd_hc2_bm), + vcov_hc2_bm_shape = dim(vcov_mpd_hc2_bm), + dof_hc2_bm = as.numeric(ct_mpd_hc2_bm$df_Satt), + reference_period = 1L, + target_period = 4L +) + output$meta <- list( source = "clubSandwich", clubSandwich_version = as.character(packageVersion("clubSandwich")), diff --git a/benchmarks/data/clubsandwich_cr2_golden.json b/benchmarks/data/clubsandwich_cr2_golden.json index e7f453e7..3524bb56 100644 --- a/benchmarks/data/clubsandwich_cr2_golden.json +++ b/benchmarks/data/clubsandwich_cr2_golden.json @@ -49,11 +49,26 @@ "vcov_cr2_shape": [12, 12], "dof_cr2": [3.000000000000001, 6.000000000000003, 1.373392274582593, 1.945786883798724, 1.867171596521997, 6.000000000000002, 6.000000000000001, 6.000000000000006, 6.000000000000002, 3.77697841726619, 3.776978417266191, 3.776978417266191] }, + "mpd_absorbed_fe_did": { + "unit": [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25], + "period": [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5], + "treated": [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0], + "y": [1.592764408921928, 1.754733008754762, 1.445348342659473, 1.373251413268618, 2.002943727920197, 0.8910220161481357, 2.215049275534196, 1.861907947387392, 1.957920128028314, 1.740338998762936, 2.241876096823999, 3.308656021852109, 2.685313932128977, 2.860108228777479, 2.324734002748834, 3.208449919760292, 2.456821239378394, 2.834211205028724, 3.660356325834778, 3.349361849633647, 3.689810962277662, 4.127892541238432, 3.17783578538435, 2.823431297385157, 2.901145241651844, 2.402548759405412, 1.359176318152682, 2.010189900649211, 2.106061746325425, 1.818844511540937, 2.40593658927693, 3.198416773173767, 3.224595168703096, 3.116222819740233, 2.527135596407027, 2.74559413963628, 2.437956710631441, 1.868974878070683, 3.683866925436484, 2.912900524323904, 3.564255416795089, 1.909820969301479, 2.669867223923734, 3.768570270091454, 3.827225860165277, 4.230364701552045, 2.893450611055403, 3.983701626712241, 4.091593826717843, 3.246600583267788, 1.929806963068948, 3.27384633236602, 2.426795135178727, 2.675831420277619, 2.264511730503031, 2.838976847403072, 3.145585636437442, 3.311897664353043, 4.072532510240145, 1.926528010818313, 3.274795990307604, 2.628734259135228, 3.67665153775195, 4.294981421347705, 3.30656020280645, 2.783811346895301, 4.24406971632601, 4.696744236178137, 4.258427335331193, 3.452164159852588, 3.727307787664688, 3.407675313404228, 3.375323590521867, 5.165255981202357, 4.801352691433715, 2.871300425345763, 2.913129143576215, 2.194229754904459, 2.938124140376827, 3.310629203740546, 3.222691537103039, 3.52157177577572, 2.94781544365317, 4.438555458481485, 3.785610336518722, 4.333549592283361, 3.836021234507394, 3.446023309581331, 3.968261858355437, 4.212435032675837, 3.418049260493099, 3.572458744450104, 5.043473471012831, 4.004090314002409, 3.809683525525007, 4.243666050278521, 3.747478241239552, 5.178859908253208, 3.900101218087583, 3.952726653364625, 2.711962703760644, 2.121888334879613, 3.0112092638935, 2.237622372075537, 3.070542156430142, 2.831976000754527, 3.044196957710926, 4.078054821300628, 3.175983354352765, 3.660561768746984, 2.984913876559905, 3.037970653922923, 4.430621137194148, 4.559615859821168, 3.959623119315983, 3.847455096864254, 4.173923203334838, 4.614430312554218, 5.490001198742418, 4.465491346740392, 4.872431801018798, 4.860726008581316, 4.461075787681401, 5.165477559867579, 4.289380153330713], + "coef": [3.169677593270309, -1.650765824192107, 0.3451191042735751, 0.751359907257881, 0.8605432486572097, 1.385795096508473, 0.04453054634629026, 0.1845443111033822, 0.4405835249932665, -0.09865132402835917, -0.01712202455468425, -0.08412217264346931, -0.03549843363573817, 0.191817613926427, 0.1364439460365977, 0.7232874316258666, 0.7831775569789207, 1.050566202068045, 0.4263582901502527, 0.8144692870216864, 0.7481911239092748, 0.9198403068448744, 1.571690206647111, 1.041573924273659, 1.062727556061063, 0.9697421224576215, 1.670068637638985, 1.934612436699907, 1.479108696282838, 0.3144751959486066, -0.03079652150595019, 0.197118275100877, 0.390930912042367, "NA"], + "coef_names": ["(Intercept)", "treated", "period_f2", "period_f3", "period_f4", "period_f5", "treated_period_2", "treated_period_3", "treated_period_4", "treated_period_5", "factor(unit)2", "factor(unit)3", "factor(unit)4", "factor(unit)5", "factor(unit)6", "factor(unit)7", "factor(unit)8", "factor(unit)9", "factor(unit)10", "factor(unit)11", "factor(unit)12", "factor(unit)13", "factor(unit)14", "factor(unit)15", "factor(unit)16", "factor(unit)17", "factor(unit)18", "factor(unit)19", "factor(unit)20", "factor(unit)21", "factor(unit)22", "factor(unit)23", "factor(unit)24", "factor(unit)25"], + "vcov_hc2": [0.1252468846815752, -0.1252468846815742, -0.09222011170385074, -0.07148136166764696, -0.1034229382377807, -0.08762470380756525, 0.09222011170385076, 0.07148136166764682, 0.1034229382377812, 0.08762470380756457, -2.284057600544696e-16, -3.935203512346296e-16, -6.002658826905006e-16, -4.634078012180163e-16, -3.007987730240803e-16, -7.888467547673029e-16, -4.431244005441231e-16, -3.962322994631656e-16, -4.863238283566331e-16, -7.570546595385548e-17, -3.305529977021927e-16, -4.837559349940205e-16, -2.713067713127314e-16, -5.320678729997817e-16, -3.130048789217813e-17, -3.936037076214141e-16, -2.049057035700561e-18, -6.958144873580223e-16, -7.478355423816202e-17, -0.05754930861900105, -0.01486518641465552, -0.07271227289046649, -0.04065221212243775, -0.1252468846815749, 0.1600702480026458, 0.0922201117038505, 0.07148136166764697, 0.1034229382377804, 0.08762470380756529, -0.1009626150969271, -0.07979996026837496, -0.1147091843204617, -0.09888698692055442, -0.03001625204693945, -0.02799172562910784, -0.02789920728561605, -0.02684424770207596, -0.02147484156378036, -0.02836182667263606, -0.02738263616848014, -0.02594397138776956, -0.03020241429904486, -0.0291638252787008, -0.01736612784171559, -0.03045806555672743, -0.02719993606590916, -0.02677176698500658, -0.02567728622482477, -0.03166785092083725, -0.02946662496338875, -0.02947104356900903, -0.02743731070006746, 0.05754930861900082, 0.01486518641465526, 0.07271227289046613, 0.04065221212243749, -0.09222011170385082, 0.09222011170385014, 0.1476032727940036, 0.1044802964758183, 0.1044802964758181, 0.1044802964758185, -0.1476032727940036, -0.1044802964758185, -0.1044802964758185, -0.1044802964758182, 1.825219968405273e-16, 1.462585520207421e-16, 5.372073176933113e-16, 4.076058449478632e-16, 2.107659825342415e-16, 6.137700969857237e-16, 4.392286830569629e-16, 3.493674864778499e-16, 2.914532185864956e-16, 9.903602989212875e-17, 3.103937531218319e-16, 4.253395924152742e-16, 2.639920630913651e-16, 4.451729164062866e-16, 7.162833476071757e-17, 3.807707439908925e-16, 6.283152369859241e-17, 4.873479979642651e-16, 1.465924208672432e-16, 0.007422539042761633, -0.03932590186857911, 0.002396408782484242, -0.03179396981650389, -0.07148136166764676, 0.07148136166764625, 0.1044802964758181, 0.2636050474961812, 0.1044802964758178, 0.1044802964758181, -0.1044802964758183, -0.2636050474961817, -0.1044802964758183, -0.1044802964758179, 1.11671041157713e-16, 1.524916887083024e-16, 6.173116109981238e-16, 3.225352696959346e-16, 2.368919524433174e-16, 5.790259762575267e-16, 3.298435240137784e-16, 3.334633946079716e-16, 1.887553444506295e-16, 1.261592856987744e-16, 3.091094689723224e-16, 3.706319108794598e-16, 2.54017076929229e-16, 4.551475907584926e-16, 4.318088223881376e-17, 2.532614531178493e-16, 5.488103860621264e-17, 4.702548460770176e-16, 6.019235227728802e-17, -0.03042926664240113, -0.09050745478694976, -0.0116719668252837, -0.03238598578622085, -0.1034229382377807, 0.10342293823778, 0.104480296475818, 0.1044802964758178, 0.1830162159474743, 0.1044802964758179, -0.104480296475818, -0.1044802964758179, -0.1830162159474747, -0.1044802964758175, 1.724374339296564e-16, 2.241879626362903e-16, 4.984744794805693e-16, 5.17141657567021e-16, 2.533462853411204e-16, 6.600221810714608e-16, 5.075087593331421e-16, 4.866575752585108e-16, 2.643601885406811e-16, 1.933816647804064e-16, 3.686186955241591e-16, 4.602409849513467e-16, 3.406595122645797e-16, 4.873228081117761e-16, 2.07393783017759e-16, 4.527275736049338e-16, 1.617517176409348e-16, 5.456376986148697e-16, 2.395435747012338e-16, 0.002126003053302255, -0.05303516158320955, 0.05848688611096394, -0.01286451877124116, -0.08762470380756487, 0.08762470380756421, 0.1044802964758181, 0.1044802964758179, 0.1044802964758176, 0.1242750090127531, -0.1044802964758181, -0.104480296475818, -0.104480296475818, -0.1242750090127527, 3.053627642121331e-16, 2.383849420826669e-16, 4.731709302629929e-16, 4.293601217712507e-16, 2.297313746972919e-16, 6.270557425758743e-16, 3.899474797622955e-16, 3.726663291743884e-16, 3.114328736248038e-16, 1.75776692540923e-16, 3.043231177021744e-16, 4.396037648520271e-16, 2.445308819657429e-16, 4.5505930323003e-16, 9.851361811939864e-17, 3.522518564065572e-16, 1.470442570338681e-16, 4.722988773945257e-16, 1.284864041214332e-16, -0.002561545524894237, -0.05399436285422081, 0.003161223217932053, -0.03088327818008186, 0.09222011170385075, -0.1009626150969265, -0.1476032727940035, -0.1044802964758185, -0.1044802964758181, -0.1044802964758185, 0.1743420771174607, 0.1153699831323211, 0.1153699831323211, 0.1153699831323207, -0.0002532314627845022, -0.002686545541977227, -0.002584427713055136, -0.003219694942485645, -0.009111893876678934, 0.001074905424384323, -0.001952542159556722, -0.00773543912673435, -0.00320551410364791, -0.003242543795346812, -0.01678329714604433, 0.008172866437558114, -0.003292598114964238, -0.003852127726411447, -0.003682257904125903, 0.009189156536837419, 0.006071308178256183, -0.003628725133903056, -0.002221063097847072, -0.007422539042761571, 0.03932590186857921, -0.002396408782484143, 0.03179396981650397, 0.07148136166764697, -0.07979996026837448, -0.1044802964758184, -0.263605047496182, -0.104480296475818, -0.1044802964758185, 0.1153699831323211, 0.2903039978499587, 0.1153699831323211, 0.1153699831323207, 0.004802502984349069, -0.003515057792523033, -0.002848901956082749, -0.002757267926859517, -0.01043450785775177, -0.003550789151425393, -0.006032606453195888, -0.004243733847770334, 0.007883848577084939, -0.003664777011072829, -0.01348978782522216, -0.003154893225227057, -0.003305640565628235, -0.00381215532233186, -0.001640905900692437, 0.00526849594828588, -0.000908132038397042, -0.003892085577326768, -0.002125366173710049, 0.03042926664240121, 0.09050745478694998, 0.0116719668252838, 0.03238598578622096, 0.1034229382377808, -0.1147091843204612, -0.1044802964758182, -0.1044802964758182, -0.1830162159474746, -0.1044802964758181, 0.1153699831323208, 0.1153699831323208, 0.2065924369306855, 0.1153699831323202, -0.0006448294778879513, 0.003896570180875877, 0.0001878794390891855, 8.819134332471639e-05, -0.00671726729153591, 3.371787528246186e-05, 0.003139694810939287, -0.0006395861908108351, 0.003310274931136943, 0.006356745526257809, -0.01329334427135956, 0.0003734837006687619, 0.001697674634352693, 0.0009253163784994255, -0.002290833551714637, 0.00064528856962497, 0.002410256538923121, 0.00688668708421375, 0.001565268293682015, -0.002126003053302193, 0.05303516158320974, -0.05848688611096382, 0.0128645187712413, 0.08762470380756485, -0.0988869869205545, -0.1044802964758182, -0.104480296475818, -0.1044802964758176, -0.1242750090127532, 0.1153699831323208, 0.1153699831323207, 0.1153699831323206, 0.1481409994577317, 0.004968022767920427, 0.001054865876064719, 0.003532691235030964, -0.001098785386697391, -0.007570918578232068, 0.003042503791840368, 0.0005498392211161055, 0.001129820681064981, 0.001814666667549845, 0.005160906250564132, -0.01081172697189723, 0.005690075447537407, -0.0003085510473141659, -0.0006109938278217947, -0.005208366942445865, 0.002027518126337383, -0.001449103284940274, 0.006780546048962117, -0.001241080944889068, 0.002561545524894277, 0.05399436285422091, -0.003161223217931945, 0.03088327818008195, 4.062509823806105e-17, -0.03001625204693958, 1.009095288742695e-16, -7.805559022406642e-17, -1.642664089675595e-17, -4.585605491446902e-17, -0.0002532314627842821, 0.004802502984349415, -0.0006448294778879744, 0.004968022767920465, 0.1073749780756504, 0.02824175908462014, 0.02824175908461987, 0.02824175908461977, 0.02824175908462029, 0.02824175908462012, 0.02824175908461981, 0.0282417590846198, 0.02824175908462045, 0.02824175908462035, 0.02824175908462035, 0.02824175908462022, 0.02824175908462011, 0.02824175908461998, 0.02824175908462055, 0.02824175908462025, 0.02824175908462031, 0.02824175908462022, 0.02824175908462037, -1.568933258096708e-16, -1.8185054474974e-16, -1.367224173832291e-16, -1.977487319512773e-16, -1.237008914912177e-16, -0.02799172562910796, 2.03217263867653e-18, 9.773449082223245e-18, 1.661782213984008e-18, -7.514628508811563e-17, -0.002686545541976977, -0.003515057792522766, 0.003896570180875884, 0.001054865876064711, 0.02824175908462014, 0.05007686871951633, 0.02824175908461982, 0.02824175908461969, 0.02824175908462021, 0.02824175908462004, 0.02824175908461973, 0.02824175908461974, 0.02824175908462037, 0.02824175908462031, 0.02824175908462024, 0.02824175908462015, 0.02824175908462004, 0.02824175908461992, 0.02824175908462046, 0.02824175908462021, 0.02824175908462025, 0.02824175908462013, 0.02824175908462028, 8.774357855143607e-17, 1.588320175417508e-16, 1.877257432480739e-16, 1.527517369301608e-16, -3.501492444743265e-16, -0.02789920728561637, 3.603631463501666e-16, 4.418226782312757e-16, 2.688893004463125e-16, 2.16192107612452e-16, -0.002584427713054802, -0.002848901956082394, 0.000187879439089251, 0.003532691235030944, 0.02824175908462009, 0.02824175908462002, 0.05139075500105788, 0.02824175908461965, 0.02824175908462018, 0.02824175908462, 0.02824175908461967, 0.02824175908461971, 0.0282417590846203, 0.02824175908462028, 0.02824175908462026, 0.02824175908462008, 0.02824175908462, 0.02824175908461988, 0.02824175908462042, 0.02824175908462012, 0.0282417590846202, 0.0282417590846201, 0.02824175908462024, 6.844147437577848e-17, 3.11963423164566e-17, 1.740140215103348e-16, 8.330662439272845e-17, -4.218462774356293e-16, -0.02684424770207623, 4.101452586848311e-16, 3.202138783063954e-16, 4.323700603542044e-16, 2.552321981813559e-16, -0.003219694942485456, -0.002757267926859305, 8.81913433246668e-05, -0.001098785386697456, 0.0282417590846201, 0.02824175908462002, 0.02824175908461977, 0.04956881971904628, 0.02824175908462019, 0.02824175908462, 0.02824175908461967, 0.02824175908461971, 0.02824175908462032, 0.02824175908462027, 0.02824175908462028, 0.02824175908462011, 0.02824175908462002, 0.02824175908461988, 0.02824175908462043, 0.02824175908462015, 0.02824175908462022, 0.02824175908462008, 0.02824175908462025, 9.628920926447364e-17, 4.041737597790488e-17, 2.53971754623282e-16, 9.197106651785481e-17, -2.783184175201069e-16, -0.02147484156378001, 2.803771715868325e-16, 2.122493698169977e-16, 2.070660234188376e-16, 1.324855969168432e-16, -0.009111893876678946, -0.01043450785775166, -0.006717267291536132, -0.007570918578232353, 0.02824175908462015, 0.02824175908462006, 0.02824175908461981, 0.02824175908461971, 0.07344827063775429, 0.02824175908462005, 0.02824175908461974, 0.02824175908461976, 0.02824175908462039, 0.02824175908462031, 0.02824175908462029, 0.02824175908462015, 0.02824175908462005, 0.02824175908461993, 0.02824175908462048, 0.02824175908462018, 0.02824175908462025, 0.02824175908462015, 0.0282417590846203, 6.363992077551269e-17, 6.21115722968449e-17, 1.682862667248403e-16, 7.127490588439686e-17, -4.060374808001935e-16, -0.02836182667263635, 4.321063070305119e-16, 3.691327953578351e-16, 3.619059126999611e-16, 2.562318610875498e-16, 0.001074905424384656, -0.003550789151425036, 3.371787528255065e-05, 0.003042503791840417, 0.02824175908462019, 0.02824175908462011, 0.02824175908461986, 0.02824175908461974, 0.02824175908462026, 0.05764352940979834, 0.02824175908461977, 0.02824175908461977, 0.02824175908462043, 0.02824175908462035, 0.02824175908462034, 0.02824175908462017, 0.02824175908462009, 0.02824175908461997, 0.02824175908462049, 0.02824175908462025, 0.02824175908462028, 0.02824175908462018, 0.02824175908462034, 8.814493025944747e-17, 5.119889396839972e-17, 2.13482742667446e-16, 8.473524154431741e-17, -3.110011644032646e-16, -0.02738263616848048, 3.266594272128619e-16, 2.521674882308006e-16, 2.89249803317383e-16, 1.944495603227771e-16, -0.001952542159556433, -0.006032606453195618, 0.003139694810939374, 0.0005498392211160444, 0.02824175908462011, 0.02824175908462003, 0.02824175908461975, 0.02824175908461964, 0.02824175908462016, 0.02824175908461999, 0.08885512135058865, 0.02824175908461969, 0.02824175908462034, 0.02824175908462026, 0.02824175908462025, 0.02824175908462007, 0.02824175908461998, 0.02824175908461986, 0.02824175908462042, 0.02824175908462005, 0.02824175908462018, 0.02824175908462005, 0.02824175908462023, 5.751548328189774e-17, 3.095112939229784e-17, 1.842450176300661e-16, 5.230595998499911e-17, -2.64336494580313e-16, -0.02594397138776993, 2.754390020278343e-16, 2.198542128853604e-16, 2.554590487453482e-16, 1.572711145121441e-16, -0.007735439126734086, -0.00424373384777001, -0.0006395861908107437, 0.001129820681064953, 0.02824175908462018, 0.0282417590846201, 0.02824175908461987, 0.02824175908461973, 0.02824175908462027, 0.02824175908462007, 0.02824175908461975, 0.1075142860920249, 0.02824175908462041, 0.02824175908462037, 0.02824175908462032, 0.02824175908462018, 0.0282417590846201, 0.02824175908461998, 0.02824175908462055, 0.02824175908462025, 0.02824175908462033, 0.02824175908462017, 0.02824175908462036, 3.710112187995967e-17, 1.513455695633787e-17, 1.609467132428958e-16, 3.743324420457153e-17, -1.842859753366553e-16, -0.03020241429904472, 1.977184275963449e-16, 6.164400699264568e-17, 1.404485576874576e-16, 6.705281640719477e-17, -0.00320551410364776, 0.00788384857708512, 0.003310274931136773, 0.001814666667549691, 0.02824175908462028, 0.02824175908462019, 0.02824175908461992, 0.02824175908461981, 0.02824175908462035, 0.02824175908462017, 0.02824175908461984, 0.02824175908461986, 0.1053130769660629, 0.02824175908462044, 0.02824175908462044, 0.02824175908462026, 0.02824175908462018, 0.02824175908462006, 0.02824175908462061, 0.02824175908462028, 0.02824175908462037, 0.02824175908462025, 0.02824175908462041, 6.200940925719407e-17, 5.182383682223791e-17, 1.73801515045023e-16, 6.532186449252204e-17, -2.60926033753841e-16, -0.02916382527870022, 2.472229282741561e-16, 1.984208675314487e-16, 1.923647111955151e-16, 1.133051229726601e-16, -0.003242543795346912, -0.003664777011072833, 0.006356745526257541, 0.005160906250563809, 0.02824175908462019, 0.02824175908462013, 0.02824175908461989, 0.02824175908461975, 0.02824175908462027, 0.0282417590846201, 0.02824175908461981, 0.02824175908461982, 0.02824175908462043, 0.07569298551405899, 0.02824175908462035, 0.02824175908462021, 0.02824175908462012, 0.02824175908461999, 0.02824175908462053, 0.02824175908462024, 0.02824175908462032, 0.02824175908462019, 0.02824175908462035, 7.718275209344477e-17, 7.36054128079999e-17, 1.863894777644187e-16, 8.508720215426157e-17, -2.036010404841931e-16, -0.01736612784171539, 2.453806060231797e-16, 1.976630328239743e-16, 1.804731387212408e-16, 1.100227266156839e-16, -0.01678329714604417, -0.01348978782522193, -0.01329334427135961, -0.01081172697189739, 0.02824175908462026, 0.02824175908462008, 0.0282417590846199, 0.02824175908461982, 0.02824175908462032, 0.02824175908462014, 0.02824175908461983, 0.0282417590846198, 0.02824175908462049, 0.02824175908462042, 0.1230179216890002, 0.02824175908462025, 0.02824175908462015, 0.02824175908462004, 0.02824175908462056, 0.02824175908462028, 0.02824175908462036, 0.02824175908462024, 0.02824175908462041, 4.761893042235767e-18, 9.573464848132841e-18, 1.293240963609233e-16, 1.723893367905454e-17, -2.408378471350186e-16, -0.0304580655567275, 2.905348053553198e-16, 2.221608763369375e-16, 2.049881232719335e-16, 1.412735706223768e-16, 0.00817286643755836, -0.003154893225226808, 0.0003734837006687314, 0.005690075447537355, 0.02824175908462027, 0.02824175908462016, 0.02824175908461989, 0.0282417590846198, 0.02824175908462033, 0.02824175908462012, 0.02824175908461981, 0.02824175908461982, 0.02824175908462046, 0.02824175908462041, 0.0282417590846204, 0.1153724203808198, 0.02824175908462017, 0.02824175908462002, 0.02824175908462057, 0.02824175908462027, 0.02824175908462035, 0.02824175908462024, 0.02824175908462042, 3.309856527454699e-17, 1.348333436514057e-17, 1.388199898633107e-16, 3.003700937075362e-17, -2.764008005792364e-16, -0.0271999360659091, 2.813946650714977e-16, 2.132072724645399e-16, 2.007541632113417e-16, 1.233270081971404e-16, -0.003292598114964092, -0.003305640565628035, 0.001697674634352655, -0.0003085510473143397, 0.02824175908462023, 0.02824175908462014, 0.02824175908461988, 0.02824175908461978, 0.02824175908462029, 0.02824175908462011, 0.02824175908461981, 0.02824175908461983, 0.02824175908462045, 0.02824175908462039, 0.02824175908462036, 0.02824175908462025, 0.0481926248924397, 0.02824175908461998, 0.02824175908462054, 0.02824175908462028, 0.02824175908462033, 0.02824175908462022, 0.02824175908462035, 8.030193591075768e-17, 7.571350779865124e-17, 1.586948757721954e-16, 8.532020130739456e-17, -3.113652668660108e-16, -0.0267717669850069, 2.889755342513486e-16, 2.182748849471877e-16, 2.142562521723722e-16, 1.290860672309059e-16, -0.003852127726411073, -0.00381215532233142, 0.0009253163784995566, -0.00061099382782171, 0.02824175908462019, 0.02824175908462012, 0.02824175908461986, 0.02824175908461975, 0.02824175908462028, 0.0282417590846201, 0.02824175908461979, 0.02824175908461981, 0.02824175908462043, 0.02824175908462037, 0.02824175908462036, 0.02824175908462022, 0.0282417590846201, 0.04027437281127354, 0.02824175908462053, 0.02824175908462024, 0.02824175908462032, 0.0282417590846202, 0.02824175908462034, 1.008504891163701e-16, 1.014861516527886e-16, 2.045726784812991e-16, 1.117047809402094e-16, -2.001524565288355e-16, -0.02567728622482404, 2.238099905597254e-16, 1.277584695452614e-16, 1.558641084172538e-16, 8.869337996001103e-17, -0.003682257904125998, -0.001640905900692484, -0.002290833551714905, -0.005208366942446251, 0.02824175908462025, 0.02824175908462016, 0.02824175908461992, 0.02824175908461983, 0.02824175908462034, 0.02824175908462013, 0.02824175908461984, 0.02824175908461986, 0.02824175908462048, 0.02824175908462044, 0.02824175908462039, 0.02824175908462027, 0.02824175908462017, 0.02824175908462005, 0.1071512475831359, 0.02824175908462031, 0.02824175908462035, 0.02824175908462024, 0.02824175908462041, 3.815992711118558e-17, 3.431722700040473e-17, 1.421551424445783e-16, 3.060184341927672e-17, -2.133354712261508e-16, -0.03166785092083726, 2.674372631136414e-16, 1.080030515721389e-16, 1.711241410168758e-16, 8.673116543519833e-17, 0.009189156536837682, 0.005268495948286192, 0.0006452885696250306, 0.002027518126337286, 0.02824175908462025, 0.0282417590846202, 0.02824175908461991, 0.02824175908461982, 0.02824175908462033, 0.02824175908462016, 0.02824175908461977, 0.02824175908461987, 0.02824175908462048, 0.02824175908462044, 0.02824175908462042, 0.02824175908462027, 0.0282417590846202, 0.02824175908462005, 0.02824175908462059, 0.1238366984018747, 0.02824175908462036, 0.02824175908462025, 0.02824175908462042, 5.627645690933536e-17, 4.968281046652755e-17, 1.578711837460467e-16, 4.454051304505052e-17, -2.355320291964027e-16, -0.02946662496338833, 2.70191427654562e-16, 1.62518780961581e-16, 1.709488567917577e-16, 1.139020940865724e-16, 0.006071308178256169, -0.0009081320383970045, 0.002410256538922886, -0.001449103284940578, 0.02824175908462027, 0.0282417590846202, 0.02824175908461994, 0.02824175908461982, 0.02824175908462034, 0.02824175908462017, 0.02824175908461987, 0.02824175908461989, 0.0282417590846205, 0.02824175908462043, 0.02824175908462042, 0.02824175908462027, 0.02824175908462017, 0.02824175908462005, 0.02824175908462056, 0.02824175908462027, 0.1126668815484987, 0.02824175908462027, 0.02824175908462042, 5.041929260163529e-17, 4.771233846556506e-17, 1.622709761226476e-16, 3.790677375916962e-17, -2.406693249409194e-16, -0.02947104356900927, 2.460365525878227e-16, 1.740505860643879e-16, 1.846577137097897e-16, 1.197737347784893e-16, -0.003628725133902694, -0.003892085577326368, 0.006886687084213896, 0.006780546048962155, 0.0282417590846202, 0.02824175908462012, 0.02824175908461987, 0.02824175908461974, 0.02824175908462029, 0.02824175908462009, 0.02824175908461977, 0.02824175908461981, 0.02824175908462041, 0.02824175908462038, 0.02824175908462035, 0.02824175908462023, 0.02824175908462013, 0.02824175908462, 0.02824175908462052, 0.02824175908462024, 0.02824175908462032, 0.09089165877157075, 0.02824175908462035, 6.595294862800735e-17, 6.09264944455641e-17, 1.681125991314724e-16, 5.615253715536687e-17, -2.627747123245314e-16, -0.0274373107000669, 2.613004045303211e-16, 1.77197089140025e-16, 1.875685272812674e-16, 1.208207598325717e-16, -0.002221063097847094, -0.002125366173710064, 0.001565268293681834, -0.001241080944889416, 0.02824175908462016, 0.0282417590846201, 0.02824175908461983, 0.02824175908461972, 0.02824175908462026, 0.02824175908462008, 0.02824175908461976, 0.0282417590846198, 0.02824175908462041, 0.02824175908462035, 0.02824175908462034, 0.02824175908462022, 0.02824175908462009, 0.02824175908461996, 0.02824175908462051, 0.02824175908462023, 0.02824175908462031, 0.02824175908462017, 0.07144364413363367, 7.883835816956167e-17, 7.091745861896508e-17, 1.87128285189014e-16, 7.212101343351941e-17, -0.05754930861900074, 0.05754930861900012, 0.007422539042761498, -0.03042926664240109, 0.002126003053302038, -0.002561545524894167, -0.00742253904276145, 0.03042926664240138, -0.002126003053302154, 0.002561545524894633, -6.496084403481151e-18, 2.393626045131996e-16, 1.863845196146732e-16, 1.258495972996973e-16, 1.058143561076547e-16, 2.971431846485409e-16, 1.064383522959465e-16, 8.561048469239574e-17, 2.753598658153981e-16, -3.577615119704615e-17, 5.50403596323341e-17, 1.429697859982588e-16, 5.482499598118724e-17, 1.663565468467466e-16, -5.841991590385133e-17, 1.143948882936918e-16, -8.118155660051131e-17, 3.023017418328888e-16, -4.106296625579747e-17, 0.09171627403279942, 0.06223776263324736, 0.06223776263324723, 0.06223776263324731, -0.01486518641465545, 0.01486518641465491, -0.03932590186857914, -0.09050745478694958, -0.05303516158320959, -0.05399436285422056, 0.0393259018685792, 0.09050745478694999, 0.05303516158320957, 0.05399436285422106, -1.419204942937526e-16, 1.899918453316857e-16, 2.779201096382317e-17, -5.23688133005627e-17, -1.461887830836616e-17, 1.344088796674636e-16, -3.858679832115659e-17, -5.792550142768003e-17, 1.436048721839882e-16, -1.618110900986449e-16, -5.772068686957214e-17, -5.320293465201895e-18, -7.211000944461016e-17, 4.797630114334946e-17, -1.838320372110861e-16, -2.398223117212109e-17, -2.054579319330353e-16, 1.717090635653412e-16, -1.734398668668729e-16, 0.06223776263324744, 0.1707482888252618, 0.06223776263324736, 0.06223776263324739, -0.07271227289046665, 0.07271227289046607, 0.002396408782484369, -0.01167196682528345, 0.05848688611096402, 0.003161223217932436, -0.002396408782484293, 0.01167196682528372, -0.05848688611096414, -0.003161223217931935, -9.679236692724132e-17, 2.188855710380094e-16, 1.706096901577015e-16, 1.611855653448159e-16, 9.155581611962979e-17, 2.966927283665105e-16, 1.147070899166128e-16, 8.788665485887916e-17, 2.655825504067739e-16, -4.902702514222516e-17, 6.202994464321924e-17, 1.200163620329693e-16, 1.087135852893464e-17, 1.510628279718605e-16, -7.59941217669121e-17, 8.420614210739884e-17, -9.089929427595236e-17, 2.7889516825125e-16, -5.722904029682359e-17, 0.06223776263324747, 0.06223776263324753, 0.1654719168952946, 0.06223776263324751, -0.04065221212243773, 0.04065221212243724, -0.03179396981650387, -0.03238598578622064, -0.01286451877124114, -0.03088327818008154, 0.03179396981650388, 0.03238598578622094, 0.01286451877124105, 0.03088327818008201, -1.578186814952884e-16, 1.839115647200969e-16, 7.990229304009638e-17, -8.151227606108364e-19, -5.455544720813103e-18, 1.679452272433827e-16, -1.723196772845384e-17, -3.562681417944449e-17, 1.571028998542737e-16, -1.50329300752382e-16, -5.005521803864926e-17, 1.123338154041231e-17, -6.250331593586586e-17, 5.819493043077164e-17, -1.875474207922135e-16, -2.912452859359707e-17, -2.152634966394298e-16, 1.669351062751451e-16, -1.722363120523177e-16, 0.06223776263324743, 0.06223776263324744, 0.06223776263324739, 0.1638354639693444], + "vcov_hc2_shape": [33, 33], + "vcov_hc2_bm": [0.1252468846815764, -0.1252468846815755, -0.09222011170385148, -0.0714813616676476, -0.1034229382377814, -0.08762470380756532, 0.09222011170385136, 0.07148136166764699, 0.1034229382377816, 0.08762470380756511, -2.322915406406573e-16, -3.515539209037987e-16, -5.06896126319525e-16, -7.387595911980773e-16, -4.752937401108498e-16, -7.864042641131272e-16, -7.584277395376675e-16, -3.535997353175598e-16, -5.796380735763776e-16, -3.977811653976134e-16, -3.765162309216743e-16, -7.919066482061346e-17, -2.709875821931517e-16, -5.32923959035176e-16, -2.768879361519443e-17, -3.516372772905834e-16, -3.887536143004506e-16, -2.841437898270143e-16, -3.495984473024078e-16, -0.057549308619002, -0.01486518641465643, -0.07271227289046743, -0.04065221212243889, -0.1252468846815758, 0.1600702480026473, 0.0922201117038516, 0.07148136166764769, 0.1034229382377812, 0.08762470380756512, -0.1009626150969274, -0.07979996026837453, -0.1147091843204624, -0.09888698692055517, -0.03001625204693956, -0.02799172562910814, -0.02789920728561669, -0.02684424770207598, -0.02147484156378025, -0.02836182667263636, -0.02738263616848017, -0.02594397138777023, -0.03020241429904494, -0.02916382527870059, -0.01736612784171587, -0.0304580655567283, -0.02719993606590962, -0.026771766985007, -0.02567728622482476, -0.03166785092083775, -0.02946662496338853, -0.02947104356900972, -0.02743731070006744, 0.05754930861900182, 0.01486518641465621, 0.07271227289046714, 0.04065221212243869, -0.09222011170385129, 0.09222011170385089, 0.1476032727940043, 0.1044802964758188, 0.1044802964758185, 0.1044802964758186, -0.1476032727940041, -0.1044802964758186, -0.1044802964758187, -0.1044802964758185, 1.828550637479149e-16, 5.810684386550462e-17, 4.047577108555305e-16, 4.964401667908976e-16, 4.301139598158673e-16, 6.113276063315485e-16, 6.124234748984875e-16, 3.067349223322441e-16, 3.776620364486393e-16, 3.180830326506722e-16, 3.563569863413138e-16, 1.202503052482814e-16, 2.69001944489986e-16, 4.415881103431803e-16, 1.17754631986941e-16, 3.121589610690579e-16, 3.322965295629259e-16, 1.609424287244692e-16, 2.810751236188693e-16, 0.007422539042762107, -0.03932590186857871, 0.002396408782484675, -0.03179396981650336, -0.07148136166764715, 0.07148136166764693, 0.1044802964758186, 0.2636050474961819, 0.1044802964758182, 0.1044802964758182, -0.1044802964758188, -0.263605047496182, -0.1044802964758186, -0.1044802964758182, 1.262149627803026e-16, 1.531578225230774e-16, 5.29270925145349e-16, 4.548903341042753e-16, 4.64677624712094e-16, 6.227687634277577e-16, 5.456708800009092e-16, 2.908308304623658e-16, 3.069385854219776e-16, 2.990210106329113e-16, 3.515199885130036e-16, 9.751704682167143e-17, 2.590269583278502e-16, 4.098183989694802e-16, 5.910911319523279e-17, 2.219531638234201e-16, 2.817134803249401e-16, 2.006926956980296e-16, 1.911223413501136e-16, -0.03042926664240078, -0.09050745478694946, -0.01167196682528334, -0.0323859857862204, -0.1034229382377811, 0.1034229382377807, 0.1044802964758185, 0.1044802964758183, 0.1830162159474749, 0.1044802964758179, -0.1044802964758183, -0.104480296475818, -0.1830162159474751, -0.1044802964758177, 2.189557786614505e-16, 1.804451754660592e-16, 3.660248726427881e-16, 6.508289896049118e-16, 4.735824410424462e-16, 6.149471262716793e-16, 6.398473438684609e-16, 3.907343059308975e-16, 4.091887821030329e-16, 3.804542444297454e-16, 4.128055719042404e-16, 1.942315482511592e-16, 2.945991345304436e-16, 5.339200827617272e-16, 1.580409001262189e-16, 3.619113301905961e-16, 3.828109996316114e-16, 2.227848430538743e-16, 3.269528212087532e-16, 0.002126003053302559, -0.05303516158320918, 0.0584868861109644, -0.01286451877124064, -0.08762470380756525, 0.08762470380756482, 0.1044802964758187, 0.1044802964758184, 0.1044802964758181, 0.1242750090127531, -0.1044802964758186, -0.1044802964758181, -0.1044802964758183, -0.124275009012753, 2.562909065237013e-16, 1.526757245816047e-16, 3.859074005274557e-16, 5.621757552624631e-16, 4.525999732734123e-16, 6.292761886251246e-16, 5.622540931841198e-16, 3.229283376711814e-16, 3.914799537002778e-16, 3.083373216811668e-16, 3.962495841411376e-16, 1.420639942524852e-16, 2.896336923411398e-16, 4.589919213500694e-16, 1.068333519101835e-16, 3.091752030511012e-16, 3.6845048371974e-16, 2.173916709405898e-16, 2.633368682499664e-16, -0.002561545524893856, -0.05399436285422038, 0.003161223217932501, -0.0308832781800813, 0.09222011170385118, -0.1009626150969271, -0.1476032727940043, -0.104480296475819, -0.1044802964758185, -0.1044802964758185, 0.1743420771174612, 0.1153699831323212, 0.1153699831323213, 0.115369983132321, -0.0002532314627846589, -0.002686545541977302, -0.002584427713055054, -0.003219694942485881, -0.009111893876679402, 0.001074905424384164, -0.001952542159557032, -0.007735439126734414, -0.003205514103648201, -0.00324254379534725, -0.01678329714604454, 0.008172866437558317, -0.003292598114964349, -0.003852127726411557, -0.00368225790412616, 0.009189156536837368, 0.006071308178255707, -0.003628725133902871, -0.002221063097847328, -0.007422539042762049, 0.03932590186857878, -0.002396408782484586, 0.03179396981650343, 0.07148136166764729, -0.07979996026837495, -0.104480296475819, -0.2636050474961827, -0.1044802964758185, -0.1044802964758185, 0.1153699831323216, 0.2903039978499589, 0.1153699831323214, 0.1153699831323211, 0.004802502984348891, -0.003515057792523197, -0.002848901956082726, -0.002757267926859805, -0.01043450785775225, -0.003550789151425618, -0.006032606453196265, -0.004243733847770398, 0.007883848577084591, -0.003664777011073223, -0.01348978782522236, -0.003154893225226922, -0.003305640565628331, -0.003812155322331938, -0.001640905900692666, 0.005268495948285772, -0.0009081320383974741, -0.00389208557732664, -0.002125366173710306, 0.03042926664240087, 0.09050745478694965, 0.01167196682528344, 0.03238598578622051, 0.1034229382377812, -0.1147091843204618, -0.1044802964758188, -0.1044802964758187, -0.1830162159474751, -0.1044802964758181, 0.1153699831323211, 0.1153699831323208, 0.206592436930686, 0.1153699831323205, -0.0006448294778880832, 0.003896570180875846, 0.0001878794390893359, 8.819134332450619e-05, -0.006717267291536295, 3.371787528243382e-05, 0.003139694810939125, -0.0006395861908107355, 0.003310274931136681, 0.006356745526257506, -0.01329334427135966, 0.0003734837006689682, 0.001697674634352719, 0.0009253163784993505, -0.002290833551714732, 0.0006452885696249918, 0.002410256538922769, 0.006886687084214013, 0.001565268293681884, -0.002126003053302492, 0.05303516158320935, -0.05848688611096432, 0.01286451877124076, 0.08762470380756518, -0.098886986920555, -0.1044802964758187, -0.1044802964758186, -0.1044802964758181, -0.1242750090127531, 0.1153699831323211, 0.1153699831323207, 0.1153699831323208, 0.148140999457732, 0.004968022767920389, 0.001054865876064726, 0.003532691235031052, -0.001098785386697607, -0.007570918578232464, 0.003042503791840269, 0.0005498392211158767, 0.00112982068106501, 0.001814666667549623, 0.005160906250563876, -0.01081172697189742, 0.005690075447537677, -0.0003085510473142359, -0.0006109938278218251, -0.005208366942446001, 0.002027518126337359, -0.001449103284940621, 0.00678054604896233, -0.001241080944889259, 0.002561545524893901, 0.05399436285422048, -0.003161223217932386, 0.03088327818008138, -2.117778451429222e-16, -0.03001625204693982, -2.146359937795525e-16, -1.472703792891827e-16, 6.645045079439408e-17, 1.319544566258688e-16, -0.0002532314627846735, 0.004802502984348843, -0.0006448294778878684, 0.00496802276792055, 0.1073749780756507, 0.02824175908462047, 0.02824175908462044, 0.02824175908462015, 0.02824175908462054, 0.02824175908462051, 0.02824175908462023, 0.02824175908462049, 0.02824175908462075, 0.02824175908462062, 0.02824175908462079, 0.02824175908462074, 0.02824175908462061, 0.02824175908462047, 0.02824175908462066, 0.02824175908462071, 0.02824175908462062, 0.02824175908462058, 0.02824175908462071, -1.903734888960203e-16, -2.26432938082343e-16, -1.813048107158317e-16, -2.423311252838785e-16, -3.71420081487064e-16, -0.02799172562910818, -2.876659702230913e-16, -7.522732361428264e-17, 8.401845686234134e-17, 9.1908940901166e-17, -0.002686545541977382, -0.003515057792523303, 0.003896570180876001, 0.001054865876064812, 0.02824175908462035, 0.05007686871951669, 0.02824175908462039, 0.02824175908462007, 0.02824175908462045, 0.02824175908462044, 0.02824175908462015, 0.02824175908462041, 0.02824175908462065, 0.02824175908462057, 0.02824175908462067, 0.02824175908462067, 0.02824175908462055, 0.02824175908462039, 0.02824175908462058, 0.02824175908462068, 0.02824175908462056, 0.02824175908462049, 0.02824175908462062, 5.790633476463828e-17, 1.067903132624483e-16, 1.356840389687718e-16, 1.007100326508596e-16, -5.989092685557592e-16, -0.02789920728561662, 6.893028001242377e-17, 3.634138547434836e-16, 3.53674587961038e-16, 3.756145503074354e-16, -0.0025844277130552, -0.002848901956082928, 0.0001878794390893609, 0.003532691235031054, 0.0282417590846203, 0.02824175908462037, 0.05139075500105845, 0.02824175908462003, 0.02824175908462041, 0.0282417590846204, 0.0282417590846201, 0.02824175908462037, 0.02824175908462059, 0.02824175908462054, 0.0282417590846207, 0.02824175908462061, 0.0282417590846205, 0.02824175908462035, 0.02824175908462054, 0.02824175908462059, 0.0282417590846205, 0.02824175908462044, 0.02824175908462057, 3.999200936976312e-17, -1.945758318206378e-17, 1.233600960118149e-16, 3.265269889420976e-17, -6.748848343402332e-16, -0.02684424770207645, 1.17439810047133e-16, 2.391446937378669e-16, 5.119362696611285e-16, 4.216504553943028e-16, -0.003219694942485868, -0.00275726792685985, 8.819134332478375e-05, -0.001098785386697344, 0.02824175908462031, 0.02824175908462036, 0.02824175908462034, 0.04956881971904663, 0.02824175908462043, 0.02824175908462041, 0.0282417590846201, 0.02824175908462037, 0.02824175908462061, 0.02824175908462052, 0.02824175908462071, 0.02824175908462063, 0.02824175908462052, 0.02824175908462036, 0.02824175908462056, 0.02824175908462061, 0.02824175908462051, 0.02824175908462042, 0.02824175908462058, 6.662543782527495e-17, -1.145085595379905e-17, 2.021035226915786e-16, 4.010283458615287e-17, -5.296696847937794e-16, -0.02147484156378024, -1.26277879010141e-17, 1.357379001310545e-16, 2.862026176149137e-16, 2.948353854775379e-16, -0.009111893876679364, -0.01043450785775222, -0.006717267291536021, -0.007570918578232252, 0.02824175908462036, 0.0282417590846204, 0.02824175908462038, 0.02824175908462009, 0.0734482706377546, 0.02824175908462044, 0.02824175908462016, 0.02824175908462043, 0.02824175908462066, 0.02824175908462057, 0.02824175908462071, 0.02824175908462068, 0.02824175908462055, 0.02824175908462041, 0.02824175908462061, 0.02824175908462064, 0.02824175908462056, 0.02824175908462049, 0.02824175908462063, 3.25883705555319e-17, 8.855561584359066e-18, 1.150302560123553e-16, 1.801889517191293e-17, -6.506341685392828e-16, -0.0283618266726366, 1.392856619119872e-16, 2.879484143084804e-16, 4.439156426531239e-16, 4.212054189056598e-16, 0.001074905424384257, -0.003550789151425578, 3.371787528266405e-05, 0.003042503791840527, 0.02824175908462041, 0.02824175908462045, 0.02824175908462043, 0.02824175908462012, 0.0282417590846205, 0.05764352940979874, 0.02824175908462019, 0.02824175908462044, 0.02824175908462071, 0.02824175908462061, 0.02824175908462078, 0.0282417590846207, 0.02824175908462059, 0.02824175908462045, 0.0282417590846206, 0.02824175908462071, 0.02824175908462059, 0.02824175908462053, 0.02824175908462067, 6.247102281499497e-17, -2.230589091683622e-18, 1.600532596073635e-16, 3.130575848423578e-17, -5.705164740357531e-16, -0.02738263616848071, 3.765517374148636e-17, 1.595339322399984e-16, 3.875659339447278e-16, 3.476269985042443e-16, -0.001952542159556831, -0.006032606453196166, 0.003139694810939496, 0.0005498392211161699, 0.02824175908462032, 0.02824175908462036, 0.02824175908462033, 0.02824175908462001, 0.0282417590846204, 0.02824175908462039, 0.08885512135058918, 0.02824175908462036, 0.02824175908462061, 0.02824175908462051, 0.02824175908462069, 0.02824175908462059, 0.02824175908462048, 0.02824175908462033, 0.02824175908462053, 0.0282417590846205, 0.02824175908462048, 0.02824175908462042, 0.02824175908462056, 3.513755044180115e-17, -1.3631263940304e-17, 1.396626242974651e-16, 7.723566652399269e-18, -5.207293019560442e-16, -0.02594397138777017, -5.585523454047818e-18, 1.296492697609256e-16, 3.534282346774978e-16, 3.208568935494722e-16, -0.007735439126734498, -0.004243733847770535, -0.0006395861908106531, 0.001129820681065051, 0.0282417590846204, 0.02824175908462044, 0.02824175908462044, 0.02824175908462012, 0.02824175908462051, 0.02824175908462047, 0.02824175908462017, 0.1075142860920257, 0.02824175908462071, 0.02824175908462064, 0.02824175908462074, 0.02824175908462071, 0.0282417590846206, 0.02824175908462046, 0.02824175908462067, 0.02824175908462071, 0.02824175908462064, 0.02824175908462052, 0.0282417590846207, 1.472318903986308e-17, -2.944783637626466e-17, 1.163643199102944e-16, -7.14914912802857e-18, -4.467503148783049e-16, -0.03020241429904497, -9.11123535274339e-17, -1.13871737196187e-17, 2.195092577314604e-16, 2.400061022147979e-16, -0.003205514103648166, 0.007883848577084584, 0.00331027493113688, 0.001814666667549773, 0.02824175908462049, 0.02824175908462053, 0.0282417590846205, 0.0282417590846202, 0.0282417590846206, 0.02824175908462057, 0.02824175908462026, 0.02824175908462053, 0.1053130769660632, 0.02824175908462069, 0.02824175908462086, 0.02824175908462079, 0.02824175908462067, 0.02824175908462053, 0.02824175908462074, 0.02824175908462075, 0.02824175908462068, 0.02824175908462059, 0.02824175908462075, 2.852924617084607e-17, -3.860786756615791e-18, 1.181168914661702e-16, 9.637240913670486e-18, -5.077778620116985e-16, -0.02916382527870045, -4.334257632559929e-17, 1.153282906585189e-16, 2.780173604482291e-16, 2.80013404266747e-16, -0.003242543795347311, -0.003664777011073359, 0.006356745526257655, 0.005160906250563923, 0.02824175908462041, 0.02824175908462047, 0.02824175908462047, 0.02824175908462014, 0.02824175908462051, 0.02824175908462051, 0.02824175908462022, 0.02824175908462049, 0.02824175908462072, 0.07569298551405927, 0.02824175908462079, 0.02824175908462074, 0.02824175908462062, 0.02824175908462047, 0.02824175908462065, 0.02824175908462072, 0.02824175908462062, 0.02824175908462054, 0.02824175908462069, 4.300869961670616e-17, 1.722689983875569e-17, 1.300109647951754e-16, 2.870868918501938e-17, -4.702287163681865e-16, -0.01736612784171564, -5.316462656606823e-17, 1.121418430846767e-16, 2.692482902307135e-16, 2.718737821770353e-16, -0.01678329714604459, -0.01348978782522247, -0.0132933442713595, -0.0108117269718973, 0.02824175908462047, 0.02824175908462042, 0.02824175908462048, 0.0282417590846202, 0.02824175908462056, 0.02824175908462054, 0.02824175908462026, 0.02824175908462047, 0.02824175908462078, 0.02824175908462069, 0.1230179216890006, 0.02824175908462079, 0.02824175908462066, 0.02824175908462051, 0.02824175908462068, 0.02824175908462077, 0.02824175908462067, 0.0282417590846206, 0.02824175908462074, -2.247326553059559e-17, -3.986615421720414e-17, 7.988447729558744e-17, -3.220068538628081e-17, -5.293230388163199e-16, -0.03045806555672772, 2.744858317127916e-18, 1.307416267793191e-16, 3.183963481402765e-16, 2.913285065470862e-16, 0.008172866437557953, -0.003154893225227345, 0.0003734837006688714, 0.005690075447537473, 0.02824175908462048, 0.0282417590846205, 0.02824175908462045, 0.02824175908462017, 0.02824175908462057, 0.02824175908462052, 0.02824175908462024, 0.02824175908462048, 0.02824175908462075, 0.02824175908462066, 0.02824175908462084, 0.1153724203808203, 0.02824175908462067, 0.0282417590846205, 0.02824175908462068, 0.02824175908462075, 0.02824175908462065, 0.02824175908462059, 0.02824175908462074, 5.863406701715203e-18, -3.595628470019616e-17, 8.938037079797444e-17, -1.940260969458184e-17, -5.310155163920901e-16, -0.02719993606590934, -1.008156935314569e-17, 1.264284082051591e-16, 2.785571887352608e-16, 2.902954980126238e-16, -0.003292598114964506, -0.003305640565628574, 0.001697674634352774, -0.0003085510473142336, 0.02824175908462044, 0.02824175908462047, 0.02824175908462046, 0.02824175908462016, 0.02824175908462052, 0.02824175908462051, 0.02824175908462023, 0.0282417590846205, 0.02824175908462073, 0.02824175908462065, 0.02824175908462079, 0.02824175908462077, 0.04819262489244023, 0.02824175908462046, 0.02824175908462065, 0.02824175908462075, 0.02824175908462063, 0.02824175908462057, 0.02824175908462069, 4.942385803837416e-17, 2.263096943376313e-17, 1.056123374073076e-16, 3.22376629425076e-17, -5.612420191851024e-16, -0.02677176698500714, -4.11616141029794e-18, 1.340438957931476e-16, 2.987704891439765e-16, 2.935175239627729e-16, -0.00385212772641149, -0.003812155322331954, 0.0009253163784996736, -0.000610993827821603, 0.02824175908462041, 0.02824175908462046, 0.02824175908462044, 0.02824175908462013, 0.02824175908462052, 0.0282417590846205, 0.0282417590846202, 0.02824175908462048, 0.02824175908462072, 0.02824175908462063, 0.02824175908462079, 0.02824175908462074, 0.0282417590846206, 0.04027437281127403, 0.02824175908462066, 0.0282417590846207, 0.02824175908462061, 0.02824175908462056, 0.02824175908462067, 7.101324532957305e-17, 4.944444737348678e-17, 1.525309742019972e-16, 5.966307666090863e-17, -4.653754129731023e-16, -0.02567728622482432, -6.364995244221687e-17, 4.499589670858168e-17, 2.319588378931503e-16, 2.409695750856002e-16, -0.003682257904126399, -0.001640905900693015, -0.002290833551714763, -0.005208366942446135, 0.02824175908462046, 0.0282417590846205, 0.02824175908462049, 0.0282417590846202, 0.02824175908462058, 0.02824175908462052, 0.02824175908462027, 0.02824175908462054, 0.02824175908462076, 0.02824175908462069, 0.02824175908462083, 0.02824175908462079, 0.02824175908462067, 0.02824175908462052, 0.1071512475831361, 0.02824175908462078, 0.02824175908462065, 0.02824175908462059, 0.02824175908462075, 9.536989757572696e-18, -2.76124010919652e-17, 8.02255143522092e-17, -3.132778467309175e-17, -4.64871052869145e-16, -0.03166785092083749, -3.371005468957221e-17, 3.653312198825358e-17, 2.406438619430057e-16, 2.572558383764341e-16, 0.009189156536837264, 0.005268495948285632, 0.0006452885696251377, 0.00202751812633739, 0.02824175908462048, 0.02824175908462054, 0.02824175908462048, 0.02824175908462021, 0.02824175908462058, 0.02824175908462056, 0.02824175908462019, 0.02824175908462055, 0.02824175908462077, 0.02824175908462071, 0.02824175908462085, 0.02824175908462079, 0.0282417590846207, 0.02824175908462053, 0.0282417590846207, 0.1238366984018752, 0.02824175908462065, 0.0282417590846206, 0.02824175908462076, 2.279629382298697e-17, -6.001813112325813e-18, 1.021865601671939e-16, -1.11441105338014e-17, -4.870676108393981e-16, -0.02946662496338853, -4.344589917568514e-17, 7.855884235066259e-17, 2.585097018680466e-16, 2.705489792199935e-16, 0.00607130817825572, -0.0009081320383975568, 0.002410256538923007, -0.00144910328494046, 0.02824175908462049, 0.02824175908462054, 0.02824175908462052, 0.02824175908462021, 0.02824175908462059, 0.02824175908462056, 0.0282417590846203, 0.02824175908462056, 0.02824175908462079, 0.02824175908462069, 0.02824175908462086, 0.02824175908462081, 0.02824175908462068, 0.02824175908462054, 0.02824175908462069, 0.02824175908462074, 0.112666881548499, 0.02824175908462061, 0.02824175908462075, 2.179635524802237e-17, -1.421728962680465e-17, 1.003413480302783e-16, -2.402285433319918e-17, -5.100725583864753e-16, -0.02947104356900953, -4.973312243986257e-17, 9.68560690097791e-17, 2.734328652192619e-16, 2.760736752167149e-16, -0.003628725133903104, -0.003892085577326923, 0.006886687084214018, 0.00678054604896225, 0.02824175908462042, 0.02824175908462046, 0.02824175908462044, 0.02824175908462012, 0.02824175908462053, 0.02824175908462049, 0.02824175908462019, 0.02824175908462048, 0.0282417590846207, 0.02824175908462064, 0.02824175908462078, 0.02824175908462076, 0.02824175908462064, 0.02824175908462048, 0.02824175908462065, 0.02824175908462071, 0.02824175908462062, 0.09089165877157111, 0.02824175908462069, 3.351361962724547e-17, 6.28270495229714e-18, 1.134688096382055e-16, 1.508747662100891e-18, -5.245993725844136e-16, -0.02743731070006713, -2.261894075209746e-17, 8.964844133817987e-17, 2.743162207281915e-16, 2.786168992688272e-16, -0.002221063097847483, -0.002125366173710592, 0.001565268293681952, -0.001241080944889302, 0.02824175908462037, 0.02824175908462043, 0.02824175908462041, 0.02824175908462011, 0.0282417590846205, 0.02824175908462048, 0.02824175908462018, 0.02824175908462047, 0.02824175908462069, 0.02824175908462061, 0.02824175908462077, 0.02824175908462074, 0.02824175908462058, 0.02824175908462044, 0.02824175908462064, 0.0282417590846207, 0.0282417590846206, 0.02824175908462051, 0.07144364413363401, 5.195014429192529e-17, 1.627366912569812e-17, 1.324844956957471e-16, 1.747722394025347e-17, -0.05754930861900152, 0.05754930861900098, 0.007422539042761807, -0.0304292666424009, 0.002126003053302426, -0.002561545524894131, -0.007422539042761708, 0.03042926664240132, -0.002126003053302353, 0.002561545524894295, -5.607905983781046e-18, 2.392515822107379e-16, 1.862734973122125e-16, 3.042624373569612e-16, 1.085899136692178e-16, 2.971570624363495e-16, 2.887369729393986e-16, 8.605457390224788e-17, 2.791346240991241e-16, 1.40527265113429e-16, 5.50403596323351e-17, -3.510998715161472e-17, 5.426988446887624e-17, 1.685769928959984e-16, -6.286080800235194e-17, 1.157271559232434e-16, 1.000068410183149e-16, 1.202251657943639e-16, 1.372388514990036e-16, 0.09171627403280007, 0.06223776263324798, 0.06223776263324783, 0.06223776263324803, -0.0148651864146562, 0.01486518641465574, -0.03932590186857886, -0.09050745478694941, -0.05303516158320923, -0.05399436285422053, 0.03932590186857896, 0.09050745478694998, 0.05303516158320941, 0.05399436285422073, -1.410323158740528e-16, 1.89880823029224e-16, 2.768098866136224e-17, 1.251558483370008e-16, -1.53960344256038e-17, 1.344227574552717e-16, 1.366063949646939e-16, -5.748141221782824e-17, 1.473796304677137e-16, 1.44923262118301e-17, -5.772068686957143e-17, -1.83400066615076e-16, -7.266512095692159e-17, 4.664403351380041e-17, -1.882729293095869e-16, -2.264996354256992e-17, -3.137496167181043e-17, -3.262085115582882e-18, 4.861950887927753e-18, 0.06223776263324811, 0.1707482888252625, 0.06223776263324798, 0.06223776263324812, -0.07271227289046744, 0.07271227289046694, 0.002396408782484629, -0.01167196682528321, 0.05848688611096449, 0.003161223217932482, -0.002396408782484543, 0.01167196682528365, -0.05848688611096442, -0.003161223217932279, -9.590418850754165e-17, 2.187745487355473e-16, 1.704986678552404e-16, 3.387102269823788e-16, 9.07786600023921e-17, 2.967066061543189e-16, 2.899002832024631e-16, 8.833074406873062e-17, 2.693573086904996e-16, 1.272763911682496e-16, 6.202994464321979e-17, -5.806341111690503e-17, 1.031624701662306e-17, 1.497305603423111e-16, -8.043501386541281e-17, 8.553840973694979e-17, 8.31836759852726e-17, 1.039240195703259e-16, 1.210727774579771e-16, 0.06223776263324814, 0.06223776263324816, 0.1654719168952953, 0.06223776263324824, -0.04065221212243855, 0.04065221212243812, -0.03179396981650356, -0.03238598578622041, -0.01286451877124073, -0.03088327818008149, 0.03179396981650363, 0.03238598578622086, 0.01286451877124084, 0.03088327818008167, -1.569305030755884e-16, 1.838005424176349e-16, 7.97912707376353e-17, 1.767095388769523e-16, -6.23270083805076e-18, 1.67959105031191e-16, 1.579612255573966e-16, -3.518272496959301e-17, 1.608776581379992e-16, 2.597411555809279e-17, -5.005521803864867e-17, -1.66846391609462e-16, -6.305842744817744e-17, 5.686266280122226e-17, -1.919883128907142e-16, -2.779226096404592e-17, -4.118052637820484e-17, -8.036042405779079e-18, 6.065505702483043e-18, 0.0622377626332481, 0.06223776263324808, 0.06223776263324801, 0.1638354639693452], + "vcov_hc2_bm_shape": [33, 33], + "dof_hc2_bm": [4.853932584269628, 9.64773087908587, 7.529411764705872, 7.529411764705888, 7.529411764705889, 7.529411764705892, 11.611917494270315, 11.611917494270454, 11.611917494270452, 11.611917494270408, 7.977900552486196, 7.977900552486197, 7.977900552486195, 7.977900552486195, 7.977900552486197, 7.977900552486192, 7.977900552486195, 7.977900552486199, 7.977900552486192, 7.977900552486197, 7.977900552486199, 7.977900552486194, 7.977900552486194, 7.977900552486194, 7.977900552486196, 7.977900552486194, 7.977900552486197, 7.977900552486196, 7.977900552486194, 7.529411764705895, 7.529411764705896, 7.529411764705891, 7.529411764705896], + "reference_period": 1, + "target_period": 4 + }, "meta": { "source": "clubSandwich", "clubSandwich_version": "0.7.0", "R_version": "R version 4.5.2 (2025-10-31)", - "generated_at": "2026-05-16 21:17:21 UTC", + "generated_at": "2026-05-17 10:36:19 UTC", "note": "CR2 Bell-McCaffrey cluster-robust parity target for diff_diff._compute_cr2_bm" } } diff --git a/diff_diff/estimators.py b/diff_diff/estimators.py index 84e712cf..d98eab32 100644 --- a/diff_diff/estimators.py +++ b/diff_diff/estimators.py @@ -1449,16 +1449,9 @@ def fit( # type: ignore[override] n_treated_raw = int(np.sum(data[treatment].values.astype(float))) n_control_raw = len(data) - n_treated_raw - # Reject multi-absorb with survey weights (single-pass demeaning is - # not the correct weighted FWL projection for N > 1 dimensions) - if absorb and len(absorb) > 1 and survey_weights is not None: - raise ValueError( - f"Multiple absorbed fixed effects (absorb={absorb}) with survey " - "weights is not supported. Single-pass sequential demeaning is not " - "the correct weighted FWL projection for multiple absorbed dimensions. " - "Use absorb with a single variable, or use fixed_effects= instead." - ) - + # Mutual-exclusion check runs ABOVE the auto-route below so that the + # `absorb=..., fixed_effects=...` combination still rejects rather + # than being silently merged. if absorb and fixed_effects: raise ValueError( "Cannot use both absorb and fixed_effects. " @@ -1468,19 +1461,50 @@ def fit( # type: ignore[override] "or fixed_effects alone (for low-dimensional FE)." ) - # Reject HC2 / HC2 + Bell-McCaffrey on absorbed-FE fits (see the - # matching guard in DifferenceInDifferences.fit / twfe.py for the - # methodology reasoning: HC2/CR2 leverage corrections depend on the - # full FE hat matrix, not the residualized design from within- - # transformation). Tracked in TODO.md. + # Auto-route absorb → fixed_effects when vcov_type needs the FULL FE + # hat matrix. Mirrors the identical pattern in + # DifferenceInDifferences.fit (PR #458). HC2 leverage and CR2 + # Bell-McCaffrey DOF both depend on the full-design hat; FWL + # preserves coefficients and residuals but not the hat matrix, so + # the demeaned design's leverage is wrong for these vcov families. + # Building the full-dummy design and routing through the existing + # fixed_effects= branch produces the algebraically correct vcov. + # Empirically matches `lm() + sandwich::vcovHC` and + # `lm() + clubSandwich::vcovCR` (singleton-cluster trick for one-way + # HC2-BM; PT2018 §3.3 unweighted CR2 algebra) at ~1e-15. + # Conley vcov is unaffected: the absorb+Conley path computes the + # panel sandwich on demeaned scores, which is FWL-correct because + # Conley's meat uses only residuals (no leverage term). + # HC1/CR1 paths remain on the demeaned design (no leverage term). + # + # Survey-replicate scope: this also short-circuits the absorb-refit + # replicate-variance branch below (search "compute_replicate_refit_variance"). + # Correct: with a fixed full-dummy design, replicate variance doesn't + # need per-replicate refit — the standard compute_replicate_vcov + # path applies directly because the design matrix does not depend + # on the replicate weights. + # + # Placement: this auto-route runs BEFORE the multi-absorb + + # survey-weights guard because that guard's rationale ("single-pass + # demeaning is not the correct weighted FWL projection for N > 1 + # dimensions") doesn't apply when we're about to swap absorb for + # fixed_effects: the fixed_effects= path builds the full-dummy + # design and solves WLS directly, with no within-transform step. if absorb and self.vcov_type in ("hc2", "hc2_bm"): - raise NotImplementedError( - f"MultiPeriodDiD(absorb=..., vcov_type={self.vcov_type!r}) " - "is not yet supported: absorbed fixed effects are handled " - "by demeaning, and the HC2 / CR2 Bell-McCaffrey leverage " - "corrections depend on the full FE hat matrix, not the " - "residualized one. Use vcov_type='hc1' with absorb=, or " - "switch to fixed_effects= dummies for a full-dummy design." + fixed_effects = list(fixed_effects or []) + list(absorb) + absorb = None + n_absorbed_effects = 0 + + # Reject multi-absorb with survey weights (single-pass demeaning is + # not the correct weighted FWL projection for N > 1 dimensions). + # Only fires when absorb is still set — i.e., the auto-route above + # didn't consume it. + if absorb and len(absorb) > 1 and survey_weights is not None: + raise ValueError( + f"Multiple absorbed fixed effects (absorb={absorb}) with survey " + "weights is not supported. Single-pass sequential demeaning is not " + "the correct weighted FWL projection for multiple absorbed dimensions. " + "Use absorb with a single variable, or use fixed_effects= instead." ) # MultiPeriodDiD is intrinsically a multi-period panel estimator; diff --git a/docs/methodology/REGISTRY.md b/docs/methodology/REGISTRY.md index 478716b8..5b7aeefe 100644 --- a/docs/methodology/REGISTRY.md +++ b/docs/methodology/REGISTRY.md @@ -2550,8 +2550,8 @@ Shipped in `diff_diff/had_pretests.py` as `stute_joint_pretest()` (residuals-in - [x] Phase 1a: HC2 + Bell-McCaffrey DOF correction in `diff_diff/linalg.py` via `vcov_type="hc2_bm"` enum (both one-way and CR2 cluster-robust with Imbens-Kolesar / Pustejovsky-Tipton Satterthwaite DOF). Weighted cluster CR2 raises `NotImplementedError` and is tracked as Phase 2+ in `TODO.md`. - **Note (scope limitation on absorbed FE):** HC2 and HC2 + Bell-McCaffrey on within-transformed designs still depend on the FULL FE hat matrix because FWL preserves coefficients and residuals but NOT the hat matrix: `h_ii = x_i' (X'X)^{-1} x_i` on the reduced design is not the diagonal of the full FE projection, and CR2's block adjustment `A_g = (I - H_gg)^{-1/2}` likewise depends on the full cluster-block hat matrix. The status across the three estimators that previously rejected this combination: - **`DifferenceInDifferences(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** When the user pairs `absorb=` with HC2 / HC2-BM, `DiD.fit()` internally promotes the absorb columns to `fixed_effects=` so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection. Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=..., type="CR2")` (singleton-cluster CR2 trick for one-way HC2-BM Satterthwaite DOF; PT2018 §3.3 unweighted CR2 algebra). **User-visible surface change**: under the auto-route, the entire `DiDResults` (coefficients, vcov, residuals, fitted_values, r_squared) reflect the full-dummy fit rather than the within-transformed fit — the FE-dummy entries are included in `result.coefficients` / `result.vcov`, `r_squared` is computed on the un-demeaned outcome, and `residuals` / `fitted_values` are on the original scale. `result.att` is unaffected (FWL-equivalent). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). **Survey-design scope**: when `survey_design=` is supplied, the existing survey variance path (Taylor-series linearization / replicate weights) takes precedence over the analytical HC2/HC2-BM sandwich; the auto-route only changes the FE handling (removing the prior reject) and does not redirect to the analytical small-sample sandwich on survey fits. - - **`TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` — still rejects.** TWFE is a standalone class with no `fixed_effects=` equivalent path, so the same auto-route surgery would require building the full-dummy design inline. Tracked as a follow-up in `TODO.md`. - - **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` — still rejects.** Same algebraic situation as DiD, but the `MultiPeriodDiD.fit()` body has many additional branches (cohort, event-study, survey, etc.) that complicate the auto-route surgery. Tracked as a follow-up in `TODO.md`. + - **`TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` — still rejects.** TWFE is a standalone class with no `fixed_effects=` equivalent path, so the same auto-route surgery used for DiD-absorb and MPD-absorb is not directly applicable; lifting requires building the full-dummy design inline or refactoring TWFE to delegate to DiD. Tracked as a follow-up in `TODO.md`. + - **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** Same auto-route pattern as `DifferenceInDifferences`: `MultiPeriodDiD.fit()` internally promotes the absorb columns to `fixed_effects=` for HC2 / HC2-BM callers, so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture; the parity target is a per-period interaction `treated:period_X` (the `treated` main-effect coefficient becomes NaN under rank deficiency because MPD requires a time-invariant ever-treated indicator that is collinear with unit dummies — expected behavior). Same `MultiPeriodDiDResults` surface change as DiD: `vcov`, `residuals`, `fitted_values`, `r_squared`, and `coefficients` reflect the full-dummy fit, with `period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` invariant by FWL. HC1/CR1 paths on `absorb=` are unchanged (no leverage term). Same survey-design scope as DiD: replicate-weight variance routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design rather than the per-replicate refit branch (which targets the demeaning path); since the auto-routed design does not depend on replicate weights, no refit is needed. - Workarounds for the still-rejecting paths: use `vcov_type="hc1"` (HC1/CR1 have no leverage term and survive FWL), or switch to `fixed_effects=` dummies so the hat matrix is computed on the full design. - [x] Phase 1a: `vcov_type` enum threaded through `DifferenceInDifferences` (`MultiPeriodDiD`, `TwoWayFixedEffects` inherit); `robust=True` <=> `vcov_type="hc1"`, `robust=False` <=> `vcov_type="classical"`. Conflict detection at `__init__`. Results summary prints the variance-family label. - **Note (deviation from the fully-symmetric enum):** `MultiPeriodDiD(cluster=..., vcov_type="hc2_bm")` is intentionally **not supported** and raises `NotImplementedError`. The scalar-coefficient `DifferenceInDifferences` path handles the cluster + CR2 Bell-McCaffrey combination (`_compute_cr2_bm` returns a per-coefficient Satterthwaite DOF that is valid for the single-ATT contrast), but `MultiPeriodDiD` also reports a post-period-average ATT constructed as a *contrast* of the event-study coefficients. The cluster-aware CR2 BM DOF for that contrast (i.e., the Pustejovsky-Tipton 2018 per-cluster adjustment matrices applied to an arbitrary aggregation contrast) is not yet implemented. Pairing CR2 cluster-robust SEs with the one-way Imbens-Kolesar (2016) contrast DOF would be a broken hybrid, so the combination fails fast with a clear workaround message (drop the cluster for one-way HC2+BM, or use `vcov_type="hc1"` with cluster for CR1 Liang-Zeger). Tracked in `TODO.md` under Methodology/Correctness. Applies only to `MultiPeriodDiD`; `DifferenceInDifferences(cluster=..., vcov_type="hc2_bm")` works. diff --git a/tests/test_estimators_vcov_type.py b/tests/test_estimators_vcov_type.py index 2f346049..dcc5fb48 100644 --- a/tests/test_estimators_vcov_type.py +++ b/tests/test_estimators_vcov_type.py @@ -744,33 +744,6 @@ def test_did_fixed_effects_dummies_still_accept_hc2_and_hc2_bm(self): assert np.isfinite(res.att) assert np.isfinite(res.se) - def test_multi_period_absorb_rejects_hc2_and_hc2_bm(self): - """MultiPeriodDiD with absorb= rejects HC2/HC2+BM for the same - methodology reason as the base class.""" - rng = np.random.default_rng(20260420) - n_units, n_time = 30, 4 - rows = [] - for i in range(n_units): - treated = int(i >= n_units // 2) - for t in range(n_time): - y = rng.normal(0.0, 1.0) + 0.3 * treated + 0.5 * treated * (t >= 2) - rows.append({"unit": i, "time": t, "treated": treated, "y": y}) - data = pd.DataFrame(rows) - - for bad in ("hc2", "hc2_bm"): - with pytest.raises( - NotImplementedError, - match="MultiPeriodDiD.*absorb.*not yet supported", - ): - MultiPeriodDiD(vcov_type=bad).fit( - data, - outcome="y", - treatment="treated", - time="time", - absorb=["unit"], - unit="unit", - ) - def test_summary_suppresses_variance_line_under_wild_bootstrap(self): """When inference_method='wild_bootstrap', the Variance label is omitted. @@ -1217,3 +1190,206 @@ def test_absorb_hc2_bm_df_sensitive_inference(self): f"HC2 CI width ({width_hc2:.6f}) — BM Satterthwaite DOF is " "smaller than n-k, so the critical value is larger." ) + + +class TestMPDAbsorbedFERParity: + """R-parity for `MultiPeriodDiD(absorb=..., vcov_type in {hc2, hc2_bm})`. + + Mirrors `TestDiDAbsorbedFERParity`. The auto-route promotes `absorb=` to + `fixed_effects=` internally; MPD's existing `fixed_effects=` code path + builds the full-dummy design that R's `lm()` produces. + + Collinearity note: MPD's `treated` is a time-invariant ever-treated + indicator, so it is perfectly collinear with the sum of treated-cohort + unit dummies post-auto-route. `solve_ols` resolves this by dropping + one column (typically a unit dummy from the never-treated cohort); the + `treated` main-effect coefficient absorbs the dropped column's effect. + Tests pin parity on a per-period interaction (`treated:period_4`) where + no collinearity exists, exposed as `result.period_effects[4]`. + """ + + def _load_golden(self): + import json + from pathlib import Path + + golden_path = ( + Path(__file__).parent.parent / "benchmarks" / "data" / "clubsandwich_cr2_golden.json" + ) + if not golden_path.exists(): + pytest.skip( + "Golden JSON not present; run `Rscript " + "benchmarks/R/generate_clubsandwich_golden.R` to generate." + ) + with open(golden_path) as f: + golden = json.load(f) + if "mpd_absorbed_fe_did" not in golden: + pytest.skip( + "Golden JSON does not yet include `mpd_absorbed_fe_did` scenario; " + "regenerate via the R script." + ) + return golden["mpd_absorbed_fe_did"] + + def _make_data(self, d): + return pd.DataFrame( + { + "unit": d["unit"], + "period": d["period"], + "treated": d["treated"], + "y": d["y"], + } + ) + + def _fit_absorb(self, d, vcov_type, absorb_cols=("unit",)): + return MultiPeriodDiD(vcov_type=vcov_type).fit( + self._make_data(d), + outcome="y", + treatment="treated", + time="period", + absorb=list(absorb_cols), + reference_period=int(d["reference_period"]), + unit="unit", + ) + + def _fit_fixed_effects(self, d, vcov_type, fe_cols=("unit",)): + return MultiPeriodDiD(vcov_type=vcov_type).fit( + self._make_data(d), + outcome="y", + treatment="treated", + time="period", + fixed_effects=list(fe_cols), + reference_period=int(d["reference_period"]), + unit="unit", + ) + + def test_absorb_hc2_matches_fixed_effects_dummies_single_absorb(self): + """`absorb=["unit"]` + hc2 produces the same per-period SE as + `fixed_effects=["unit"]` + hc2 (auto-route invariant).""" + d = self._load_golden() + target_period = int(d["target_period"]) + res_a = self._fit_absorb(d, "hc2", absorb_cols=("unit",)) + res_f = self._fit_fixed_effects(d, "hc2", fe_cols=("unit",)) + pe_a = res_a.period_effects[target_period] + pe_f = res_f.period_effects[target_period] + np.testing.assert_allclose(pe_a.effect, pe_f.effect, atol=1e-12) + np.testing.assert_allclose(pe_a.se, pe_f.se, atol=1e-12) + + def test_absorb_hc2_matches_fixed_effects_dummies_multi_absorb(self): + """`absorb=["unit","time"]` invariant: with both unit and time + FE auto-routed, the period dummies collide with time-FE dummies; + `solve_ols` handles rank deficiency, slope SE on the + target per-period interaction stays well-defined and matches the + explicit `fixed_effects=` path.""" + d = self._load_golden() + target_period = int(d["target_period"]) + # Use "period" as the time/FE column name to match the data column. + res_a = MultiPeriodDiD(vcov_type="hc2").fit( + self._make_data(d), + outcome="y", + treatment="treated", + time="period", + absorb=["unit", "period"], + reference_period=int(d["reference_period"]), + unit="unit", + ) + res_f = MultiPeriodDiD(vcov_type="hc2").fit( + self._make_data(d), + outcome="y", + treatment="treated", + time="period", + fixed_effects=["unit", "period"], + reference_period=int(d["reference_period"]), + unit="unit", + ) + pe_a = res_a.period_effects[target_period] + pe_f = res_f.period_effects[target_period] + assert np.isfinite(pe_a.se) + np.testing.assert_allclose(pe_a.effect, pe_f.effect, atol=1e-12) + np.testing.assert_allclose(pe_a.se, pe_f.se, atol=1e-12) + + def test_absorb_hc2_bm_matches_fixed_effects_dummies(self): + """`absorb=` + hc2_bm equals `fixed_effects=` + hc2_bm bit-for-bit + on both per-period SE and inference (DOF transfers identically).""" + d = self._load_golden() + target_period = int(d["target_period"]) + res_a = self._fit_absorb(d, "hc2_bm", absorb_cols=("unit",)) + res_f = self._fit_fixed_effects(d, "hc2_bm", fe_cols=("unit",)) + pe_a = res_a.period_effects[target_period] + pe_f = res_f.period_effects[target_period] + np.testing.assert_allclose(pe_a.effect, pe_f.effect, atol=1e-12) + np.testing.assert_allclose(pe_a.se, pe_f.se, atol=1e-12) + np.testing.assert_allclose(pe_a.p_value, pe_f.p_value, atol=1e-12) + + def test_absorb_hc2_matches_sandwich_vcovhc(self): + """`absorb=` + hc2 matches `lm() + sandwich::vcovHC(type="HC2")` + at 1e-10 on the target per-period interaction.""" + d = self._load_golden() + target_period = int(d["target_period"]) + target_coef = f"treated_period_{target_period}" + coef_names = d["coef_names"] + idx = coef_names.index(target_coef) + expected_vcov = np.asarray(d["vcov_hc2"]).reshape(d["vcov_hc2_shape"]) + expected_se = float(np.sqrt(expected_vcov[idx, idx])) + expected_coef = float(d["coef"][idx]) + + res = self._fit_absorb(d, "hc2", absorb_cols=("unit",)) + pe = res.period_effects[target_period] + np.testing.assert_allclose(pe.effect, expected_coef, atol=1e-10) + np.testing.assert_allclose(pe.se, expected_se, atol=1e-10) + + def test_absorb_hc2_bm_matches_clubsandwich_singleton_cluster(self): + """`absorb=` + hc2_bm matches `clubSandwich::vcovCR(cluster=1:n, type="CR2")` + at 1e-10 on the target per-period interaction (singleton-cluster + CR2 = one-way HC2-BM by PT2018 §3.3).""" + d = self._load_golden() + target_period = int(d["target_period"]) + target_coef = f"treated_period_{target_period}" + coef_names = d["coef_names"] + idx = coef_names.index(target_coef) + expected_vcov = np.asarray(d["vcov_hc2_bm"]).reshape(d["vcov_hc2_bm_shape"]) + expected_se = float(np.sqrt(expected_vcov[idx, idx])) + + res = self._fit_absorb(d, "hc2_bm", absorb_cols=("unit",)) + pe = res.period_effects[target_period] + np.testing.assert_allclose(pe.se, expected_se, atol=1e-10) + + def test_absorb_plus_fixed_effects_still_rejected_under_hc2_bm(self): + """Mutual-exclusion of `absorb=` and `fixed_effects=` is preserved + on MPD across all vcov_types (the auto-route does NOT silently merge).""" + d = self._load_golden() + for vcov in ("hc1", "hc2", "hc2_bm"): + with pytest.raises(ValueError, match="Cannot use both absorb and fixed_effects"): + MultiPeriodDiD(vcov_type=vcov).fit( + self._make_data(d), + outcome="y", + treatment="treated", + time="period", + absorb=["unit"], + fixed_effects=["period"], + reference_period=int(d["reference_period"]), + unit="unit", + ) + + def test_absorb_hc2_bm_df_sensitive_inference(self): + """HC2 vs HC2-BM produce the same SE on the target per-period + interaction but different `p_value` / `conf_int` because the BM + Satterthwaite DOF differs from n-k. Guards against an unwired DOF + path (R1 review on PR #458 caught the analogous gap on DiD).""" + d = self._load_golden() + target_period = int(d["target_period"]) + res_hc2 = self._fit_absorb(d, "hc2", absorb_cols=("unit",)) + res_hc2_bm = self._fit_absorb(d, "hc2_bm", absorb_cols=("unit",)) + pe_hc2 = res_hc2.period_effects[target_period] + pe_hc2_bm = res_hc2_bm.period_effects[target_period] + np.testing.assert_allclose(pe_hc2.effect, pe_hc2_bm.effect, atol=1e-12) + np.testing.assert_allclose(pe_hc2.se, pe_hc2_bm.se, atol=1e-12) + assert pe_hc2.p_value != pe_hc2_bm.p_value, ( + "HC2 and HC2-BM should have different p_values " + "because the BM Satterthwaite DOF differs from n-k." + ) + width_hc2 = float(pe_hc2.conf_int[1] - pe_hc2.conf_int[0]) + width_hc2_bm = float(pe_hc2_bm.conf_int[1] - pe_hc2_bm.conf_int[0]) + assert width_hc2_bm > width_hc2, ( + f"HC2-BM CI width ({width_hc2_bm:.6f}) should exceed " + f"HC2 CI width ({width_hc2:.6f}) — BM Satterthwaite DOF is " + "smaller than n-k, so the critical value is larger." + ) From 226710826e110b68c1bb7aa9f7f5aa628dbc423a Mon Sep 17 00:00:00 2001 From: igerber Date: Sun, 17 May 2026 06:54:11 -0400 Subject: [PATCH 2/4] R1 polish: correct rank-deficiency claim + add survey/avg_att MPD tests MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Codex local review surfaced two findings on the MPD-absorb gate lift: P3 (methodology accuracy): REGISTRY/CHANGELOG/test-class docstring claimed the `treated` main-effect coefficient becomes NaN under rank deficiency. Empirically false — in the shipped parity fixture solve_ols drops a never-treated unit dummy (`unit_25`) and keeps `treated` finite. The pivot order is data-dependent. Rewrite to say one column from the collinear set is dropped under R-style rank-deficiency handling; per-period interactions and avg_att are identified and invariant to that choice. P2 (test coverage): the new MPD test class missed two surfaces that the DiD analogue covers: 1. Survey-weighted multi-absorb auto-route bypass of the `len(absorb) > 1 + survey_weights` reject — exercised by new `test_absorb_hc2_bm_survey_multi_absorb_auto_routes` with parity assertion against the explicit `fixed_effects=` path on both `period_effects[target]` and `avg_att`. 2. The MPD-specific `avg_att` (post-period-average) contrast did not have a direct HC2-vs-HC2-BM inference pin. Added `test_absorb_hc2_bm_avg_att_df_sensitive_inference` asserting same avg_se but different avg_p_value / wider avg_conf_int under HC2-BM (the contract guard the per-period df_sensitive test cannot reach). Co-Authored-By: Claude Opus 4.7 (1M context) --- CHANGELOG.md | 2 +- docs/methodology/REGISTRY.md | 2 +- tests/test_estimators_vcov_type.py | 98 ++++++++++++++++++++++++++++-- 3 files changed, 96 insertions(+), 6 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 2ab0c506..76810979 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,7 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] ### Added -- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that is collinear with unit dummies, so the `treated` main-effect coefficient becomes NaN under solve_ols's rank-deficiency handling — this is expected; per-period interaction parity targets (`treated:period_X`) are unaffected. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. +- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that is collinear with the sum of treated-cohort unit dummies post-auto-route, so `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling. Which specific column is dropped is pivot-order dependent (in the shipped parity fixture it is a never-treated unit dummy, not the `treated` main effect itself). The per-period interaction coefficients (`treated:period_X`) and `avg_att` are identified and invariant to that choice; parity tests target those rather than the `treated` main effect. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. - **`DifferenceInDifferences(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:382`). Previously raised `NotImplementedError` because the HC2 leverage correction and CR2 Bell-McCaffrey DOF depend on the FULL FE hat matrix, while within-transformation (FWL) preserves coefficients and residuals but not the hat. Lift via internal auto-route: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, the fit promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov. Empirically matches `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=..., type="CR2")` at ~1e-10 (verified via new `tests/test_estimators_vcov_type.py::TestDiDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `absorbed_fe_did`, with the R generator using the singleton-cluster CR2 trick for one-way HC2-BM Satterthwaite DOF). HC1/CR1 paths unchanged. `MultiPeriodDiD(absorb=...)` and `TwoWayFixedEffects` rejections remain as follow-ups (different fit-path structure). **Behavioral note (full `DiDResults` surface change under auto-route):** under the auto-route, the entire returned `DiDResults` reflects the full-dummy fit rather than the within-transformed fit. Specifically, `result.coefficients` and `result.vcov` include the FE-dummy entries (matching the `fixed_effects=` path), `result.residuals` and `result.fitted_values` are on the un-demeaned outcome scale, and `result.r_squared` is computed on the un-demeaned outcome (so it absorbs the FE variance and will typically be higher than the within-R²). `result.att` is invariant to this routing (FWL guarantee). Downstream consumers reading `result.att` are unaffected; consumers reading the broader result surface should expect the full-dummy values. **Survey-design scope:** the auto-route changes the FE handling (and removes the prior absorbed-FE rejection), but `survey_design=` continues to drive its own variance path (Taylor-series linearization or replicate-weight variance, per the existing survey contract) rather than the analytical HC2/HC2-BM sandwich. The auto-route is therefore methodologically meaningful for non-survey fits and for the FE-handling side of survey fits; analytical small-sample inference under `vcov_type in {"hc2","hc2_bm"}` is bypassed when a survey design is supplied. - **BaconDecomposition R parity goldens.** Closes the PR-B deferral row in `TODO.md`. JSON goldens at `benchmarks/data/r_bacondecomp_golden.json` generated from the committed `benchmarks/R/generate_bacon_golden.R` script (3 fixtures: `uniform_3groups_with_never_treated`, `two_groups_no_never_treated`, `always_treated_remapped`) against `bacondecomp 0.1.1` on R 4.5.2. `tests/test_methodology_bacon.py::TestBaconParityR` now active (4 tests, no skips): TWFE coefficient parity at `atol=1e-6` across all 3 fixtures; weights-sum parity at `atol=1e-6` across all 3 fixtures; per-component estimate + weight parity at `atol=1e-6` on the 2 non-remap fixtures **and on the 6 timing-vs-timing rows of `always_treated_remapped`** (carve-out narrowed to U-bucket rows only); plus a dedicated fold-back test (`test_always_treated_remapped_fold_back_matches_r`) that pins the **documented convention divergence** on `always_treated_remapped` (R keeps `first_treat=1` as a distinct timing cohort and emits `Later vs Always Treated` comparisons; Python's paper-footnote-11 convention remaps those units to `U` and folds them into a single `treated_vs_never` cell per treated cohort) by aggregating R's split rows per cohort and asserting they match Python's single fold at `atol=1e-6`. The aggregate is invariant per Theorem 1; the per-component breakdown differs structurally between conventions but the fold-back is now directly asserted. New `**Note (R parity convention divergence on always-treated)**` and `**Deviation (first-period boundary extension on always-treated remap)**` in `docs/methodology/REGISTRY.md`. **First-period boundary deviation:** the paper uses strict `t_i < 1` for the always-treated bucket; the library uses the inclusive `first_treat <= min(time)` rule and folds `first_treat == min(time)` cohorts into `U`. R does NOT apply this fold (it keeps such cohorts as their own bucket). When `min(time) > 1` the rules coincide. Explicitly labeled in REGISTRY's Deviations block and mirrored in `METHODOLOGY_REVIEW.md` and `bacon.py`. METHODOLOGY_REVIEW.md tracker row promoted `**Complete** (R parity goldens pending)` → `**Complete**`. - **`generate_ddd_panel_data` — panel-structured DGP for Triple-Difference power analysis** (`diff_diff/prep_dgp.py`). New public function exported from `diff_diff` and `diff_diff.prep` for panel DDD simulations. Cross-sectional `generate_ddd_data` remains available unchanged. Produces a balanced panel of `n_units × n_periods` with two unit-level binary dimensions (`group`, `partition`) and a derived `post = 1[period >= treatment_period]` indicator; columns: `unit, period, outcome, group, partition, post, treated, true_effect` (+ `x1, x2` when `add_covariates=True`). DDD-CPT identification holds because the `group * partition` interaction enters as a unit-level (time-invariant) term, leaving the triple-interaction `treatment_effect * group * partition * post` as the sole source of differential group × partition trend. Compatible with `TripleDifference(cluster="unit").fit(..., time="post")` (the cluster kwarg is required because `TripleDifference` is the repeated-cross-section `panel=FALSE` estimator and unclustered SE on panel-generated rows understates variance under within-unit serial correlation; the point estimate `att` is invariant to clustering — see the new `TripleDifference` REGISTRY note on panel-shaped input). Users get panel-realistic unit fixed effects and within-unit serial correlation while the binary 2×2×2 estimator surface is unchanged. **Stratified allocation:** the partition split is drawn stratified-by-group at the requested `partition_frac` so every `(group, partition)` cell receives at least one unit; a targeted `ValueError` is raised at fit-time when the rounded cell counts (`n_units`, `group_frac`, `partition_frac`) would leave any cell empty. This guarantees the 2x2x2 DDD surface is populated for any valid input — independent marginal sampling (the cross-sectional `generate_ddd_data` convention) could collapse cells when marginals are small (e.g., `n_units=4, group_frac=partition_frac=0.25`). Validates `1 <= treatment_period < n_periods`, `group_frac` and `partition_frac` strictly in `(0, 1)`, and `n_units >= 4`. Deterministic recovery (`noise_sd=0`) matches `treatment_effect` to ~1e-15 (covered by `tests/test_prep.py::TestGenerateDddPanelData`, 16 tests including infeasible-config rejection and smallest-feasible-config round-trip through `TripleDifference.fit`). `power.simulate_power` is NOT yet auto-routed to the panel DGP for `TripleDifference` (the existing `_ddd_dgp_kwargs` registry entry still ignores `n_periods` and the existing `_check_ddd_dgp_compat` warning still fires on non-default kwargs) — that wiring is tracked as a follow-up in TODO.md. diff --git a/docs/methodology/REGISTRY.md b/docs/methodology/REGISTRY.md index 5b7aeefe..a77961ca 100644 --- a/docs/methodology/REGISTRY.md +++ b/docs/methodology/REGISTRY.md @@ -2551,7 +2551,7 @@ Shipped in `diff_diff/had_pretests.py` as `stute_joint_pretest()` (residuals-in - **Note (scope limitation on absorbed FE):** HC2 and HC2 + Bell-McCaffrey on within-transformed designs still depend on the FULL FE hat matrix because FWL preserves coefficients and residuals but NOT the hat matrix: `h_ii = x_i' (X'X)^{-1} x_i` on the reduced design is not the diagonal of the full FE projection, and CR2's block adjustment `A_g = (I - H_gg)^{-1/2}` likewise depends on the full cluster-block hat matrix. The status across the three estimators that previously rejected this combination: - **`DifferenceInDifferences(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** When the user pairs `absorb=` with HC2 / HC2-BM, `DiD.fit()` internally promotes the absorb columns to `fixed_effects=` so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection. Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=..., type="CR2")` (singleton-cluster CR2 trick for one-way HC2-BM Satterthwaite DOF; PT2018 §3.3 unweighted CR2 algebra). **User-visible surface change**: under the auto-route, the entire `DiDResults` (coefficients, vcov, residuals, fitted_values, r_squared) reflect the full-dummy fit rather than the within-transformed fit — the FE-dummy entries are included in `result.coefficients` / `result.vcov`, `r_squared` is computed on the un-demeaned outcome, and `residuals` / `fitted_values` are on the original scale. `result.att` is unaffected (FWL-equivalent). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). **Survey-design scope**: when `survey_design=` is supplied, the existing survey variance path (Taylor-series linearization / replicate weights) takes precedence over the analytical HC2/HC2-BM sandwich; the auto-route only changes the FE handling (removing the prior reject) and does not redirect to the analytical small-sample sandwich on survey fits. - **`TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` — still rejects.** TWFE is a standalone class with no `fixed_effects=` equivalent path, so the same auto-route surgery used for DiD-absorb and MPD-absorb is not directly applicable; lifting requires building the full-dummy design inline or refactoring TWFE to delegate to DiD. Tracked as a follow-up in `TODO.md`. - - **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** Same auto-route pattern as `DifferenceInDifferences`: `MultiPeriodDiD.fit()` internally promotes the absorb columns to `fixed_effects=` for HC2 / HC2-BM callers, so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture; the parity target is a per-period interaction `treated:period_X` (the `treated` main-effect coefficient becomes NaN under rank deficiency because MPD requires a time-invariant ever-treated indicator that is collinear with unit dummies — expected behavior). Same `MultiPeriodDiDResults` surface change as DiD: `vcov`, `residuals`, `fitted_values`, `r_squared`, and `coefficients` reflect the full-dummy fit, with `period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` invariant by FWL. HC1/CR1 paths on `absorb=` are unchanged (no leverage term). Same survey-design scope as DiD: replicate-weight variance routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design rather than the per-replicate refit branch (which targets the demeaning path); since the auto-routed design does not depend on replicate weights, no refit is needed. + - **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** Same auto-route pattern as `DifferenceInDifferences`: `MultiPeriodDiD.fit()` internally promotes the absorb columns to `fixed_effects=` for HC2 / HC2-BM callers, so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture; the parity target is a per-period interaction `treated:period_X` because MPD requires the `treated` column to be a time-invariant ever-treated indicator, which is exactly collinear with the sum of treated-cohort unit dummies post-auto-route. `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling; in the shipped parity fixture (4 ever-treated cohorts of 5 units + 1 never-treated cohort of 5 units) it drops a unit dummy from the never-treated cohort (`unit_25`) and the `treated` main effect remains finite, but the specific column that gets NaN'd is pivot-order dependent and users should not rely on `treated` being either kept or dropped on other fixtures. Either way, the slope coefficients (`treated:period_X`) and the post-period-average `avg_att` are identified and invariant to which column was dropped. Same `MultiPeriodDiDResults` surface change as DiD: `vcov`, `residuals`, `fitted_values`, `r_squared`, and `coefficients` reflect the full-dummy fit, with `period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` invariant by FWL. HC1/CR1 paths on `absorb=` are unchanged (no leverage term). Same survey-design scope as DiD: replicate-weight variance routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design rather than the per-replicate refit branch (which targets the demeaning path); since the auto-routed design does not depend on replicate weights, no refit is needed. - Workarounds for the still-rejecting paths: use `vcov_type="hc1"` (HC1/CR1 have no leverage term and survive FWL), or switch to `fixed_effects=` dummies so the hat matrix is computed on the full design. - [x] Phase 1a: `vcov_type` enum threaded through `DifferenceInDifferences` (`MultiPeriodDiD`, `TwoWayFixedEffects` inherit); `robust=True` <=> `vcov_type="hc1"`, `robust=False` <=> `vcov_type="classical"`. Conflict detection at `__init__`. Results summary prints the variance-family label. - **Note (deviation from the fully-symmetric enum):** `MultiPeriodDiD(cluster=..., vcov_type="hc2_bm")` is intentionally **not supported** and raises `NotImplementedError`. The scalar-coefficient `DifferenceInDifferences` path handles the cluster + CR2 Bell-McCaffrey combination (`_compute_cr2_bm` returns a per-coefficient Satterthwaite DOF that is valid for the single-ATT contrast), but `MultiPeriodDiD` also reports a post-period-average ATT constructed as a *contrast* of the event-study coefficients. The cluster-aware CR2 BM DOF for that contrast (i.e., the Pustejovsky-Tipton 2018 per-cluster adjustment matrices applied to an arbitrary aggregation contrast) is not yet implemented. Pairing CR2 cluster-robust SEs with the one-way Imbens-Kolesar (2016) contrast DOF would be a broken hybrid, so the combination fails fast with a clear workaround message (drop the cluster for one-way HC2+BM, or use `vcov_type="hc1"` with cluster for CR1 Liang-Zeger). Tracked in `TODO.md` under Methodology/Correctness. Applies only to `MultiPeriodDiD`; `DifferenceInDifferences(cluster=..., vcov_type="hc2_bm")` works. diff --git a/tests/test_estimators_vcov_type.py b/tests/test_estimators_vcov_type.py index dcc5fb48..4c740fec 100644 --- a/tests/test_estimators_vcov_type.py +++ b/tests/test_estimators_vcov_type.py @@ -1202,10 +1202,14 @@ class TestMPDAbsorbedFERParity: Collinearity note: MPD's `treated` is a time-invariant ever-treated indicator, so it is perfectly collinear with the sum of treated-cohort unit dummies post-auto-route. `solve_ols` resolves this by dropping - one column (typically a unit dummy from the never-treated cohort); the - `treated` main-effect coefficient absorbs the dropped column's effect. - Tests pin parity on a per-period interaction (`treated:period_4`) where - no collinearity exists, exposed as `result.period_effects[4]`. + one column from that collinear set under R-style rank-deficiency + handling. In the shipped parity fixture the dropped column is a unit + dummy from the never-treated cohort (`unit_25`) — the `treated` main + effect remains finite there — but the specific column dropped is + pivot-order dependent and not guaranteed across fixtures. Tests + therefore pin parity on a per-period interaction (`treated:period_4`) + which is identified independent of that choice, exposed as + `result.period_effects[4]`. """ def _load_golden(self): @@ -1393,3 +1397,89 @@ def test_absorb_hc2_bm_df_sensitive_inference(self): f"HC2 CI width ({width_hc2:.6f}) — BM Satterthwaite DOF is " "smaller than n-k, so the critical value is larger." ) + + def test_absorb_hc2_bm_avg_att_df_sensitive_inference(self): + """The post-period-average ATT (`avg_att`, the MPD-specific + contrast that does NOT have a DiD analogue) must also reflect + the Satterthwaite DOF under HC2-BM: HC2 and HC2-BM share + `avg_se` but differ in `avg_p_value` and `avg_conf_int`. This is + the MPD-specific inference pin that the DiD test class cannot + cover.""" + d = self._load_golden() + res_hc2 = self._fit_absorb(d, "hc2", absorb_cols=("unit",)) + res_hc2_bm = self._fit_absorb(d, "hc2_bm", absorb_cols=("unit",)) + np.testing.assert_allclose(res_hc2.avg_att, res_hc2_bm.avg_att, atol=1e-12) + np.testing.assert_allclose(res_hc2.avg_se, res_hc2_bm.avg_se, atol=1e-12) + assert res_hc2.avg_p_value != res_hc2_bm.avg_p_value, ( + "HC2 and HC2-BM should produce different `avg_p_value` because " + "the BM Satterthwaite DOF on the post-period-average contrast " + "differs from n-k. Same p_value indicates the DOF was not " + "propagated to the avg_att inference path." + ) + width_hc2 = float(res_hc2.avg_conf_int[1] - res_hc2.avg_conf_int[0]) + width_hc2_bm = float(res_hc2_bm.avg_conf_int[1] - res_hc2_bm.avg_conf_int[0]) + assert width_hc2_bm > width_hc2, ( + f"HC2-BM avg_att CI width ({width_hc2_bm:.6f}) should exceed " + f"HC2 avg_att CI width ({width_hc2:.6f}) — BM Satterthwaite " + "DOF is smaller than n-k, so the critical value is larger." + ) + + def test_absorb_hc2_bm_survey_multi_absorb_auto_routes(self): + """Survey-weighted multi-absorb + HC2-BM should auto-route, not reject. + + Mirrors the DiD-class test of the same name: the legacy guard at + `estimators.py:1505-1512` rejects `survey_design + len(absorb) > 1` + because single-pass demeaning is not the correct weighted FWL + projection for multiple absorbed dimensions. But when the auto-route + fires (hc2/hc2_bm), absorb is swapped for fixed_effects= BEFORE the + survey guard sees it, so the demeaning rationale doesn't apply. The + auto-route placement is precisely tuned for this case; this test + pins it on the MPD path.""" + from diff_diff import SurveyDesign + + d = self._load_golden() + rng = np.random.default_rng(20260420) + n = len(d["y"]) + data = pd.DataFrame( + { + "unit": d["unit"], + "period": d["period"], + "treated": d["treated"], + "y": d["y"], + "weight": rng.uniform(0.5, 2.0, size=n), + } + ) + sd = SurveyDesign(weights="weight", weight_type="aweight") + # Multi-absorb (unit + period) + survey + hc2_bm: auto-route fires + # and the multi-absorb-survey guard is bypassed cleanly. + res = MultiPeriodDiD(vcov_type="hc2_bm").fit( + data, + outcome="y", + treatment="treated", + time="period", + absorb=["unit", "period"], + reference_period=int(d["reference_period"]), + unit="unit", + survey_design=sd, + ) + # Parity invariant: the explicit fixed_effects= path on the same + # data must produce the same per-period SE. + res_fe = MultiPeriodDiD(vcov_type="hc2_bm").fit( + data, + outcome="y", + treatment="treated", + time="period", + fixed_effects=["unit", "period"], + reference_period=int(d["reference_period"]), + unit="unit", + survey_design=sd, + ) + target_period = int(d["target_period"]) + pe_a = res.period_effects[target_period] + pe_f = res_fe.period_effects[target_period] + assert np.isfinite(pe_a.effect) + assert np.isfinite(pe_a.se) + np.testing.assert_allclose(pe_a.effect, pe_f.effect, atol=1e-12) + np.testing.assert_allclose(pe_a.se, pe_f.se, atol=1e-12) + np.testing.assert_allclose(res.avg_att, res_fe.avg_att, atol=1e-12) + np.testing.assert_allclose(res.avg_se, res_fe.avg_se, atol=1e-12) From f25373a971050199d9f8582963de0950ed01986d Mon Sep 17 00:00:00 2001 From: igerber Date: Sun, 17 May 2026 07:01:58 -0400 Subject: [PATCH 3/4] R2 polish: replicate-weight regression + corrected collinearity wording MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Codex R2 returned ✅ with two P3s. P3 (test coverage upgraded to actionable per feedback_test_coverage_gap_treat_as_actionable.md): the survey test added in R1 used aweight (generic survey-vcov path), but the CHANGELOG/REGISTRY claim specifically that the auto-route short-circuits the absorb-refit replicate-variance branch at estimators.py:1693. Added test_absorb_hc2_bm_replicate_weights_auto_routes using SurveyDesign(replicate_method="JK1", replicate_weights=[...]) that exercises the replicate path and pins SE parity vs the explicit fixed_effects= form on both period_effects[target_period] and avg_att. Passing at atol=1e-12 confirms the documented short-circuit works as claimed. P3 (doc accuracy): REGISTRY/CHANGELOG/test-class docstring described the `treated` alias as "exactly collinear with the sum of treated-cohort unit dummies". Under pd.get_dummies(drop_first=True) the exact alias depends on the omitted FE reference category (and the intercept), not just on the cohort-dummy sum. Rewrite to say `treated` lies in the span of the intercept plus the post-auto-route unit FE dummies; which specific nuisance column gets dropped is pivot-order and dummy-coding dependent. Co-Authored-By: Claude Opus 4.7 (1M context) --- CHANGELOG.md | 2 +- docs/methodology/REGISTRY.md | 2 +- tests/test_estimators_vcov_type.py | 92 +++++++++++++++++++++++++++--- 3 files changed, 85 insertions(+), 11 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 76810979..6aee904a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,7 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] ### Added -- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that is collinear with the sum of treated-cohort unit dummies post-auto-route, so `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling. Which specific column is dropped is pivot-order dependent (in the shipped parity fixture it is a never-treated unit dummy, not the `treated` main effect itself). The per-period interaction coefficients (`treated:period_X`) and `avg_att` are identified and invariant to that choice; parity tests target those rather than the `treated` main effect. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. +- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that lies in the span of the intercept and the post-auto-route unit FE dummies (the exact alias depends on the omitted FE reference category under `pd.get_dummies(drop_first=True)`, not just on "the sum of treated-cohort unit dummies"), so `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling. Which specific column is dropped is pivot-order and dummy-coding dependent (in the shipped parity fixture it is a never-treated unit dummy, not the `treated` main effect itself). The per-period interaction coefficients (`treated:period_X`) and `avg_att` are identified and invariant to that choice; parity tests target those rather than the `treated` main effect. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. - **`DifferenceInDifferences(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:382`). Previously raised `NotImplementedError` because the HC2 leverage correction and CR2 Bell-McCaffrey DOF depend on the FULL FE hat matrix, while within-transformation (FWL) preserves coefficients and residuals but not the hat. Lift via internal auto-route: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, the fit promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov. Empirically matches `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=..., type="CR2")` at ~1e-10 (verified via new `tests/test_estimators_vcov_type.py::TestDiDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `absorbed_fe_did`, with the R generator using the singleton-cluster CR2 trick for one-way HC2-BM Satterthwaite DOF). HC1/CR1 paths unchanged. `MultiPeriodDiD(absorb=...)` and `TwoWayFixedEffects` rejections remain as follow-ups (different fit-path structure). **Behavioral note (full `DiDResults` surface change under auto-route):** under the auto-route, the entire returned `DiDResults` reflects the full-dummy fit rather than the within-transformed fit. Specifically, `result.coefficients` and `result.vcov` include the FE-dummy entries (matching the `fixed_effects=` path), `result.residuals` and `result.fitted_values` are on the un-demeaned outcome scale, and `result.r_squared` is computed on the un-demeaned outcome (so it absorbs the FE variance and will typically be higher than the within-R²). `result.att` is invariant to this routing (FWL guarantee). Downstream consumers reading `result.att` are unaffected; consumers reading the broader result surface should expect the full-dummy values. **Survey-design scope:** the auto-route changes the FE handling (and removes the prior absorbed-FE rejection), but `survey_design=` continues to drive its own variance path (Taylor-series linearization or replicate-weight variance, per the existing survey contract) rather than the analytical HC2/HC2-BM sandwich. The auto-route is therefore methodologically meaningful for non-survey fits and for the FE-handling side of survey fits; analytical small-sample inference under `vcov_type in {"hc2","hc2_bm"}` is bypassed when a survey design is supplied. - **BaconDecomposition R parity goldens.** Closes the PR-B deferral row in `TODO.md`. JSON goldens at `benchmarks/data/r_bacondecomp_golden.json` generated from the committed `benchmarks/R/generate_bacon_golden.R` script (3 fixtures: `uniform_3groups_with_never_treated`, `two_groups_no_never_treated`, `always_treated_remapped`) against `bacondecomp 0.1.1` on R 4.5.2. `tests/test_methodology_bacon.py::TestBaconParityR` now active (4 tests, no skips): TWFE coefficient parity at `atol=1e-6` across all 3 fixtures; weights-sum parity at `atol=1e-6` across all 3 fixtures; per-component estimate + weight parity at `atol=1e-6` on the 2 non-remap fixtures **and on the 6 timing-vs-timing rows of `always_treated_remapped`** (carve-out narrowed to U-bucket rows only); plus a dedicated fold-back test (`test_always_treated_remapped_fold_back_matches_r`) that pins the **documented convention divergence** on `always_treated_remapped` (R keeps `first_treat=1` as a distinct timing cohort and emits `Later vs Always Treated` comparisons; Python's paper-footnote-11 convention remaps those units to `U` and folds them into a single `treated_vs_never` cell per treated cohort) by aggregating R's split rows per cohort and asserting they match Python's single fold at `atol=1e-6`. The aggregate is invariant per Theorem 1; the per-component breakdown differs structurally between conventions but the fold-back is now directly asserted. New `**Note (R parity convention divergence on always-treated)**` and `**Deviation (first-period boundary extension on always-treated remap)**` in `docs/methodology/REGISTRY.md`. **First-period boundary deviation:** the paper uses strict `t_i < 1` for the always-treated bucket; the library uses the inclusive `first_treat <= min(time)` rule and folds `first_treat == min(time)` cohorts into `U`. R does NOT apply this fold (it keeps such cohorts as their own bucket). When `min(time) > 1` the rules coincide. Explicitly labeled in REGISTRY's Deviations block and mirrored in `METHODOLOGY_REVIEW.md` and `bacon.py`. METHODOLOGY_REVIEW.md tracker row promoted `**Complete** (R parity goldens pending)` → `**Complete**`. - **`generate_ddd_panel_data` — panel-structured DGP for Triple-Difference power analysis** (`diff_diff/prep_dgp.py`). New public function exported from `diff_diff` and `diff_diff.prep` for panel DDD simulations. Cross-sectional `generate_ddd_data` remains available unchanged. Produces a balanced panel of `n_units × n_periods` with two unit-level binary dimensions (`group`, `partition`) and a derived `post = 1[period >= treatment_period]` indicator; columns: `unit, period, outcome, group, partition, post, treated, true_effect` (+ `x1, x2` when `add_covariates=True`). DDD-CPT identification holds because the `group * partition` interaction enters as a unit-level (time-invariant) term, leaving the triple-interaction `treatment_effect * group * partition * post` as the sole source of differential group × partition trend. Compatible with `TripleDifference(cluster="unit").fit(..., time="post")` (the cluster kwarg is required because `TripleDifference` is the repeated-cross-section `panel=FALSE` estimator and unclustered SE on panel-generated rows understates variance under within-unit serial correlation; the point estimate `att` is invariant to clustering — see the new `TripleDifference` REGISTRY note on panel-shaped input). Users get panel-realistic unit fixed effects and within-unit serial correlation while the binary 2×2×2 estimator surface is unchanged. **Stratified allocation:** the partition split is drawn stratified-by-group at the requested `partition_frac` so every `(group, partition)` cell receives at least one unit; a targeted `ValueError` is raised at fit-time when the rounded cell counts (`n_units`, `group_frac`, `partition_frac`) would leave any cell empty. This guarantees the 2x2x2 DDD surface is populated for any valid input — independent marginal sampling (the cross-sectional `generate_ddd_data` convention) could collapse cells when marginals are small (e.g., `n_units=4, group_frac=partition_frac=0.25`). Validates `1 <= treatment_period < n_periods`, `group_frac` and `partition_frac` strictly in `(0, 1)`, and `n_units >= 4`. Deterministic recovery (`noise_sd=0`) matches `treatment_effect` to ~1e-15 (covered by `tests/test_prep.py::TestGenerateDddPanelData`, 16 tests including infeasible-config rejection and smallest-feasible-config round-trip through `TripleDifference.fit`). `power.simulate_power` is NOT yet auto-routed to the panel DGP for `TripleDifference` (the existing `_ddd_dgp_kwargs` registry entry still ignores `n_periods` and the existing `_check_ddd_dgp_compat` warning still fires on non-default kwargs) — that wiring is tracked as a follow-up in TODO.md. diff --git a/docs/methodology/REGISTRY.md b/docs/methodology/REGISTRY.md index a77961ca..774cbacd 100644 --- a/docs/methodology/REGISTRY.md +++ b/docs/methodology/REGISTRY.md @@ -2551,7 +2551,7 @@ Shipped in `diff_diff/had_pretests.py` as `stute_joint_pretest()` (residuals-in - **Note (scope limitation on absorbed FE):** HC2 and HC2 + Bell-McCaffrey on within-transformed designs still depend on the FULL FE hat matrix because FWL preserves coefficients and residuals but NOT the hat matrix: `h_ii = x_i' (X'X)^{-1} x_i` on the reduced design is not the diagonal of the full FE projection, and CR2's block adjustment `A_g = (I - H_gg)^{-1/2}` likewise depends on the full cluster-block hat matrix. The status across the three estimators that previously rejected this combination: - **`DifferenceInDifferences(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** When the user pairs `absorb=` with HC2 / HC2-BM, `DiD.fit()` internally promotes the absorb columns to `fixed_effects=` so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection. Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=..., type="CR2")` (singleton-cluster CR2 trick for one-way HC2-BM Satterthwaite DOF; PT2018 §3.3 unweighted CR2 algebra). **User-visible surface change**: under the auto-route, the entire `DiDResults` (coefficients, vcov, residuals, fitted_values, r_squared) reflect the full-dummy fit rather than the within-transformed fit — the FE-dummy entries are included in `result.coefficients` / `result.vcov`, `r_squared` is computed on the un-demeaned outcome, and `residuals` / `fitted_values` are on the original scale. `result.att` is unaffected (FWL-equivalent). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). **Survey-design scope**: when `survey_design=` is supplied, the existing survey variance path (Taylor-series linearization / replicate weights) takes precedence over the analytical HC2/HC2-BM sandwich; the auto-route only changes the FE handling (removing the prior reject) and does not redirect to the analytical small-sample sandwich on survey fits. - **`TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` — still rejects.** TWFE is a standalone class with no `fixed_effects=` equivalent path, so the same auto-route surgery used for DiD-absorb and MPD-absorb is not directly applicable; lifting requires building the full-dummy design inline or refactoring TWFE to delegate to DiD. Tracked as a follow-up in `TODO.md`. - - **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** Same auto-route pattern as `DifferenceInDifferences`: `MultiPeriodDiD.fit()` internally promotes the absorb columns to `fixed_effects=` for HC2 / HC2-BM callers, so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture; the parity target is a per-period interaction `treated:period_X` because MPD requires the `treated` column to be a time-invariant ever-treated indicator, which is exactly collinear with the sum of treated-cohort unit dummies post-auto-route. `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling; in the shipped parity fixture (4 ever-treated cohorts of 5 units + 1 never-treated cohort of 5 units) it drops a unit dummy from the never-treated cohort (`unit_25`) and the `treated` main effect remains finite, but the specific column that gets NaN'd is pivot-order dependent and users should not rely on `treated` being either kept or dropped on other fixtures. Either way, the slope coefficients (`treated:period_X`) and the post-period-average `avg_att` are identified and invariant to which column was dropped. Same `MultiPeriodDiDResults` surface change as DiD: `vcov`, `residuals`, `fitted_values`, `r_squared`, and `coefficients` reflect the full-dummy fit, with `period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` invariant by FWL. HC1/CR1 paths on `absorb=` are unchanged (no leverage term). Same survey-design scope as DiD: replicate-weight variance routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design rather than the per-replicate refit branch (which targets the demeaning path); since the auto-routed design does not depend on replicate weights, no refit is needed. + - **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2","hc2_bm"})` — SUPPORTED (auto-route).** Same auto-route pattern as `DifferenceInDifferences`: `MultiPeriodDiD.fit()` internally promotes the absorb columns to `fixed_effects=` for HC2 / HC2-BM callers, so the existing full-dummy code path computes the algebraically correct vcov from the full FE projection on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture; the parity target is a per-period interaction `treated:period_X` because MPD requires the `treated` column to be a time-invariant ever-treated indicator, which lies in the span of the intercept and the post-auto-route unit FE dummies (under `pd.get_dummies(..., drop_first=True)` the dropped reference unit is implicit in the intercept, so the exact alias relation depends on the omitted FE category — it is NOT simply "the sum of treated-cohort unit dummies"). `solve_ols` drops one column from the collinear set under R-style rank-deficiency handling; in the shipped parity fixture (4 ever-treated cohorts of 5 units + 1 never-treated cohort of 5 units) it drops a unit dummy from the never-treated cohort (`unit_25`) and the `treated` main effect remains finite, but the specific column that gets NaN'd is pivot-order and dummy-coding dependent. Either way, the slope coefficients (`treated:period_X`) and the post-period-average `avg_att` are identified and invariant to which column was dropped. Same `MultiPeriodDiDResults` surface change as DiD: `vcov`, `residuals`, `fitted_values`, `r_squared`, and `coefficients` reflect the full-dummy fit, with `period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` invariant by FWL. HC1/CR1 paths on `absorb=` are unchanged (no leverage term). Same survey-design scope as DiD: replicate-weight variance routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design rather than the per-replicate refit branch (which targets the demeaning path); since the auto-routed design does not depend on replicate weights, no refit is needed. - Workarounds for the still-rejecting paths: use `vcov_type="hc1"` (HC1/CR1 have no leverage term and survive FWL), or switch to `fixed_effects=` dummies so the hat matrix is computed on the full design. - [x] Phase 1a: `vcov_type` enum threaded through `DifferenceInDifferences` (`MultiPeriodDiD`, `TwoWayFixedEffects` inherit); `robust=True` <=> `vcov_type="hc1"`, `robust=False` <=> `vcov_type="classical"`. Conflict detection at `__init__`. Results summary prints the variance-family label. - **Note (deviation from the fully-symmetric enum):** `MultiPeriodDiD(cluster=..., vcov_type="hc2_bm")` is intentionally **not supported** and raises `NotImplementedError`. The scalar-coefficient `DifferenceInDifferences` path handles the cluster + CR2 Bell-McCaffrey combination (`_compute_cr2_bm` returns a per-coefficient Satterthwaite DOF that is valid for the single-ATT contrast), but `MultiPeriodDiD` also reports a post-period-average ATT constructed as a *contrast* of the event-study coefficients. The cluster-aware CR2 BM DOF for that contrast (i.e., the Pustejovsky-Tipton 2018 per-cluster adjustment matrices applied to an arbitrary aggregation contrast) is not yet implemented. Pairing CR2 cluster-robust SEs with the one-way Imbens-Kolesar (2016) contrast DOF would be a broken hybrid, so the combination fails fast with a clear workaround message (drop the cluster for one-way HC2+BM, or use `vcov_type="hc1"` with cluster for CR1 Liang-Zeger). Tracked in `TODO.md` under Methodology/Correctness. Applies only to `MultiPeriodDiD`; `DifferenceInDifferences(cluster=..., vcov_type="hc2_bm")` works. diff --git a/tests/test_estimators_vcov_type.py b/tests/test_estimators_vcov_type.py index 4c740fec..58542b8a 100644 --- a/tests/test_estimators_vcov_type.py +++ b/tests/test_estimators_vcov_type.py @@ -1200,15 +1200,18 @@ class TestMPDAbsorbedFERParity: builds the full-dummy design that R's `lm()` produces. Collinearity note: MPD's `treated` is a time-invariant ever-treated - indicator, so it is perfectly collinear with the sum of treated-cohort - unit dummies post-auto-route. `solve_ols` resolves this by dropping - one column from that collinear set under R-style rank-deficiency - handling. In the shipped parity fixture the dropped column is a unit - dummy from the never-treated cohort (`unit_25`) — the `treated` main - effect remains finite there — but the specific column dropped is - pivot-order dependent and not guaranteed across fixtures. Tests - therefore pin parity on a per-period interaction (`treated:period_4`) - which is identified independent of that choice, exposed as + indicator, so it lies in the span of the intercept and the + post-auto-route unit FE dummies (under `pd.get_dummies(drop_first=True)` + the dropped reference unit is folded into the intercept; the exact + alias relation depends on the omitted category and is NOT simply + "the sum of treated-cohort unit dummies"). `solve_ols` resolves this + by dropping one column from the collinear set under R-style + rank-deficiency handling. In the shipped parity fixture the dropped + column is a unit dummy from the never-treated cohort (`unit_25`) and + the `treated` main effect remains finite, but the specific column + dropped is pivot-order and dummy-coding dependent. Tests therefore + pin parity on a per-period interaction (`treated:period_4`) which is + identified independent of that choice, exposed as `result.period_effects[4]`. """ @@ -1483,3 +1486,74 @@ def test_absorb_hc2_bm_survey_multi_absorb_auto_routes(self): np.testing.assert_allclose(pe_a.se, pe_f.se, atol=1e-12) np.testing.assert_allclose(res.avg_att, res_fe.avg_att, atol=1e-12) np.testing.assert_allclose(res.avg_se, res_fe.avg_se, atol=1e-12) + + def test_absorb_hc2_bm_replicate_weights_auto_routes(self): + """Replicate-weight survey design + absorb + HC2-BM auto-routes + through `compute_replicate_vcov` on the full-dummy design. + + The CHANGELOG/REGISTRY claim that, under the auto-route, the + survey-replicate absorb-refit branch at `estimators.py:1693` is + short-circuited (no per-replicate refit needed because the + full-dummy design does not depend on replicate weights — the + standard `compute_replicate_vcov` path applies directly). This + test pins the parity invariant on a JK1 fixture: `absorb=` + + replicate weights must produce the same `period_effects` + and `avg_att` SEs as the explicit `fixed_effects=` form.""" + from diff_diff import SurveyDesign + + d = self._load_golden() + rng = np.random.default_rng(20260420) + n = len(d["y"]) + data = pd.DataFrame( + { + "unit": d["unit"], + "period": d["period"], + "treated": d["treated"], + "y": d["y"], + "weight": rng.uniform(0.5, 2.0, size=n), + } + ) + # 10 JK1 jackknife replicate-weight columns; weights drawn from + # {0.5, 1.5} match the BRR pattern of the existing replicate + # tests in this file. + rep_cols = [f"rep{r}" for r in range(10)] + for col in rep_cols: + data[col] = rng.choice([0.5, 1.5], size=n) + sd = SurveyDesign( + weights="weight", + replicate_weights=rep_cols, + replicate_method="JK1", + replicate_scale=1.0, + ) + res_absorb = MultiPeriodDiD(vcov_type="hc2_bm").fit( + data, + outcome="y", + treatment="treated", + time="period", + absorb=["unit", "period"], + reference_period=int(d["reference_period"]), + unit="unit", + survey_design=sd, + ) + res_fe = MultiPeriodDiD(vcov_type="hc2_bm").fit( + data, + outcome="y", + treatment="treated", + time="period", + fixed_effects=["unit", "period"], + reference_period=int(d["reference_period"]), + unit="unit", + survey_design=sd, + ) + target_period = int(d["target_period"]) + pe_a = res_absorb.period_effects[target_period] + pe_f = res_fe.period_effects[target_period] + assert np.isfinite(pe_a.effect) + assert np.isfinite(pe_a.se) + # The auto-route short-circuits the absorb-refit branch and routes + # both calls through the standard replicate-vcov path; SE parity + # is therefore exact (bit-identical, not just to within 1e-10). + np.testing.assert_allclose(pe_a.effect, pe_f.effect, atol=1e-12) + np.testing.assert_allclose(pe_a.se, pe_f.se, atol=1e-12) + np.testing.assert_allclose(res_absorb.avg_att, res_fe.avg_att, atol=1e-12) + np.testing.assert_allclose(res_absorb.avg_se, res_fe.avg_se, atol=1e-12) From 1e25f306f8576501eb861507007bbb686a06b9ce Mon Sep 17 00:00:00 2001 From: igerber Date: Sun, 17 May 2026 07:19:49 -0400 Subject: [PATCH 4/4] R1 CI: fix duplicate period_X dummy names on MPD multi-absorb path MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit CI Codex review on PR #459 surfaced a P1 newly exposed by the auto-route: when MPD(absorb=["unit","period"]) auto-routes to fixed_effects=["unit", "period"], the existing fixed_effects= expansion loop adds `period_X` dummies via `pd.get_dummies(prefix="period")` that collide on name with the event-study period dummies MPD already builds for non-reference periods. The duplicate `var_names` entries silently collapse in `coef_dict = {name: coef for name, coef in zip(var_names, coefficients)}`, overwriting the real event-study coefficients with the rank-deficient NaN drops on the redundant FE block. Result: `len(coefficients) < vcov.shape[0]` and `coefficients["period_X"] = NaN` even though `period_effects[X]` (read by position) was correct. Bug was pre-existing on MPD's `fixed_effects=[]` path; the auto-route just made it newly reachable via `absorb=`. Fix: in MPD's fixed_effects expansion at estimators.py:1604, skip entries where `fe == time` — MPD's design already absorbs the time dimension via non-reference period dummies, so the FE-block dummies would be perfectly redundant anyway (NaN'd by solve_ols, dropping nothing useful while corrupting the result surface). Empirical evidence: - Pre-fix: `MPD(absorb=["unit","period"])` -> len(coefs)=34, vcov.shape=(38,38), coefs["period_2"]=NaN. - Post-fix: same call -> len(coefs)=34, vcov.shape=(34,34), coefs["period_2"]=0.345 (finite, matches MPD's event-study fit). Tests: new `test_absorb_hc2_result_surface_invariants_multi_absorb` asserts `len(coefficients) == vcov.shape[0]`, no duplicate names, and finite event-study `period_X` on BOTH the auto-route and the explicit `fixed_effects=` paths (Codex P2: regression coverage for the result- surface contract on the newly reachable path). 11/11 MPD tests pass; 249/249 in the broader sweep (test_estimators.py / test_linalg_hc2_bm.py unchanged). REGISTRY/CHANGELOG: documented the time-FE skip rule for both auto-route and pre-existing `fixed_effects=[]` invocations. Co-Authored-By: Claude Opus 4.7 (1M context) --- CHANGELOG.md | 2 +- diff_diff/estimators.py | 15 +++++++- docs/methodology/REGISTRY.md | 2 +- tests/test_estimators_vcov_type.py | 60 ++++++++++++++++++++++++++++++ 4 files changed, 76 insertions(+), 3 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 6aee904a..896aecb9 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,7 +8,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] ### Added -- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that lies in the span of the intercept and the post-auto-route unit FE dummies (the exact alias depends on the omitted FE reference category under `pd.get_dummies(drop_first=True)`, not just on "the sum of treated-cohort unit dummies"), so `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling. Which specific column is dropped is pivot-order and dummy-coding dependent (in the shipped parity fixture it is a never-treated unit dummy, not the `treated` main effect itself). The per-period interaction coefficients (`treated:period_X`) and `avg_att` are identified and invariant to that choice; parity tests target those rather than the `treated` main effect. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. +- **`MultiPeriodDiD(absorb=..., vcov_type in {"hc2", "hc2_bm"})` now supported** (`diff_diff/estimators.py:1476`). Mirrors the DiD-absorb auto-route shipped earlier in this release: when `absorb=` is paired with `vcov_type in {"hc2","hc2_bm"}`, `MultiPeriodDiD.fit()` promotes the absorb columns to `fixed_effects=` internally so the existing full-dummy-design code path computes the algebraically correct vcov on the event-study design (`treated + period_X dummies + treated:period_X interactions + factor(unit)`). Verified at ~1e-10 vs `lm() + sandwich::vcovHC(type="HC2")` and `lm() + clubSandwich::vcovCR(cluster=1:n, type="CR2")` on a 5-cohort × 5-period event-study fixture (new `tests/test_estimators_vcov_type.py::TestMPDAbsorbedFERParity` against `benchmarks/data/clubsandwich_cr2_golden.json` scenario `mpd_absorbed_fe_did`). HC1/CR1 paths on `absorb=` are unchanged (no leverage term). `TwoWayFixedEffects(vcov_type in {"hc2","hc2_bm"})` rejection remains as a follow-up (different fit-path structure — no `fixed_effects=` equivalent inside TWFE). **Behavioral note (full `MultiPeriodDiDResults` surface change under auto-route):** under the auto-route, the entire returned `MultiPeriodDiDResults` reflects the full-dummy fit rather than the within-transformed fit — `result.coefficients`, `result.vcov`, `result.residuals`, `result.fitted_values`, `result.r_squared` all include the FE-dummy entries / un-demeaned values. `result.period_effects[t].effect` / `.se` / `.p_value` / `.conf_int` and `result.avg_att` / `.avg_se` are invariant to this routing (FWL guarantee). MPD requires a time-invariant ever-treated indicator that lies in the span of the intercept and the post-auto-route unit FE dummies (the exact alias depends on the omitted FE reference category under `pd.get_dummies(drop_first=True)`, not just on "the sum of treated-cohort unit dummies"), so `solve_ols` drops one column from that collinear set under R-style rank-deficiency handling. Which specific column is dropped is pivot-order and dummy-coding dependent (in the shipped parity fixture it is a never-treated unit dummy, not the `treated` main effect itself). The per-period interaction coefficients (`treated:period_X`) and `avg_att` are identified and invariant to that choice; parity tests target those rather than the `treated` main effect. **Survey-design scope (replicate weights):** when `survey_design=` uses replicate weights, the auto-route short-circuits the absorb-refit branch at `estimators.py:1693` and routes through the standard `compute_replicate_vcov` path on the fixed full-dummy design — correct because the design does not depend on replicate weights so no per-replicate refit is needed. **Redundant time-FE skip:** when the routed (or directly-supplied) `fixed_effects` list contains the `time` column, MPD silently skips emitting `