从lme4 mer模型对象中提取随机效应方差

dyn*_*amo 41 random r effects lme4

我有一个具有固定和随机效果的mer对象.如何提取随机效应的方差估计?这是我的问题的简化版本.

study <- lmer(Reaction ~ Days + (1|Subject), data = sleepstudy)
study
Run Code Online (Sandbox Code Playgroud)

这样可以提供长输出 - 在这种情况下不会太长.无论如何,我如何明确选择

Random effects:
Groups   Name        Variance Std.Dev.
Subject  (Intercept) 1378.18  37.124  
Residual              960.46  30.991  
Run Code Online (Sandbox Code Playgroud)

部分输出?我想要自己的价值观.

我长期看看

str(study)
Run Code Online (Sandbox Code Playgroud)

那里什么都没有!还检查了lme4包中的任何提取器功能都无济于事.请帮忙!

Ben*_*ker 72

其他一些答案是可行的,但我声称最好的答案是使用专为此设计的访问器方法 - VarCorr(这与lme4前一个nlme程序包相同).

更新版本lme4(版本1.1-7,但下面的所有内容可能适用于版本> = 1.0),VarCorr比以前更灵活,应该做你想要的一切,而不必在拟合的模型对象内钓鱼.

library(lme4)
study <- lmer(Reaction ~ Days + (1|Subject), data = sleepstudy)
VarCorr(study)
##  Groups   Name        Std.Dev.
##  Subject  (Intercept) 37.124  
##  Residual             30.991
Run Code Online (Sandbox Code Playgroud)

默认情况下VarCorr()打印标准偏差,但如果您愿意,可以改为:

print(VarCorr(study),comp="Variance")
##  Groups   Name        Variance
##  Subject  (Intercept) 1378.18 
##  Residual              960.46 
Run Code Online (Sandbox Code Playgroud)

(comp=c("Variance","Std.Dev.")将打印两者).

为了更加灵活,您可以使用该as.data.frame方法转换VarCorr对象,该对象提供分组变量,效果变量,方差/协方差或标准差/相关性:

as.data.frame(VarCorr(study))
##        grp        var1 var2      vcov    sdcor
## 1  Subject (Intercept) <NA> 1378.1785 37.12383
## 2 Residual        <NA> <NA>  960.4566 30.99123
Run Code Online (Sandbox Code Playgroud)

最后,VarCorr对象的原始形式(如果你不需要,你可能不应该把它弄乱)是一个方差 - 协方差矩阵列表,其中附加(冗余)信息编码标准偏差和相关性,以及attributes("sc")给出残差标准差并指定模型是否具有估计的尺度参数("useSc").

unclass(VarCorr(fm1))
## $Subject
##             (Intercept)      Days
## (Intercept)  612.089748  9.604335
## Days           9.604335 35.071662
## attr(,"stddev")
## (Intercept)        Days 
##   24.740448    5.922133 
## attr(,"correlation")
##             (Intercept)       Days
## (Intercept)  1.00000000 0.06555134
## Days         0.06555134 1.00000000
## 
## attr(,"sc")
## [1] 25.59182
## attr(,"useSc")
## [1] TRUE
## 
Run Code Online (Sandbox Code Playgroud)

  • (1)很容易平衡标准偏差; (2)`print(VarCorr(fitted_model),comp ="Variance")`或`as.data.frame(VarCorr(fitted_model))`将很容易检索方差; (3)报告差异与标准偏差是依赖于上下文的 - 如果试图考虑var分解/比例解释和std devs,如果试图与固定效应的大小进行比较,我通常更喜欢差异 (3认同)

Tom*_*mmy 14

lmer 返回一个S4对象,所以这应该工作:

remat <- summary(study)@REmat
print(remat, quote=FALSE)
Run Code Online (Sandbox Code Playgroud)

哪个印刷品:

 Groups   Name        Variance Std.Dev.
 Subject  (Intercept) 1378.18  37.124  
 Residual              960.46  30.991  
Run Code Online (Sandbox Code Playgroud)

...通常,您可以查看"mer"对象的来源printsummary方法:

class(study) # mer
selectMethod("print", "mer")
selectMethod("summary", "mer")
Run Code Online (Sandbox Code Playgroud)

  • 如果你想要这些值,那么VarCorr()会更有效率.看看Ben Bolker的帖子 (6认同)
  • 现在这有点过时了(尽管最初的问题确实指的是“mer 对象”,根据定义,它们与 1.0 之前的 `lme4` 相关 - 该类现在称为 `merMod`。 (2认同)

Joh*_*lby 5

> attributes(summary(study))$REmat
 Groups     Name          Variance  Std.Dev.
 "Subject"  "(Intercept)" "1378.18" "37.124"
 "Residual" ""            " 960.46" "30.991"
Run Code Online (Sandbox Code Playgroud)

  • 我可能错了,因为在 `attributes(summary(study))` 中似乎没有 `REmat` (2认同)