为什么summary()显示的标准误差与coeftest()不同？

Question

我使用鲁棒标准误差运行glm()。在后续的模型比较中，我计算了两个回归模型（系数和se）之间的差异。为此计算，我使用了summary()函数。然而，模型的summary函数显示的标准误差与我从coeftest()获取的标准误差不同。系数的值保持不变。

输入：

mod.01 <- glm(dep ~ indep1 + indep2 + indep3,
          family = binomial (link = "logit"), data = data)
coeftest(mod.01, vcov. = vcovHC, type= "HC3", df = NULL)

summary(mod.01, robust=T)

输出：

coeftest()      
                           Estimate  Std. Error  t value  Pr(>|t|)    
    (Intercept)         -2.72917626  0.16367787 -16.6741 < 2.2e-16 ***
    indep1               0.00427870  0.41928906   0.0102  0.991859    
    indep2               2.00243724  0.19757861  10.1349 < 2.2e-16 ***
    indep3               0.36385098  0.32783817   1.1098  0.267159    


summary()
                           Estimate Std. Error z value Pr(>|z|)    
    (Intercept)          -2.7291763  0.1758744 -15.518  < 2e-16 ***
    indep1                0.0042787  0.3389472   0.013  0.98993    
    indep2                2.0024372  0.1746829  11.463  < 2e-16 ***
    indep3                0.3638510  0.2604196   1.397  0.16236  

如何将模型 01 的稳健标准误 (robust se) 加入到摘要函数中？这样最终我的后续计算中就可以包括正确的稳健标准误，并显示在回归表中。

提前感谢！

- P_MA1

Answer 1

lmtest::coeftest的优点在于可以使用与lm()计算出的不同协方差矩阵。

fit <- glm(am ~ hp + cyl, family=binomial(link="logit"), mtcars)

summary(fit)
#             Estimate Std. Error z value Pr(>|z|)   
# (Intercept)  5.83220    2.06595   2.823  0.00476 **
# hp           0.02775    0.01366   2.031  0.04228 * 
# cyl         -1.70306    0.60286  -2.825  0.00473 **

如果你将 lm() 的初始协方差矩阵放入 coeftest，标准误差就会相同。

library(sandwich); library(lmtest)

coeftest(fit, vcov.=vcov(fit))
#              Estimate Std. Error z value Pr(>|z|)   
# (Intercept)  5.832204   2.065949  2.8230 0.004757 **
# hp           0.027748   0.013664  2.0307 0.042282 * 
# cyl         -1.703064   0.602862 -2.8250 0.004729 **

计算稳健标准误差

现在，我们可以使用coeftest来指定一个不同的协方差矩阵。（旁注： 实际上，带有“稳健”标准误差的协方差矩阵也被称为白色标准误差，并使用type='HC0'进行指定，其他的'HC*'是改进，例如Stata使用'HC1'）。

(ct <- coeftest(fit, vcov.=vcovHC(fit, type='HC0')))
#              Estimate Std. Error z value Pr(>|z|)   
# (Intercept)  5.832204   2.139307  2.7262 0.006407 **
# hp           0.027748   0.012254  2.2643 0.023553 * 
# cyl         -1.703064   0.572045 -2.9772 0.002909 **

这就是你可以从 coeftest 中获得不同标准误的方法。

使用 texreg::screenreg 可以将纠正后的标准误放入表格中。您还可以使用 texreg::texreg 生成 LaTeX 或使用 texreg::htmlreg 生成 HTML，但基本结构保持不变。

这里是使用传统标准误的表格：

texreg::screenreg(fit)
# =========================
#                 Model 1  
# -------------------------
# (Intercept)       5.83 **
#                  (2.07)  
# hp                0.03 * 
#                  (0.01)  
# cyl              -1.70 **
#                  (0.60)  
# -------------------------
# AIC              34.63   
# BIC              39.02   
# Log Likelihood  -14.31   
# Deviance         28.63   
# Num. obs.        32      
# =========================
# *** p < 0.001; ** p < 0.01; * p < 0.05

在使用鲁棒性标准误时，需要覆盖p值和标准误。

texreg::screenreg(fit, 
                  override.pvalues=ct[, 4],
                  override.se=ct[, 3])
# =========================
#                 Model 1  
# -------------------------
# (Intercept)       5.83 **
#                  (2.73)  
# hp                0.03 * 
#                  (2.26)  
# cyl              -1.70 **
#                 (-2.98)  
# -------------------------
# AIC              34.63   
# BIC              39.02   
# Log Likelihood  -14.31   
# Deviance         28.63   
# Num. obs.        32      
# =========================
# *** p < 0.001; ** p < 0.01; * p < 0.05

- jay.sf