Skip to content

Commit 90b524b

Browse files
igerberclaude
andauthored
feat: RegressionDiscontinuity - sharp RD estimator with rdrobust 4.0.0 parity (RDD PR-2) (#684)
* feat: RegressionDiscontinuity - sharp RD estimator with rdrobust parity Estimation-path port + public estimator + docs sweep (PR-2 of 2). Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01QGca52n6H8oDDXALjjrsp4 * fix: run mass-point detection before the bandwidth branch + reject bools on integer knobs rdrobust.R:365-380 runs the mass-point check in rdrobust() itself, ahead of the manual-vs-data-driven bandwidth split, so the warning must fire on manual-h fits too (verified on installed 4.0.0). fit() now performs the detection and the port's rdbwselect-level warning is suppressed on the estimator path to mirror R's single warning per fit. Integer knobs (p/q/nnmatch/bwcheck) now reject bool (int subclass) inputs. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01QGca52n6H8oDDXALjjrsp4 * chore: drop stale mypy ignores after rebase onto #683 (blocking mypy job) The type: ignore[assignment] comments on the beta_p_* dataclass fields are unused under the pinned CI mypy config, which now fails on unused ignores (warn_unused_ignores = true). Co-Authored-By: Claude Fable 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01QGca52n6H8oDDXALjjrsp4 --------- Co-authored-by: Claude Fable 5 <noreply@anthropic.com>
1 parent 8dc3421 commit 90b524b

23 files changed

Lines changed: 2888 additions & 133 deletions

CHANGELOG.md

Lines changed: 33 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -8,6 +8,31 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
88
## [Unreleased]
99

1010
### Added
11+
- **`RegressionDiscontinuity` - sharp regression discontinuity estimation with robust
12+
bias-corrected inference (alias `RDD`).** Local-polynomial sharp RD per Calonico,
13+
Cattaneo & Titiunik (2014), parity-targeting R `rdrobust` 4.0.0 end-to-end: all 10
14+
data-driven bandwidth selectors (`mserd` default, `msetwo`/`msesum`/comb and the
15+
CER-optimal variants), triangular/epanechnikov/uniform kernels, `masspoints`
16+
adjust/check/off, manual `h`/`b`/`rho` with R-exact resolution semantics (including
17+
`rho`-without-`h` applying to selected bandwidths and the unconditional N<20
18+
full-range fallback), and the three-row Conventional / Bias-Corrected / Robust
19+
output. **Canonical binding:** `att`/`se`/`t_stat`/`p_value`/`conf_int` are ONE
20+
coherent row - the robust bias-corrected row (`att = tau_bc`, CI centered on it,
21+
`t_stat == att/se`), preserving the library-wide field identities; rdrobust's
22+
printed headline coefficient is exposed as `att_conventional` with a full inference
23+
row, and `summary()` prints the familiar three-row table. Estimation-path port
24+
(`rdrobust_fit_sharp`: Q_q bias-correction score matrix, conventional/robust NN
25+
sandwiches) validated against a new estimates golden
26+
(`benchmarks/data/rdrobust_estimates_golden.json`, 16 configurations incl. the
27+
Senate anchors) at rtol=1e-9 in `tests/test_rdd_parity.py`; R-free methodology
28+
anchors (CCT 2014 Remark 7 bias-corrected == local-quadratic equivalence at rel
29+
1e-10 across all kernels, invariances, joint-NaN degenerate contracts) in
30+
`tests/test_rdd_methodology.py`; API/validation suite in `tests/test_rdd.py`.
31+
Deviations from R (each labeled in the REGISTRY section): warn-instead-of-silent
32+
NaN drops, warn-and-ignore `b`-without-`h`, fail-closed targeted errors on
33+
degenerate designs, and the canonical-binding note above. Sharp designs only:
34+
fuzzy RD, covariates (CCFT 2019 - review on file), cluster-robust variance,
35+
weights, kink estimands, and rdplot/density diagnostics are documented follow-ups.
1136
- **Internal: mypy enforced at zero errors.** Triaged the 184 pre-existing
1237
`mypy diff_diff` errors to an enforceable zero and added a blocking Mypy job to
1338
the Lint CI workflow (pinned `mypy==2.1.0` + pinned numpy/pandas/scipy for stub
@@ -24,7 +49,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
2449
already-stale ignores were removed. Tightening tracked in TODO.md. No public
2550
API or numerical behavior change.
2651
- **Internal: rdrobust sharp-RD bandwidth-selection port (`diff_diff/_rdrobust_port.py`).**
27-
Step-1 machinery for the upcoming `RegressionDiscontinuity` estimator: a faithful
52+
Step-1 machinery for the `RegressionDiscontinuity` estimator above: a faithful
2853
pure-Python NumPy/SciPy port of R `rdrobust` 4.0.0's `rdbwselect` sharp-RD path (all 10 data-driven
2954
selectors - mserd default, msetwo/msesum/comb and the CER-optimal variants -
3055
triangular/epanechnikov/uniform kernels, `masspoints` adjust/check/off, IK-style
@@ -36,10 +61,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
3661
Cattaneo-Frandsen-Titiunik 2015) validated at rtol=1e-9 across 17 configurations x 10
3762
selectors in `tests/test_rdrobust_port.py` (kernels, p/q orders, deriv, masspoints
3863
modes, stdvars, bwrestrict, scaleregul, nnmatch), including the 2017 Stata Journal
39-
published Senate anchors under `masspoints="off"`. **No public API change** - the module is private;
40-
the estimator, robust bias-corrected inference, and docs surfaces land in the follow-up PR.
41-
Port-level deviations documented in `docs/methodology/REGISTRY.md` (new
42-
RegressionDiscontinuity section stub).
64+
published Senate anchors under `masspoints="off"`. The module is private; the public
65+
estimator, robust bias-corrected inference, and docs surfaces are the
66+
`RegressionDiscontinuity` entry above. Port-level deviations documented in
67+
`docs/methodology/REGISTRY.md` (RegressionDiscontinuity section).
4368
- **Internal: ungated `Lint` CI workflow.** Check-only `ruff check` + `black --check`
4469
(at the versions pinned in the `dev` extra) run on every PR push and on pushes to
4570
main — no `ready-for-ci` label needed; the aggregate `Lint Gate` job is the single
@@ -65,6 +90,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
6590
`benchmarks/data/qte_golden.json`.
6691

6792
### Changed
93+
- `diff_diff/guides/llms-autonomous.txt` no longer lists regression discontinuity as
94+
out of scope: sharp RD routes to `RegressionDiscontinuity`; only fuzzy RD is
95+
referred to external tooling.
6896
- **Internal: repo-wide lint normalization + pinned tooling.** black/ruff/mypy are now
6997
pinned exactly in the `dev` extra (`black==26.3.1`, `ruff==0.15.13`, `mypy==2.1.0`;
7098
the tools require Python >= 3.10 — the library floor stays 3.9); full `black` +

README.md

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -112,6 +112,7 @@ Full guide: `diff_diff.get_llm_guide("practitioner")`.
112112
- [TripleDifference](https://diff-diff.readthedocs.io/en/stable/api/triple_diff.html) - triple difference (DDD) estimator for designs requiring two criteria for treatment eligibility
113113
- [ContinuousDiD](https://diff-diff.readthedocs.io/en/stable/api/continuous_did.html) - Callaway, Goodman-Bacon & Sant'Anna (2024) continuous treatment DiD with dose-response curves
114114
- [HeterogeneousAdoptionDiD](https://diff-diff.readthedocs.io/en/stable/api/had.html) - de Chaisemartin, Ciccia, D'Haultfœuille & Knau (2026) for designs where **no unit remains untreated**; local-linear estimator at the dose support boundary returning Weighted Average Slope (WAS) on Design 1' (`d̲ = 0` / QUG) or `WAS_{d̲}` on Design 1 (`d̲ > 0`, continuous-near-d̲ or mass-point), with a multi-period event-study extension (last-treatment cohort, pointwise CIs). **Panel-only** in this release - repeated cross-sections rejected by the validator. Alias `HAD`.
115+
- [RegressionDiscontinuity](https://diff-diff.readthedocs.io/en/stable/api/regression_discontinuity.html) - Calonico, Cattaneo & Titiunik (2014) sharp regression discontinuity with robust bias-corrected inference and rdrobust-parity bandwidth selection (all 10 selectors, mass-point handling). **Sharp designs only** in this release; canonical `att` is the bias-corrected estimate with a coherent robust CI (rdrobust's printed headline is `att_conventional`). Alias `RDD`.
115116
- [StackedDiD](https://diff-diff.readthedocs.io/en/stable/api/stacked_did.html) - Wing, Freedman & Hollingsworth (2024) stacked DiD with Q-weights and sub-experiments; optional covariate balancing (Ustyuzhanin 2026)
116117
- [EfficientDiD](https://diff-diff.readthedocs.io/en/stable/api/efficient_did.html) - Chen, Sant'Anna & Xie (2025) efficient DiD with optimal weighting for tighter SEs
117118
- [TROP](https://diff-diff.readthedocs.io/en/stable/api/trop.html) - Triply Robust Panel estimator (Athey et al. 2025) with nuclear norm factor adjustment
Lines changed: 136 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,136 @@
1+
# Golden-value generator for the diff-diff sharp-RD ESTIMATION port
2+
# (diff_diff/_rdrobust_port.py::rdrobust_fit_sharp and the public
3+
# RegressionDiscontinuity estimator).
4+
#
5+
# Deliberately a SEPARATE file/JSON from generate_rdrobust_golden.R so the
6+
# bandwidth fixtures reviewed in the machinery PR are never regenerated.
7+
# Same parity target: CRAN rdrobust 4.0.0 (tarball sha256
8+
# 78f0d6b4bdec4091cc8f42f6f1598704747f95926446d3aaee381ea1d613a36f);
9+
# the GitHub 4.1.0-dev tree must NOT be installed when regenerating.
10+
#
11+
# Outputs benchmarks/data/rdrobust_estimates_golden.json: per config the
12+
# embedded inputs (17 significant digits), selected/echoed bandwidths, the
13+
# full three-row (Conventional / Bias-Corrected / Robust) coef/se/z/pv/ci
14+
# block, effective counts, and per-side order-p coefficient vectors.
15+
#
16+
# Run from the repo root: Rscript benchmarks/R/generate_rdrobust_estimates_golden.R
17+
18+
library(rdrobust)
19+
library(jsonlite)
20+
21+
stopifnot(packageVersion("rdrobust") == "4.0.0")
22+
23+
TARBALL_SHA256 <- "78f0d6b4bdec4091cc8f42f6f1598704747f95926446d3aaee381ea1d613a36f"
24+
25+
run_estimate <- function(y, x, c = 0, masspoints = "adjust", kernel = "tri",
26+
p = 1, q = 2, h = NULL, b = NULL, rho = NULL,
27+
level = 95, bwselect = "mserd") {
28+
args <- list(y = y, x = x, c = c, masspoints = masspoints, kernel = kernel,
29+
p = p, q = q, level = level, bwselect = bwselect)
30+
if (!is.null(h)) args$h <- h
31+
if (!is.null(b)) args$b <- b
32+
if (!is.null(rho)) args$rho <- rho
33+
r <- suppressWarnings(do.call(rdrobust, args))
34+
list(
35+
c = c, masspoints = masspoints, kernel = kernel, p = p, q = q,
36+
bwselect = bwselect,
37+
h_in = if (is.null(h)) NA else h,
38+
b_in = if (is.null(b)) NA else b,
39+
rho_in = if (is.null(rho)) NA else rho,
40+
level = level,
41+
h_l = r$bws["h", "left"], h_r = r$bws["h", "right"],
42+
b_l = r$bws["b", "left"], b_r = r$bws["b", "right"],
43+
tau_cl = unname(r$coef[1, 1]), tau_bc = unname(r$coef[2, 1]),
44+
se_cl = unname(r$se[1, 1]), se_rb = unname(r$se[3, 1]),
45+
z = unname(r$z[, 1]), pv = unname(r$pv[, 1]),
46+
ci_lower = unname(r$ci[, 1]), ci_upper = unname(r$ci[, 2]),
47+
N = unname(r$N), N_h = unname(r$N_h), N_b = unname(r$N_b),
48+
bias = unname(r$bias),
49+
beta_p_l = unname(as.vector(r$beta_Y_p_l)),
50+
beta_p_r = unname(as.vector(r$beta_Y_p_r))
51+
)
52+
}
53+
54+
golden <- list()
55+
56+
golden$metadata <- list(
57+
rdrobust_version = as.character(packageVersion("rdrobust")),
58+
rdrobust_tarball_sha256 = TARBALL_SHA256,
59+
seeds = list(dgp_lee_smooth = 42L, dgp_ties_moderate = 123L,
60+
dgp_asymmetric_scaled = 777L),
61+
generator = "benchmarks/R/generate_rdrobust_estimates_golden.R",
62+
algorithm = paste(
63+
"rdrobust() sharp-RD estimation blocks (three-row coef/se/z/pv/ci,",
64+
"counts, per-side beta_p) for the vce='nn' no-covariate path,",
65+
"complementing the bandwidth fixtures in rdrobust_golden.json."
66+
),
67+
r_version = R.version.string
68+
)
69+
70+
# Same seeded DGPs as generate_rdrobust_golden.R (kept in sync by seed).
71+
set.seed(42)
72+
n <- 1000
73+
x1 <- 2 * rbeta(n, 2, 4) - 1
74+
y1 <- 0.48 + 1.27 * x1 + 7.18 * x1^2 + 20.21 * x1^3 + 21.54 * x1^4 +
75+
7.33 * x1^5 + 0.04 * (x1 >= 0) + rnorm(n, sd = 0.1295)
76+
77+
golden$dgp_lee_smooth <- list(
78+
x = x1, y = y1,
79+
configs = list(
80+
default = run_estimate(y1, x1),
81+
manual_h = run_estimate(y1, x1, h = 0.15),
82+
h_rho2 = run_estimate(y1, x1, h = 0.15, rho = 2),
83+
rho_only = run_estimate(y1, x1, rho = 2),
84+
p2q3 = run_estimate(y1, x1, p = 2, q = 3),
85+
p0q1 = run_estimate(y1, x1, p = 0, q = 1),
86+
epa = run_estimate(y1, x1, kernel = "epa"),
87+
uni = run_estimate(y1, x1, kernel = "uni"),
88+
level90 = run_estimate(y1, x1, level = 90),
89+
msetwo = run_estimate(y1, x1, bwselect = "msetwo"),
90+
cercomb2 = run_estimate(y1, x1, bwselect = "cercomb2")
91+
)
92+
)
93+
94+
set.seed(123)
95+
n2 <- 800
96+
x2 <- round(2 * rbeta(n2, 2, 4) - 1, 2)
97+
y2 <- 0.5 + 0.8 * x2 + (x2 >= 0) * 1.0 + rnorm(n2, sd = 0.3)
98+
99+
golden$dgp_ties_moderate <- list(
100+
x = x2, y = y2,
101+
configs = list(
102+
adjust = run_estimate(y2, x2, masspoints = "adjust"),
103+
off = run_estimate(y2, x2, masspoints = "off")
104+
)
105+
)
106+
107+
set.seed(777)
108+
n3 <- 300
109+
x3 <- 40 * rbeta(n3, 5, 2)
110+
y3 <- 2 + 0.05 * x3 - 0.001 * x3^2 + 0.8 * (x3 >= 28) + rnorm(n3, sd = 0.5)
111+
112+
golden$dgp_asymmetric_scaled <- list(
113+
x = x3, y = y3,
114+
configs = list(
115+
default = run_estimate(y3, x3, c = 28)
116+
)
117+
)
118+
119+
senate_path <- "benchmarks/data/rdrobust_senate.csv"
120+
stopifnot(file.exists(senate_path))
121+
senate <- read.csv(senate_path)
122+
ok <- complete.cases(senate$vote, senate$margin)
123+
sv <- senate$vote[ok]
124+
sm <- senate$margin[ok]
125+
126+
golden$senate <- list(
127+
csv = "benchmarks/data/rdrobust_senate.csv",
128+
configs = list(
129+
adjust = run_estimate(sv, sm, masspoints = "adjust"),
130+
off = run_estimate(sv, sm, masspoints = "off")
131+
)
132+
)
133+
134+
out_path <- "benchmarks/data/rdrobust_estimates_golden.json"
135+
write_json(golden, out_path, auto_unbox = TRUE, pretty = TRUE, digits = I(17))
136+
cat("Wrote", out_path, "\n")

benchmarks/data/rdrobust_estimates_golden.json

Lines changed: 517 additions & 0 deletions
Large diffs are not rendered by default.

diff_diff/__init__.py

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -206,6 +206,10 @@
206206
TreatmentDoseShape,
207207
profile_panel,
208208
)
209+
from diff_diff.rdd import (
210+
RegressionDiscontinuity,
211+
RegressionDiscontinuityResults,
212+
)
209213
from diff_diff.results import (
210214
DiDResults,
211215
MultiPeriodDiDResults,
@@ -308,6 +312,7 @@
308312
DCDH = ChaisemartinDHaultfoeuille
309313
HAD = HeterogeneousAdoptionDiD
310314
CiC = ChangesInChanges
315+
RDD = RegressionDiscontinuity
311316

312317
__version__ = "3.7.0"
313318
__all__ = [
@@ -508,6 +513,10 @@
508513
"HeterogeneousAdoptionDiDResults",
509514
"HeterogeneousAdoptionDiDEventStudyResults",
510515
"HAD",
516+
# RegressionDiscontinuity (sharp RD, rdrobust parity)
517+
"RegressionDiscontinuity",
518+
"RegressionDiscontinuityResults",
519+
"RDD",
511520
# HeterogeneousAdoptionDiD pre-test diagnostics (Phase 3)
512521
"qug_test",
513522
"stute_test",

0 commit comments

Comments
 (0)