FE 回归中的 Driscoll 和 Kraay 标准误差：在 R 中再现 Stata xtscc 输出答案

【问题标题】：Driscoll and Kraay standard errors in FE regression: reproducing Stata xtscc output in RFE 回归中的 Driscoll 和 Kraay 标准误差：在 R 中再现 Stata xtscc 输出
【发布时间】：2021-03-21 05:58:43
【问题描述】：

我正在尝试使用包plm 在 R 中复制 Stata 命令xtscc 提供的结果，但我在看到相同的标准错误时遇到了一些麻烦我也在 Stata 中使用 plm 包中的数据集进行复制。

# code to obtain dataset
library(lmtest)
library(car)
library(tidyverse)
data("Produc", package="plm")
write.dta(Produc,"test.dta")

我的目标是使用 Driscoll 和 Kraay 标准误差运行双向固定效应面板模型估计。 Stata中的套路如下

use "test.dta", clear \\ to import data
** i declare the panel 
xtset state year
* create the dummies for the time fixed effects
quietly tab year, gen(yeardum)
* run a two way fixed effect regression model with Driscoll and Kraay standard errors
xi: xtscc gsp pcap emp unemp yeardum*,fe 
* results are the following
                    Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
        pcap |  -.1769881    .265713    -0.67   0.515    -.7402745    .3862983
         emp |   40.61522   2.238392    18.14   0.000     35.87004     45.3604
       unemp |   23.59849   85.10647     0.28   0.785    -156.8192    204.0161

在 R 中，我使用以下例程：

# I declare the panel
Produc <- pdata.frame(Produc, index = c("state","year"), drop.index = FALSE)
# run a two way fixed effect model
femodel <- plm(gsp~pcap+emp+unemp, data=Produc,effect = "twoway", 
               index = c("iso3c","year"), model="within")
# compute Driscoll and Kraay standard errors using vcovSCC
coeftest(femodel, vcovSCC(femodel))

pcap  -0.17699    0.25476 -0.6947   0.4874    
emp   40.61522    2.14610 18.9252   <2e-16 ***
unemp 23.59849   81.59730  0.2892   0.7725

虽然点估计与 Stata 中的相同，但标准误不同。

为了检查我是否对标准误差使用了“错误的”小样本调整，我还尝试使用所有可用调整运行 coeftest，但没有一个产生与 xtscc 相同的值。

library(purrr)
results <- map(c("HC0", "sss", "HC1", "HC2", "HC3", "HC4"),~coeftest(femodel, vcovSCC(femodel,type = .x)))
walk(results,print)
# none of the estimated standard errors is the same as xtscc

有人知道我如何在 R 中复制 Stata 的结果吗？

【问题讨论】：

可能由于 Stata 在其内部模型中包含截距而 plm 没有，因此小样本调整略有不同（有限情况；在大样本中不相关）。在此处或 Stackoverflow 上的小样本调整中讨论了 Stata 与 plm 的差异，也许这会有所帮助。您还可以查看plm::within_intercept 的编码（或文档），以了解如何计算具有截距的 FE 模型（文档是参考）。我怀疑xtssc 以某种方式使用了Stata 的xtreg。不知道xtssc是否使用了和xtreg一样的小样本调整。
感谢您的评论，实际数据的问题在于，虽然 R 系数很重要，但 Stata 系数不重要。调整为 0.03。还有你知道Stata中使用什么类型的标准错误吗？从 xtscc 的文档看不清楚..
xtscc 是用户贡献的命令。如果其文档没有显示所采用的小样本调整，您可能需要阅读函数源代码或联系其作者。在 plm 中，type = sss 模仿 Stata 的小样本调整（因此得名）。
请注意，您对plm 的调用中的index 参数将被忽略，因为您的数据已经是pdata.frame。此外，索引参数包含一个不在数据中的变量 (iso3c)。

标签： r panel stata panel-data plm

【解决方案1】：

从plm 2.4 版开始，它的函数within_intercept(., return.model = TRUE) 可以返回带有截距的内部模型的完整模型，就像在Stata 中一样。这样，就可以准确复制 Stata 用户贡献的命令 xtscc 的结果。

xtscc 的工作方式似乎是将双向 FE 模型估计为单向 FE 模型 + 时间维度的假人。所以让我们用plm复制它：

data("Produc", package="plm")
Produc <- pdata.frame(Produc, index = c("state","year"), drop.index = FALSE)

femodel <- plm(gsp ~ pcap + emp + unemp + factor(year), data = Produc, model="within")
femodelint  <- within_intercept(femodel,  return.model = TRUE)

lmtest::coeftest(femodelint, vcov. = function(x) vcovSCC(x, type = "sss"))
#                     Estimate  Std. Error t value              Pr(>|t|)    
# (Intercept)      -6547.68816  3427.47163 -1.9104             0.0564466 .  
# pcap                -0.17699     0.26571 -0.6661             0.5055481    
# emp                 40.61522     2.23839 18.1448 < 0.00000000000000022 ***
# unemp               23.59849    85.10647  0.2773             0.7816356    
# [...]

【讨论】：