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
lmeto 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
predicton 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,
Is this the wanted behaviour for
predictin this situation? Am I missing something really obvious?
If you look at the help for
predict.lme you will see that it has a
level argument that determines which level to make the predictions at. The default is the highest or innermost which means that if you don’t specify the level then it is trying to predict at the subject level. If you specify
level=0 as part of your first
predict call (without
subject) then it will give the prediction at the population level and not need a subject number.