我们可以使用sandwich::vcovCL
来获得与lfe::felm
基本相同的标准误差.让我们估计一些模型.
> est1 <- lfe::felm(y ~ x1 + x2 | id + firm | 0 | firm, data=d) ## id + firm FE, clustered by firm
> est2 <- lfe::felm(y ~ x1 | id + firm | 0 | firm, data=d) ## same, restricted model
> est3 <- lm(y ~ x1 + x2 + as.factor(id) + as.factor(firm), data=d) ## as est1 w/o clustering
> est4 <- lm(y ~ x1 + as.factor(id) + as.factor(firm), data=d) ## as restricted est3
比较est3和est1的标准误差,
> lmtest::coeftest(est3, vcov.=\(x)
+ sandwich::vcovCL(x, cluster=d$firm, type='HC0'))[1:3, ]
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.7416642 0.10554929 35.44945 5.454680e-177
x1 1.0432612 0.02965723 35.17730 3.631800e-175
x2 0.4904104 0.03186679 15.38939 5.822822e-48
> coef(summary(est1))
Estimate Cluster s.e. t value Pr(>|t|)
x1 1.0432612 0.02968696 35.14207 1.799173e-13
x2 0.4904104 0.03189874 15.37398 2.931632e-09
yield 率基本相同.
因此,根据Cross validated [1]上的一篇帖子,我们可以使用lmtest::waldtest
来比较est1和est2型号(还有一种lfe::waldtest
,它的工作方式不同).
> lmtest::waldtest(est3, est4, vcov=\(x)
+ sandwich::vcovCL(x, cluster=d$firm, type='HC0'))
Wald test
Model 1: y ~ x1 + x2 + as.factor(id) + as.factor(firm)
Model 2: y ~ x1 + as.factor(id) + as.factor(firm)
Res.Df Df F Pr(>F)
1 966
2 967 -1 236.83 < 2.2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
希望这能让你走得更远.
顺便说一句:你绝对可以把这件事作为作者GitHub分的一个问题来提交.
Data:个
> set.seed(42)
> n <- 1e3
> d <- data.frame(x1=rnorm(n), x2=rnorm(n), id=factor(sample(20, n, replace=TRUE)),
+ firm=factor(sample(13, n, replace=TRUE)), u=rnorm(n))
> id.eff <- rnorm(nlevels(d$id))
> firm.eff <- rnorm(nlevels(d$firm))
> d$y <- d$x1 + 0.5 * d$x2 + id.eff[d$id] + firm.eff[d$firm] + d$u