【发布时间】：2017-03-02 09:14:16
【问题描述】：
我正在使用多层次模型来尝试描述纵向变化的不同模式。 Dingemanse et al (2010) 描述了随机效应完全相关时的“扇出”模式。然而，我发现当随机效应之间的关系是非线性的但在观察到的时间间隔内单调增加时，会出现类似的模式。在这种情况下，随机效应不是完全相关的，而是由函数描述的。 请参阅下面的示例以了解此情况。该示例仍然具有较高的截距-斜率相关性 (>.9)，但有可能获得低于 0.7 的相关性，同时仍保持完美的截距-斜率关系。
我的问题：有没有办法使用 nlme 或其他 R 包在多级模型中包含完美（非线性）随机效应之间的关系？ MLwiN 有一种方法来约束斜率截距协方差，这将是一个开始，但不足以包括非线性关系。到目前为止，我一直无法找到 nlme 的解决方案，但也许您知道其他一些可以将其包含在模型中的包。
为草率的编码道歉。我希望我的问题足够清楚，但如果有什么需要澄清的，请告诉我。非常感谢任何帮助或替代解决方案。
set.seed(123456)
# Change function, quadratic
# Yit = B0ij + B1ij*time + B2ij*time^2
chn <- function(int, slp, slp2, time){
score<-int + slp * time+ slp2 * time^2
return(score)
}
# Set N, random intercept, time and ID
N<-100
start<-rnorm(N,100,15) # Random intercept
time<- matrix(1:15,ncol = 15, nrow = 100,byrow = T) # Time, balanced panel data
ID<-1:N # ID variable
# Random intercept, linear slope: exp(intercept/25)/75, quadratic slope: exp(intercept/25)/250
score3<-matrix(NA,ncol = ncol(time), nrow = N)
for(x in ID){
score3[x,]<-chn(start[x],exp(start[x]/25)/75,exp(start[x]/25)/250,time[x,])
}
#Create dataframe
df<- data.frame(ID,score3)
df2<- melt(df,id = 'ID')
df2$variable<-as.vector(time)
# plot
ggplot(df2, aes(x= variable, y= value)) + geom_line(aes(group = ID)) +
geom_smooth(method = "lm", formula = y ~ x + I(x^2), se =F, size = 2, col ="red" )
# Add noise and estimate model
df2$value2<-df2$value + rnorm(N*ncol(time),0,2)
# Random intercept
mod1<-lme(value2 ~ variable + I(variable^2),
random= list(ID = ~1),
data=df2,method="ML",na.action=na.exclude)
summary(mod1)
# Random slopes
mod2<-update(mod1,.~.,random= list(ID = ~variable + I(variable^2)))
summary(mod2)
pairs(ranef(mod2))
【问题讨论】：
-
总是可以去贝叶斯...
-
我正在考虑使用贝叶斯算法，但对 R2WinBUGS 和 R2jags 的经验有限。如果你能给我一个例子，将不胜感激。
-
@MattTyers 花了一些时间，但我尝试使用 rjags 包。欢迎任何反馈
标签： r mixed-models multi-level nlme r2jags