提取nlme中的隨機效應設計矩陣

線性混合效應模型傳統上按照以下方式制定。 Ri = Xi×β+ Zi×bi +εi 其中β代表估計的固定效應，Z代表隨機效應。 X因此是經典的設計矩陣。使用R，我希望能夠從nlme包中使用lme擬合模型後提取這兩個矩陣。例如，在nlme軟件包中也可以找到數據集「Rails」，其中包含對隨機選擇的6條軌道上的超聲波傳播時間的三次單獨測量。我可以用一個簡單的模型來擬合每個鐵路的攔截固定效應和隨機效應，並給出以下幾點。提取nlme中的隨機效應設計矩陣

library(nlme) 
lmemodel<-lme(travel ~ 1, random = ~ 1 | Rail, data=Rail)

的X設計矩陣僅僅是一個18×1矩陣（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([email protected]) ##takes the transpose of [email protected] 
[email protected] ## extracts the X matrix

不過，我很茫然，以如何提取這馬特里x從一個lme擬合的模型。

來源

2013-10-31 iantist

如果內部表示在未來發生變化，最好使用getME（lmermodel，「Z」）或getME（lmermodel，「X」）「，而不是直接訪問對象中的槽。 –

感謝您的提示。 – iantist

據我所見，Z矩陣不存儲在lme對象的任何地方。最好的選擇是在modelStruct$reStruct組件（嘗試names(modelfit); str(modelfit); sapply(modelfit,class)等）來探索，但它不在我所知的範圍內。事實上，一些暗示lme.default的內容暗示矩陣可能永遠不會被明確地構造;內部lme似乎與分組結構一起使用。當然，你可以做

Z <- model.matrix(~Rail-1,data=Rail)

但是這可能不是你腦子裏想什麼......

來源

2013-10-31 16:32:42

model.matrix(formula(lmemodel$modelStruct$reStr)[[1]],data=lmemodel$data)

1是有點兒具體到這個例子，因爲只有一個隨機效應。當你有多個隨機效果時，你可以做更多的自動編程來將不同的Z_i堆疊在一起。

來源

2013-11-23 03:13:07 user3023981

提取nlme中的隨機效應設計矩陣

回答

相關問題