grantmcdermott / etwfe

Extended two-way fixed effects
https://grantmcdermott.com/etwfe/
Other
50 stars 11 forks source link

Can I wild cluster bootstrap this? =) #10

Closed s3alfisc closed 1 year ago

s3alfisc commented 1 year ago

Hi Grant,

Thanks for sharing this package - it looks really cool! Naturally, I am wondering if I can wild (cluster) bootstrap this. =)

After taking a very brief (really very, very brief) glance at both code and the Wooldridge paper, my very topline understanding of the estimation procedure is that etwfe implements a linear regression / glm with a set of interacted covariates and a few estimated properties (e.g. equation 5.7 in the Wooldridge paper attached below, where $\bar{x}_{1}$ is an estimated mean and therefore has a sampling variance for which one has to correct inference for).

image

After estimating the 'full' model, one then computes marginal effects, e.g. for the effect on a specific cohort, you compute the marginal effect for a given cohort. If we assume that there are $C$ cohorts, for linear models, inference for all $C$ 'cohort effects' should then boil down to C linear multi-parameter tests of the form $R{c}'\beta{c} = r$ vs $R{c}'\beta{c} \neq r$, where $R{c}$ is a vector of dimension $k{c}$, with $k_{c}$ denoting all regression parameters that depend on cohort $c$?

Here's a quick formulaic example to (hopefully) clarify:

We estimate a linear model

$y = \beta0 + \beta{1} X{1} + \beta{2} X{1} X{2} + \epsilon_{i}$

via OLS.

We then compute the average marginal effect at the mean for $X_{1}$ as

$E[\frac{\partial E(y|X)}{\partial X{1}}] = \beta{1} + \beta{2} \bar{X}{2}.$

Inference for the marginal effect then follows a standard $R\beta = r$ schema, where $R = [1, \bar{X}_{2}]$, $\beta = [\beta_1, \beta_2]$ and $r = 0$.

Do I get this right? If inference can be squeezed into the $R\beta = r$ schema, it will be fairly straightforward to make fwildclusterboot work with etwfe - but only under the assumption that $\bar{x}_{1}$ is a noiseless estimate of $x_1$.

Best, Alex

grantmcdermott commented 1 year ago

Hi Alex. Sorry for the super slow response. I've been on vacation and then straight into work travel...

I think you're intuition is correct and everything should pass through correctly. However, I must confess, that I haven't thought about this passed reading your comment and the time it's taken to type up this response. I'll let you test (simulate) and decide when is best to close/update ;-)

s3alfisc commented 1 year ago

Hi Grant, thanks for your feedback! No worries about the delay, I hope you had a relaxing vacation =) I will close this for now and will update you once I manage to implement support for the wcb and etwfe + have some simulations at hand that look promising!

jtorcasso commented 1 year ago

Hi all! Checking into this issue to see if there was any progress. I tried doing this myself but got an error.

s3alfisc commented 1 year ago

Hi @jtorcasso - I have never looked into this in more detail. But it should definitely be possible, and if I understand everything correctly, implementing support for the wild cluster bootstrap / etwfe should be around one Saturday of work for me (I am alwazys optimistic with these things, though). But I am committed to too many projects at the moment, so I cannot promise to take a closer look in the next weeks.

In principle, there are two options how to add support for wild cluster bootstrapping to etwfe - one of them is sketched out above and requires to loop over boottest() multiple times, similar as done for the boot_aggregate() implementation of the wcb for fixest::sunab() here. The other option is to simply pass a bootstrapped covariance matrix to the feols call / marginaleffects(). You could do this via sandwich::vcovBS (potentially slow, and note that this way you would not compute a percentile-t bootstrap and would not impose the null on the bootstrap dgp). Alternatively, I have a rough implementation of the classical pairs bootstrap in the dev branch of fwildclusterboot, which computes a bootstrapped vcov matrix (but I have not really tested it). Beyond, I unfortunately don't think it is currently possible to get a full bootstrap vcov matrix out of fwildclusterboot::boottest(). If you post your code as an issue / PR to fwildclusterboot, I can take a look tomorrow evening, in case you think it might help. You can also send me an email if you'd prefer that (you could find it e.g. here). =)

s3alfisc commented 1 year ago

Also, just to point this out again - if the Sun-Abrams method is also appropriate for your use case, boot_aggregate does support the wcb for fixest::sunab(). =)

s3alfisc commented 1 year ago

Ok, I took a quick look at this, there are (at least) two issues:

s3alfisc commented 1 year ago

Actually, here is a hacky preliminary solution, in five steps, which works without the PR to fixest mentioned above:

library(etwfe)
library(fixest)
library(fwildclusterboot)

data("mpdta", package="did")
head(mpdta)

mod = etwfe(
  fml = lemp ~ lpop,
  tvar  = year,
  gvar = first.treat,
  data = mpdta,
  se = "hetero"
)

# get it from assigning a global var data_etwfe <<- data in `etwfe()`
data_etwfe

# get formulas, write one part formula, varying slopes as var1 / var2
mod$fml_all
# new formula
fml <- lemp ~ .Dtreat:i(first.treat, i.year, ref = 0, ref2 = 2003)/lpop_dm + first.treat + first.treat/lpop + year + year/lpop
# fit model
fit <- feols(fml, data_etwfe, cluster ~ countyreal, ssc = ssc(adj = FALSE, fixef.K = "none"))

#aggregate
agg <- aggregate(
  fit, 
  c("ATT"="^.Dtreat:first.treat::[0-9]{4}:year::[0-9]{4}$")
)
agg
#       Estimate Std. Error t value    Pr(>|t|)
# ATT 0.09483085 0.02676169 3.54353 0.000431911

boot_aggregate(
  fit, 
  c("ATT"="^.Dtreat:first.treat::[0-9]{4}:year::[0-9]{4}$"), 
  B = 999,
  boot_ssc = ssc(adj = FALSE),
  clustid = ~countyreal 
)
# Estimate  t value Pr(>|t|)    [0.025%   0.975%]
# [1,] 0.09483085 3.547078        0 0.04099321 0.1487666

t-stats, pvalue and confidence interval of boot_aggregate() are computed via boottest(). Note that the reported non-bootstrapped t-stat of aggregate() and boot_aggregate() should match exactly, differences are likely / potentially due to incorrect handling of small sample adjustments in boot_aggregate (partially fixed in the dev branch).

Keep in mind that none of this is properly tested, that fixed effect varying slopes syntax does not work (wrong results), and that there is a bug for boot_aggregate() and the HC bootstrap (fixed in dev). Anyways, I hope this helps!