[R-sig-ME] how to compute Confident intervals for BLUPS for lme function in nlme library
Reinhold Kliegl
reinhold.kliegl at gmail.com
Wed Jan 14 07:50:22 CET 2009
Well, intervals() is a pretty good nlme response, I think, even if it
was not exactly answering the question. The lme4 equivalent is
doplot() which produces a so-called "caterpillar plot" of conditional
means.
For example:
?ranef
fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
rr1 <- ranef(fm1, postVar = TRUE)
dotplot(rr1,scales = list(x = list(relation = 'free')))[["Subject"]]
It is important to be clear about the distinction between random
effect estimates of variances and correlations, provided as part of
the model output and the conditional means for groups/units based on
("predicted with") them. At least, correlations computed from
conditional means can differ strongly from the correlation estimates.
The short story: You should only trust the model estimates. As Greg
Lee said, conditional means are to be handled with care; you can not
apply any of the usual inferential statistics to them. Nevertheless,
the caterpillar plots are highly informative and diagnostic about, for
example, whether you need a random effect to account for unit-by-unit
variability. They may also lead you to discover distinct subgroups
suggestive of a fixed effect (for a future model). There is a
manuscript at the top of my publications page for download (and
constructive feedback) walking through an example.
Reinhold Kliegl
On Wed, Jan 14, 2009 at 2:03 AM, Greg Lee <sp8ial at gmail.com> wrote:
> Hello again Maria,
>
> I was on autopilot when I answered the first time. The intervals()
> function in nlme provides confidence intervals for the estimated model
> components (fixed effects, random effect variances and correlation
> coefficients, if present), which is not what you were asking.
>
> Actually, the more I think about your question the less certain I am. The
> coefficients in a mixed effects model are (potentially) comprised of both
> fixed- and random-effects. But the random effects are not "estimated" so
> much as "predicted". The only (assumed) requirement on them is that they be
> normally distributed, with standard deviation estimates as reported by
> intervals().
>
> Does it even make sense to seek confidence intervals in this context? I am
> not sure, and leave it for wiser heads to respond.
>
> Regards,
> Greg
>
>
>
> 2009/1/14 Greg Lee <sp8ial at gmail.com>
>
>> Hi Maria,
>>
>> Try:
>>
>> intervals(ssb_model2.1)
>>
>> Cheers,
>> Greg
>>
>>
>>
>> 2009/1/14 Maria Jose Juan Jorda <mjuanjorda at gmail.com>
>>
>>> Dear all,
>>>
>>>
>>> How do you compute the confident intervals for the estimations of the
>>> random
>>> coefficients? I have not seem a lot of discussion about this in the forum.
>>> I
>>> wonder why. Is there any theory background about this?
>>>
>>> Description of my data:
>>>
>>> I have 29 time series of data, biomass over time for 29 populations. I
>>> want
>>> a hierarchical linear model with varying intercepts and slopes. So it
>>> computes 29 intercepts, 29 slopes and the overall mean of the intercept
>>> and
>>> slope and the variation among intercepts and slopes.
>>>
>>> *Here you have my linear mixed model:
>>> *
>>> ssb_model2.1<-lme(log(SSB_tonnes)~year2,random=~year2|id, data=ssb)
>>>
>>> *Data:*
>>> SSB_tonnes: spawning stock biomass (ssb), continuous variable
>>> year2: time in years, continuous variable
>>> id: there are 29 id, in other words 29 populations of fish
>>>
>>> *Here, in the summary of the ouput* I find the oveall mean of the
>>> interceps
>>> and slopes and their variation.
>>> *> summary(ssb_model2.1)*
>>> Linear mixed-effects model fit by maximum likelihood
>>> Data: ssb
>>> AIC BIC logLik
>>> 1575.647 1605.880 -781.8235
>>>
>>> Random effects:
>>> Formula: ~year2 | id
>>> Structure: General positive-definite, Log-Cholesky parametrization
>>> StdDev Corr
>>> (Intercept) 3.30025584 (Intr)
>>> year2 0.02960239 -0.674
>>> Residual 0.42259598
>>>
>>> Fixed effects: log(SSB_tonnes) ~ year2
>>> Value Std.Error DF t-value p-value
>>> (Intercept) 13.190739 0.6271420 1110 21.033098 0.0000
>>> year2 -0.019053 0.0058418 1110 -3.261529 0.0011
>>> Correlation:
>>> (Intr)
>>> year2 -0.69
>>>
>>> Standardized Within-Group Residuals:
>>> Min Q1 Med Q3 Max
>>> -5.34384192 -0.39116381 0.04824214 0.41871509 4.60300635
>>>
>>> Number of Observations: 1140
>>> Number of Groups: 29
>>>
>>> *Here I have the estimations for the random coefficient. 29 interceps and
>>> 29
>>> slopes. How do you compute the confidence intervals for these estimations?
>>> please any advise.
>>> *
>>> *> coef(ssb_model2.1)*
>>> (Intercept) year2
>>> 1 11.490207 -0.017008685
>>> 2 11.169138 0.008690365
>>> 3 13.373769 -0.030298582
>>> 4 13.680485 -0.024891089
>>> 5 15.498426 -0.033080606
>>> 6 14.869951 -0.012568273
>>> 7 13.016336 -0.025133918
>>> 8 13.968938 -0.030511346
>>> 9 13.079165 -0.017713727
>>> 10 12.935920 -0.054212415
>>> 11 9.124597 -0.026398183
>>> 12 11.493282 -0.005128911
>>> 13 10.371121 0.012636602
>>> 14 18.273836 -0.081414152
>>> 15 15.679493 -0.015946371
>>> 16 5.397496 0.027787751
>>> 17 14.479558 -0.024691465
>>> 18 16.388451 -0.030843356
>>> 19 18.326202 -0.098614391
>>> 20 12.031710 -0.010605874
>>> 21 10.352634 -0.008145436
>>> 22 14.501495 -0.038909443
>>> 23 14.415786 0.009461663
>>> 24 6.915476 0.025257749
>>> 25 7.247051 0.016988203
>>> 26 14.519202 -0.026469286
>>> 27 16.403880 -0.035931279
>>> 28 18.120565 0.014417736
>>> 29 15.407256 -0.019264309
>>>
>>> Thanks
>>>
>>> --
>>> >))):) >))):) >))):) >))):) >))):) >))):) >))):) >))):)
>>>
>>> Maria Jose Juan Jorda
>>>
>>> AZTI - Tecnalia / Unidad de Investigación Marina
>>> Herrera Kaia Portualde z/g
>>> 20110 Pasaia, Gipuzkoa, Spain
>>>
>>> Recursos Marinos y Pesquerias
>>> Depart. Biologia Animal, Vegetal y Ecologia
>>> Universidade A Coruña
>>> Campus A Zapateira s/n
>>> 15071, A Coruña, Spain
>>>
>>> Tel. Oficina +34-981167000 ext. 2204
>>> Tel. Mobil + 34-671072900
>>> Fax. +34981167065
>>> mjuan at pas.azti.es
>>> mjuanjorda at gmail.com
>>>
>>> [[alternative HTML version deleted]]
>>>
>>>
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>>
>>
>>
>>
>
>
> --
> --------------
> Greg Lee
> Biometrician
> Tasmanian Institute of Agricultural Research
> New Town Research Laboratories
> University of Tasmania
> 13 St Johns Avenue,
> New Town, 7008
> Australia
> Ph: +613 6233 6858
> Fax: +613 6233 6145
>
> [[alternative HTML version deleted]]
>
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>
More information about the R-sig-mixed-models
mailing list