I'm taking a GLMM course and working through some homework that conceptually isn't that difficult, but is very challenging with respect to implementation in R.

The dataset is composed of measurements taken on 30 subjects at 5 different time points, so the overall design is RCBD. The time points are treated as a continuous measurement and the individual subjects are treated as a factor. For a linear model, we have individual acting as a blocking factor, while in the GLMM model, individual is a random effect.

I can compute the models easily enough, but then I have questions like:

  • The mean weight at 15 days for individual 3.
  • The mean weight at 15 days for the average individual.
  • The difference in mean weight at 15 days between individual 3 and the overall mean.

While I can calculate these using predict() using the fit models and artificial data, this is not very satisfying as it doesn't provide a true contrast or inference for the model, especially for the last question above in which I'm only subtracting the two previous answers. I've tried to use contrasts, but I cannot find examples of creating contrasts when I have a two variable model. I assume that if I can create a proper contrast that is acceptable to emmeans(), I can have a better answer to the question that provides some confidence in the form of standard errors, etc.

Below are a couple of examples of my current code, but I have not provided the data as it is a bit long to include here.

Fixed Effects example for first question, above:

fit.fixed <- lm(weight ~ day + id, data=dat)

nd.1b <- data.frame( id = factor('3', levels(dat$id)), 
              day = 15)
predict(fit.fixed, newdata = nd.1b)

#alternatively, using emmeans
mns <-  emmeans(fit.fixed, specs='day', by='id', at=list(day=15))

#how can I create a proper contrast for this?

Random Effects Example for first question above

fit.random.reml <- lmer(weight ~ day + (day | id), data=dat)

nd.3b <- data.frame( id = factor('3', levels(dat$id)), 
                     day = 15)

predict(fit.random.reml, newdata = nd.3b)

Your Answer

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

Browse other questions tagged or ask your own question.