Skip to content

Residuals for a binomial lmer model

3 messages · Andy Fugard, Douglas Bates

Original
Andy Fugard
# 
Hello,

I was wondering how to get residuals out of binomial lmers as the 
"residuals" function isn't implemented.  I've noticed a few other people 
ask the question too but get no response.  (Or at least I haven't found 
a response.)

I guess the answer is just to use the "fitted" function, which is 
implemented for binomial GLMMs.

Take the sleepstudy data, dichotomize Reaction (just to give us a 
dataset), and fit a multilevel logistic model:


med = median(sleepstudy$Reaction)
sleepstudy$bin = (sleepstudy$Reaction > med) + 0
M2 = lmer(bin ~ Days + (1|Subject) + (0+Days|Subject),
           data = sleepstudy, family = binomial)


We can pull out the fitted values and, say, plot fitted (post-inverse
logit) against data using a boxplot:


ilog = function(x) { 1/(1 + exp(-x)) }
boxplot(ilog(fitted(fm3)) ~ bin, data = sleepstudy)


Not sure now how useful this is, but I had some reason for wanting to peek!

Andy
Andy Fugard
# 
Sorry, that works modulo variable names...  Second try.

med = median(sleepstudy$Reaction)
sleepstudy$bin = (sleepstudy$Reaction > med) + 0
mod = lmer(bin ~ Days + (1|Subject) + (0+Days|Subject),
            data = sleepstudy, family = binomial)

ilog = function(x) { 1/(1 + exp(-x)) }
boxplot(ilog(fitted(mod)) ~ bin, data = sleepstudy)
Andy Fugard wrote:

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Douglas Bates
    

    

      
      #
    

  


  

      
The problem with defining residuals for a generalized linear model or
a generalized linear mixed model is that you have so many possible
definitions.  Take a look at

?residuals.glm

Well, actually that isn't good advice because that help page refers
you to "the reference", of which there are two, and presumably you
need to sift through the descriptions in the reference to decide what
all these different types of residuals are.

In the development version of the lme4 package, the fitted model is an
object of class "mer" and has slots named "y", "eta", "mu", "resid",
"var", "pWt" and "sqrtrWt" which, for a generalized linear mixed
model, correspond to

y -> observed responses (remember that in the case of a binomial
response this is the fraction of positive responses, not the count)

eta -> linear predictor, depending on the fixed-effects parameters and
the random effects

mu -> conditional means of the responses, given the parameters and the
values of the random effects

resid -> sqrtrWt * (y - mu), i.e. the residuals whose squared length
is the weighted residual sum of squares

var -> the values of the variance function for the glm family.  This
is somewhat of a misnomer because it is the variance of the response
up to the prior weights.

pWt -> the prior weights.  For example, in a binomial response where
each "observation" is the number of successes in several trials, this
is the number of trials for each observation.

sqrtrWt -> the square root of the weights used to define the weighted
sum of squares of the residuals, sqrt(pWt/var)

Using these you can get several of the types of residuals defined in
residuals.glm.  The deviance residuals can also be calculated but I
never use them externally so I don't have an extractor for them.

So, what type of residuals do you want?

        
On Fri, Mar 14, 2008 at 7:41 AM, Andy Fugard <a.fugard at ed.ac.uk> wrote: