run_es() — diagnostic warning for ambiguous timing = NA (classic TWFE,
staggered = TRUE only): a warning() is now emitted when any unit has
treatment = 1 in at least one row but timing = NA. Such units are silently
absorbed as never-treated controls — correct when NA is intentional, but
dangerous when NA represents genuinely missing treatment timing. The warning
names up to 5 affected unit IDs (when unit is supplied) for quick diagnosis.
@param timing documentation clarified: all staggered estimators (cs, sa,
bjs, twm, flex, classic with staggered = TRUE) share the same NA = never treated convention as did::att_gt() and fixest::sunab(). The roxygen
documentation now states this explicitly and notes the risk of silent
misclassification when NA represents missing data rather than "never treated".
honest_sensitivity() — honest robust inference (Rambachan & Roth 2023)
(R/honest_did.R): sensitivity analysis for event-study / DiD
designs when parallel trends may be violated. Instead of assuming parallel
trends holds exactly, it reports confidence sets for a post-treatment effect
under progressively weaker restrictions on the difference in trends, and a
breakdown value (the largest restriction at which the effect is still
significant):
res <- run_es(df, outcome = y, treatment = treat, time = year, timing = 6,
fe = ~ id + year)
h <- honest_sensitivity(res, type = "relative_magnitude",
Mvec = c(0, 0.5, 1, 1.5, 2))
plot_honest(h)
"relative_magnitude" (\eqn{\Delta^{RM}(\bar M)}) and
"smoothness" (\eqn{\Delta^{SD}(M)}).run_es() es_result (estimator "twfe": classic or
method = "sunab", which now carry the event-study coefficient covariance
via the new es_vcov attribute), or from a raw betahat / sigma pair for
any estimator.plot_honest() / autoplot() draw the "top-down" sensitivity plot;
print.honest_result() reports the breakdown value."Conditional"): \eqn{\Delta^{SD}} (single post period) matches to machine
precision; \eqn{\Delta^{RM}} matches to grid resolution. Heavy numeric
helpers (lpSolveAPI, Rglpk, TruncatedNormal, Matrix, pracma) are in
Suggests and gated with requireNamespace().Estimator-core consolidation (behaviour-preserving): factored the
boilerplate duplicated across run_es(), calc_att(), run_did(), and the
five estimator backends into shared helpers in
R/utils-internal.R — .resolve_col() (NSE column
resolution), .add_ci_columns() (CI construction), .validate_panel_cols(),
.compute_cohort_sizes(), .model_vcov_full(), and .aggregate_iw(). The
calc_att() aggregation loops were vectorised. All existing tests pass
unchanged.
run_did() — basic TWFE DiD estimator (R/run_did.R):
New function for classic two-way fixed-effects difference-in-differences.
Accepts a pre-built binary treatment indicator D_it and returns a
did_result S3 object with full modelsummary / tinytable compatibility:
df$D <- as.integer(df$treated & df$year >= 2006)
result <- run_did(df, outcome = y, treatment = D, fe = ~ id + year)
print(result)
modelsummary::modelsummary(result) # tinytable output (modelsummary >= 2.0)
run_es() and calc_att() — specify variables
by name, not by formula construction.fe auto-inferred as ~ unit + time when both unit and time are
supplied and fe = NULL.tidy.did_result() delegates to broom::tidy.fixest() for full
regression table (all regressors); glance.did_result() delegates to
broom::glance.fixest() for model-level stats (within.r.squared,
nobs, AIC, etc.).print.did_result() shows a clean header: estimator, sample sizes,
FE spec, VCOV type.log(y) supported.modelsummary and tinytable added to Suggests.calc_att() — ATT aggregation function (R/calc_att.R): New
function that separates "calculating ATT" from "running an event study".
While run_es() produces a full dynamic event-study curve (effects by
relative time), calc_att() yields aggregated ATT estimates with a clear,
minimal interface:
calc_att(data, outcome = y, time = year, timing = g, unit = id,
estimator = "cs",
aggregation = "simple") # overall ATT
calc_att(..., aggregation = "by_cohort") # ATT per cohort
calc_att(..., aggregation = "by_time") # ATT per calendar period
"cs" (Callaway-Sant'Anna 2021) and "bjs" (Borusyak
et al. 2024). SA, TWM, FLEX are event-study-only estimators and are not
supported (use run_es()).did::aggte() for point estimates:
"simple": exposure-weighted average over all post-treatment (g,t) pairs."by_cohort": mean ATT(g,t) over post-treatment periods per cohort."by_time": cohort-size-weighted ATT per calendar period.did::aggte SEs, which use influence-function variance
accounting for within-cohort time correlation. Cluster-robust SE for BJS
aggregations planned for v0.12.att_result S3 object (data.frame subclass) with columns
group, estimate, std.error, statistic, p.value,
conf_low_XX/conf_high_XX. Attributes: aggregation, estimator,
N, N_units, N_treated, control_group (CS only), att_gt (raw
CS ATT(g,t) table), tau_it (BJS unit-time effects table).print.att_result() S3 method provides a concise summary.tests/testthat/test-att.R:
did::aggte() at tolerance 1e-4mean(tau_it$tau_hat) at tolerance 1e-12conf.level support verifiedrun_es.R refactoring: Extracted .make_tidy() and .es_finalize() private
helpers that centralise the repeated tidy-construction / baseline-row / range-filter /
attr-stamping pattern that was previously duplicated across the cs, sa, bjs,
twm, and flex branches (~360 lines removed). No API change; all 593 tests pass.inst/profile/benchmark.R with profile_estimator()
and benchmark_all() utilities. Run interactively to generate flamegraph HTML and
wall-clock summary across estimators before starting Phase 3 Rcpp work.src/compute_att_gt.cpp (shipped in v0.8.1 but
missing from the Rcpp Acceleration Roadmap), and updated Package Structure to list all
13 R source files, 5 Rcpp files, and 16 test files that now exist.test-twm.R Test 16 — trends=TRUE with exactly 2 pre-treatment periods succeeds
without warning (boundary condition of the ≥ 2 requirement).test-twm.R Test 17 — explicit lead_range/lag_range correctly trims TWM output
and preserves estimates in the shared window (exercises .es_finalize() filter).test-flex.R Test 12 — explicit lead_range/lag_range trims FLEX output correctly.test-flex.R Test 13 — FLEX succeeds with an all-treated panel (no never-treated
groups; N_nevertreated == 0).test-sa.R — explicit lead_range/lag_range trims SA output and preserves
estimates within the shared window.run_es(estimator = "twm", trends = TRUE): Cohort-specific linear trend
detrending (Wooldridge 2025, Section 8). Adds d_g * t regressors so each
cohort's counterfactual trend can deviate linearly from the common time
trend. Implemented as post-treatment-only cells + trend columns; requires
= 2 pre-treatment periods per cohort. Output shows
relative_time >= 0only (pre-trend testing is deferred to the no-trend model).
run_es(estimator = "twm", covariates = ~ x1 + x2): Full Wooldridge (2025)
Procedure 5.1 covariate interactions. Adds treatment-cell × cohort-demeaned
covariate terms (ẋ_{ig} = x_i - x̄_g) and i(time, x_j) conditional
parallel-trends controls (Eqs. 5.2-5.3).run_es(estimator = "flex", covariates = ~ x1 + x2): Full Deb et al. (2024)
Eq. 3.1 covariate interactions for RCS data. Cell-level centering
X_{i,t} - X̄_{g,t} (Eq. 2.11) with i(time, x_j) and i(group, x_j)
conditional PT controls.src/indicator_matrix.cpp — build_indicator_matrix_cpp(): Replaces the
pure-R nested for-loop that fills the 0/1 indicator matrix in the SA, TWM,
and FLEX estimators. Single shared Rcpp utility, O(N×K) integer fill with
no R allocation overhead per column.src/iw_aggregation.cpp — aggregate_iw_cpp(): Replaces the R event-time
aggregation loop and the t(w) %*% V_sub %*% w quadratic-form VCOV step in
SA, TWM, and FLEX. Uses RcppArmadillo arma::mat::submat() for O(1) VCOV
block extraction and arma::as_scalar(w.t() * V_sub * w) for the quadratic
form.src/cov_demeaning.cpp — build_cov_interactions_cpp(): Replaces the
pure-R covariate demeaning loops in TWM (cohort-level centering) and FLEX
(group×time cell-level centering) with a single std::unordered_map-based
O(N) grouping pass followed by column-wise mean subtraction, then builds the
N × (K×p) treatment-cell × centred-covariate interaction matrix in C++.src/bootstrap_cs.cpp — bootstrap_cs_cpp(): Replaces the R-level
apply/sweep calls in the CS multiplier bootstrap with RcppArmadillo DGEMM
for both the pilot and main stages. Uses R::runif() for Mammen weight
generation (column-major fill order matches R's matrix()) so results are
numerically identical to the pure-R path under the same set.seed().src/Makevars and src/Makevars.win linking $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS) for cross-platform RcppArmadillo support.RcppArmadillo added to LinkingTo in DESCRIPTION.run_es(estimator = "twm"): Wooldridge (2025) Two-Way
Mundlak (TWM) estimator for panel data. Implements
Procedure 5.1 via POLS on cohort x calendar-time
treatment cells with two-way FE. Algebraically
identical to estimator = "sa" in the no-covariate
base case (verified numerically at tolerance 1e-10).
The trends = TRUE extension (cohort-specific linear
trend detrending, Wooldridge 2025 Section 8) is
reserved for v0.9.1.run_es(estimator = "flex", group = <group_id>):
Deb, Norton, Wooldridge & Zabel (2024) FLEX estimator
for repeated cross-section (RCS) data. Uses group
x calendar-time OLS with group+time FE (no unit FE).
Algebraically equivalent to the multi-step imputation
estimator (Proposition 2.1). New group argument
identifies the treatment group (R_{ig}) each
individual belongs to. Default clustering is by group.
Covariate interactions deferred to v0.9.1.papers/Deb et al. (2024).txt,
papers/Woodridge (2025).txt.run_es(estimator = "cs"): Callaway-Sant'Anna (2021)
estimator. Computes group-time ATT(g,t) via the
unconditional DiD estimand (eq. 2.8) with
never-treated or not-yet-treated control groups.
Supports control_group = "nevertreated" (default)
or "notyettreated". The full ATT(g,t) matrix is
stored as the att_gt attribute of the result.run_es(estimator = "sa"): Sun-Abraham (2021)
interaction-weighted estimator. Aggregates cohort x
relative-time interactions by cohort share weights.
Numerically identical to fixest::sunab() to
machine precision.run_es(estimator = "bjs"): Borusyak, Jaravel &
Spiess (2024) imputation estimator. Fits TWFE on
untreated observations only, imputes counterfactuals
for treated observations, and averages treatment
effects by horizon. Handles singleton unit fixed
effects via closed-form recovery.run_es(bootstrap = TRUE): Multiplier bootstrap for
simultaneous confidence bands (Algorithm 1, Callaway
& Sant'Anna 2021). Adds conf_low_sim /
conf_high_sim columns to the result and stores the
(g,t)-level bootstrap object as
attr(result, "bootstrap"). Controlled by B
(draws, default 999) and boot_seed arguments.
Available only when estimator = "cs".plot_att_gt(): Visualize the full ATT(g,t) matrix
from CS results as a heatmap (type = "heatmap") or
cohort-faceted time series (type = "facet").
Requires estimator = "cs" result as input.
When bootstrap = TRUE was used, the heatmap adds
a subtitle with the simultaneous critical value and
open-diamond markers for simultaneously significant
cells; the facet adds a lighter simultaneous CI
ribbon.plot_es(show_simultaneous = TRUE): Overlays the
simultaneous bootstrap CI (lighter band, alpha 0.15)
alongside the pointwise CI (alpha 0.3) with a
two-entry legend. Errors informatively if bootstrap
was not run.plot_es_interactive(show_simultaneous = TRUE):
Adds a second lighter ribbon trace for the
simultaneous CI and extends the hover tooltip with
simultaneous CI bounds.run_es() bootstrap path: base::merge() silently
dropped all custom attributes (e.g. att_gt,
N_units, lead_range) from the es_result object.
Attributes are now saved and restored around the
merge, so plot_att_gt() works correctly after a
bootstrap = TRUE run.did and didimputation to Suggests for
numerical agreement tests against reference
implementations.cluster in run_es() did not produce clustered standard errorsvcov = "HC1" was being applied after model estimation, overriding the clustered SE from feols()cluster is specified and vcov is left at its default ("HC1"), the model's clustered standard errors are usedcluster is set, pass vcov = "HC1" together with cluster = NULLclassic and sunab methodsrun_es() results and equivalent direct fixest::feols() callsLead/lag range filtering now enforced:
lead_range and lag_range parameters were ignored, causing all estimated coefficients to be returned regardless of specified ranges[-lead_range, lag_range]classic and sunab methodssunab() namespace error resolved:
'sunab' is not an exported object from 'namespace:fixest' error when using method = "sunab"sunab function from fixest packageBaseline row now included in sunab results:
baseline parameter (default: -1) now applies to both classic and sunab methodsestimate = 0, std.error = 0, and is_baseline = TRUEbaseline is used for both methodsImproved term column formatting:
"-9", "0", "3")"year::-9") for better readabilityModernized code with native pipe operator:
%>%) with base R pipe (|>) throughout the package|>)relative_time was NA for all coefficients except baseline when using non-staggered treatment timingfixest::i() now correctly corresponds to the baseline parametertiming = 5 and baseline = -1, period 4 (not period 5) is now used as the referenceplot_es_interactive() for creating interactive event study plotses_result objects from run_es()sapply() approach for significant performance gainsplotly package (optional, suggested dependency)run_es() and plot_es() now support returning and visualizing multiple confidence levels (e.g., 90%, 95%, 99%) in a single analysis.plot_es() adds ci_level selection, more theme options (theme_style), and improved ribbon/error bar display.NA in timing):
timing are now retained as never-treated controls (staggered only).weights input:
~ popwt), bare names (popwt), or character strings ("popwt").cluster input handling:
unit is supplied without time_transform = TRUE.Added support for staggered treatment timing:
staggered = TRUE option allows timing to vary by unit (e.g., treatment year column).NA in timing are safely retained as untreated.Added support for observation weights:
weights argument (e.g., ~ popwt) to run weighted regressions.Automatic lead/lag range detection:
lead_range or lag_range is NULL, the function computes the maximum feasible range from the data.unit is specified without time_transform = TRUE.plot_es():
ggplot2::scale_x_continuous() with integer breaks spaced by 1, aligned to the relative_time range.Date:
time_transform = TRUE to automatically convert the time variable into a unit-level sequential index (1, 2, 3, ...) for event study estimation.unit argument to specify the panel unit identifier required when time_transform = TRUE.Date class in the time variable and converts it automatically to numeric if time_transform = FALSE.unit is missing or time is of unsupported type.@examples in the function documentation to include Date-based examples.time_transform usage.README.md to describe irregular time handling and demonstrate new use cases.time_transform, unit handling, and Date conversion edge cases.lead1, lag0) already exist in the dataset to prevent accidental overwriting.lead_range, lag_range, and interval) has fewer than 10 rows, helping users identify overly narrow estimation windows.treatment variable: it is now coerced to logical using as.logical() to support both binary numeric (0/1) and logical (TRUE/FALSE) formats.fe argument (e.g., ~ id + year) were combined using model_formula | fe_text, which caused evaluation errors during tests.as.formula() to ensure compatibility with fixest::feols().run_es():
~ x1 + x2).fe and cluster arguments must now be specified using a one-sided formula (e.g., ~ id + year).cluster is still accepted.fe_var argument now supports additive notation (firm_id + year) instead of character vectors.plot_es() efficiency and documentation.fe notation.This version introduced several enhancements and refinements to improve usability and maintainability.
outcome_var, treated_var, and time_var are now processed using rlang::ensym() for better robustness.fe_var and cluster_var handling improved for more reliable column referencing.plot_es() function:
relative_time, estimate, etc.) are present.baseline handling could lead to incorrect sorting of lead/lag terms.This is the first release of the fixes package, providing tools for estimating and visualizing event study models with fixed effects.
run_es(): A function to estimate event study models using fixest::feols(), generating lead and lag variables automatically.
fe_var as character vector).cluster_var.interval argument.plot_es(): A function to visualize event study results with ggplot2.
type = "ribbon", default).type = "errorbar").fixest::feols().interval argument).c("firm_id", "year"))."state_id").fe_var.