Dans quelle mesure la somme des effets aléatoires doit-elle être proche de zéro dans GLMM (avec lme4)

8

J'utilise le lme4package dans R pour faire de la modélisation logistique à effets mixtes.
J'ai cru comprendre que la somme de chaque effet aléatoire devait être nulle.

Lorsque je crée des modèles mixtes linéaires de jouets à l'aide de lmer, les effets aléatoires sont généralement < confirmant ma conviction que le Mais dans les modèles binomiaux jouets (et dans les modèles de mes données binomiales réelles) une partie de la somme des effets aléatoires à ~ 0.9.1010colSums(ranef(model)$groups) ~ 0

Dois-je m'inquiéter? Comment interpréter cela?

Voici un exemple de jouet linéaire

toylin<-function(n=30,gn=10,doplot=FALSE){
 require(lme4)
 x=runif(n,0,1000)
 y1=matrix(0,gn,n)
 y2=y1
 for (gx in 1:gn)
 {
   y1[gx,]=2*x*(1+(gx-5.5)/10) + gx-5.5  + rnorm(n,sd=10)
   y2[gx,]=3*x*(1+(gx-5.5)/10) * runif(1,1,10)  + rnorm(n,sd=20)
 }
 c1=y1*0;
 c2=y2*0+1;
 y=c(t(y1[c(1:gn),]),t(y2[c(1:gn),]))
 g=rep(1:gn,each=n,times=2)
 x=rep(x,times=gn*2)
 c=c(c1,c2)
 df=data.frame(list(x=x,y=y,c=factor(c),g=factor(g)))
 (m=lmer(y~x*c + (x*c|g),data=df))
 if (doplot==TRUE)
  {require(lattice)
   df$fit=fitted(m)
   plot1=xyplot(fit ~ x|g,data=df,group=c,pch=19,cex=.1)
   plot2=xyplot(y ~ x|g,data=df,group=c)
   print(plot1+plot2)
  }
 print(colMeans(ranef(m)$g))
 m
}

Dans ce cas, les colMeans sortent toujours<106

Voici un exemple de jouet binomial (je partagerais mes données réelles, mais elles sont soumises pour publication et je ne suis pas sûr de ce que la politique de la revue est en train de publier au préalable):


toybin<-function(n=100,gn=4,doplot=FALSE){
  require(lme4)
x=runif(n,-16,16) y1=matrix(0,gn,n) y2=y1 for (gx in 1:gn) { com=runif(1,1,5) ucom=runif(1,1,5) y1[gx,]=tanh(x/(com+ucom) + rnorm(1)) > runif(x,-1,1) y2[gx,]=tanh(2*(x+2)/com + rnorm(1)) > runif(x,-1,1) } c1=y1*0; c2=y2*0+1; y=c(t(y1[c(1:gn),]),t(y2[c(1:gn),])) g=rep(1:gn,each=n,times=2) x=rep(x,times=gn*2) c=c(c1,c2) df=data.frame(list(x=x,y=y,c=factor(c),g=factor(g))) (m=lmer(y~x*c + (x*c|g),data=df,family=binomial)) if (doplot==TRUE) {require(lattice) df$fit=fitted(m) print(xyplot(fit ~ x|g,data=df,group=c,pch=19,cex=.1)) } print(colMeans(ranef(m)$g)) m }

Maintenant, les colMeans dépassent parfois 0,3, et nettement plus, en moyenne que l'exemple linéaire.

jerlich
la source
3
Pouvez-vous inclure du code pour reproduire ces exemples de jouets ici? Cela aiderait à explorer ce comportement intéressant.
Aaron a quitté Stack Overflow le
J'ai également vu ce même comportement avec mes expériences. Dans le cas gaussien, il existe une contrainte somme à zéro, mais pas dans les cas non gausiens. Je ne sais pas si la somme à zéro est une condition nécessaire, tant que la valeur attendue des effets aléatoires est zéro. Cela peut être utile dans certains cas, et apparemment il est facile de coder en cas gaussien donc c'est là ... J'espère que quelqu'un avec une meilleure compréhension des carillons.
Jouni

Réponses:

3

Depuis que le code de @ Hemmo a été légèrement altéré dans la case "Bounty", j'ajoute cette version reformatée en "wiki communautaire". Si ce n'est pas une utilisation appropriée du wiki, je m'excuse à l'avance. N'hésitez pas à le retirer.

library(mvabund)
library(lme4) 
data(spider) 
Y <- as.matrix(spider$abund)
X <- spider$x 
X <- X[ ,c(1, 4, 5, 6)] 
X <- rbind(X, X, X, X, X, X, X, X, X, X, X, X) 
site <- rep(seq(1, 28), 12) 
dataspider <- data.frame(c(Y), X, site) 
names(dataspider) <- c("Y","soil.dry", "moss", "herb.layer", "reflection", "site") 
fit <- glmer(
  Y ~ soil.dry + moss + herb.layer + reflection + (1|site), 
  family = poisson(link = log), 
  data = dataspider,
  control = glmerControl(optimizer = "bobyqa")
) 
David J. Harris
la source
1
Eh bien, il semble que cette question n'ait toujours pas reçu suffisamment d'attention. Ma propre conclusion est qu'il n'y a vraiment rien de mal ici, il n'y a vraiment pas de condition de somme à zéro, mais cela se produit uniquement dans les cas gaussiens où tout est linéaire. L'attente des effets aléatoires doit être 0, c'est l'hypothèse réelle, et non pas que les sommes réelles des effets estimés soient nulles. Je vais devoir récompenser quelqu'un, alors je vous en prie. :)
Jouni
2
@Hemmo Yikes, maintenant je sens que je devrais réellement contribuer quelque chose. Vous avez raison, rien n'est vraiment faux. La réponse courte (que j'espérais écrire mais que je n'ai pas trouvé le temps), est que la moyenne sera nulle si la surface de probabilité est gaussienne. De manière informelle, nous pouvons «prouver» cela en notant que les erreurs gaussiennes * Les effets aléatoires gaussiens conduisent à un autre gaussien. Lorsque vous avez un glmm avec une fonction d'erreur non gaussienne (par exemple Poisson dans votre cas), la surface de probabilité peut devenir non gaussienne et tous les paris sont désactivés. J'espère que cela t'aides.
David J. Harris