Standard errors in fetwfe: the variance decomposition and se_type

library(fetwfe)

This vignette documents how fetwfe(), betwfe(), etwfe(), and twfeCovs() compute their standard errors: the tight Gaussian variance used by default under (Ψ-IF), the two variance components \(V_1 + V_2\) the package exposes, the three values se_type can take (including an experimental cluster-robust option), and how the bridge penalty lambda is selected. The four estimators share the same inferential machinery, so the discussion applies to all of them; we use fetwfe() in the running example.

Two companion vignettes cover the rest of the package’s inference toolkit:

  • vignette("simultaneous_cis_vignette") — family-wise simultaneous confidence bands (simultaneousCIs()) and the selection-consistency zero-effect test.
  • vignette("high_dimensional_vignette") — high-dimensional (\(p \geq NT\)) debiased inference (debiasedATT()).

1. The headline: tight Gaussian CIs by default under (Ψ-IF)

By default, the package’s standard errors — the att_se slot on the returned object and the entries of catt_ses — are the tight Gaussian variance from the paper’s Theorem te.asym.norm.thm(c\('\)), valid under the influence-function condition (Ψ-IF) (Faletto 2025, arXiv:2312.05985, Assumption (Ψ-IF) at paper Eq. psi.if.assum). The combined ATT variance has two pieces,

\[ V_1 \;=\; \sigma_\varepsilon^2 \, \boldsymbol{\alpha}^\top \boldsymbol{\Sigma}^{-1} \boldsymbol{\alpha}, \qquad V_2 \;=\; v_\psi(\boldsymbol{\beta}_0, \mathcal{S}), \]

where \(V_1\) is the Kock-2013 regression-coefficient variance contribution and \(V_2\) is the propensity-score variance contribution; the names \(V_1\) and \(V_2\) are the paper’s, and the package exposes both of them by name on the fitted object under result$internal$variance_components. Theorem (c\('\)) shows that under (Ψ-IF) the same-data confidence interval is asymptotically Gaussian with the tight variance \(V_1 + V_2\), identical to the independent-data variance in part (b) of the theorem. Concretely, the package returns

\[ \widehat{\text{SE}}(\hat{\tau}_{\text{ATT}}) \;=\; \sqrt{\widehat{\text{Var}}_1 + \widehat{\text{Var}}_2} \;=\; \sqrt{\tilde{v}_N \,/\, N}, \]

where \(\tilde{v}_N := \hat{v}_N / T\) is the unit-scaled paper-notation variance estimator (paper line 2004), and the Wald confidence interval is

\[ \text{CI}_N(\tilde{v}_N;\,\alpha) \;=\; \left[\, \hat{T}_N \;\pm\; z_{1 - \alpha/2}\,\sqrt{\tilde{v}_N \,/\, N} \,\right] \]

per paper Eq. conf.int.form (line 1259). The same combined variance backs the per-cohort entries of catt_ses and the catt_df confidence-interval columns.

1.1 Which propensity-score estimators satisfy (Ψ-IF)?

Assumption (Ψ-IF) requires the cohort-weight estimator \(\widehat{\boldsymbol{\psi}}_N\) to admit an i.i.d.-across-\(i\) unit-level influence-function representation. It is satisfied by

  • Cohort sample proportions \(\widehat{\pi}_g = N_g / N\) — the package’s default propensity-score estimator. So out of the box the tight Gaussian CI is asymptotically exact.
  • Multinomial-logit MLE for \(P(W = g \mid \boldsymbol{X})\).
  • Any GLM fit by maximum likelihood on \(W \mid \boldsymbol{X}\).
  • Kernel or series regression of \(\mathbf{1}\{W = g\}\) on \(\boldsymbol{X}\).

(See the discussion at paper line 1268: “When (Ψ-IF) holds — which is the case for every standard parametric or semi-parametric propensity-score estimator listed in Remark psi.if.coverage (cohort sample proportions, multinomial logit, any GLM, kernel/series regression) — Theorem te.asym.norm.thm(c\('\)) gives asymptotic Gaussianity with the tight variance \(V_1 + V_2\).”)

The (Ψ-IF) condition fails for Robins-Rotnitzky-style augmented doubly-robust estimators that augment the propensity score with outcome residuals. The package does not currently implement these; when it does, users will be able to opt out of the tight CI via se_type = "conservative" (see the conservative-fallback note in Section 2.3).

1.2 What changed in v1.12.0

Versions \(\le\) 1.11.7 of the package returned the conservative Cauchy-Schwarz upper bound

\[ \widehat{\text{SE}}(\hat{\tau}_{\text{ATT}}) \;=\; \sqrt{\widehat{\text{Var}}_1 + \widehat{\text{Var}}_2 + 2\sqrt{\widehat{\text{Var}}_1\,\widehat{\text{Var}}_2}} \]

as the same-data default (from paper Theorem (c)). Theorem (c\('\)) was added in the paper’s same-data Gaussianity upgrade, and v1.12.0 switches the package default to match: the tight Gaussian variance is now the headline number. The conservative formula remains available as an opt-in fallback (Section 2.3) for any users whose propensity-score estimator violates (Ψ-IF).

2. The two variance pieces and how the package exposes them

The combined ATT variance decomposes as \(V_1 + V_2\), with both pieces estimated by the package and exposed as named slots.

2.1 \(V_1\): regression-coefficient variance

The regression-coefficient piece is

\[ \widehat{\text{Var}}_1 \;=\; \frac{\sigma_\varepsilon^2}{NT}\;\psi_{\text{att}}^\top \widehat G^{-1} \psi_{\text{att}}, \]

where \(\widehat G^{-1}\) is the Gram inverse on the bridge-selected support and \(\psi_{\text{att}}\) encodes the cohort-weighted ATT contrast. The per-cohort SEs (catt_ses) have the same form with a cohort selector \(\psi_g\) in place of \(\psi_{\text{att}}\). This piece is stored on the fitted object as result$internal$variance_components$att_var_1; the paper-notation per-unit-scaled version is result$internal$variance_components$V_1 (with V_1 = N * att_var_1).

2.2 \(V_2\): propensity-score variance

The propensity-score piece accounts for variability from estimating the cohort-membership probabilities \(\widehat{\pi}_g\). It scales like \(1/N\) and is unrelated to the regression residuals. This piece is stored as result$internal$variance_components$att_var_2, with paper-notation per-unit-scaled V_2 = N * att_var_2.

2.3 The combined variance: tilde_v_N and the paper’s catalogue

The unit-scaled total variance is

\[ \tilde{v}_N \;=\; N \cdot \widehat{\text{Var}}(\hat{T}_N), \]

stored as result$internal$variance_components$tilde_v_N. The original (un-scaled) paper notation hat_v_N = T * tilde_v_N is also exposed.

The paper catalogues six variance estimators at line 2006; the package exposes the four that apply to the ATT (the other two apply to the CATT(x), which connects to the on-hold predict() work in issue #33):

Slot Paper Eq. Meaning
tilde_v_N_C v.n.r.t.att.const Fixed-π exact (no propensity-score noise)
tilde_v_N_C_pi_hat v.n.r.t.att.rand Random-π exact (tight Gaussian under (Ψ-IF))
tilde_v_N_C_pi_hat_cons v.n.r.t.att.rand.cons Random-π conservative (Cauchy-Schwarz)
tilde_v_N_cons var.est.kock.wooldridge.subgauss.cons Subgaussian conservative bound

The default att_se value is computed from tilde_v_N_C_pi_hat (via sqrt(tilde_v_N_C_pi_hat / N)); se_type = "conservative" switches the headline to tilde_v_N_C_pi_hat_cons. When indep_counts is supplied (two-sample regime), the two-sample-exact formula tilde_v_N_C_pi_hat is in force regardless of se_type.

Conservative fallback (se_type = "conservative"). If a cohort-weight (propensity-score) estimator violates (Ψ-IF) — the canonical case being a Robins-Rotnitzky-style augmented doubly-robust estimator that injects outcome-residual information into the weights — Theorem (c\('\)) no longer applies, and se_type = "conservative" switches the headline SE to the Cauchy-Schwarz upper bound \(\sqrt{\widehat{\text{Var}}_1 + \widehat{\text{Var}}_2 + 2\sqrt{\widehat{\text{Var}}_1\,\widehat{\text{Var}}_2}}\) (tilde_v_N_C_pi_hat_cons, paper Theorem (c)). The package ships no non-(Ψ-IF) estimator — its default cohort-sample-proportions weights \(\widehat{\pi}_g = N_g/N\) satisfy (Ψ-IF) — so under the default estimator there is no practical reason to set it, and default users never need it. It exists for users who supply their own propensity-score estimates (issue #32) or who want to reproduce the pre-v1.12.0 conservative default; see ?fetwfe.

2.4 A running example

The rest of this vignette uses the following small simulated panel and its default fit; the cluster-robust comparison in Section 3.3 reuses them.

set.seed(2026)
sim_coefs <- genCoefs(G = 3, T = 6, d = 2, density = 0.5, eff_size = 2)
# sim_coefs was built without a seed, so draw from the ambient generator
# (seeded just above); seed = NA selects the ambient generator silently.
sim_data <- simulateData(
  sim_coefs,
  N = 120,
  sig_eps_sq = 1,
  sig_eps_c_sq = 0.5,
  seed = NA
)

res_default <- fetwfe(
  pdata = sim_data$pdata,
  time_var = sim_data$time_var,
  unit_var = sim_data$unit_var,
  treatment = sim_data$treatment,
  covs = sim_data$covs,
  response = sim_data$response,
  sig_eps_sq = sim_data$sig_eps_sq,
  sig_eps_c_sq = sim_data$sig_eps_c_sq
)

3. Experimental: cluster-robust standard errors via se_type = "cluster"

The package’s default standard errors rely on the model-based covariance structure (paper Assumption F1: i.i.d. units, idiosyncratic noise plus a unit random effect). In applied DiD work it is not unusual to suspect one or more violations of F1:

  • Within-unit serial correlation beyond the random effect. The unit random effect \(c_i\) in F1 absorbs a single, time-constant deviation per unit. It does not model serial correlation in the idiosyncratic shocks. Bertrand, Duflo, and Mullainathan (2004) is the classic warning that ignoring within-unit serial correlation can drastically understate DiD standard errors.
  • Heteroskedasticity across units. F1 imposes a single \(\sigma_\varepsilon^2\) across all units.
  • Higher-level clustering. F1 treats units as i.i.d.; if observations within a state-year or industry-year share unobserved shocks, the unit-level F1 variance will understate the true sampling variability.

Starting in version 1.6.0, all four estimators accept an experimental se_type argument with value "cluster":

fetwfe(..., se_type = "cluster")

Setting se_type = "cluster" swaps the model-based regression-coefficient variance \(V_1\) for a unit-clustered Liang-Zeger CR1 sandwich computed on the bridge-selected support. The \(V_2\) piece and the (Ψ-IF)-tight combination \(V_1 + V_2\) are unchanged — only the \(V_1\) formula switches.

3.1 The formula

Let \(\widehat{S}\) be the support selected by the bridge regression, \(X_{\widehat{S}}\) the corresponding design matrix in the coordinate system the regression was solved in (GLS-transformed for ETWFE/twfeCovs, fusion-then-GLS-transformed for FETWFE/BETWFE), and \(\widehat\varepsilon\) the residuals from OLS on that selected support. The cluster-robust variance is

\[ V_{\text{CR}} \;=\; \frac{N}{N-1}\; (X_{\widehat{S}}^\top X_{\widehat{S}})^{-1}\; \left(\sum_{i=1}^N X_{i\cdot\widehat{S}}^\top \widehat\varepsilon_{i\cdot} \widehat\varepsilon_{i\cdot}^\top X_{i\cdot\widehat{S}}\right)\; (X_{\widehat{S}}^\top X_{\widehat{S}})^{-1}, \]

with units \(i = 1, \dots, N\) as clusters and an \(N/(N-1)\) small-sample adjustment (matching sandwich::vcovCL(cadjust = TRUE, type = "HC0")). The CATT SE for cohort \(g\) is \(\sqrt{\psi_g^\top V_{\text{CR}} \psi_g}\) (using a zero-padded \(\psi_g\) on the full selected support); the ATT regression-coefficient variance is \(\psi_{\text{att}}^\top V_{\text{CR}} \psi_{\text{att}}\), replacing \(V_1\) above. The propensity-score piece \(V_2\) is unchanged because it depends on empirical cohort proportions, not regression residuals.

For FETWFE and BETWFE, se_type = "cluster" is only meaningful when q < 1 (the bridge oracle property is required); for q >= 1 the cluster path returns NA just like the default. ETWFE and twfeCovs have no q argument, so the cluster path always runs when the Gram matrix is invertible.

3.2 Why we call this experimental

The CR1 sandwich is a textbook estimator and the package’s implementation matches sandwich::vcovCL() to numerical precision on a clean panel without selection. What is not yet covered by the paper’s theory is:

  • Verification that the bridge oracle property of Theorem 6.2 still holds under the relaxed covariance structure that motivates cluster-robust SEs in the first place.
  • Sandwich consistency after model selection, i.e., that the CR1 sandwich evaluated at the bridge-selected support is a consistent estimator of the true asymptotic variance.

These extensions are mechanically routine but conceptually non-trivial, and they are explicitly outside the package’s current scope. Until they land, se_type = "cluster" is exposed as an opt-in, clearly-labelled experimental feature.

Recommendation. Until the theory lands, treat se_type = "cluster" as a sensitivity check: report both se_type = "default" and se_type = "cluster" in applied work, comment on the gap, and lean on the default for headline numbers.

3.3 Worked cluster comparison

res_cluster <- fetwfeWithSimulatedData(sim_data, se_type = "cluster")

c(
  default = res_default$att_se,
  cluster = res_cluster$att_se
)
#>    default    cluster 
#> 0.10486600 0.09988543

On this F1-conforming simulated panel the two SEs are similar by construction: the data-generating process satisfies F1, so the model-based SE is already valid and the cluster-robust SE estimates the same underlying variance. Under a deliberately serially-correlated DGP (or under heteroskedasticity, or higher-level clustering) the cluster-robust SE will typically be larger.

The print() and summary() methods label the SE so it is clear which one was used:

print(res_cluster)
#> Fused Extended Two-Way Fixed Effects Results
#> ===========================================
#> 
#> Overall Average Treatment Effect (ATT):
#>   Estimate:   -0.7777
#>   Std. Error (cluster-robust): 0.0999
#>   P-value:    6.912e-15
#>   Selected:   TRUE
#>   95% CI:    [-0.9735, -0.5819]
#> 
#> Cohort Average Treatment Effects (CATT) [simultaneous 95% CI]:
#>  cohort   estimate         se     ci_low    ci_high      p_value selected
#>       2 -0.4115949 0.06057451 -0.5468233 -0.2763665 2.168443e-11     TRUE
#>       3 -1.6779650 0.13059071 -1.9694998 -1.3864302 0.000000e+00     TRUE
#>       4  0.0000000 0.00000000  0.0000000  0.0000000           NA    FALSE
#> 
#> Event-Study Average Treatment Effects (per event time) [simultaneous 95% CI]:
#>  event_time n_cohorts   estimate        se     ci_low    ci_high      p_value
#>           0         3  0.0000000 0.0000000  0.0000000  0.0000000           NA
#>           1         3 -0.6342543 0.1102732 -0.8954173 -0.3730913 2.246783e-08
#>           2         3 -0.6342543 0.1102732 -0.8954173 -0.3730913 2.628733e-08
#>           3         2 -3.0242895 0.1920647 -3.4791616 -2.5694175 0.000000e+00
#>           4         1  0.1683984 0.2513507 -0.4268825  0.7636793 8.582874e-01
#> 
#> Model Details:
#>   Units (N)           : 120
#>   Time periods (T)    : 6
#>   Treated cohorts (G) : 3
#>   Covariates (d)      : 2
#>   Features (p)        : 62
#>   Selected size       : 26
#>   Lambda*             : 0.0122

The CATT SEs and confidence intervals in catt_df are recomputed from the same cluster-robust sandwich; the CATT p-values follow accordingly.

4. Choosing the bridge penalty parameter lambda

As of v1.13.0, fetwfe() and betwfe() select the bridge penalty lambda via 10-fold cross-validation by default (using grpreg::cv.grpreg). The prior default was BIC over the same grpreg lambda grid. The change addresses a finite-sample bias issue documented in simulation studies (see issue #164 in the package repository): under the prior BIC default, the overall ATT estimator was biased toward zero at moderate sample sizes, producing 95% confidence intervals whose empirical coverage was as low as 0.00 in some regimes. Cross-validation restores near-nominal coverage in every regime tested and additionally has the lowest MSE in the high-dim regime FETWFE was originally designed for.

The lambda_star and lambda_star_model_size result slots are populated identically under both selection paths; only the underlying selection criterion differs. The CV path stores its provenance in three new top-level slots: lambda_selection, cv_folds, and cv_seed. The seed defaults to as.integer(N * T) so consecutive calls on the same dataset are reproducible without the user having to specify a seed.

To select lambda by BIC instead, pass lambda_selection = "bic":

fit_bic <- fetwfe(..., lambda_selection = "bic")

To control the CV fold assignment explicitly, pass cv_folds and/or cv_seed:

fit_cv <- fetwfe(..., cv_folds = 20L, cv_seed = 42L)

References

Bertrand, M., Duflo, E., & Mullainathan, S. (2004). “How much should we trust differences-in-differences estimates?” Quarterly Journal of Economics 119(1), 249–275.

Faletto, G. (2025). “Fused Extended Two-Way Fixed Effects for Difference-in-Differences with Staggered Adoptions.” arXiv preprint arXiv:2312.05985.

Kock, A. B. (2013). “Oracle Efficient Variable Selection in Random and Fixed Effects Panel Data Models.” Econometric Theory 29(1), 115–152.

Liang, K.-Y., & Zeger, S. L. (1986). “Longitudinal data analysis using generalized linear models.” Biometrika 73(1), 13–22.