我对 R 中的模型设置有疑问,经过长时间彻底的搜索后,我没有找到回答我的两个问题中的任何一个的线程:
我先描述一下设置:
它是一个重复测量数据集,具有两种不同的干预措施(食物和培训),每个干预措施都有两个级别。所有参与者都经历了各种条件组合。在一天中收集激素样品。年龄和 BMI 作为协变量包括在内。
一个可重现的数据集是这个:
ID <- rep(rep(c("A","B","C"),each=5),4)
TIME <- rep(paste(sprintf("%02i",9:13),"00",sep=":"),12)
training <- rep(rep(rep(c("T1","T2"),each=5),each=3),2)
food <- c(rep(rep(c("F1","F2"),each=5),each=3),rep(rep(c("F2","F1"),each=5),each=3))
hormone <- rnorm(n = 60,mean=5,sd=2)
BMI <- as.numeric(rep(rep(c(22,27,25),each=5),4))
age <- as.numeric(rep(rep(c(26,27,23),each=5),4))
DF <- as.data.frame(cbind(ID,TIME,training,food,hormone,BMI,age))
DF$BMI <- scale(as.numeric(DF$BMI))
DF$age <- scale(as.numeric(DF$age))
DF$hormone <- as.numeric(as.character(DF$hormone))
现在我有两个主要问题:
1)首先,我想评估两种干预措施对基线激素水平的相互作用。为此,我设置了以下模型:
model.df <- DF[DF$TIME=="09:00",]
m1 <- lmer(hormone~training*food+age+BMI+(1+training*food|ID),model.df)
但是,无法设置此模型,因为
Error: number of observations (=16) <= number of random effects (=16) for term (1 + training * food | ID); the random-effects parameters and the residual variance (or scale parameter) are probably unidentifiable
当我从随机效应中排除交互项时,它确实有效:
m1 <- lmer(激素~训练*食物+年龄+BMI+(1+训练+食物|ID),model.df)
我现在想知道这是否仍然是测试交互的有效模型?因此,我的空模型将是:
m0 <- lmer(hormone~training+food+age+BMI+(1+training+food|ID),model.df)
现在到第二点:
2)
我们还想监测激素随时间的变化。
但是,我不确定如何在模型中包含一天中的时间。
正如该线程中指出的那样https://stats.stackexchange.com/questions/245866/is-hour-of-day-a-categorical-variable
但是,它可以作为循环变量包含在内,因为我的采样没有涵盖一整天,所以我不确定如何在我的情况下实现它。任何人都可以帮忙吗?
另外,我不确定如何设置包括时间变量在内的模型。
我们仍然对这两种干预措施的相互作用感兴趣,所以我会设置类似以下模型的东西。(现在假设时间为数字)
model.df <- DF
model.df$TIME <- as.numeric(model.df$TIME)
m2 <- lmer(hormone~training*food*TIME+age+BMI+(1+training*food*TIME||ID),model.df)
然而:1)这个模型不收敛
Warning messages:
1: In commonArgs(par, fn, control, environment()) :
maxfun < 10 * length(par)^2 is not recommended.
2: In optwrap(optimizer, devfun, getStart(start, rho$lower, rho$pp),:
convergence code 1 from bobyqa: bobyqa -- maximum number of function evaluations exceeded
3: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0035331 (tol = 0.002, component 1)
2)你认为这是一个合适的模型吗?
然后将针对以下空模型对其进行测试:
m02 <- lmer(hormone~training*TIME+food*TIME+age+BMI+(1+training*TIME+food*TIME||ID),model.df)
我希望这些不是一个线程的太多问题,我会对任何指针感到非常高兴。非常感谢。