[R-sig-ME] lmer predicted and fitted values differ

[ONKELINX, Thierry]

Dear Amy,

Here is a solution. It should work with nested random effects. Untested
with random slopes and crossed random effects.

Best regards,

Thierry

library(lme4)
cbpp$Obs <- factor(seq_len(nrow(cbpp)))
(gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 |
herd/Obs), family = binomial, data = cbpp))
str(gm1)

predict.lmerBin <- function(object, X, Zt){
	if(missing(X)){
		X <- object at X
	}
	b <- fixef(object)
	if(missing(Zt)){
		Zt <- as.matrix(object at Zt)
	}
	z <- unlist(ranef(object))
	plogis(X %*% b + t(Zt) %*% z)
}

all.equal(as.vector(predict.lmerBin(gm1)), fitted(gm1))

PS Replying to all should keep the thread intact.

------------------------------------------------------------------------
----
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek
team Biometrie & Kwaliteitszorg
Gaverstraat 4
9500 Geraardsbergen
Belgium

Research Institute for Nature and Forest
team Biometrics & Quality Assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium

tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to
say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of
data.
~ John Tukey
  

> -----Oorspronkelijk bericht-----
> Van: Amy Wade [mailto:a.s.wade at reading.ac.uk] 
> Verzonden: woensdag 24 november 2010 14:44
> Aan: ONKELINX, Thierry; r-sig-mixed-models at r-project.org
> Onderwerp: lmer predicted and fitted values differ
> 
> Theirry,
> 
> Thanks for that clarification. I would like to predict data 
> for the existing levels of the random effects.
> 
> 
> Many thanks,
> Amy
> 
> P.S. not sure how to continue a thread - hope this works....
> 
> 
> 
> 
> 
> .......................................
> Amy S. I. Wade
> Centre for Agri-Environmental Research
> University of Reading
> 
> 
> 
> 
> 
> -----Original Message-----
> From: ONKELINX, Thierry [mailto:Thierry.ONKELINX at inbo.be]
> Sent: 24 November 2010 12:47
> To: Amy Wade; r-sig-mixed-models at r-project.org
> Subject: RE: [R-sig-ME] lmer predicted and fitted values differ
> 
> Dear Amy,
> 
> The functions that you used take only the fixed effects into account.
> While fitted() takes both the fixed and the random effects 
> into account.
> 
> What do you want to predict? The 'average' data (all random 
> effect levels = 0)? Data for existing levels of the random 
> effects? Data for new levels of the random effects?
> 
> Best regards,
> 
> Thierry
> 
> --------------------------------------------------------------
> ----------
> ----
> ir. Thierry Onkelinx
> Instituut voor natuur- en bosonderzoek
> team Biometrie & Kwaliteitszorg
> Gaverstraat 4
> 9500 Geraardsbergen
> Belgium
> 
> Research Institute for Nature and Forest team Biometrics & 
> Quality Assurance Gaverstraat 4 9500 Geraardsbergen Belgium
> 
> tel. + 32 54/436 185
> Thierry.Onkelinx at inbo.be
> www.inbo.be
> 
> To call in the statistician after the experiment is done may 
> be no more than asking him to perform a post-mortem 
> examination: he may be able to say what the experiment died of.
> ~ Sir Ronald Aylmer Fisher
> 
> The plural of anecdote is not data.
> ~ Roger Brinner
> 
> The combination of some data and an aching desire for an 
> answer does not ensure that a reasonable answer can be 
> extracted from a given body of data.
> ~ John Tukey
>   
> 
> > -----Oorspronkelijk bericht-----
> > Van: r-sig-mixed-models-bounces at r-project.org
> > [mailto:r-sig-mixed-models-bounces at r-project.org] Namens Amy Wade
> > Verzonden: woensdag 24 november 2010 13:29
> > Aan: r-sig-mixed-models at r-project.org
> > Onderwerp: [R-sig-ME] lmer predicted and fitted values differ
> > 
> > Hello,
> > 
> >  
> > 
> > I am trying to predict response values for a model with new, 
> > unmeasured parameter values. It is an lmer model 
> constructed using the 
> > lme4 package. I realise there is not a 'predict.lmer()' function or 
> > similar for lmer models.
> > Having searched on the internet I found 3 possible methods for 
> > calculating predicted values given unmeasured values for a given 
> > parameter in the model.
> > 
> > 
> >  
> > 
> > This is my model:
> > 
> > data.m2<-data.frame(pods.bind,FUNGICIDE,N.TENTS,Day,HUMIDITY,P
> > LOT,TREE)
> > #where 'pods.bind' is cbind(Black.pods,healthy.pods)
> > 
> > m2<-lmer(pods.bind~FUNGICIDE+N.TENTS*Day+HUMIDITY*Day+(1|PLOT/
> > TREE),family="
> > quasibinomial")
> > 
> >  
> > 
> >  
> > 
> > Here are the possible methods:
> > 
> > #1
> > 
> > mm = model.matrix(terms(m2),data.m2)
> > 
> > data.m2$pods.bind = mm %*% fixef(m2)
> > 
> > predicted.1<-exp(data.m2$pods.bind)/(1+exp(data.m2$pods.bind))
> > 
> >  
> > 
> > #2
> > 
> > predict.lmerBin <- function(object, X){
> > 
> > if(missing(X))
> > 
> > X <- object at X
> > 
> > b <- fixef(object)
> > 
> > plogis(X %*% b)
> > 
> > }
> > 
> > predicted.2<-predict.lmerBin(m2)
> > 
> >  
> > 
> > #3
> > 
> > predicted.3<-exp(model.matrix(terms(m2),data.m2)%*%fixef(m2))
> > 
> >  
> > 
> >  
> > 
> > All three methods result in exactly the same predicted values. 
> > However, when I plug in the real data I get different 
> values than I do 
> > for fitted(). This makes me doubt the validity of my 
> predicted values. 
> > Could anybody explain why the predicted values, using these 
> methods, 
> > differ from the fitted values?
> > 
> >  
> > 
> > Thanks for your help.
> > 
> >  
> > 
> > Amy
> > 
> >  
> > 
> > ..............................................................
> > ..............
> > ...............
> > 
> > Amy S. I. Wade
> > 
> > Centre for Agri-Environmental Research
> > 
> > University of Reading
> > 
> >  
> > 
> > 
> > 	[[alternative HTML version deleted]]
> > 
> > _______________________________________________
> > R-sig-mixed-models at r-project.org mailing list 
> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> > 
> 
>