Preserve pyfixest as the IW event-study engine

Implement eventreg as a thin interaction-weighted layer over pyfixest. The first-stage interacted regression and second-stage cohort-share regressions run through pyfixest so native solver, vcov, demeaning backend, fixed-effect, and split behavior stay visible. IW aggregation adds the linear-combination and delta-method pieces needed for Stata/Julia eventstudyinteract parity.

Constraint: pyfixest event_study saturated only exposes unit/time fixed effects, so this package keeps explicit rel_varlist plus arbitrary fixest formula fixed effects.

Constraint: pyfixest does not expose multivariate cross-equation covariance for share regressions; cross blocks reuse pyfixest bread, scores, and CRV1 meat internals.

Rejected: Reimplement high-dimensional demeaning in this package | would lose pyfixest performance and backend choices.

Rejected: Claim R/JL vcov equality for the single-cohort saturated fixture | their saturated/small-sample paths differ or are collinear, so the script reports explicit skips.

Confidence: high

Scope-risk: broad

Directive: Keep first-stage and share regressions pyfixest-backed; do not add a separate regression engine without a benchmark and parity reason.

Tested: .venv/bin/python -m ruff check .

Tested: .venv/bin/python -m ruff format --check .

Tested: .venv/bin/python -m pytest (13 passed)

Tested: .venv/bin/python -m compileall src tests scripts

Tested: .venv/bin/python scripts/check_external_parity.py (Python/R/Stata/Julia coefficient parity; Stata V_iw parity)

Not-tested: Large nlswork benchmark parity and GPU/JAX/CuPy backend runtime.

Author: xiaobaaaa

.github/workflows/ci.yml

@@ -13,7 +13,7 @@ jobs:
     strategy:
       fail-fast: false
       matrix:
-        python-version: ["3.9", "3.10", "3.11", "3.12"]
+        python-version: ["3.10", "3.11", "3.12"]
 
     steps:
       - uses: actions/checkout@v4

README.md

@@ -1,47 +1,92 @@
 # EventStudyInteracts.py
 
-Python project scaffold for a `pyfixest`-backed implementation of the
-Sun and Abraham interaction-weighted event-study estimator.
-
-This repository is currently in the planning/scaffold stage. The intended
-package import name is `eventstudyinteracts`, and the primary public function
-will be `eventreg(...)`, mirroring the workflow of
-[`EventStudyInteracts.jl`](https://github.com/xiaobaaaa/EventStudyInteracts.jl)
-while preserving the native features of
+`EventStudyInteracts.py` is a Python implementation of the Sun and Abraham
+interaction-weighted event-study estimator, built as a thin layer over
 [`pyfixest`](https://github.com/py-econometrics/pyfixest).
 
-## Design Goal
+The estimator keeps `pyfixest.feols()` as the regression engine. First-stage
+coefficients and covariance matrices come from pyfixest, including its fast
+fixed-effect demeaning backends and solver choices. The second-stage cohort-share
+regressions are also fit with pyfixest; the package only adds the
+interaction-weighted aggregation and the cross-equation delta-method covariance
+that pyfixest does not expose directly.
 
-The package should feel like a thin, Python-native layer over `pyfixest`, not a
-separate regression engine. In particular, the implementation plan keeps
-`pyfixest.feols()` as the estimation backend and forwards important backend
-options such as:
+## Installation
 
-- `vcov`, `vcov_kwargs`, `ssc`, and cluster-robust inference settings
-- `weights` and `weights_type`
-- `fixef_rm`, `fixef_tol`, `fixef_maxiter`, and `collin_tol`
+```bash
+python -m pip install -e ".[dev]"
+```
+
+The package requires Python 3.10+ and `pyfixest>=0.50.1`.
+
+## Basic Usage
+
+```python
+from eventstudyinteracts import eventreg
+
+fit = eventreg(
+    "ln_wage ~ age + tenure | idcode + year",
+    data=df,
+    rel_varlist=["g_4", "g_3", "g_2", "g0", "g1", "g2", "g3"],
+    cohort="first_treat",
+    control_cohort="never_treat",
+    vcov={"CRV1": "idcode"},
+    demeaner_backend="rust-cg",
+)
+
+fit.coef()
+fit.se()
+fit.vcov()
+fit.tidy()
+```
+
+`control_cohort` is a boolean column identifying never-treated/control rows.
+`rel_varlist` should contain the relative-time indicator columns to aggregate,
+in the output order you want.
+
+## pyfixest Passthrough
+
+`eventreg(...)` mirrors pyfixest naming where possible and forwards the important
+performance and inference options:
+
+- `vcov`, `vcov_kwargs`, `ssc`
+- `weights`, `weights_type`
+- `fixef_rm`, `fixef_tol`, `fixef_maxiter`, `collin_tol`
 - `solver`
-- `demeaner_backend` with `numba`, `rust`, `jax`, `cupy`, `cupy32`, `cupy64`,
-  and `scipy`
-- `copy_data`, `store_data`, `lean`, and `use_compression`
+- `demeaner_backend`, including `numba`, `rust`, `rust-cg`, `jax`, `cupy`,
+  `cupy32`, `cupy64`, and `scipy`
+- `copy_data`, `store_data`, `lean`, `use_compression`
 - `split` and `fsplit`
 
-## Current Layout
+The IW share covariance currently supports the EventStudyInteracts.jl parity
+modes: `iid`/`simple`, `hetero`/`HC1`, and one- or two-way `{"CRV1": ...}`.
 
-- `src/eventstudyinteracts/`: future package source
+## Project Layout
+
+- `src/eventstudyinteracts/`: package source
+- `tests/`: unit and estimator regression tests
+- `scripts/check_external_parity.py`: optional local parity harness for Python,
+  R fixest, Stata `eventstudyinteract`, and Julia `EventStudyInteracts.jl`
 - `docs/requirements.md`: detailed requirements
-- `docs/design/api.md`: proposed Python API and result object design
-- `docs/design/upstream-analysis.md`: notes from `EventStudyInteracts.jl` and
-  `pyfixest`
+- `docs/design/api.md`: API and result object design
+- `docs/design/upstream-analysis.md`: upstream implementation notes
 - `docs/implementation-plan.md`: milestone plan
-- `docs/workflows/`: development, maintenance, GitHub, and release workflows
-- `tests/`: test suite placeholder
-- `references/`: external parity baselines
-- `benchmarks/`: performance scripts and outputs
-- `.github/workflows/ci.yml`: initial CI skeleton
+- `.github/workflows/ci.yml`: lint, format, test, and compile checks
+
+## Verification
+
+Current local checks:
 
-## Status
+```bash
+python -m ruff check src tests scripts
+python -m ruff format --check src tests scripts
+python -m pytest
+python -m compileall src tests
+python scripts/check_external_parity.py
+```
 
-No estimator implementation is claimed yet. The next implementation step is to
-add `eventreg(...)` with tests that first lock the single-cohort equivalence
-case before implementing the full multi-cohort IW covariance.
+The external parity script compares coefficients against R `fixest`, Stata
+`eventstudyinteract`, and Julia `EventStudyInteracts.jl` when those tools are
+available. It also compares a single-cohort IW covariance matrix against Stata;
+Julia/R saturated paths use different small-sample conventions in that fixture
+and are reported explicitly as skipped for vcov.

docs/design/api.md

@@ -39,8 +39,8 @@ def eventreg(
     weights: str | None = None,
     ssc: dict[str, str | bool] | None = None,
     fixef_rm: Literal["singleton", "none"] = "singleton",
-    fixef_tol: float = 1e-8,
-    fixef_maxiter: int = 100_000,
+    fixef_tol: float = 1e-6,
+    fixef_maxiter: int = 10_000,
     collin_tol: float = 1e-9,
     drop_intercept: bool = False,
     copy_data: bool = True,
@@ -95,7 +95,7 @@ fit.interacted_model_
 fit.stage2_models_
 fit.esample_
 fit.demeaner_backend_
-fit.vcov_type_
+fit.vcov_spec
 ```
 
 ### 3.3 `tidy()` 输出
@@ -140,7 +140,22 @@ y ~ __esi_n_g_2_x_1980 + __esi_n_g_2_x_1981 + ... + x1 + x2 | id + year
 - 公式拆分应使用稳健 parser 或严格的 formula section splitter，不能只做
   脆弱字符串拼接。
 
-## 5. Split / Fsplit
+## 5. Stage 2 与 pyfixest
+
+第二阶段 cohort-share 回归使用 `pyfixest.feols` 的多因变量公式：
+
+```text
+share_1980 + share_1981 + ... ~ 0 + n_rel0 + n_rel1 + ...
+```
+
+对角协方差块直接来自每个 pyfixest share 模型的 `_vcov`。由于 pyfixest 当前不
+公开多方程交叉协方差，包装层只使用 `_bread`、`_scores` 和 CRV1 meat 函数拼接
+delta-method 需要的 cross-equation blocks。
+如果 pyfixest 模型暴露了额外系数（例如旧版本或公式差异导致的 intercept），
+share covariance 会先按 `rel_varlist` 对应的内部列名做子集和重排，避免聚合块
+错位。
+
+## 6. Split / Fsplit
 
 `split` 或 `fsplit` 存在时，`eventreg(...)` 应返回
 `EventStudyInteractMulti`：
@@ -152,10 +167,9 @@ fits["east"].tidy()
 fits.summary()
 ```
 
-初版可以先实现单模型路径，但必须在检测到 `split`/`fsplit` 时明确报错，并在
-计划中将其列为保留 pyfixest 特性的必做项。
+当前实现会逐组调用 `eventreg(...)`，并用 `"all"` 键保存 `fsplit` 的全样本结果。
 
-## 6. 与 pyfixest 保持一致的错误策略
+## 7. 与 pyfixest 保持一致的错误策略
 
 - 后端参数非法时，让 `pyfixest` 自身报错。
 - 本项目只在 IW 语义层报错，如缺少 `rel_varlist` 列、没有 treated cohort、

docs/design/upstream-analysis.md

@@ -73,7 +73,7 @@ Julia 版 `eventfit.jl` 的关键路径：
 
 - GitHub: <https://github.com/py-econometrics/pyfixest>
 - Docs: <https://pyfixest.org/>
-- PyPI stable package inspected: `pyfixest==0.40.1`
+- Local package inspected: `pyfixest==0.50.1`
 
 ### 2.1 `feols` 是本项目核心后端
 
@@ -88,8 +88,8 @@ feols(
     weights=None,
     ssc=None,
     fixef_rm="singleton",
-    fixef_tol=1e-08,
-    fixef_maxiter=100_000,
+    fixef_tol=1e-06,
+    fixef_maxiter=10_000,
     collin_tol=1e-09,
     drop_intercept=False,
     copy_data=True,
@@ -113,6 +113,7 @@ feols(
 
 - `numba`：默认 CPU alternating projections。
 - `rust`：Rust 实现的 CPU alternating projections。
+- `rust-cg`：Rust conjugate-gradient-schwarz 后端，适合困难稀疏固定效应。
 - `jax`：JAX CPU/GPU alternating projections，需要 `jax`/`jaxlib`。
 - `cupy` / `cupy64`：CuPy GPU float64 FWL sparse backend。
 - `cupy32`：CuPy GPU float32 FWL sparse backend。
@@ -129,12 +130,20 @@ feols(
 - `estimator="did2s"`
 - `estimator="saturated"`
 
+官方 reference 中 `event_study` 的 `estimator` 参数列出 `"did2s"`、
+`"twfe"` 和 `"saturated"`，并示例调用 `fit_twfe_saturated.aggregate()`。
+本地 `pyfixest.did.saturated_twfe.SaturatedEventStudy.aggregate()` 的核心做法是
+构造线性组合向量 `R`，再用 `R @ model._vcov @ R` 计算聚合标准误。
+
 但该接口按 `yname/idname/tname/gname` 组织，不是
 `rel_varlist + cohort + control_cohort + formula` 的 Stata/Julia 工作流。因此：
 
 - 本项目应提供显式 `rel_varlist` 工作流。
 - 后续可以用 `pyfixest.event_study(..., estimator="saturated")` 做 sanity check，
   但不能直接替代 `eventreg(...)`。
+- 当前 Python 实现已保存 `aggregation_matrix_`，即与 saturated `aggregate()` 同类
+  的线性组合矩阵；第一阶段聚合协方差走 `R @ V @ R.T`，再叠加 Stata/Julia
+  `eventstudyinteract` 要求的 cohort-share 不确定性。
 
 ## 3. 初步风险
 

docs/implementation-plan.md

@@ -2,7 +2,7 @@
 
 ## Phase 0: 项目骨架
 
-状态：已开始。
+状态：已完成。
 
 交付物：
 
@@ -16,10 +16,12 @@
 
 - `python -m compileall src` 通过。
 - `pytest` 最小导入测试通过。
-- 文档清楚标注尚未实现 estimator。
+- 文档清楚标注已实现范围和剩余外部 parity 缺口。
 
 ## Phase 1: API 与结果对象骨架
 
+状态：已完成。
+
 任务：
 
 - 新增 `eventstudyinteracts.api.eventreg`。
@@ -34,6 +36,8 @@
 
 ## Phase 2: 数据准备与第一阶段 pyfixest 回归
 
+状态：已完成。
+
 任务：
 
 - 实现安全临时列前缀。
@@ -50,6 +54,8 @@
 
 ## Phase 3: 单 cohort 等价测试
 
+状态：已完成基础路径。
+
 任务：
 
 - 构造 deterministic panel fixture。
@@ -65,10 +71,12 @@
 
 ## Phase 4: 多 cohort IW 权重与协方差
 
+状态：已完成 MVP。
+
 任务：
 
 - 实现 treated sample cohort-share 估计。
-- 实现 share 协方差：
+- 实现 pyfixest-backed share 协方差：
   - iid
   - hetero/HC1
   - one-way/two-way CRV1
@@ -83,13 +91,15 @@
 
 ## Phase 5: pyfixest 特性完整化
 
+状态：部分完成。
+
 任务：
 
-- `split`/`fsplit` 返回 `EventStudyInteractMulti`。
+- `split`/`fsplit` 返回 `EventStudyInteractMulti`。已完成。
 - `use_compression=True` 路径验证。
 - `lean`、`store_data`、`copy_data` 行为验证。
 - `demeaner_backend` 覆盖测试：
-  - 必测：`numba`、`rust`、`scipy`
+  - 必测：`numba`、`rust-cg`/`rust`、`scipy`
   - 可选环境：`jax`、`cupy*`
 
 验收标准：
@@ -99,11 +109,16 @@
 
 ## Phase 6: 外部基准与文档示例
 
+状态：进行中。
+
 任务：
 
 - 将 `EventStudyInteracts.jl` 的 `nlswork` 数据/结果转成 Python reference。
 - 增加 `references/` 中的 baseline 文件。
-- 编写 README usage 示例。
+- 编写 README usage 示例。已完成基础示例。
+- 增加 `scripts/check_external_parity.py`。已完成 Python fixture、R fixest
+  Sun-Abraham parity、Stata `eventstudyinteract` parity、Julia
+  `EventStudyInteracts.jl` parity。
 - 增加 `examples/` 脚本或 notebook。
 
 验收标准：

docs/requirements.md

@@ -21,8 +21,9 @@ interaction-weighted (IW) 估计流程，但底层回归、固定效应吸收、
   `method`、`nthreads`、`tol`、`maxiter`、`weights`、`drop_singletons` 等
   后端参数透传。
 - Python 版应以 `pyfixest.feols(...)` 为唯一高维固定效应回归后端。
-- `pyfixest` 当前稳定版为 `0.40.1`，`feols` 暴露了 `demeaner_backend`：
-  `numba`、`rust`、`jax`、`cupy`、`cupy32`、`cupy64`、`scipy`。
+- 本项目当前以 `pyfixest>=0.50.1` 为基线，`feols` 暴露了
+  `demeaner_backend`：`numba`、`rust`、`rust-cg`、`jax`、`cupy`、
+  `cupy32`、`cupy64`、`scipy`。
 - `pyfixest` 已有 `event_study(..., estimator="saturated")`，但它不是
   Stata/Julia `eventstudyinteract` 风格的显式 `rel_varlist` 工作流。本项目
   不替代它，而是补足更贴近原 Stata/Jl 包的接口。
@@ -46,8 +47,8 @@ eventreg(
     weights=None,
     ssc=None,
     fixef_rm="singleton",
-    fixef_tol=1e-8,
-    fixef_maxiter=100_000,
+    fixef_tol=1e-6,
+    fixef_maxiter=10_000,
     collin_tol=1e-9,
     drop_intercept=False,
     copy_data=True,
@@ -71,7 +72,8 @@ eventreg(
   `"ln_wage ~ controls | idcode + year"`。
 - `rel_varlist` 是已在数据中生成好的相对时间 dummy 列名列表，和 Julia 版
   `rel_varlist::Vector{Symbol}` 对齐。
-- `cohort` 是首次处理时间 cohort 列；从未处理样本应为缺失值。
+- `cohort` 是首次处理时间 cohort 列；控制组的 cohort 值不会进入 treated
+  cohort 列表。
 - `control_cohort` 是控制组二元列，可以表示 never-treated 或 last-treated
   控制组。
 - 所有 `pyfixest.feols` 后端参数必须透传，不能被本包装层吞掉。
@@ -112,7 +114,8 @@ eventreg(
 - 收敛控制：`fixef_tol`、`fixef_maxiter` 必须透传。
 - 固定效应样本处理：`fixef_rm` 必须透传。
 - 推断：`vcov`、`vcov_kwargs`、`ssc` 必须透传至第一阶段；第二阶段
-  cohort-share 协方差应实现与主 `vcov` 兼容的路径。
+  cohort-share 回归也由 `pyfixest.feols` 估计，IW 包装层只补 pyfixest 未暴露的
+  cross-equation delta-method 协方差。
 - 权重：`weights`、`weights_type` 必须透传；第二阶段 share 估计必须明确处理
   加权样本。
 - 性能：`solver`、`use_compression`、`copy_data`、`store_data`、`lean` 必须可用。