Skip to contents

This vignette shows the shortest path from a data frame to an ATT estimate. The default workflow uses cross-fitted DML with built-in nuisance learners.

Simulate Example Data

set.seed(123)
dat <- simulate_comp_did(n = 1000, delta = 0.25)

head(dat)
#>          y d post         x1          x2 x3 x4 x5 x6
#> 1 179.8029 0    1 -0.4248450 -0.45275453  0  0  1  1
#> 2 263.3234 0    1  0.5766103  0.18773387  0  1  2  3
#> 3 179.6477 0    1 -0.1820462 -0.67963037  0  0  2  2
#> 4 336.6057 0    1  0.7660348  0.70686048  1  1  2  2
#> 5 348.8698 1    0  0.8809346  0.69547832  0  0  1  1
#> 6 145.9946 1    1 -0.9088870 -0.04422637  1  1  1  0

The analysis data need four ingredients:

  • an outcome column;
  • a binary treatment-group column coded as 0/1;
  • a binary post-period column coded as 0/1;
  • optional pre-treatment covariates.

Estimate the ATT

The formula supplies the outcome and covariates. The treatment and post-period columns are supplied separately because they define the DiD design.

fit <- comp_did(
  y ~ x1 + x2 + x3 + x4 + x5 + x6,
  data = dat,
  tname = "post",
  dname = "d",
  K = 2,
  seed = 123
)

fit
#> Compositional-change DiD
#> Nuisance method: dml 
#> Nuisance backend: custom 
#> ATT: 1.782 
#> SE : 2.676 
#> 95% CI: [-3.462, 7.026]

Work With the Result

summary() returns a compact table object, as.data.frame() returns a one-row data frame, and confint() returns the confidence interval.

summary(fit)
#> Compositional-change DiD summary
#> Observations: 1000 
#> Nuisance method: dml 
#> Nuisance backend: custom 
#> 
#>  term estimate std.error conf.low conf.high
#>   ATT    1.782     2.676   -3.462     7.026
as.data.frame(fit)
#>   term estimate std.error  conf.low conf.high method backend
#> 1  ATT 1.781903   2.67559 -3.462252  7.026059    dml  custom
confint(fit)
#>          2.5%    97.5%
#> ATT -3.462156 7.025963

Add the Stationarity Test

Set stationary = TRUE to compute the stationary estimator as well, and stationarity_test = TRUE to compare it with the non-stationary estimator.

fit_test <- comp_did(
  y ~ x1 + x2 + x3 + x4 + x5 + x6,
  data = dat,
  tname = "post",
  dname = "d",
  K = 2,
  seed = 123,
  stationary = TRUE,
  stationarity_test = TRUE
)

summary(fit_test)
#> Compositional-change DiD summary
#> Observations: 1000 
#> Nuisance method: dml 
#> Nuisance backend: custom 
#> 
#>  term estimate std.error conf.low conf.high
#>   ATT    1.782     2.676   -3.462     7.026
#> 
#> Stationarity test p-value: 0.04051
fit_test$stationarity_test
#> Hausman-type test based on the difference in estimated influence functions
#> Null: ATT(nonstationary) = ATT(stationary) 
#> Statistic: 4.1962 
#> p-value: 0.0405147 
#>  alpha critical_value reject
#>   0.10       2.705543   TRUE
#>   0.05       3.841459   TRUE
#>   0.01       6.634897  FALSE

For clustered inference, pass a cluster identifier to drdid_cluster_bootstrap(). The example below uses artificial clusters only to show the workflow.

set.seed(321)
cluster_id <- sample(1:50, nrow(dat), replace = TRUE)

clustered <- drdid_cluster_bootstrap(
  nonstationary_fit = fit_test$att,
  stationary_fit = fit_test$stationary,
  group = cluster_id,
  nboot = 99,
  seed = 123,
  keep_boots = FALSE
)

clustered
#> Mammen cluster bootstrap for DR DiD estimators
#> Clusters: 50   Observations: 1000   Draws: 99 
#> Non-stationary cluster SE: 2.33 
#> Stationary cluster SE    : 2.003 
#> Stationarity clustered p-value: 0.0404

Local-Polynomial Nuisances

The local-polynomial path is available through nuisance_method = "loo". If you use a formula, every right-hand-side variable is treated as a continuous local-polynomial covariate, and the package requires those formula covariates to be numeric.

For mixed continuous and discrete kernels, use covariate block selectors instead of a formula. This lets the nonparametric engine use the right kernel for each covariate type.

fit_lp <- comp_did(
  yname = "y",
  tname = "post",
  dname = "d",
  data = dat,
  continuous = c("x1", "x2"),
  unordered = c("x3", "x4"),
  ordered = c("x5", "x6"),
  nuisance_method = "loo",
  list_control_ps = list(n_start = 1, ps_min = 1e-5, lp_order = 1),
  list_control_or = list(n_start = 1, bw_constrained = FALSE, or_order = 1)
)