I am just dipping my toe into the ocean that is linear effects models and am working through Barr et al.'s 'Keeping it Maximal' paper, trying to figure out the best way to fit a lmem for my experiment. Say you have three groups given three different types of drug over three days: 100mg on first day, 50mg on second day, 10mg on last day. The outcome measure is how they feel that next day on some scale (e.g. mood), before they are given their daily dose (i.e. so we are measuring the effects of the previous day's dose). However participants don't come in at exactly the same time each day, thus as each time of measurement the drug will have had less time to take effect.
I would like to know how best to include that random effect of 'time elapsed since dose' into this model, and just how best to fit the model really.
This is a toy dataset. I have not built any trends into it.
dose100mg <- c(6,2,9,4,6,5,2,4,6,7,3,2)
dose50mg <- c(1,2,4,3,6,1,3,3,2,1,4,1)
dose10mg <- c(8,9,7,9,6,7,8,9,8,7,1,3)
timeD1 <- c(24.2,20.5,26,30,22,26,19,23,29,30,24,16)
timeD2 <- c(24,16,28,20,19,28,30,20,18,15,27,32)
timeD3 <- c(21,28,29,30,29,17,23,18,24,16,28,21)
subject <- c(1,2,3,4,5,6,7,8,9,10,11,12)
group <- factor(c(0,1,0,2,1,2,0,2,1,2,1,0))
df <- data.frame(subject, group, dose100mg, dose50mg, dose10mg, cov)
Turn it from wide to long
require(tidyr)
df <- gather(df, dose, score, dose100mg:dose50mg:dose10mg)
Now add the 'hours elapsed since last dose' variable to the dataframe (btw: if anyone knows how to build this into the gather function above I'd appreciate it)
df$hrsElapsed <- c(timeD1, timeD2, timeD3)
Now fit a model. First group*dose plus with random intercepts for subject.
require(lme4)
# random intercepts
anDf_randomintercepts <- lmer(score ~ group*dose + (1|subject), data = df)
anova(anDf_randomintercepts)
Next random slopes, and my first question. Is it better to include hrsElasped as a covariate, like this?
anDf_randomSlopes <- lmer(score ~ group*dose + hrsElapsed + (1|subject) + (1+hrsElapsed|subject), data = df)
Or to include it as a random effect? like this
nDf_randomSlopes <- lmer(score ~ group*dose + (1|subject) + (1+hrsElapsed|subject), data = df)
I know it's not the latter because I get an error message.
Warning messages:
1: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0450795 (tol = 0.002, component 1)
2: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model is nearly unidentifiable: large eigenvalue ratio
- Rescale variables?
But I don't know WHY this doesn't work. I would have thought time elapsed would be exactly the sort of variable you'd want to assign to random effects.
What am I doing wrong?
An ancillary question pertains to fitting random slopes for the group-by-subject effect
anDf_randomSlopes <- lmer(score ~ group*dose + (1|subject) + (1+group|subject), data = df)
When I run this i get the error message
Error: number of observations (=36) <= number of random effects (=36) for term (1 + group | subject); the random-effects parameters and the residual variance (or scale parameter) are probably unidentifiable
Why doesn't this work? What does it mean?