[R] Hierarchical multi-level model with lmer: why are the highest-level random adjustments 0?

Stefan Th. Gries stgries at gmail.com
Sun Jul 7 22:14:08 CEST 2013


Hi all

I have a hopefully not too stupid question about multi-level /
mixed-effects modeling. I was trying to test a strategy from Crawley's
2013 R Book on a data set with the following structure:

- dependent variable: CONSTRUCTION (a factor with 2 levels)
- independent fixed effect: LENGTH (an integer in the interval [1, 61])
- random effects with the following hierarchical structure: MODE >
REGISTER > SUBREGISTER > FILE. Specifically:

MODE: S
  REGISTER: monolog
    SUBREGISTER: scripted
    SUBREGISTER: unscripted
  REGISTER: dialog
    SUBREGISTER: private
    SUBREGISTER: public
  REGISTER: mix
    SUBREGISTER: broadcast
MODE: W
  REGISTER: printed
    SUBREGISTER: academic
    SUBREGISTER: creative
    SUBREGISTER: instructional
    SUBREGISTER: nonacademic
    SUBREGISTER: persuasive
    SUBREGISTER: reportage
  REGISTER: nonprinted
    SUBREGISTER: letters
    SUBREGISTER: nonprofessional

with various levels of FILE in each level of SUBREGISTER. Here's the
head of the relevant data frame (best viewed with a non-proportional
font):

  CASE MODE   REGISTER     SUBREGISTER    FILE CONSTRUCTION LENGTH
1    1    W    printed       reportage W2C-002 V_P_DO       11
2    2    W    printed     nonacademic W2B-035 V_P_DO        8
3    3    W nonprinted nonprofessional W1A-014 V_P_DO        8
4    4    W    printed       reportage W2C-005 V_P_DO        6
5    5    S     dialog         private S1A-073 V_DO_P        2
6    6    S     dialog         private S1A-073 V_DO_P        2

And here's the unique-types distribution of FILE in the design:
tapply(FILE, list(SUBREGISTER, REGISTER, MODE), function (qwe)
length(unique(qwe)))

, , S
                dialog mix monolog nonprinted printed
academic             .   .       .          .       .
broadcast            .  20       .          .       .
creative             .   .       .          .       .
instructional        .   .       .          .       .
letters              .   .       .          .       .
nonacademic          .   .       .          .       .
nonprofessional      .   .       .          .       .
persuasive           .   .       .          .       .
private             96   .       .          .       .
public              77   .       .          .       .
reportage            .   .       .          .       .
scripted             .   .      25          .       .
unscripted           .   .      66          .       .

, , W
                dialog mix monolog nonprinted printed
academic             .   .       .          .      26
broadcast            .   .       .          .       .
creative             .   .       .          .      20
instructional        .   .       .          .      19
letters              .   .       .         28       .
nonacademic          .   .       .          .      37
nonprofessional      .   .       .         17       .
persuasive           .   .       .          .       9
private              .   .       .          .       .
public               .   .       .          .       .
reportage            .   .       .          .      20
scripted             .   .       .          .       .
unscripted           .   .       .          .       .

# I would usually have done this (using lme4)
model.1.tog <- lmer(CONSTRUCTION ~ LENGTH +
(1|MODE/REGISTER/SUBREGISTER), family=binomial)

# but Crawley (2013:692ff.) suggests this:
LEVEL2 <- MODE:REGISTER
LEVEL3 <- MODE:REGISTER:SUBREGISTER
model.1.sep <- lmer(CONSTRUCTION ~ LENGTH + (1|MODE) + (1|LEVEL2) +
(1|LEVEL3), family=binomial)

The results are the same for fixed and random effects, ok, but what I
don't understand in this case is why the random adjustments to
intercepts at the highest level of hierarchical organization (MODE)
are 0:

ranef(model.1.sep)
$LEVEL3
                             (Intercept)
S:dialog:private              -0.9482442
S:dialog:public                0.1216021
[...]

$LEVEL2
             (Intercept)
S:dialog      -0.4746389
[...]

$MODE
  (Intercept)
S           0
W           0

I am guessing that's got something to do with the fact that what the
random adjustments to intercepts at the level of MODE would do is
already taken care of by the random adjustments to intercepts on the
lower levels of the hierarchical organization of the data - but when I
run the code from Crawley on his data (a linear mixed-effects model,
not a generalized linear mixed-effects model as mine), the highest
hierarchical level of organization does not return 0 as random
adjustment. What am I missing?

Thanks for any input you may have!
STG



More information about the R-help mailing list