线性混合效应模型传统上以下列方式配制 . Ri = Xi×βZi×biεi其中β表示估计的固定效应，Z表示随机效应 . 因此，X是经典的设计矩阵 . 使用R，我希望能够在使用来自nlme包的lme拟合模型后提取这两个矩阵 . 例如，也可以在nlme包中找到的数据集“Rails”包含在6个随机选择的铁路轨道上的三个单独的超声波传播时间测量值 . 我可以使用截距固定效果和每个轨道的随机效果拟合一个简单的模型，具体如下 .
library(nlme)
lmemodel<-lme(travel ~ 1, random = ~ 1 | Rail, data=Rail)
X设计矩阵只是一个18x1矩阵（6个轨道* 3个测量）的统一，可以通过以下方式轻松提取：
model.matrix(lmemodel, data=Rail)
(Intercept)
1 1
2 1
3 1
4 1
5 1
6 1
7 1
8 1
9 1
10 1
11 1
12 1
13 1
14 1
15 1
16 1
17 1
18 1
attr(,"assign")
[1] 0
我想要做的是提取随机效果设计矩阵Z.我意识到如果我使用lme4包适合相同的模型，这可以通过以下方式完成：
library(lme4)
lmermodel<-lmer(travel ~ 1 + (1|Rail),data=Rail)
t(lmermodel@Zt) ##takes the transpose of lmermodel@Zt
lmermodel@X ## extracts the X matrix
但是，我不知道如何从lme拟合模型中提取这个矩阵 .
2 回答
1对于这个例子是特定的，因为只有一个随机效应 . 当你有多个随机效果时，你可以做一些更自动的编程来将不同的Z_i叠加在一起 .
据我所知，
Z矩阵未存储在lme对象中的任何位置 . 最好的选择是在modelStruct$reStruct组件中（尝试names(modelfit); str(modelfit); sapply(modelfit,class)等进行探索），但据我所知，它并不存在 . 事实上，一些人认为lme.default可能永远不会明确构建Z矩阵 . 内部lme似乎适用于分组结构 . 你当然可以但那可能不是你想到的......