在 R 中替换 Reg 模型中的标准错误

Question

我正在寻找一种方法来直接用我自己的标准误差替换回归模型中的标准误差，以便在另一个 R 包中使用稳健模型，该包没有自己的稳健选项，只能被喂养特定类型的模型而不是 coeftest 格式。

假设我有一个线性模型：

model <- lm(data=dat, Y ~ X1 + X2 + X3)

然后我想要稳健的标准误差：

robust <- coeftest(model, vcov=sandwich)

接下来，我需要在一个特定的包中使用这个模型，这个包不能提供 coeftest 并且没有自己强大的标准错误选项。我想在将原始模型的标准误差（以及与此相关的 p 值、t 统计量等）输入包之前替换它们，以便将它们计算在内。

要访问 原始模型 中的标准误差，我使用：

summary(model)$coefficients[,2]

要从 coeftest 中提取标准误差，我使用：

coeftest.se <- robust[, 2]

然而，以下方法在尝试替换模型的标准错误时returns出错，因为它将"summary"本身视为命令：

summary(model)$coefficients[,2] <- coeftest.se


Error in summary(M3)$coefficients[, 2] <- seM3 : could not find function "summary<-"

细节

我正在尝试运行使用 Mediation R 包进行中介分析。该包将使用 "mediate" 函数执行一种方式的集群标准错误，但我想要两种方式的集群标准错误。

为了获得两种方式的聚类标准误差，我正在使用 Mahmood Arai 的 mclx 函数（可以找到代码 here（第 4 页）。

我的想法是将已经报告正确标准错误的模型提供给包的中介函数。根据文档，中介包接受以下类模型：lm、polr、bayespolr、glm、bayesglm、gam、rq、survreg 和 merMod。

Answer 1

对我有用的是：

    model <- lm(data=dat, Y ~ X1 + X2 + X3)

    library(sandwich)

    SE_robust <- sqrt(diag(vcovHC(model, type="HC2")))

    model2 <- summary(model)

    model2$coefficients[,2] <- SE_robust

Answer 2

@MySeppo 的回答很好，但我认为值得在这里提供一个函数，它可以将所有内容更改为健壮的（SE、t-stats、p-值）以获得更多的通用性（我也使用 HC1 来与 Stata 一致):

model <- lm(mpg~cyl,data=mtcars)

robustsummary <- function(model) { 
  library(sandwich)
  library(lmtest)
  coeftest <- coeftest(model, vcov = vcovHC(model, type="HC1"))
  summ <- summary(model)
  summ$coefficients[,2] <- coeftest[,2]
  summ$coefficients[,3] <- coeftest[,3]
  summ$coefficients[,4] <- coeftest[,4]
  summ
}
summary(model)
robustsummary(model)

或者如果您不想使用 coeftest

robustsummary <- function(model) { 
  library(sandwich)
  SE_robust <- sqrt(diag(vcovHC(model, type="HC1")))
  summ <- summary(model)
  summ$coefficients[,2] <- SE_robust
  summ$coefficients[,3] <- summ$coefficients[,1]/SE_robust # get t-stats
  summ$coefficients[,4] <- 2*(1-pt(abs(summ$coefficients[,3]),summ$df[2])) # apply t distribution with the right degrees of freedom to get p-values
  summ
}

在 R 中替换 Reg 模型中的标准错误

Replacing Standard Errors in a Reg Model in R

r

regression

standard-error