0

I am using LME model defined like: 

mod4.lme <- lme(pRNFL ~  Init.Age + Status + I(Time^2), random= ~1|Patient/EyeID,data = long1, na.action = na.omit)

The output is: 

> summary(mod4.lme)
Linear mixed-effects model fit by REML
 Data: long1 
       AIC      BIC    logLik
  2055.295 2089.432 -1018.647

Random effects:
 Formula: ~1 | Patient
        (Intercept)
StdDev:    7.949465

 Formula: ~1 | EyeID %in% Patient
        (Intercept) Residual
StdDev:    12.10405 2.279917

Fixed effects: pRNFL ~ Init.Age + Status + I(Time^2) 
                Value Std.Error  DF   t-value p-value
(Intercept)  97.27827  6.156093 212 15.801950  0.0000
Init.Age      0.02114  0.131122  57  0.161261  0.8725
StatusA     -27.32643  3.762155 212 -7.263504  0.0000
StatusF     -23.31652  3.984353 212 -5.852023  0.0000
StatusN      -0.28814  3.744980  57 -0.076940  0.9389
I(Time^2)    -0.06498  0.030223 212 -2.149921  0.0327
 Correlation: 
          (Intr) Int.Ag StatsA StatsF StatsN
Init.Age  -0.921                            
StatusA   -0.317  0.076                     
StatusF   -0.314  0.088  0.834              
StatusN   -0.049 -0.216  0.390  0.365       
I(Time^2) -0.006 -0.004  0.001 -0.038 -0.007

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-2.3565641 -0.4765840  0.0100608  0.4670792  2.7775392 

Number of Observations: 334
Number of Groups: 
           Patient EyeID %in% Patient 
                60                119

I wanted to get comparisons between my 'Status' factors (named A, N, F and H). So I did a emmeans model using this code: 

emmeans(mod4.lme, pairwise ~ Status, adjust="bonferroni")

The output for this, is: 

> emmeans(mod4.lme, pairwise ~ Status, adjust="bonferroni")
$emmeans
 Status   emmean       SE df lower.CL  upper.CL
 H      98.13515 2.402248 57 93.32473 102.94557
 A      70.80872 2.930072 57 64.94135  76.67609
 F      74.81863 3.215350 57 68.38000  81.25726
 N      97.84701 2.829706 57 92.18062 103.51340

Degrees-of-freedom method: containment 
Confidence level used: 0.95 

$contrasts
 contrast    estimate       SE  df t.ratio p.value
 H - A     27.3264289 3.762155 212   7.264  <.0001
 H - F     23.3165220 3.984353 212   5.852  <.0001
 H - N      0.2881375 3.744980  57   0.077  1.0000
 A - F     -4.0099069 2.242793 212  -1.788  0.4513
 A - N    -27.0382913 4.145370  57  -6.523  <.0001
 F - N    -23.0283844 4.359019  57  -5.283  <.0001
Lili
  • 547
  • 6
  • 19
  • 2
    Not sure what your question is. If it is ["Does `emmeans` support `lme` models?"](https://cran.r-project.org/web/packages/emmeans/vignettes/models.html) then the answer is "yes" (as you demonstrated yourself). – Maurits Evers Sep 28 '18 at 10:24
  • Sorry for the late reply! I am wondering if the multiple comparisons in this case is giving me the pure group differences (taking into account both the variables affecting the group and the nesting from the random effects) or is giving me a differences for interaction effects (Status*Time). Because what I want are the pure effects. Does this make sense?? Thanks! :) – Lili Nov 08 '18 at 15:40
  • It is giving you the differences between Status based on your model that takes into account the interactions. You can add time in the pairwise comparisons/contrasts by specifying this in your emmeans: `emmeans(mod4.lme, pairwise ~ Status | Time, adjust="bonferroni")` and then it should return the differences between `Status` for each `Time`. I am not sure if this answers your question? – Nakx Nov 08 '18 at 18:49
  • I think you reply to my question, yes. Because what I want to report are the group differences for the interaction effects and not only of the pure group (status) differences. So I think yes, you answered my question thanks. Just another clarification, if my model is: 'mod4.lme <- lme(pRNFL ~ Init.Age + Status* Time, random= ~1|Patient/EyeID,data = long1, na.action = na.omit)' rather than the one I explained before how could I define the interaction effect in the emmeans model? Thanks!! – Lili Nov 09 '18 at 11:38
  • Have you seen this link: https://cran.r-project.org/web/packages/emmeans/vignettes/interactions.html ? I think it may help you a lot with building and testing the effects from your model, as well as making figures to illustrate your results. – Nakx Nov 09 '18 at 21:24

1 Answers1

1

The answer is yes, emmeans does the calculation based on the model

Lili
  • 547
  • 6
  • 19