Prédiction sur les modèles à effets mixtes: que faire des effets aléatoires?

13

Prenons cet ensemble de données hypothétique:

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)

nous pouvons utiliser lmepour modéliser la réponse avec un modèle à effet aléatoire:

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

Je voudrais utiliser predictsur le résultat de ce modèle pour obtenir, par exemple, la réponse d'un sujet générique du groupe 1 à une dose de 10:

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

Cependant, avec ce code, j'obtiens l'erreur suivante:

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

Pour m'en débarrasser, je dois faire, par exemple

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

Cependant, cela n'a pas vraiment de sens pour moi ... le sujet est un facteur de nuisance dans mon modèle, alors dans quel sens faut-il l'inclure predict? Si je mets un numéro de sujet non présent dans l'ensemble de données, predictretourne NA.

Est-ce le comportement recherché predictdans cette situation? Suis-je en train de manquer quelque chose de vraiment évident?

Nico
la source
modelXβ+ZγyN(Xβ+Zγ,σ2je)Z ). C'est pourquoi l'ajustement () vous donne des résultats "avec des nuisances" en premier lieu. (Et en fait je ne pense pas que ce soit gênant mais plutôt un extra info mais OK ...)
@ user11852: juste pour clarifier, je pense à cela comme un modèle à utiliser, par exemple, en cas de mesures répétées pour le même sujet.
nico
Étant donné que vous recherchez des estimations pour le même sujet, pourquoi l'exclure alors? Si vous voulez des estimations au niveau de la population (nonZl'information) alors c'est une question de différence. Comme le dit Greg dans sa réponse, vous pouvez obtenir des estimations pour la population si vous le souhaitez, mais cela ne sera pas spécifique au sujet.
usεr11852
2
@ user11852: Je ne recherche pas d'estimation sur le même sujet. Je fais des mesures répétées sur divers sujets afin d'obtenir des estimations de population. Je ne me soucie pas des sujets que j'ai déjà testés car j'ai déjà la réponse expérimentale ... Je veux pouvoir prédire comment un nouveau sujet d'un groupe spécifique répondra au stimulus. La réponse de Greg résout en effet le problème.
nico

Réponses:

17

Si vous regardez l'aide, predict.lmevous verrez qu'elle a un levelargument qui détermine à quel niveau faire les prédictions. La valeur par défaut est la plus élevée ou la plus intérieure, ce qui signifie que si vous ne spécifiez pas le niveau, elle essaie de prédire au niveau du sujet. Si vous spécifiez level=0dans le cadre de votre premier predictappel (sans subject), il donnera la prédiction au niveau de la population et n'aura pas besoin d'un numéro de sujet.

Greg Snow
la source