[R] Extracting random parameters from summary lme and lmer

Chuck Cleland ccleland at optonline.net
Tue Jul 31 12:58:21 CEST 2007


Rense Nieuwenhuis wrote:
> LS,
> 
> I'm estimating multilevel regression models, using the lme-function  
> from the nlme-package. Let's say that I estimated a model and stored  
> it inside the object named 'model'. The summary of that model is  
> shown below:
> 
> Using summary(model)$tTable , I receive the following output:
> 
>  > summary(model)$tTable
>                      Value  Std.Error   DF     t-value       p-value
> (Intercept)    0.23268607 0.09350662 3990   2.4884449  1.287080e-02
> sexM          -0.15338225 0.03169762 3990  -4.8389206  1.354802e-06
> standLRT       0.38593558 0.01677195 3990  23.0107762 4.005182e-110
> vrmid 50%      0.07606394 0.09389376   61   0.8101064  4.210281e-01
> vrtop 25%      0.24561327 0.10483374   61   2.3428838  2.241317e-02
> intakemid 50% -0.41469716 0.03177240 3990 -13.0521199  3.698344e-38
> intaketop 25% -0.75920783 0.05357980 3990 -14.1696648  1.666780e-44
> typeSngl       0.15680532 0.07173835   61   2.1857949  3.267903e-02
> 
> 
> All looks fine to me. The output above is simply  a section from the  
> full summary shown below. Now, I want to extract from the summary (or  
> the full model) the part stating the random parameters. More  
> specifically, I want to extract from the summary the following:
> 
> (Intercept) 0.2869401 (Intr)
> typeSngl    0.2791040 -0.617
> Residual    0.7302233
> 
> How could this be done, and how can the same be done using the lmer- 
> function from the lme4-package?

?VarCorr

> library(nlme)

> fm1 <- lme(distance ~ age, data = Orthodont) # random is ~ age

> summary(fm1)
Linear mixed-effects model fit by REML
 Data: Orthodont
       AIC      BIC    logLik
  454.6367 470.6173 -221.3183

Random effects:
 Formula: ~age | Subject
 Structure: General positive-definite
            StdDev    Corr
(Intercept) 2.3270339 (Intr)
age         0.2264276 -0.609
Residual    1.3100399

Fixed effects: distance ~ age
                Value Std.Error DF   t-value p-value
(Intercept) 16.761111 0.7752461 80 21.620375       0
age          0.660185 0.0712533 80  9.265334       0
 Correlation:
    (Intr)
age -0.848

Standardized Within-Group Residuals:
         Min           Q1          Med           Q3          Max
-3.223106017 -0.493760867  0.007316632  0.472151090  3.916032743

Number of Observations: 108
Number of Groups: 27

> VarCorr(fm1)
Subject = pdSymm(age)
            Variance   StdDev    Corr
(Intercept) 5.41508688 2.3270339 (Intr)
age         0.05126947 0.2264276 -0.609
Residual    1.71620459 1.3100399

> library(lme4)

> (fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy))
Linear mixed-effects model fit by REML
Formula: Reaction ~ Days + (Days | Subject)
   Data: sleepstudy
  AIC  BIC logLik MLdeviance REMLdeviance
 1754 1770 -871.8       1752         1744
Random effects:
 Groups   Name        Variance Std.Dev. Corr
 Subject  (Intercept) 610.835  24.7151
          Days         35.056   5.9208  0.067
 Residual             655.066  25.5943
number of obs: 180, groups: Subject, 18

Fixed effects:
            Estimate Std. Error t value
(Intercept)  251.405      6.820   36.86
Days          10.467      1.546    6.77

Correlation of Fixed Effects:
     (Intr)
Days -0.137

> VarCorr(fm1)
$Subject
2 x 2 Matrix of class "dpoMatrix"
            (Intercept)      Days
(Intercept)  610.834546  9.738707
Days           9.738707 35.056337

attr(,"sc")
   scale
25.59426

> Thanks for the effort,
> 
> Rense Nieuwenhuis
> 
> Linear mixed-effects model fit by REML
>   Data: Exam
>        AIC      BIC   logLik
>    9158.56 9234.241 -4567.28
> 
> Random effects:
>   Formula: ~type | school
>   Structure: General positive-definite, Log-Cholesky parametrization
>              StdDev    Corr
> (Intercept) 0.2869401 (Intr)
> typeSngl    0.2791040 -0.617
> Residual    0.7302233
> 
> Fixed effects: normexam ~ sex + standLRT + vr + intake + type
>                     Value  Std.Error   DF    t-value p-value
> (Intercept)    0.2326861 0.09350662 3990   2.488445  0.0129
> sexM          -0.1533822 0.03169762 3990  -4.838921  0.0000
> standLRT       0.3859356 0.01677195 3990  23.010776  0.0000
> vrmid 50%      0.0760639 0.09389376   61   0.810106  0.4210
> vrtop 25%      0.2456133 0.10483374   61   2.342884  0.0224
> intakemid 50% -0.4146972 0.03177240 3990 -13.052120  0.0000
> intaketop 25% -0.7592078 0.05357980 3990 -14.169665  0.0000
> typeSngl       0.1568053 0.07173835   61   2.185795  0.0327
>   Correlation:
>                (Intr) sexM   stnLRT vrm50% vrt25% int50% int25%
> sexM          -0.201
> standLRT      -0.125  0.028
> vrmid 50%     -0.742  0.028 -0.035
> vrtop 25%     -0.652  0.051 -0.065  0.649
> intakemid 50% -0.246 -0.011  0.541 -0.002  0.007
> intaketop 25% -0.218 -0.018  0.676  0.014  0.013  0.660
> typeSngl      -0.421  0.080  0.007  0.033 -0.027 -0.001  0.001
> 
> Standardized Within-Group Residuals:
>          Min          Q1         Med          Q3         Max
> -3.59074329 -0.63776965  0.03829878  0.67303837  3.33952680
> 
> Number of Observations: 4059
> Number of Groups: 65
> 
> 
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 


-- 
Chuck Cleland, Ph.D.
NDRI, Inc.
71 West 23rd Street, 8th floor
New York, NY 10010
tel: (212) 845-4495 (Tu, Th)
tel: (732) 512-0171 (M, W, F)
fax: (917) 438-0894



More information about the R-help mailing list