Why and how does the inclusion of random effects in mixed models influence the fixed-effect intercept term?

Question

The question is best illustrated by this example which uses a dataset (in library faraway) and lme4 library (both in R). This intercept-only model is only for illustrative purposes.

    library(lme4)
    library(faraway)
    data(epilepsy)
    log(mean(epilepsy$seizures))
          #expected intercept in intercept only model = 2.5544
    (Ep1a=glm(seizures~1, family=poisson, data=epilepsy))
          #intercept term =2.554 as expected.
    (Ep1b=glmer(seizures ~ 1 + (1|id), family=poisson, 
           data=epilepsy))#intercept term =2.214.

-------results----
1. simple glm 

    Call:  glm(formula = seizures ~ 1, family = poisson, 
             data = epilepsy)
    Coefficients:(Intercept)  2.554  
    
    Generalized linear mixed model fit by maximum likelihood  
     ['glmerMod']
    Family: poisson ( log )
    Formula: seizures ~ 1 + (1 | id) 
    Random effects:
     Groups Name        Std.Dev.
     id     (Intercept) 0.7795  
    Number of obs: 295, groups: id, 59
    Fixed Effects:
    (Intercept)  
          2.214  

My understanding is that the inclusion of the random term (id) tells the model that there is a repeated measures across subject (in this case). I can understand that this allows for the non-independence of the data: there are fewer than n=295 independent data points. But why does the fixed-effect intercept value decrease? Are there circumstances where it would increase?

I note from the following website: http://www.danielezrajohnson.com/glasgow_workshop.R the following in relation to a model unrelated to the one i've specified above (suggest you search for "average speaker"):

"...this model has random effects for speaker and word. The fixed effects reported are for a sort of average speaker and word. However, word, especially, tends to be a very skewed variable. There will always be a few very common words, that may favor or disfavor the response. The mixed model largely counteracts this weighting."

In my real example, all the coefficients (from a two categorical factor, crossed design with random-intercept-and-slope random effect - such that there are DistClaz dependencies within ID but Treat is independent) are considerably less than the corresponding mean values for those factor combinations as apparent in the raw data.

Results----------------------------------- 

    Generalized linear mixed model fit by maximum likelihood 
     ['glmerMod']
    Family: poisson ( log )
    Formula: RCP ~ Treated * DistClaz + (1 + DistClaz | ID) 
          AIC       BIC    logLik  deviance 
     69246.20  69307.43 -34611.10  69222.20 
    Random effects:
     Groups    Name        Variance Std.Dev. Corr       
     ID        (Intercept) 12.945   3.598               
               DistClazAZE 11.418   3.379    -0.63      
               DistClazRef  8.769   2.961    -0.80  0.71
    Number of obs: 1215, groups: ID, 344
    Fixed effects:
                               Estimate Std. Error z value Pr(>|z|)    
    (Intercept)                  2.2528     0.4240   5.313 1.08e-07 ***
    TreatedTreated              -2.1244     0.4916  -4.321 1.55e-05 ***
    DistClazAZE                  1.5925     0.5506   2.892  0.00382 ** 
    DistClazRef                  2.6302     0.3614   7.277 3.41e-13 ***
    TreatedTreated:DistClazAZE   0.3613     0.6110   0.591  0.55428    
    TreatedTreated:DistClazRef   0.7068     0.4201   1.682  0.09248 .  
    ---
    Correlation of Fixed Effects:
                (Intr) TrtdTr DsCAZE DstClR TT:DCA
    TreatedTrtd -0.862                            
    DistClazAZE -0.500  0.431                     
    DistClazRef -0.828  0.714  0.556              
    TrtdT:DCAZE  0.450 -0.525 -0.901 -0.501       
    TrtdTrt:DCR  0.712 -0.830 -0.478 -0.860  0.583
    ------------------------------------------END

These are the log(means) for the combinations:

    Results-------------------------------------
                 CE      AZE      Ref
    No EMB  5.454940 5.808012 6.273650
    Treated 3.626005 4.387088 5.253717
    -------------------------------------------END 

The fixed-effect intercept value (2.2528) relates to No EMB:CE (mean =5.4549)- but is this in relation to a 'sort-of-average' ID?

Any pointers would be much appreciated.

Answer 1

You have to consider that you are not predicting a single intercept, but a distribution (random intercepts) across a nonlinear link function. Remember http://en.wikipedia.org/wiki/Jensen%27s_inequality , which basically tells you that in general f(mean(x) != mean(f(x))

If I look at the values in your example and roughly estimate that we have an intercept of 2, random effect variance of 3, and a log link, I expect a log mean on the response of

    log(mean(exp(rnorm(100, 2, 3))))

which gives 5.036582, roughly fitting to your raw values.

So, to answer your general question: if you have a nonlinear link, a random effect can change the intercept because the random term is applied on a nonlinear link function.

Model misspecification could be another reason for changes in the intercept (e.g. that the random effects are not normal), but I suspect the former reason explains most of your puzzle.