如何从多个组的 lme() 函数中提取结果,然后在 R 中组合?
How to extract result from lme() function from multiple groups and then combine in R?
首先,根据sl
变量将以下数据随机分为两组,然后运行两组的模型使用数据集下方所示的for循环
mydata
y x sl
1 5.297967 1 1
2 3.322833 2 1
3 4.969813 3 1
4 4.276666 4 1
5 5.972807 1 2
6 6.619440 2 2
7 8.045588 3 2
8 7.377759 4 2
9 6.907755 5 2
10 8.672486 6 2
11 8.283999 7 2
12 8.455318 8 2
13 7.414573 9 2
14 8.634087 10 2
15 7.356355 1 3
16 6.606247 2 3
17 6.396930 9 3
18 6.579251 10 3
19 5.521110 1 4
20 2.224221 2 4
21 6.742881 3 4
22 6.709304 4 4
23 6.875232 5 4
24 8.476371 6 4
25 7.360104 7 4
Runnign 模型对两组都使用 lme() 函数,然后将 beta
系数存储为矩阵,将 theta
[ 随机截距项 ] 存储为向量
sl.no=unique(mydata$sl)
m=length(unique(mydata$sl))
ngrp=2
set.seed(125)
idx=sample(1:ngrp, size=m, replace = T)
beta=matrix(NA, nrow = ngrp, ncol=3, byrow=T) #null matrix to store coefficients from both groups
theta=rep(0,m) #null vector to store intercepts from both groups
library(nlme)
for ( g in 1:ngrp){
rg=sl.no[idx==g]
mydata_rG=mydata[mydata$sl %in% rg,] #Data set belongs to group-g
lme_mod=lme(y~x+I(x^2),random = ~ 1|sl,
data = mydata_rG, method = "ML") #mixed effect model for each group
beta[g,]=c(unlist(lme_mod$coefficients[1])[[1]],
unlist(lme_mod$coefficients[1])[[2]],
unlist(lme_mod$coefficients[1])[[3]])
theta=c(unname(lme_mod$coefficients$random$sl))
}
我期待一个 theta
长度为 m 的向量。不幸的是,theta
是一个的大小。
任何帮助表示赞赏。
beta
和 theta
的结果
beta
[,1] [,2] [,3]
[1,] 4.895805 0.7954474 -0.05602771
[2,] 6.423533 -1.7441753 0.32049662
theta
[1] 4.264366e-21 #it should be length of m.
这只是关于你如何存储 theta
值
sl.no=unique(mydata$sl)
m=length(unique(mydata$sl))
ngrp=2
set.seed(125)
idx=sample(1:ngrp, size=m, replace = T)
beta=matrix(NA, nrow = ngrp, ncol=3, byrow=T)
theta=numeric() #Change here
library(nlme)
for ( g in 1:ngrp){
rg=sl.no[idx==g]
mydata_rG=mydata[mydata$sl %in% rg,]
lme_mod=lme(y~x+I(x^2),random = ~ 1|sl,
data = mydata_rG, method = "ML")
beta[g,]=c(unlist(lme_mod$coefficients[1])[[1]],
unlist(lme_mod$coefficients[1])[[2]],
unlist(lme_mod$coefficients[1])[[3]])
theta=c(theta, unname(lme_mod$coefficients$random$sl)) #Change here
}
首先,根据sl
变量将以下数据随机分为两组,然后运行两组的模型使用数据集下方所示的for循环
mydata
y x sl
1 5.297967 1 1
2 3.322833 2 1
3 4.969813 3 1
4 4.276666 4 1
5 5.972807 1 2
6 6.619440 2 2
7 8.045588 3 2
8 7.377759 4 2
9 6.907755 5 2
10 8.672486 6 2
11 8.283999 7 2
12 8.455318 8 2
13 7.414573 9 2
14 8.634087 10 2
15 7.356355 1 3
16 6.606247 2 3
17 6.396930 9 3
18 6.579251 10 3
19 5.521110 1 4
20 2.224221 2 4
21 6.742881 3 4
22 6.709304 4 4
23 6.875232 5 4
24 8.476371 6 4
25 7.360104 7 4
Runnign 模型对两组都使用 lme() 函数,然后将 beta
系数存储为矩阵,将 theta
[ 随机截距项 ] 存储为向量
sl.no=unique(mydata$sl)
m=length(unique(mydata$sl))
ngrp=2
set.seed(125)
idx=sample(1:ngrp, size=m, replace = T)
beta=matrix(NA, nrow = ngrp, ncol=3, byrow=T) #null matrix to store coefficients from both groups
theta=rep(0,m) #null vector to store intercepts from both groups
library(nlme)
for ( g in 1:ngrp){
rg=sl.no[idx==g]
mydata_rG=mydata[mydata$sl %in% rg,] #Data set belongs to group-g
lme_mod=lme(y~x+I(x^2),random = ~ 1|sl,
data = mydata_rG, method = "ML") #mixed effect model for each group
beta[g,]=c(unlist(lme_mod$coefficients[1])[[1]],
unlist(lme_mod$coefficients[1])[[2]],
unlist(lme_mod$coefficients[1])[[3]])
theta=c(unname(lme_mod$coefficients$random$sl))
}
我期待一个 theta
长度为 m 的向量。不幸的是,theta
是一个的大小。
任何帮助表示赞赏。
beta
和 theta
beta
[,1] [,2] [,3]
[1,] 4.895805 0.7954474 -0.05602771
[2,] 6.423533 -1.7441753 0.32049662
theta
[1] 4.264366e-21 #it should be length of m.
这只是关于你如何存储 theta
值
sl.no=unique(mydata$sl)
m=length(unique(mydata$sl))
ngrp=2
set.seed(125)
idx=sample(1:ngrp, size=m, replace = T)
beta=matrix(NA, nrow = ngrp, ncol=3, byrow=T)
theta=numeric() #Change here
library(nlme)
for ( g in 1:ngrp){
rg=sl.no[idx==g]
mydata_rG=mydata[mydata$sl %in% rg,]
lme_mod=lme(y~x+I(x^2),random = ~ 1|sl,
data = mydata_rG, method = "ML")
beta[g,]=c(unlist(lme_mod$coefficients[1])[[1]],
unlist(lme_mod$coefficients[1])[[2]],
unlist(lme_mod$coefficients[1])[[3]])
theta=c(theta, unname(lme_mod$coefficients$random$sl)) #Change here
}