Getting started with staggeredGMM

Overview

staggeredGMM estimates cohort-by-time average treatment effects (CATTs) under staggered treatment adoption by the generalized method of moments, as proposed in

Arora, P. and Bijani, R. (2026). “Estimating Treatment Effects under Staggered Timing and Non-Spherical Errors.” Available at SSRN: https://doi.org/10.2139/ssrn.6558759

Only clean two-by-two difference-in-differences comparisons enter the moment system: a treated cohort against a never-treated group, or against a cohort not yet treated at the relevant period.

This vignette covers how to call the estimator and what it returns. For the methodology, the choice among the three weighting schemes, and the simulation and empirical evidence, see the paper.

Data coding conventions

gmm_staggered() has the following input requirements:

one row per unit-period, with no duplicated unit-period pair;
a numeric outcome (factor, character and logical outcomes are errors);
integer-valued time and cohort columns;
the cohort constant within unit, with 0 for never-treated units;
no period absent from the whole panel;
no unit treated at or before the first observed period.

A missing outcome is fine, and is what an unbalanced panel looks like. A missing unit id, period or cohort is an error, since the row cannot be placed in the panel at all.

Two conventions are worth stating explicitly. First, a unit whose cohort falls after the last observed period is untreated throughout the window, so it is pooled with the never-treated group as a clean control rather than discarded. Second, a unit treated at or before the first observed period has no pre-treatment period at all, so none of its effects is identified; rather than return zeros for those cells and contaminate the aggregate, the estimator stops and names the offending cohort.

A first estimate

sim_panel is a balanced synthetic panel: 60 units over 33 periods, five treated cohorts and a never-treated group.

fit <- gmm_staggered(sim_panel, yname = "y", tname = "year",
                     idname = "unit_id", gname = "cohort")
fit
#> Staggered-adoption GMM estimator
#> Weighting: pooled Toeplitz (GMM-T)
#> 
#>   Units                       60
#>   Periods                     33
#>   Treated cohorts              5
#>   CATTs (identified)     90 (90)
#>   Clean comparisons         1935
#>   Moment-space rank          160
#>   Iterations                   7
#>   Converged                  yes
#> 
#> Treated-observation ATT       -16.6321  (se 0.2522) 
#> Cohort-equal ATT              -15.5701  (se 0.2521)

The never-treated group was detected automatically. never_treated = TRUE requires one and errors if none exists; never_treated = FALSE excludes them.

What is returned

head(fit$catt)
#>    g  t event_time  estimate std_error identified n_comparisons
#> 1 10 10          0 -16.01541 0.3200093       TRUE            45
#> 2 10 11          1 -15.67617 0.3786223       TRUE            45
#> 3 10 12          2 -15.98558 0.4040038       TRUE            45
#> 4 10 13          3 -16.42134 0.4216487       TRUE            36
#> 5 10 14          4 -16.48369 0.4310857       TRUE            36
#> 6 10 15          5 -16.92851 0.4391031       TRUE            36

catt gives one row per cohort-by-time cell, with the event time, the estimate and its standard error, whether the cell is identified, and how many clean comparisons anchor it. coef(), vcov() and confint() work as usual:

head(coef(fit))
#>   g10:t10   g10:t11   g10:t12   g10:t13   g10:t14   g10:t15 
#> -16.01541 -15.67617 -15.98558 -16.42134 -16.48369 -16.92851
head(confint(fit))
#>             2.5 %    97.5 %
#> g10:t10 -16.64262 -15.38821
#> g10:t11 -16.41826 -14.93408
#> g10:t12 -16.77741 -15.19375
#> g10:t13 -17.24775 -15.59492
#> g10:t14 -17.32860 -15.63877
#> g10:t15 -17.78914 -16.06789

The two aggregates

Two scalar summaries are reported, and they answer different questions.

CW, treated-observation weighting. Every treated unit-period gets equal weight, so a cohort observed for longer contributes more total weight. This is the aggregation used throughout the paper, and it coincides with the “simple” aggregation of Callaway and Sant’Anna (2021).
EW, cohort-equal weighting. Each treated cohort gets total weight 1/K, spread equally over its own post-treatment window, so exposure length does not affect a cohort’s influence.

c(CW = fit$aggregate$CW$estimate, EW = fit$aggregate$EW$estimate)
#>        CW        EW 
#> -16.63214 -15.57012

On this design the true values are -16.79375 and -15.82020 respectively; they differ because earlier cohorts are observed for more post-treatment periods. Neither is more correct than the other. For any other target, weight fit$catt$estimate directly.

Choosing a weighting

fit_ht <- gmm_staggered(sim_panel, yname = "y", tname = "year",
                        idname = "unit_id", gname = "cohort",
                        weighting = "cohort_toeplitz")
c(pooled = fit$aggregate$CW$estimate, cohort = fit_ht$aggregate$CW$estimate)
#>    pooled    cohort 
#> -16.63214 -16.61543

All three weightings target the same CATT vector using the same moment conditions. They do not generally return the same numbers: the system is over-identified, so a different weighting matrix gives a different estimate. Which to prefer is a methodological question the paper addresses.

Convergence

unlist(fit$convergence[c("converged", "solve_ok", "termination", "n_iter")])
#>   converged    solve_ok termination      n_iter 
#>      "TRUE"      "TRUE" "converged"         "7"

converged is TRUE only when the largest change in any CATT fell below tol. It is deliberately distinct from solve_ok, which records merely that a weighted solve succeeded; a run can solve successfully at every iteration and still exhaust max_iter without converging. termination says why the iteration stopped.

max_iter defaults to 100 because the iteration can converge linearly at a slow rate.

On this design the unrestricted weighting cannot complete a weighted step at all:

fit_u <- gmm_staggered(sim_panel, yname = "y", tname = "year",
                       idname = "unit_id", gname = "cohort",
                       weighting = "unrestricted")
#> Warning: No weighted GMM step could be completed (termination:
#> singular_normal_equations). The reported estimates are the identity-weighted
#> seed, not a unrestricted-weighted estimator, and standard errors use the
#> sandwich fallback. The estimated weighting matrix supports 60 independent
#> moment direction(s), but 90 effect(s) must be identified. A cohort-specific
#> covariance carries more parameters than a short panel with few units per cohort
#> can support; a more parsimonious `weighting`, or more units, is needed for this
#> design.

The estimate returned in that case is the identity-weighted seed, not an unrestricted-weighted GMM estimator, which is why it always warns. Failing to converge after a weighted step did run is a different matter, and is common enough under this weighting that it does not warn.

Unbalanced panels

Missing outcomes are supported. Group-period means are taken over the units observed at that period; autocovariances are divided by the number of residual pairs actually observed at each lag; and the variance of a group-period mean accounts for the overlap between the unit sets contributing at each pair of periods. All three reduce exactly to the complete-data formulas when nothing is missing, so there is one code path rather than two.

gappy <- sim_panel
set.seed(1)
gappy$y[sample(nrow(gappy), 100)] <- NA
fit_gap <- gmm_staggered(gappy, yname = "y", tname = "year",
                         idname = "unit_id", gname = "cohort")
c(balanced = fit$aggregate$CW$estimate,
  unbalanced = fit_gap$aggregate$CW$estimate)
#>   balanced unbalanced 
#>  -16.63214  -16.63035

This is valid when missingness is independent of the outcome given group and period. Selection on outcomes is not addressed by any weighting scheme.

A comparison referencing a group-period at which no unit is observed cannot be formed and is dropped, with a warning giving the count.

Partial identification

A cohort-by-time effect with no clean comparison is not estimable. Such cells are reported as NA with identified = FALSE, never as an estimated zero, and are excluded from the parameter vector so that one unestimable cell does not disable the efficient weighting everywhere else.

Where an unidentified cell carries positive aggregation weight, the aggregate itself is NA, and an identified-subset aggregate renormalised over the estimable cells is reported alongside it:

names(fit$aggregate$CW)
#> [1] "estimate"                "std_error"              
#> [3] "estimate_identified"     "std_error_identified"   
#> [5] "identified_weight_share" "weights"

Baseline covariates

covar takes a character vector of baseline, time-invariant covariate columns and applies the outcome-regression adjustment of Section 4.5 of the paper, so that parallel trends need hold only conditional on them.

fit_cov <- gmm_staggered(sim_panel, yname = "y", tname = "year",
                         idname = "unit_id", gname = "cohort",
                         covar = c("x1", "x2"))
fit_cov$aggregate$CW$estimate
#> [1] -16.66011

If a supplied covariate is not in fact constant within unit, each affected unit’s value at its earliest observed period is used and a warning names the units. Earliest by time, not first by row: the answer does not depend on how the data were sorted.

Where the adjustment cannot be identified for some control cohort and pair of periods – too few complete-case control units, or a rank-deficient fit – those comparisons alone revert to the unconditional form, and the count is reported in fit$controls$cov_n_fallback.

Testing the identifying assumptions

The over-identifying restrictions are the pre-treatment placebo moments, whose expectation is zero only under parallel trends and no anticipation. gmm_j_test() tests them jointly, weighting them by the estimated serial-correlation model.

gmm_j_test(fit)
#> Specification test for parallel trends and no anticipation
#> Hansen J on the pre-treatment placebo restrictions
#> 
#>   Restriction set   local pre-window (3 period(s) before adoption)
#>   Base period       9
#>   Weighting         pooled Toeplitz (GMM-T)
#>   Moments           14
#> 
#>   J = 9.0396, df = 14, p = 0.8285

By default the test uses a local pre-window of three periods before each cohort’s adoption. Using every placebo restriction is also available, but a large number of restrictions degrades the chi-squared approximation when the number of units is small:

gmm_j_test(fit, type = "full")
#> Specification test for parallel trends and no anticipation
#> Hansen J on the pre-treatment placebo restrictions
#> 
#>   Restriction set   all pre-treatment placebo moments
#>   Base period       9
#>   Weighting         pooled Toeplitz (GMM-T)
#>   Moments           70
#> 
#>   J = 73.1373, df = 70, p = 0.3754

The test requires a never-treated group, since the event-study basis it uses differences every cohort against a common never-treated reference.

A real design, and a coding pitfall

beck_banks is a panel of 49 US states over 1976–2006 on bank branch deregulation and income inequality. It is also a worked example of the most common way a staggered design goes wrong.

gmm_staggered(beck_banks, yname = "ln_gini", tname = "wrkyr",
              idname = "state", gname = "branch_reform")
#> Error:
#> ! 13 unit(s) in cohort(s) 1960, 1970, 1975, 1976 are treated at or before the first observed period (1976). Such units have no pre-treatment period, so none of their cohort-by-time effects is identified, and including them would silently contaminate the aggregate.

Thirteen states deregulated at or before 1976, the first observed year, ten of them coded 1960. They are treated throughout the window and have no pre-treatment period, so none of their effects is identified. Silently including them would return a plausible-looking aggregate contaminated by cells that were never estimated. Restrict to states adopting inside the window:

dat <- beck_banks[beck_banks$branch_reform > 1976, ]
fit_beck <- gmm_staggered(dat, yname = "ln_gini", tname = "wrkyr",
                          idname = "state", gname = "branch_reform")
#> Warning: 144 of 380 cohort-by-time effect(s) have no clean comparison and are
#> reported as NA, not as an estimated zero. Affected: (g=1977, t=1999), (g=1977,
#> t=2000), (g=1977, t=2001), (g=1977, t=2002), (g=1977, t=2003), (g=1977,
#> t=2004), (g=1977, t=2005), (g=1977, t=2006), ... (144 total). A never-treated
#> control group, if available, would identify more.
#> Warning: The efficient GMM reweighting did not converge (termination:
#> singular_normal_equations) after 13 of at most 100 iteration(s), with a final
#> change of 5.13e+01. Estimates are the last valid iterate. Consider a larger
#> `max_iter`, a looser `tol`, or a more parsimonious `weighting`.
fit_beck
#> Staggered-adoption GMM estimator
#> Weighting: pooled Toeplitz (GMM-T)
#> 
#>   Units                       36
#>   Periods                     31
#>   Treated cohorts             18
#>   CATTs (identified)    380 (236)
#>   Clean comparisons         5829
#>   Moment-space rank          374
#>   Iterations                  13
#>   Converged                   no
#>   Termination           singular_normal_equations
#> 
#> Treated-observation ATT           NA  (identified subset: 0.0390, se 2.0518, 60.8% of weight) 
#> Cohort-equal ATT                  NA  (identified subset: 0.0364, se 2.8117, 57.8% of weight)

This design has no never-treated group either, so identification rests entirely on not-yet-treated controls. The consequence is visible in the output: the last cohort to adopt has no later cohort to serve as its control, so none of its effects is identified, and neither is any cell at a period after the final adoption year. Of the 380 cohort-by-time cells, 144 are unestimable.

table(fit_beck$catt$identified)
#> 
#> FALSE  TRUE 
#>   144   236
fit_beck$aggregate$CW$identified_weight_share
#> [1] 0.6081633

Because unidentified cells carry positive weight, the aggregate itself is NA and the identified-subset figure is reported instead. That is the intended behaviour: an aggregate over a set of cells more than a third of which were never estimated would not mean what its name suggests. If a never-treated group is available in your own application, supplying it identifies far more of the surface.