Я думаю, что вы на правильном пути с coeftest
в пакете lmtest. Взгляните на сэндвич-пакет , который включает в себя эту функциональность и предназначен для работы рука об руку с lmtest пакетом, который вы уже нашли.
> # generate linear regression relationship
> # with Homoskedastic variances
> x <- sin(1:100)
> y <- 1 + x + rnorm(100)
> ## model fit and HC3 covariance
> fm <- lm(y ~ x)
> vcovHC(fm)
(Intercept) x
(Intercept) 0.010809366 0.001209603
x 0.001209603 0.018353076
> coeftest(fm, vcov. = vcovHC)
t test of coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.01973 0.10397 9.8081 3.159e-16 ***
x 0.93992 0.13547 6.9381 4.313e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Чтобы получить F-тест, посмотрите на функцию waldtest()
:
> waldtest(fm, vcov = vcovHC)
Wald test
Model 1: y ~ x
Model 2: y ~ 1
Res.Df Df F Pr(>F)
1 98
2 99 -1 48.137 4.313e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Вы всегда можете придумать простую функцию, чтобы объединить эти два для вас, если бы вы хотели однострочник ...
Есть множество примеров в эконометрических вычислениях с оценщиками ковариантных матриц HC и HAC виньетка, которая поставляется с пакетом сэндвичей, связывающих lmtest и бутерброд, для того, что вы хотите.
Редактировать: Однострочник может быть простым:
mySummary <- function(model, VCOV) {
print(coeftest(model, vcov. = VCOV))
print(waldtest(model, vcov = VCOV))
}
Что мы можем использовать следующим образом (на примерах выше):
> mySummary(fm, vcovHC)
t test of coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1.01973 0.10397 9.8081 3.159e-16 ***
x 0.93992 0.13547 6.9381 4.313e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Wald test
Model 1: y ~ x
Model 2: y ~ 1
Res.Df Df F Pr(>F)
1 98
2 99 -1 48.137 4.313e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1