Prediction on mixed effect models: what to do with random effects?

Let’s consider this hypothetical dataset:

set.seed(12345)

num.subjects <- 10

dose <- rep(c(1,10,50,100), num.subjects)
subject <- rep(1:num.subjects, each=4)
group <- rep(1:2, each=num.subjects/2*4)

response <- dose*dose/10 * group + rnorm(length(dose), 50, 30)

df <- data.frame(dose=dose, response=response, 
                 subject=subject, group=group)

we can use lme to model the response with a random effect model:

require(nlme)
model <- lme(response ~ dose + group + dose*group, 
             random = ~1|subject, df)

I would like to use predict on the result of this model to get, for instance, the response of a generic subject of group 1 to a dose of 10:

pred <- predict(model, newdata=list(dose=10, group=1))

However, with this code I get the following error:

Error in predict.lme(model, newdata = list(dose = 10, group = 1)) : 
cannot evaluate groups for desired levels on 'newdata'

To get rid of it I need to do, for instance

pred <- predict(model, newdata=list(dose=10, group=1, subject=5))

This, however, does not really make much sense to me… the subject is a nuisance factor in my model, so what sense does it have to include it in predict? If I put a subject number not present in the dataset, predict returns NA.

Is this the wanted behaviour for predict in this situation? Am I missing something really obvious?

Answer

Attribution
Source : Link , Question Author : nico , Answer Author : Greg Snow