Стандартные ошибки Ньюи-Уэста со средними группами / оценка Фама-Макбета
Я пытаюсь заставить стандартные ошибки Ньюи-Уэста работать с выводом pmg()
(Средняя группа / оценка Фама-Макбета) из plm
пакет.
Следуя примеру отсюда:
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 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
пакет. Вы знаете, как сделать эту работу? Должен ли я кодировать свой собственный estfun
за pmg
объекты? Код оценки Ньюи-Уэста с нуля? Или я должен обойти plm
пакет в целом?
2 ответа
В настоящее время это невозможно с 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-статистики, потому что это вычислительно быстрее. Более общим решением является вычисление скорректированной t-статистики Ньюи-Уэста и ее p-значений с помощью coeftest()
функция lmtest
пакет.
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 ***
Я думаю, что мы можем использовать tidyverse для настройки newey west, поскольку базовые методы R слишком многословны. 3 или 4 года назад тидиверс был незрелым, но сегодня мы можем относительно совершить тидиверс. вот метод 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 adjustment:
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)
regco %>% ungroup() %>% select(x1,x2,`(Intercept)`) %>% map(~coeftest(lm(.x ~
1), vcov = NeweyWest(lm(.x ~ 1), lag = lag, prewhite = FALSE)))