提取Rme中lme4模型随机效应的后验估计和可靠区间

Question

我需要从我的模型中提取随机效应的后验估计和区间。

出于说明目的，与我正在使用的数据集相似的数据集将是基数R中的ChickWeight数据集。

我为固定效果提取后验估计和间隔的方法如下：

#load package
library(lme4)

#model
m.surv<-lmer(weight ~ Time + Diet + (1|Chick), data=ChickWeight)

#load packages
library(MCMCglmm)
library(arm)

#set up for fixed effects
sm.surv<-sim(m.surv)
smfixef.surv=sm.surv@fixef
smfixef.surv=as.mcmc(smfixef.surv)

#which gives
> posterior.mode(smfixef.surv)
(Intercept)        Time       Diet2  ... 
  8.5963329   8.7034260   5.1220436  ...
> HPDinterval(smfixef.surv)
                   lower      upper
(Intercept) -0.90309142 21.3617805
Time         8.42279728  9.0306337
Diet2       -6.84371527 35.1745980
...
attr(,"Probability")
[1] 0.95
>

当我为随机效果（Chick）尝试此操作时，在第二行代码中收到以下错误：

smranef.surv=sm.surv@ranef
smranef.surv=as.mcmc(smranef.surv)

mcmc.list（x）中的错误：参数必须是mcmc对象

关于如何修改代码以提取这些值以产生随机效果的任何建议？

对其他用户的注意：如果该模型将是MCMCglmm模型，则可以像这样提取MCMC输出的随机效应的后验模式值：

posterior.mode(sm.surv$VCV[,1])
HPDinterval(sm.surv$VCV[,1])

Answer 1

要提取随机效应的估计值和95％CI，请使用以下代码：

sm.surv <-sim(m.surv)

#between Chick variance
bChick <-sm@ranef$Chick[,,1]
bvar<-as.vector(apply(bChick, 1, var)) #between ind variance posterior distribution
bvar<-as.mcmc(bvar)
posterior.mode(bvar) #mode of the distribution
HPDinterval(bvar)

这将给您：

>     posterior.mode(bvar)
     var1 
     501.24353 
>     HPDinterval(bvar)
      lower   upper
var1 412.36042 630.201
attr(,"Probability")
[1] 0.95

这意味着估计值为501，较低的95％间隔为412，较高的95％间隔为630。