在 lme 模型拟合中提取每个单元的系数及其标准误差 [英] Extracting coefficients and their standard error for each unit in an lme model fit
问题描述
如何在线性混合模型中提取系数(b0 和 b1)及其各自的标准误差(图):
How could I extract coefficients (b0 and b1) with their respectively standard errors for each experimental unit (plot )in a linear mixed model such as this one:
使用相同的数据集(df)，以及拟合模型(fitL1):我怎么能得到这样的数据框...
with this same dataset(df), and for the fitted model (fitL1): how could I get a data frame as this one...
plot b0 b0_se b1 b1_se
1 2898.69 53.85 -7.5 4.3
... ... ... ... ...
推荐答案
第一个评论是，这实际上是一个不平凡的理论问题:有一个相当r-sig-mixed-models 上的长线程，其中涉及一些技术细节；你绝对应该看看，即使它有点吓人.基本问题是每个组的估计系数值是该组的固定效应参数和 BLUP/条件模式的总和，它们是不同类别的对象(一个是参数，一个是条件均值)随机变量)，这带来了一些技术困难.
The first comment is that this is actually a non-trivial theoretical question: there is a rather long thread on r-sig-mixed-models that goes into some of the technical details; you should definitely have a look, even though it gets a bit scary. The basic issue is that the estimated coefficient values for each group are the sum of the fixed-effect parameter and the BLUP/conditional mode for that group, which are different classes of objects (one is a parameter, one is a conditional mean of a random variable), which creates some technical difficulties.
第二点是(不幸的是)我不知道在 lme 中有什么简单的方法可以做到这一点，所以我的答案使用了 lmer(来自 lmer>lme4 包).
The second point is that (unfortunately) I don't know of an easy way to do this in lme, so my answer uses lmer (from the lme4 package).
如果您愿意做最简单的事情并且忽略固定效应参数和 BLUP 之间的(可能定义不明确的)协方差，您可以使用下面的代码.
If you are comfortable doing the easiest thing and ignoring the (possibly ill-defined) covariance between the fixed-effect parameters and the BLUPs, you can use the code below.
两种选择是 (1) 使用贝叶斯分层方法(例如 MCMCglmm 包)拟合您的模型并计算每个级别的后验预测的标准偏差 (2) 使用参数自举计算 BLUP/条件模式，然后取自举分布的标准偏差.
Two alternatives would be (1) to fit your model with a Bayesian hierarchical approach (e.g. the MCMCglmm package) and compute the standard deviations of the posterior predictions for each level (2) use parametric bootstrapping to compute the BLUPs/conditional modes, then take the standard deviations of the bootstrap distributions.
请记住，像往常一样，此建议不提供任何保证.
library(lme4)
fm1 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
cc <- coef(fm1)$Subject
## variances of fixed effects
fixed.vars <- diag(vcov(fm1))
## extract variances of conditional modes
r1 <- ranef(fm1,condVar=TRUE)
cmode.vars <- t(apply(cv <- attr(r1[[1]],"postVar"),3,diag))
seVals <- sqrt(sweep(cmode.vars,2,fixed.vars,"+"))
res <- cbind(cc,seVals)
res2 <- setNames(res[,c(1,3,2,4)],
c("int","int_se","slope","slope_se"))
## int int_se slope slope_se
## 308 253.6637 13.86649 19.666258 2.7752
## 309 211.0065 13.86649 1.847583 2.7752
## 310 212.4449 13.86649 5.018406 2.7752
## 330 275.0956 13.86649 5.652955 2.7752
## 331 273.6653 13.86649 7.397391 2.7752
## 332 260.4446 13.86649 10.195115 2.7752
这篇关于在 lme 模型拟合中提取每个单元的系数及其标准误差的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
