[R-sig-ME] Question on mixed effects model syntax and lack of p-values in lme4 summary

Dennis Murphy djmuser at gmail.com
Tue May 31 01:19:01 CEST 2011


Hi:

See inline.

On Mon, May 30, 2011 at 2:06 PM, Richard Friedman
<friedman at cancercenter.columbia.edu> wrote:
> Dear R-sig-mixed-model-list,
>
>        I am a beginner in mixed-effects models. I wish to do a mixed effects
> anova analysis with
> treatment as the fixed effect and mouse and field as the random effect, I am
> not sure of the syntax of
> the commands involved neither of which gives p-values. Here is a record of
> my session:
>
> #############################################################################
>
> R version 2.13.0 (2011-04-13)
> Copyright (C) 2011 The R Foundation for Statistical Computing
> ISBN 3-900051-07-0
> Platform: i386-apple-darwin9.8.0/i386 (32-bit)
>
> R is free software and comes with ABSOLUTELY NO WARRANTY.
> You are welcome to redistribute it under certain conditions.
> Type 'license()' or 'licence()' for distribution details.
>
>  Natural language support but running in an English locale
>
> R is a collaborative project with many contributors.
> Type 'contributors()' for more information and
> 'citation()' on how to cite R or R packages in publications.
>
> Type 'demo()' for some demos, 'help()' for on-line help, or
> 'help.start()' for an HTML browser interface to help.
> Type 'q()' to quit R.
>
> [R.app GUI 1.40 (5751) i386-apple-darwin9.8.0]
>
> [Workspace restored from /Users/friedman/.RData]
> [History restored from /Users/friedman/.Rhistory]
>
>> library(lme4)
> Loading required package: Matrix
> Loading required package: lattice
>
> Attaching package: 'Matrix'
>
> The following object(s) are masked from 'package:base':
>
>    det
>
>
> Attaching package: 'lme4'
>
> The following object(s) are masked from 'package:stats':
>
>    AIC, BIC
>
>> data<-read.table("mixinp.txt",header=T,sep="\t")
>> data
>   treatment mouse field         y
> 1          A     1     1  88.70659
> 2          A     1     2  83.85923
> 3          A     2     1  89.22233
> 4          A     2     2  73.78000
> 5          A     3     1  73.60307
> 6          A     3     2  76.53985
> 7          A     4     1  68.97187
> 8          A     4     2  70.39560
> 9          A     5     1  83.85923
> 10         A     5     2  88.70659
> 11         B     1     1 111.60902
> 12         B     1     2 113.25638
> 13         B     2     1 131.72738
> 14         B     2     2 124.26105
> 15         B     3     1 126.82843
> 16         B     3     2 134.61614
> 17         B     4     1 123.26297
> 18         B     4     2 112.42667
> 19         B     5     1 127.35469
> 20         B     5     2 140.14831
> 21         C     1     1  97.12810
> 22         C     1     2 116.70144
> 23         C     2     1 112.01635
> 24         C     2     2 108.07211
> 25         C     3     1 103.34168
> 26         C     3     2 109.22591
> 27         C     4     1 121.31415
> 28         C     4     2 109.22591
> 29         D     1     1 113.85043
> 30         D     1     2 120.03858
> 31         D     2     1 120.42591
> 32         D     2     2 123.36211
> 33         D     3     1 137.63444
> 34         D     3     2 127.88533
> 35         D     4     1 136.41102
> 36         D     4     2 133.44557
>> dim(data)
> [1] 36  4
>> attach(data)
>> mouse<-factor(mouse)
>> mouse
>  [1] 1 1 2 2 3 3 4 4 5 5 1 1 2 2 3 3 4 4 5 5 1 1 2 2 3 3 4 4 1 1 2 2 3 3 4 4
> Levels: 1 2 3 4 5
>> treatment<-factor(treatment)
>> treatment
>  [1] A A A A A A A A A A B B B B B B B B B B C C C C C C C C D D D D D D D D
> Levels: A B C D
>> model1<-lmer(y~treatment+(1|treatment/mouse))
>> summary(model1)
> Linear mixed model fit by REML
> Formula: y ~ treatment + (1 | treatment/mouse)
>   AIC   BIC logLik deviance REMLdev
>  246.3 257.4 -116.2    252.3   232.3
> Random effects:
>  Groups          Name        Variance Std.Dev.
>  mouse:treatment (Intercept) 38.136   6.1755
>  treatment       (Intercept) 11.737   3.4260
>  Residual                    39.625   6.2948
> Number of obs: 36, groups: mouse:treatment, 18; treatment, 4
>
> Fixed effects:
>            Estimate Std. Error t value
> (Intercept)   79.764      4.829  16.519
> treatmentB    44.785      6.829   6.558
> treatmentC    29.864      7.038   4.243
> treatmentD    46.867      7.038   6.659
>
> Correlation of Fixed Effects:
>           (Intr) trtmnB trtmnC
> treatmentB -0.707
> treatmentC -0.686  0.485
> treatmentD -0.686  0.485  0.471
>
>> model2<-lmer(y~treatment+(1|mouse:treatment))
>> summary(model2)
> Linear mixed model fit by REML
> Formula: y ~ treatment + (1 | mouse:treatment)
>   AIC   BIC logLik deviance REMLdev
>  244.3 253.8 -116.2    249.7   232.3
> Random effects:
>  Groups          Name        Variance Std.Dev.
>  mouse:treatment (Intercept) 38.136   6.1755
>  Residual                    39.625   6.2948
> Number of obs: 36, groups: mouse:treatment, 18
>
> Fixed effects:
>            Estimate Std. Error t value
> (Intercept)   79.764      3.404  23.431
> treatmentB    44.785      4.814   9.303
> treatmentC    29.864      5.106   5.848
> treatmentD    46.867      5.106   9.178
>
> Correlation of Fixed Effects:
>           (Intr) trtmnB trtmnC
> treatmentB -0.707
> treatmentC -0.667  0.471
> treatmentD -0.667  0.471  0.444
>>
>
>> sessionInfo()
> R version 2.13.0 (2011-04-13)
> Platform: i386-apple-darwin9.8.0/i386 (32-bit)
>
> locale:
> [1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
>
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base
>
> other attached packages:
> [1] lme4_0.999375-39   Matrix_0.999375-50 lattice_0.19-23
>
> loaded via a namespace (and not attached):
> [1] grid_2.13.0   nlme_3.1-100  stats4_2.13.0
>>
>
> ##################################################
>
> I have the following questions:
>
> 1. Is model1 or model 2 correct? or is neither correct/

It's impossible to tell without a thorough description of the
statistical design.

> 2. I do not get p-values. Is there a way to get p-values for this analysis
> with a mixed effects model?

See R-FAQ 7.35 for why p-values are not reported in lmer().

> 3. Does it seem as if I am doing anything incorrect/

See the response to (1). It may be a good idea to consult with someone
locally who knows how to use the lme4 package; that would be optimal.
If no such person exists at Columbia (which I would find rather
difficult to believe, but availability might be an issue), then please
repost with more details about the design, the relationships between
the factors, etc. and perhaps someone will be able to provide a
suitable response.

HTH,
Dennis
>
> Thanks and best wishes,
> Rich
> -----------------------------------------------------------
> Richard A. Friedman, PhD
> Associate Research Scientist,
> Biomedical Informatics Shared Resource
> Herbert Irving Comprehensive Cancer Center (HICCC)
> Lecturer,
> Department of Biomedical Informatics (DBMI)
> Educational Coordinator,
> Center for Computational Biology and Bioinformatics (C2B2)/
> National Center for Multiscale Analysis of Genomic Networks (MAGNet)
> Room 824
> Irving Cancer Research Center
> Columbia University
> 1130 St. Nicholas Ave
> New York, NY 10032
> (212)851-4765 (voice)
> friedman at cancercenter.columbia.edu
> http://cancercenter.columbia.edu/~friedman/
>
> I am a Bayesian. When I see a multiple-choice question on a test and I don't
> know the answer I say "eeney-meaney-miney-moe".
>
> Rose Friedman, Age 14
>
> _______________________________________________
> 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