带均值 Groups/Fama-MacBeth 估计器的 Newey-West 标准误

Question

我试图让 Newey-West 标准误差与来自 plm 包的 pmg()（平均 Groups/Fama-MacBeth 估计器）的输出一起工作。

效仿 here 中的示例：

require(foreign)
require(plm)
require(lmtest)
test <- read.dta("http://www.kellogg.northwestern.edu/faculty/petersen/htm/papers/se/test_data.dta")

fpmg <- pmg(y~x, test, index=c("firmid", "year")) # Time index in second position, unlike the example

我可以直接使用coeftest就可以得到Fama-MacBeth标准错误：

# Regular “Fama-MacBeth” standard errors
coeftest(fpmg)

# t test of coefficients:
#   
#             Estimate Std. Error t value Pr(>|t|)    
# (Intercept) 0.032470   0.071671   0.453   0.6505    
# x           0.969212   0.034782  27.866   <2e-16 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

但是，尝试使用 Newey-West 估计器失败了：

# Newey-West standard-errors
coeftest(fpmg, vcov = NeweyWest(fpmg, lag=3))

# Error in UseMethod("estfun") : 
#   no applicable method for 'estfun' applied to an object of class "c('pmg', 'panelmodel')"

这似乎是 plm 包中的一个缺点。你知道做这个工作的方法吗？我应该为 pmg 对象编写自己的 estfun 代码吗？从头开始编写 Newey-West 估算器？或者我应该完全绕过 plm 包吗？

Answer 1

目前 plm 软件包无法做到这一点。

但是，您可以自己创建它们。

假设你有：

fpmg <- pmg(y~x, test, index = c('year', 'firmid'))
fpmg.coefficients <- fpmg$coefficients
# (Intercept)            x 
#  0.03127797   1.03558610 

coeftest(fpmg)
#             Estimate Std. Error t value Pr(>|t|)    
# (Intercept) 0.031278   0.023356  1.3392   0.1806    
# x           1.035586   0.033342 31.0599   <2e-16 ***

然后您可以自己简单地创建估算器，例如：

the.years <- unique(test$year)
a.formula <- y ~ x


first.step <-  lapply(the.years, function(a.year) {
                temp.data <- test[test$year == a.year, ]
                an.lm <- lm(a.formula, data = temp.data)
                the.coefficients <- an.lm$coef
                the.results <- as.data.frame(cbind(a.year, t(the.coefficients)))
                the.results
                }) 

first.step.df <- do.call('rbind', first.step)

second.step.coefficients <- apply(first.step.df[, -1], 2, mean)
second.step.coefficients
# (Intercept)           x 
#  0.03127797  1.03558610 

identical(fpmg.coefficients, second.step.coefficients)
# [1] TRUE

检查它们是否完全相同以防万一。最后，您可以获得 Newey-West (1987) 和一个滞后调整的 t 统计量：

library(sandwich)
second.step.NW.sigma.sq <- apply(first.step.df[, -1], 2, 
                             function(x) sqrt(NeweyWest(lm(x ~ 1), 
                               lag = 1, prewhite = FALSE)['(Intercept)',       
                                 '(Intercept)']))
second.step.NW.sigma.sq
#  (Intercept)            x 
#   0.02438398   0.02859447
t.statistics.NW.lag.1 <- second.step.coefficients / second.step.NW.sigma.sq

t.statistics.NW.lag.1
# (Intercept)           x 
#    1.282726   36.216301

更新

在我的回答中，我只包含了 t 统计量的 "manual" 计算，因为它的计算速度更快。更通用的解决方案是使用 lmtest 包的 coeftest() 函数计算 Newey-West 校正的 t 统计量及其 p 值。

coeftest(lm(first.step.df$'(Intercept)' ~ 1), vcov = NeweyWest(lm(first.step.df$'(Intercept)' ~ 1), lag = 1, prewhite = FALSE))
# t test of coefficients:
#             Estimate Std. Error t value Pr(>|t|)
# (Intercept) 0.031278   0.024384  1.2827   0.2316
coeftest(lm(first.step.df$x ~ 1), vcov = NeweyWest(lm(first.step.df$x ~ 1), lag = 1, prewhite = FALSE))
# t test of coefficients:
#             Estimate Std. Error t value  Pr(>|t|)    
# (Intercept) 1.035586   0.028594  36.216 4.619e-11 ***

Answer 2

我认为我们可以使用 tidyverse 进行 Newey–West 调整，因为基础 R 方法过于冗长。三四年前，tidyverse还很不成熟，但今天我们可以比较完美的使用tidyverse。这是 tidyverse 方法。此外，broom包在我的回答中使用。

# test data:
set.seed(1234)
a <- tibble(
            year = rep(1:20, each = 100),
            firmid = rep(1:100, times = 20),
            x1 = rnorm(2000),
            x2 = rnorm(2000),
            y = 2 + 3 * x1 - 2 * x2 + rnorm(1)
           )

Newey-West 调整：

a.formula <- formula("y ~ x1+x2")
lag <- 1

regco <- a %>% 
   group_by(year) %>%
   group_modify(~tidy(lm(a.formula, data = .x))) %>% 
   select(year,term,estimate) %>% 
   pivot_wider(names_from = term,
               values_from = estimate) %>% 
   ungroup() %>% 
   select(x1, x2, `(Intercept)`) %>% 
   map(
       ~coeftest(
                 lm(.x ~ 1), 
                 vcov = NeweyWest(
                                  lm(.x ~ 1), 
                                  lag = lag, 
                                  prewhite = FALSE
                                 )
                )
      )

带均值 Groups/Fama-MacBeth 估计器的 Newey-West 标准误

Newey-West standard errors with Mean Groups/Fama-MacBeth estimator

regression

r

standard-error

更新