r - 来自 lmer 输出的方差分量的标准误差

Question

我需要提取 standard error来自 lmer 输出的方差分量.

library(lme4)
model <- lmer(Reaction ~ Days + (1|Subject), sleepstudy)

以下产生方差分量的估计值:

s2 <- VarCorr(model)$Subject[1]

是 不是 方差的标准误。我想要标准错误。我怎么能拥有它？
编辑 :
也许我无法让您理解我所说的“方差分量的标准误差”是什么意思。所以我正在编辑我的帖子。
在 Design and Analysis of Experiments 一书中的第 12 章“随机因素实验”中由 Douglas C. Montgomery 撰写，在本章末尾，示例 12-2 由 SAS 完成。在例 12-2 中，该模型是一个二因子因子随机效应模型。输出如表 12-17 所示

我正在尝试通过 lmer 在 R 中拟合模型.

library(lme4)
fit <- lmer(y~(1|operator)+(1|part),data=dat)

用于提取 Estimate 的 R 代码，在表12-17中用4注释:

est_ope=VarCorr(fit)$operator[1]
est_part = VarCorr(fit)$part[1]
sig = summary(fit)$sigma
est_res = sig^2

现在我想提取 Std Errors的结果，由 lmer 输出的表 12-17 中的 5 注释。
非常感谢 !

Answer 1

我认为您正在寻找方差估计的 Wald 标准误差。请注意这些(Doug Bates 经常指出的)Wald 标准误差通常是 很差方差不确定性的估计，因为似然分布在方差尺度上通常远非二次......我假设你知道你在做什么并且对这些数字有很好的利用......
这可以(现在)使用 merDeriv 完成包裹。

library(lme4)
library(merDeriv)
m1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
sqrt(diag(vcov(m1, full = TRUE)))
vv <- vcov(m1, full = TRUE)
colnames(vv)
## [1] "(Intercept)"                  "Days"                        
## [3] "cov_Subject.(Intercept)"      "cov_Subject.Days.(Intercept)"
## [5] "cov_Subject.Days"             "residual"

如果我们想要方差分量的标准误差，我们取对角线的平方根，只保留最后三个元素:

sqrt(diag(vv)[3:5])
## [1] 288.78602  46.67876  14.78208

旧答案

library("lme4")
model <- lmer(Reaction ~ Days + (1|Subject), sleepstudy, REML=FALSE)

(目前对于 REML 估计要做到这一点有点困难......)
提取根据标准偏差和相关性而不是 Cholesky 因子参数化的偏差函数(注意这是一个内部函数，所以不能保证它将来会以相同的方式工作......)

 dd.ML <- lme4:::devfun2(model,useSc=TRUE,signames=FALSE)

提取参数作为原始尺度的标准偏差:

 vv <- as.data.frame(VarCorr(model)) ## need ML estimates!
 pars <- vv[,"sdcor"]
 ## will need to be careful about order if using this for
 ## a random-slopes model ...

现在计算二阶导数(Hessian)矩阵:

library("numDeriv")
hh1 <- hessian(dd.ML,pars)
vv2 <- 2*solve(hh1)  ## 2* converts from log-likelihood to deviance scale
sqrt(diag(vv2))  ## get standard errors

这些是标准差的标准误:将它们加倍以获得方差的标准误(当您变换一个值时，其标准误根据变换的导数进行缩放)。
我认为这应该这样做，但您可能需要仔细检查它......