Previsão em modelos de efeitos mistos: o que fazer com efeitos aleatórios?

13

Vamos considerar este conjunto de dados hipotéticos:

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)

podemos usar lmepara modelar a resposta com um modelo de efeito aleatório:

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

Eu gostaria de usar predicto resultado desse modelo para obter, por exemplo, a resposta de um sujeito genérico do grupo 1 a uma dose de 10:

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

No entanto, com este código, recebo o seguinte erro:

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

Para me livrar disso, preciso fazer, por exemplo

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

Isso, no entanto, realmente não faz muito sentido para mim ... o assunto é um fator incômodo no meu modelo, então em que sentido ele tem que incluí-lo predict? Se eu colocar um número de assunto não presente no conjunto de dados, predictretornará NA.

Esse é o comportamento desejado predictnessa situação? Estou perdendo algo realmente óbvio?

nico
fonte
modelXβ+ZγyN(Xβ+Zγ,σ2Eu)Z
usεr11852 diz Reinstate Monic
@ user11852: apenas para esclarecer, estou pensando nisso como um modelo a ser usado, por exemplo, no caso de medidas repetidas para o mesmo assunto.
Nico 03/03
Z
2
@ user11852: Não estou procurando uma estimativa sobre o mesmo assunto. Eu faço repetidas medidas em vários assuntos para obter estimativas populacionais. Não ligo para os assuntos que já testei, pois já tenho a resposta experimental ... Quero ser capaz de prever como um novo sujeito de um grupo específico responderá ao estímulo. A resposta de Greg resolve o problema de fato.
Nico 04/04

Respostas:

17

Se você olhar para a ajuda predict.lme, verá que ele possui um levelargumento que determina em qual nível fazer as previsões. O padrão é o mais alto ou o mais interno, o que significa que, se você não especificar o nível, ele estará tentando prever no nível do assunto. Se você especificar level=0como parte de sua primeira predictchamada (sem subject), ela fornecerá a previsão no nível da população e não precisará de um número de assunto.

Greg Snow
fonte