5

私は、回帰モデルの標準誤差を独自の標準誤差に直接置き換える方法を探しています。これは、独自の堅牢なオプションが付属しておらず、特定の種類のモデルのみを供給することができる別の 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 関数を使用しています (コードはここ(p. 4) にあります)。

私の考えは、正しい標準誤差を既に報告しているモデルをパッケージの mediate 関数に与えることです。ドキュメントによると、メディエーション パッケージは次のクラスのモデルを受け入れます: lm、polr、bayespolr、glm、bayesglm、gam、rq、survreg、および merMod。

4

2 に答える 2

0

@MySeppoの答えは素晴らしいですが、一般性を高めるためにすべてを堅牢に変更する関数(SE、t-stats、p-values)をここに与える価値があると思います(Stataと一致させるためにHC1も使用しています):

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 に答える