5

我正在寻找一种方法,用我自己的标准误差直接替换回归模型中的标准误差,以便在另一个 R 包中使用鲁棒模型,该 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]

但是,以下方法在尝试替换模型的标准错误时会返回错误,因为它将“摘要”本身视为命令:

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


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

细节

我正在尝试使用 Mediation R 包运行调解分析。该包将使用“调解”功能执行一种方式聚集标准错误,但我想要两种方式聚集标准错误。

为了获得两种聚集标准错误,我使用了 Mahmood Arai 的 mclx 函数(代码可以在这里找到(第 4 页)。

我的想法是为包的调解功能提供已经报告正确标准错误的模型。根据文档,中介包接受以下模型类别:lm、polr、bayespolr、glm、bayesglm、gam、rq、survreg 和 merMod。

4

2 回答 2

3

对我有用的是:

    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
于 2017-02-22T20:28:34.867 回答
0

@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
}
于 2021-07-10T02:52:37.087 回答