[R-sig-ME] between-blocks variation in intercepts estimated as zero by lmer in RB design???
Douglas Bates
bates at stat.wisc.edu
Tue Feb 8 17:06:20 CET 2011
On Mon, Feb 7, 2011 at 5:44 PM, Duncan Mackay
<duncan.mackay at flinders.edu.au> wrote:
> Hello,
> I am having difficulty understanding why lmer is estimating the variation among intercepts of random groups (blocks) as zero in the following case:-
>
> I have samples of the no. of exotic species taken from two positions in each of 38 sites. (see data frame "rb.data" below). The data are coded as follows:-
>> str(rb.data)
> 'data.frame': 76 obs. of 3 variables:
> $ site : Factor w/ 38 levels "AMD1","BR1","DA1-1",..: 1 2 3 4 5 6 7 8 9 10 ...
> $ exotic : num 2 8 9 12 9 5 9 11 6 6 ...
> $ position: Factor w/ 2 levels "above","below": 2 2 2 2 2 2 2 2 2 2 ...
>>
>
> I have run a standard random-blocks analysis on these data as follows:-
>
>> summary(lmer(exotic ~ position + (1|site), data=rb.data))
> Linear mixed model fit by REML
> Formula: exotic ~ position + (1 | site)
> Data: rb.data
> AIC BIC logLik deviance REMLdev
> 519.8 529.1 -255.9 516.1 511.8
> Random effects:
> Groups Name Variance Std.Dev.
> site (Intercept) 0.000 0.0000
> Residual 53.526 7.3162
> Number of obs: 76, groups: site, 38
>
> Fixed effects:
> Estimate Std. Error t value
> (Intercept) 12.526 1.187 10.554
> positionbelow -4.053 1.678 -2.415
>
> Correlation of Fixed Effects:
> (Intr)
> positionblw -0.707
>
>
> So why, oh why, does lmer estimate the variance among the site intercepts as zero??? A plot of the data obtained by:-
>> with(rb.data, interaction.plot(position,site,exotic))
> suggests to me that the intercepts are actually quite variable, however a look at the fitted coefficients shows that lmer has fit the same intercept for each site.
I would plot the data using lattice as
dotplot(reorder(site, exotic) ~ sqrt(exotic), rb.data,
group=position, type=c("p","a"), auto.key=list(columns=2,lines=TRUE))
assuming that exotic is actually a count and thus a square root
transformation would help to stabilize the variance. From the
enclosed plot you can see that in some ways the average is being
driven by the above counts. The below counts are more-or-less
constant whereas the above counts go from very low to very high. So
the above counts are what are driving the variation.
The model that you are fitting only considers the mean count for each
location as a random effect so it is modeling the data as in the
second plot where there is not a great deal of variability in the
averages compared to the residual variability, which leads to an
estimate of zero for the between-site variance.
If you fit this as a generalized linear mixed model with the Poisson
family you do get a non-zero variance estimate for the site term but
even this model is not getting at the pattern in the data.
> (gm1 <- glmer(exotic ~ position + (1|site), rb.data, poisson))
Generalized linear mixed model fit by maximum likelihood ['merMod']
Family: poisson
Formula: exotic ~ position + (1 | site)
Data: rb.data
AIC BIC logLik deviance
302.9469 309.9391 -148.4734 296.9469
Random effects:
Groups Name Variance Std.Dev.
site (Intercept) 0.09967 0.3157
Number of obs: 76, groups: site, 38
Fixed effects:
Estimate Std. Error z value
(Intercept) 2.47856 0.06922 35.81
positionbelow -0.39087 0.07234 -5.40
Correlation of Fixed Effects:
(Intr)
positionblw -0.422
The problem with even this fit is that the pattern in the data is not
that of an additive model.
>> coef(rb.exotic.lme)
> $site
> (Intercept) positionbelow
> AMD1 12.52632 -4.052632
> BR1 12.52632 -4.052632
> DA1-1 12.52632 -4.052632
> DA1-2 12.52632 -4.052632
> DA1-3 12.52632 -4.052632
> .......
>
>
> ..................................... I think I'm missing something pretty basic here!!
>
> Many thanks for any help,
> Regards,
> Duncan
>
>
>
>> rb.data
> site exotic position
> 1 AMD1 2 below
> 2 BR1 8 below
> 3 DA1-1 9 below
> 4 DA1-2 12 below
> 5 DA1-3 9 below
> 6 DA1-4 5 below
> 7 DA1-5 9 below
> 8 MAY1 11 below
> 9 PH1 6 below
> 10 PH2 6 below
> 11 RA1 9 below
> 12 RA2 11 below
> 13 RB1 14 below
> 14 RB2 20 below
> 15 RB3 13 below
> 16 RB4-1 9 below
> 17 RB6-1 9 below
> 18 RB6-2 5 below
> 19 RB7-2 5 below
> 20 RB8-1 10 below
> 21 RLCL1 9 below
> 22 RLCL2 12 below
> 23 ROW1 12 below
> 24 ROW2 9 below
> 25 RS1 5 below
> 26 RS2/RS3 5 below
> 27 SC1 5 below
> 28 SC2-2 5 below
> 29 SC5-1 10 below
> 30 SC7 8 below
> 31 SC8 5 below
> 32 TH5-1 12 below
> 33 TH5-2 11 below
> 34 WI5 1 below
> 35 WI6-1 6 below
> 36 WI6-2 13 below
> 37 WR2-1 6 below
> 38 WR2-2 6 below
> 39 AMD1 27 above
> 40 BR1 16 above
> 41 DA1-1 4 above
> 42 DA1-2 4 above
> 43 DA1-3 5 above
> 44 DA1-4 5 above
> 45 DA1-5 6 above
> 46 MAY1 13 above
> 47 PH1 22 above
> 48 PH2 19 above
> 49 RA1 16 above
> 50 RA2 13 above
> 51 RB1 7 above
> 52 RB2 5 above
> 53 RB3 4 above
> 54 RB4-1 3 above
> 55 RB6-1 1 above
> 56 RB6-2 4 above
> 57 RB7-2 13 above
> 58 RB8-1 8 above
> 59 RLCL1 6 above
> 60 RLCL2 15 above
> 61 ROW1 6 above
> 62 ROW2 7 above
> 63 RS1 27 above
> 64 RS2/RS3 29 above
> 65 SC1 15 above
> 66 SC2-2 18 above
> 67 SC5-1 7 above
> 68 SC7 25 above
> 69 SC8 33 above
> 70 TH5-1 4 above
> 71 TH5-2 3 above
> 72 WI5 20 above
> 73 WI6-1 2 above
> 74 WI6-2 3 above
> 75 WR2-1 30 above
> 76 WR2-2 31 above
>
>
> ______________________________________________________________________________
> Dr. Duncan Mackay
> School of Biological Sciences
> Flinders University
> GPO Box 2100
> Adelaide
> S.A. 5001
> AUSTRALIA
>
> Phone 61-8-82012627
> FAX 61-8-82013015
>
> http://www.flinders.edu.au/science_engineering/biology/our-school/staff-postgrads/academic-staff/mackay-duncan.cfm
>
>
>
> ______________________________________________________________________________
> Dr. Duncan Mackay
> School of Biological Sciences
> Flinders University
> GPO Box 2100
> Adelaide
> S.A. 5001
> AUSTRALIA
>
> Phone 61-8-82012627
> FAX 61-8-82013015
>
> http://www.flinders.edu.au/science_engineering/biology/our-school/staff-postgrads/academic-staff/mackay-duncan.cfm
>
>
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: rb.pdf
Type: application/pdf
Size: 12871 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20110208/5035bf5f/attachment-0004.pdf>
More information about the R-sig-mixed-models
mailing list