[R-sig-ME] How do I interpret linear mixed model contrast estimates from multcomp::glht()?
ONKELINX, Thierry
Thierry.ONKELINX at inbo.be
Wed Dec 17 10:11:15 CET 2014
Matthew,
After scanning the helpfile of contrast, I get the feeling that you can vary only one variable and not two like you did. Make sure that you understand what contrast does. I prefer to use glht() directly. You have to specify more complex contrasts yourself, but then you know exactly what contrasts you are comparing.
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium
+ 32 2 525 02 51
+ 32 54 43 61 85
Thierry.Onkelinx op inbo.be<mailto:Thierry.Onkelinx op inbo.be>
www.inbo.be<http://www.inbo.be/>
To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of.
~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data.
~ Roger Brinner
The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data.
~ John Tukey
Van: Matthew Van Scoyoc [mailto:scoyoc op gmail.com]
Verzonden: woensdag 17 december 2014 6:25
Aan: Ken Beath
CC: ONKELINX, Thierry; r-sig-mixed-models op r-project.org
Onderwerp: Re: [R-sig-ME] How do I interpret linear mixed model contrast estimates from multcomp::glht()?
Thanks Ken and Thierry,
The rows in summary.glht() correspond to the row names of the contrast matrix given in the linfct argument as Thierry stated in my post on r-sig-ecology<http://r-sig-ecology.471788.n2.nabble.com/How-do-I-interpret-linear-mixed-model-contrast-estimates-from-multcomp-glht-td7579236.html#a7579237> (thanks again Thierry). I'm using the contrast function in the contrast package to help compute the contrast matrix. The contrast function from the contrast package produces this...
> # Calculate the contrasts
> cc<-contrast(math.lm, a = list(year = c(.5, 1.5, 2.5), size = 380, female = levels(egsingle$female)),
+ b = list(year = c(.5, 1.5, 2.5), size = 800, female = levels(egsingle$female)))
> cc
lm model parameter contrast
Contrast S.E. Lower Upper t df Pr(>|t|)
0.1671791 0.01739038 0.1330889 0.2012694 9.61 7225 0
0.2148280 0.02235506 0.1710055 0.2586504 9.61 7225 0
0.2624768 0.03167964 0.2003754 0.3245781 8.29 7225 0
0.1671791 0.01739038 0.1330889 0.2012694 9.61 7225 0
0.2148280 0.02235506 0.1710055 0.2586504 9.61 7225 0
0.2624768 0.03167964 0.2003754 0.3245781 8.29 7225 0
So, how do the rows of the contrast object relate to the coefficients of the mixed or linear model? The summary of math.lm produces 5 rows of coefficients...
> summary(math.lm)
Call:
lm(formula = math ~ year * size + female, data = egsingle)
Residuals:
Min 1Q Median 3Q Max
-3.5631 -0.7775 -0.0523 0.7183 5.0195
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -5.590e-01 3.687e-02 -15.163 < 2e-16 ***
year 8.521e-01 2.391e-02 35.644 < 2e-16 ***
size -3.413e-04 4.251e-05 -8.030 1.13e-15 ***
femaleFemale -2.232e-02 2.576e-02 -0.867 0.38624
year:size -1.135e-04 2.948e-05 -3.849 0.00012 ***
I keep looking at cc$X (what I'm passing to glht() as the contrast matrix in the linfct argument), but my brain starts to melt at this point...
> cc$X
(Intercept) year size femaleFemale year:size
1 0 0 -420 0 -210
2 0 0 -420 0 -630
3 0 0 -420 0 -1050
4 0 0 -420 0 -210
5 0 0 -420 0 -630
6 0 0 -420 0 -1050
attr(,"assign")
[1] 0 1 2 3 4
attr(,"contrasts")
attr(,"contrasts")$female
[1] "contr.SAS"
This is my first try at mixed models, and it's been several years since my linear regression course that was taught in SAS, so any help is appreciated.
Cheers,
MVS
=====
Matthew Van Scoyoc
https://sites.google.com/site/scoyoc/
=====
Think SNOW!
On Tue, Dec 16, 2014 at 2:24 PM, Ken Beath <ken.beath op mq.edu.au<mailto:ken.beath op mq.edu.au>> wrote:
You first need to work out what contrasts you want to calculate (and possibly find out what a contrast is). Then find out what the appropriate contrast matrix is, or the routines that will give you one a predefine one. Then if you use glht it will probably be fairly easy to interpret what results it is producing.
On 17 December 2014 at 03:28, Matthew Van Scoyoc <scoyoc op gmail.com<mailto:scoyoc op gmail.com>> wrote:
Sorry, I meant no disrespect. I just thought that by posting the question
on a more appropriate list, someone who has come across the same problem
would answer. The problem is actually with contrast::contrast(). The
function does not produce row names that correspond to the fixed effects or
the interactions in the model. So how do I figure out which row corresponds
to what, and how do I pull that out of the contrast() object *cc?*
MVS
=====
Matthew Van Scoyoc
<mvanscoyoc op aggiemail.usu.edu<mailto:mvanscoyoc op aggiemail.usu.edu>>https://sites.google.com/site/scoyoc/
=====
Think SNOW!
On Tue, Dec 16, 2014 at 1:28 AM, ONKELINX, Thierry <Thierry.ONKELINX op inbo.be<mailto:Thierry.ONKELINX op inbo.be>
> wrote:
>
> Matthew,
>
> Just reposting exactly the same question on a other list is not very
> polite. I answered it on r-sig-ecology
> http://r-sig-ecology.471788.n2.nabble.com/How-do-I-interpret-linear-mixed-model-contrast-estimates-from-multcomp-glht-td7579236.html#a7579237
> If the answer is not clear enough, then let us know what you don't
> understand about it. Have you look at the examples in ?glht
>
> ir. Thierry Onkelinx
> Instituut voor natuur- en bosonderzoek / Research Institute for Nature and
> Forest
> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
> Kliniekstraat 25
> 1070 Anderlecht
> Belgium
> + 32 2 525 02 51<tel:%2B%2032%202%20525%2002%2051>
> + 32 54 43 61 85<tel:%2B%2032%2054%2043%2061%2085>
> Thierry.Onkelinx op inbo.be<mailto:Thierry.Onkelinx op inbo.be>
> www.inbo.be<http://www.inbo.be>
>
> To call in the statistician after the experiment is done may be no more
> than asking him to perform a post-mortem examination: he may be able to say
> what the experiment died of.
> ~ Sir Ronald Aylmer Fisher
>
> The plural of anecdote is not data.
> ~ Roger Brinner
>
> The combination of some data and an aching desire for an answer does not
> ensure that a reasonable answer can be extracted from a given body of data.
> ~ John Tukey
>
> -----Oorspronkelijk bericht-----
> Van: R-sig-mixed-models [mailto:r-sig-mixed-models-bounces op r-project.org<mailto:r-sig-mixed-models-bounces op r-project.org>]
> Namens Matthew Van Scoyoc
> Verzonden: vrijdag 12 december 2014 18:23
> Aan: r-sig-mixed-models op r-project.org<mailto:r-sig-mixed-models op r-project.org>
> Onderwerp: [R-sig-ME] How do I interpret linear mixed model contrast
> estimates from multcomp::glht()?
>
> How do the rows in the summary (e.g. "1 == 0") correspond to the model?
> The answer is buried *contrast::contrast()*, but I can't figure it out.
> Consider this modified example I stole from here <
> https://stat.ethz.ch/pipermail/r-sig-mixed-models/2009q4/003061.html>...
>
> > options(contrasts = c(factor = "contr.SAS", ordered = "contr.poly"))
> > library("mlmRev")
> > library("lme4")
> > library("lmerTest")
> > library("contrast")
> > library("multcomp")
> >
> > data("egsingle")
> > # Linear mixed model
> > math.lmm <- lmer(math ~ year * size + female + (1|childid) +
> (1|schoolid), egsingle)
> > # Linear model
> > math.lm <- lm(math ~ year * size + female, data = egsingle) #
> > Calculate contrast matrix cc<-contrast(math.lm, a = list(year = c(0.5,
> > 1.5, 2.5), size = 380,
> female = levels(egsingle$female)), +
> b = list(year = c(0.5,
> 1.5, 2.5), size = 800, female = levels(egsingle$female)))
> > # Calculate estimates
> > summary(glht(math.lmm, linfct = cc$X))
>
> Simultaneous Tests for General Linear Hypotheses
>
> Fit: lme4::lmer(formula = math ~ year * size + female + (1 | childid) +
> (1 | schoolid), data = egsingle)
>
> Linear Hypotheses:
> Estimate Std. Error z value Pr(>|z|)
> 1 == 0 0.12774 0.08020 1.593 0.1272
> 2 == 0 0.15322 0.08066 1.900 0.0669 .
> 3 == 0 0.17870 0.08178 2.185 0.0341 *
> 4 == 0 0.12774 0.08020 1.593 0.1273
> 5 == 0 0.15322 0.08066 1.900 0.0669 .
> 6 == 0 0.17870 0.08178 2.185 0.0342 *
> ---
> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 (Adjusted p
> values reported -- single-step method)
>
> The row names correspond to the levels of *year* and *female,* and are
> probably Female:0.5, Female:1.5, Female:2.5, and so on. But how do I pull
> that out of the contrast() object *cc?* It might be simple with 3 main
> effects, but my current project has 5 main effects, four 2-way
> interactions, and one 3-way interaction, and the summary table has 24 rows.
> Ultimately I would like to create a dataframe so I can plot the contrasts,
> something like this...
>
> > x = summary(glht(math.lmm, linfct = cc$X)) # Contrast data frame
> > math.contr = data.frame(Effect.Interaction =
> > reference.something.in.cc<http://reference.something.in.cc>,
> Estimate =
> x[["test"]]$coefficients, Std.Error = x[["test"]]$sigma)
>
> Thanks for the help!
> Cheers,
> MVS
> =====
> Matthew Van Scoyoc
>
> <
> https://mail.google.com/mail/?view=cm&fs=1&tf=1&to=mvanscoyoc@aggiemail.usu.edu
> >
> https://sites.google.com/site/scoyoc/
> =====
> Think SNOW!
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models op r-project.org<mailto:R-sig-mixed-models op r-project.org> mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> Disclaimer Bezoek onze website / Visit our website<
> https://drupal.inbo.be/nl/disclaimer-mailberichten-van-het-inbo>
>
[[alternative HTML version deleted]]
_______________________________________________
R-sig-mixed-models op r-project.org<mailto:R-sig-mixed-models op r-project.org> mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
--
Ken Beath
Lecturer
Statistics Department
MACQUARIE UNIVERSITY NSW 2109, Australia
Phone: +61 (0)2 9850 8516
Building E4A, room 526
http://stat.mq.edu.au/our_staff/staff_-_alphabetical/staff/beath,_ken/
CRICOS Provider No 00002J
This message is intended for the addressee named and may contain confidential information. If you are not the intended recipient, please delete it and notify the sender. Views expressed in this message are those of the individual sender, and are not necessarily the views of the Faculty of Science, Department of Statistics or Macquarie University.
Disclaimer Bezoek onze website / Visit our website<https://drupal.inbo.be/nl/disclaimer-mailberichten-van-het-inbo>
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list