有没有办法提取 coeftest() 函数的 $R^2$ ？答案

【问题标题】：Is there a way to extract $R^2$ of the coeftest() function?有没有办法提取 coeftest() 函数的 $R^2$ ？
【发布时间】：2019-11-25 23:11:15
【问题描述】：

我有一个关于函数 coeftest() 的问题。我试图弄清楚从哪里可以得到这个函数的 R 平方的任何结果。我正在拟合一个标准的多元线性回归，如下所示：

Wetterstation.lm <- lm(temp~t+temp_auto+dum.jan+dum.feb+dum.mar+dum.apr+dum.may+dum.jun+dum.aug+dum.sep+dum.oct+dum.nov+dum.dec+
                      dum.jan*t+dum.feb*t+dum.mar*t+dum.apr*t+dum.may*t+dum.jun*t+dum.aug*t+dum.sep*t+dum.oct*t+dum.nov*t+dum.dec*t)

预先我分别定义了这些变量中的每一个，我的结果如下：


> summary(Wetterstation.lm)

Call:
lm(formula = temp ~ t + temp_auto + dum.jan + dum.feb + dum.mar + 
    dum.apr + dum.may + dum.jun + dum.aug + dum.sep + dum.oct + 
    dum.nov + dum.dec + dum.jan * t + dum.feb * t + dum.mar * 
    t + dum.apr * t + dum.may * t + dum.jun * t + dum.aug * t + 
    dum.sep * t + dum.oct * t + dum.nov * t + dum.dec * t)

Residuals:
     Min       1Q   Median       3Q      Max 
-10.9564  -1.3214   0.0731   1.3621   9.9312 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  3.236e+00  9.597e-02  33.714  < 2e-16 ***
t            1.206e-05  3.744e-06   3.221  0.00128 ** 
temp_auto    8.333e-01  2.929e-03 284.503  < 2e-16 ***
dum.jan     -3.550e+00  1.252e-01 -28.360  < 2e-16 ***
dum.feb     -3.191e+00  1.258e-01 -25.367  < 2e-16 ***
dum.mar     -2.374e+00  1.181e-01 -20.105  < 2e-16 ***
dum.apr     -1.582e+00  1.142e-01 -13.851  < 2e-16 ***
dum.may     -7.528e-01  1.106e-01  -6.809 9.99e-12 ***
dum.jun     -3.283e-01  1.106e-01  -2.968  0.00300 ** 
dum.aug     -2.144e-01  1.094e-01  -1.960  0.05005 .  
dum.sep     -8.009e-01  1.103e-01  -7.260 3.96e-13 ***
dum.oct     -1.752e+00  1.123e-01 -15.596  < 2e-16 ***
dum.nov     -2.622e+00  1.181e-01 -22.198  < 2e-16 ***
dum.dec     -3.287e+00  1.226e-01 -26.808  < 2e-16 ***
t:dum.jan    2.626e-06  5.277e-06   0.498  0.61877    
t:dum.feb    2.479e-06  5.404e-06   0.459  0.64642    
t:dum.mar    1.671e-06  5.277e-06   0.317  0.75145    
t:dum.apr    1.357e-06  5.320e-06   0.255  0.79872    
t:dum.may   -3.173e-06  5.276e-06  -0.601  0.54756    
t:dum.jun    2.481e-06  5.320e-06   0.466  0.64098    
t:dum.aug    5.998e-07  5.298e-06   0.113  0.90985    
t:dum.sep   -5.997e-06  5.321e-06  -1.127  0.25968    
t:dum.oct   -5.860e-06  5.277e-06  -1.110  0.26683    
t:dum.nov   -4.215e-06  5.320e-06  -0.792  0.42820    
t:dum.dec   -2.526e-06  5.277e-06  -0.479  0.63217    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.12 on 35744 degrees of freedom
Multiple R-squared:  0.9348,    Adjusted R-squared:  0.9348 
F-statistic: 2.136e+04 on 24 and 35744 DF,  p-value: < 2.2e-16

现在我尝试使用函数 coeftest() 和 vcovHAC 来调整异方差和自相关，如下所示：

library("lmtest")
library("sandwich")
Wetterstation.lm.HAC <- coeftest(Wetterstation.lm, vcov = vcovHAC)

这些结果是：

> Wetterstation.lm.HAC

t test of coefficients:

               Estimate  Std. Error  t value  Pr(>|t|)    
(Intercept)  3.2356e+00  7.8816e-02  41.0529 < 2.2e-16 ***
t            1.2059e-05  2.7864e-06   4.3280 1.509e-05 ***
temp_auto    8.3334e-01  2.9798e-03 279.6659 < 2.2e-16 ***
dum.jan     -3.5505e+00  1.1843e-01 -29.9789 < 2.2e-16 ***
dum.feb     -3.1909e+00  1.2296e-01 -25.9507 < 2.2e-16 ***
dum.mar     -2.3741e+00  1.0890e-01 -21.8002 < 2.2e-16 ***
dum.apr     -1.5821e+00  9.5952e-02 -16.4881 < 2.2e-16 ***
dum.may     -7.5282e-01  8.8987e-02  -8.4599 < 2.2e-16 ***
dum.jun     -3.2826e-01  8.2271e-02  -3.9899 6.622e-05 ***
dum.aug     -2.1440e-01  7.7966e-02  -2.7499  0.005964 ** 
dum.sep     -8.0094e-01  8.4456e-02  -9.4835 < 2.2e-16 ***
dum.oct     -1.7519e+00  9.2919e-02 -18.8538 < 2.2e-16 ***
dum.nov     -2.6224e+00  1.0028e-01 -26.1504 < 2.2e-16 ***
dum.dec     -3.2873e+00  1.1393e-01 -28.8546 < 2.2e-16 ***
t:dum.jan    2.6256e-06  5.2429e-06   0.5008  0.616517    
t:dum.feb    2.4790e-06  5.5284e-06   0.4484  0.653850    
t:dum.mar    1.6713e-06  4.8632e-06   0.3437  0.731107    
t:dum.apr    1.3567e-06  4.5670e-06   0.2971  0.766423    
t:dum.may   -3.1734e-06  4.2970e-06  -0.7385  0.460209    
t:dum.jun    2.4809e-06  4.1490e-06   0.5979  0.549880    
t:dum.aug    5.9983e-07  4.0379e-06   0.1485  0.881910    
t:dum.sep   -5.9975e-06  4.1675e-06  -1.4391  0.150125    
t:dum.oct   -5.8595e-06  4.4635e-06  -1.3128  0.189265    
t:dum.nov   -4.2151e-06  4.6555e-06  -0.9054  0.365263    
t:dum.dec   -2.5257e-06  4.9871e-06  -0.5065  0.612539    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

但由于我想将 R 平方添加到总结我的结果的表格中，我无法弄清楚如何获得它。现在我想知道是否有人可以帮助解决这个问题并告诉我可以从哪里获得信息。也许我只是愚蠢，但由于我对 R 很陌生，我很乐意得到任何帮助。

提前非常感谢您。

【问题讨论】：

coeftest() 来自 lmtest 包？
当我进一步处理该主题时，我发现 HAC 只会调整标准误差，而不会改变估计系数的方式。这是正确的还是我只是被误导了？如果这是真的，R^2 不会改变，因为残差保持不变。那不是真的吗？如果是这样，coeftests 显示不同的参数估计值不是很奇怪吗？我现在完全糊涂了。
@mischva11 是的，它来自包“lmtest”

标签： r regression linear-regression lm

【解决方案1】：

简单的回答：不，没有。而且也没有理由这样做。 coeftest() 函数正在使用给定模型的值。使用stats4::coef，coeftest 函数获取模型的系数。

如果函数打算这样做，则可以提取r^2 值。另外，imtest coeftest() 只返回一个粘贴，所以你不能提取值。

总结一下：lmtest::coeftest() 正在使用lm 的值，因此r^2 没有改变

关于标准误差的背景：lm 使用稍微不同的方法来计算标准误差。在source code可以提取：

  R <- chol2inv(Qr$qr[p1, p1, drop = FALSE])
  se <- sqrt(diag(R))

所以这让我们关注：lm 使用 Cholesky 组合（我也认为 R 默认使用它）。

coeftest() 同时使用残差的方差-协方差矩阵的sqrt()（见here）。所以自相关。 vcov 提取给定模型的方差-协方差矩阵（如lm）

se <- vcov.
se <- sqrt(diag(se))

我个人认为lm 的用户出于速度原因使用 chalesky 组合，因为您不必反转矩阵。但是lmtest 包的作者并没有意识到这一点。但这只是猜测。

由于 t 值是在两个包中使用估计值和标准误差计算的，如下所示：

tval <- as.vector(est)/se

而p-value 是基于tval：

2*pt(abs(tval), rdf, lower.tail = FALSE)

所有差异均基于不同的估计误差。请注意，估计值是相同的，因为 coeftest() 只是继承了它们。

【讨论】：

当你提到coef()函数时，你的意思是coeftest()吗？所以总结一下，这是否意味着唯一改变的是标准误差估计？如果是这样，为什么 coeftest() 函数会返回其他参数估计值？
@Leon 更新了答案，很抱歉给您带来不便，是的，我的意思是 coeftest。