3
$\begingroup$

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?

$\endgroup$
  • $\begingroup$ 1. cov is undefined. 2. Please provide your toy set without the requirement to install packages. 3 In the first model that includes hrsElapsed, you have it both as a fixed and as a random term. $\endgroup$ – Hans Ekbrand Nov 14 '15 at 15:31
  • $\begingroup$ Thank you Hans. Yes I must have been underslept when i posted this. I think I might start a new post with a new dataset. $\endgroup$ – llewmills Jan 2 '16 at 11:37
2
$\begingroup$

I the following model might be good, but the toy data set may be too small for this kind of model. There is a general (linear) effect of hrsElapsed, and this effect is allowed to be different between individuals.

anDf_randomSlopes <- lmer(score ~ group*dose + hrsElapsed + (1|subject) + (1+hrsElapsed|subject), data = df)

However, I would simplify it to this:

anDf_randomSlopes <- lmer(score ~ group*dose + hrsElapsed + (1+hrsElapsed|subject), data = df)

...which I believe is exactly the same model. In the first model, the random effect of subject is included twice.

The model below probably fails to converge because you include the slope of hrsElapsed as a random effect but not as a fixed effect:

nDf_randomSlopes <- lmer(score ~ group*dose + (1|subject) + (1+hrsElapsed|subject), data = df)

I think there must be a fixed effect whose slope can vary between individuals.

As for the model below, the problem is that you use group|subject as a random slope, which gives 36 different combinations (1/group 0, subject 1/group 1, subject 1/group 2, subject 2/group 0 etc.) plus 12 random intercepts for the subjects. The model will not run because you need to have more observations than random effects to be estimated even though only 12 of the 36 random slope combinations are present (subject 1 always belongs to group 0, subject 2 always belongs to group 1 etc).

anDf_randomSlopes <- lmer(score ~ group*dose + (1|subject) + (1+group|subject), data = df)
$\endgroup$
1
$\begingroup$

About the error message for hrsElapsed random slope, you need to rescale your continuous variables before using them to estimate random coefficients (for instance, using the scale function).

Check this guide to troubleshoot convergence warnings for lme4: https://rstudio-pubs-static.s3.amazonaws.com/33653_57fc7b8e5d484c909b615d8633c01d51.html

$\endgroup$
  • 1
    $\begingroup$ Can you expand on exactly why that is necessary? $\endgroup$ – mdewey Nov 8 '16 at 21:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.