require("emmeans") knitr::opts_chunk$set(fig.width = 4.5, class.output = "ro")
Much of what you do with the emmeans package involves these three basic steps:
EMM <- emmeans(...)
(see scenarios below) to obtain estimates of
means or marginal meanscontrast(EMM, ...)
or pairs(EMM)
one or more times to obtain
estimates of contrasts or pairwise comparisons among the means.Note: A lot of users have developed the habit of running something like
emmeans(model, pairwise ~ factor(s))
, which conflates steps 1 and 2. We recommend
against doing this because it often yields output you don't want or need -- especially
when there is more than one factor. You are better off keeping steps 1 and 2 separate.
What you do in step 2 depends on how many factors you have, and how they relate.
If one-factor model fits well and the factor is named treatment
, do
EMM <- emmeans(model, "treatment") # or emmeans(model, ~ treatment) EMM # display the means ### pairwise comparisons contrast(EMM, "pairwise") # or pairs(EMM)
You may specify other contrasts in the second argument of the contrast()
call, e.g. "trt.vs.ctrl", ref = 1
(compare each mean to the first),
or "consec"
(compare 2 vs 1, 3 vs 2, etc.), or "poly", max.degree = 3
(polynomial contrasts)
If the model fits well and factors are named treat
and dose
, and they don't interact, follow
the same steps as for one factor at a time.
That is, something like
(EMM1 <- emmeans(model, ~ treat)) pairs(EMM1) (EMM2 <- emmeans(model, ~ dose)) pairs(EMM2)
These analyses will yield the estimated marginal means for each factor, and comparisons/contrasts thereof.
In this case, unless the interaction effect is negligible, we usually want to do "simple comparisons" of the cell means. That is, compare or contrast the means separately, holding one factor fixed at each level.
EMM <- emmeans(model, ~ treat * dose) EMM # display the cell means ### Simple pairwise comparisons... pairs(EMM, simple = "treat") # compare treats for each dose -- "simple effects" pairs(EMM, simple = "dose") # compare doses for each treat
The default is to apply a separate Tukey adjustment to the P values in each
by
group (so if each group has just 2 means, no adjustment at all is applied).
If you want to adjust the whole family combined, you need to undo the by
variable
and specify the desired adjustment (which can't be Tukey because that method
is invalid when you have more than one set of pairwise comparisons.) For example
test(pairs(EMM, by = "dose"), by = NULL, adjust = "mvt")
If the "diagonal" comparisons (where both factors differ) are of
interest, you would do pairs(EMM)
without a by
variable. But you get a lot more
comparisons this way.
Sometimes you may want to examine interaction contrasts, which are contrasts of contrasts. The thing to know here is that contrast()
or (pairs()
) creates the same kind of object as emmeans()
, so you can run them multiple times. For example,
CON <- pairs(EMM, by = "dose") contrast(CON, "consec", by = NULL) # by = NULL is essential here!
Or equivalently, the named argument interaction
can be used
contrast(EMM, interaction = c("pairwise", "consec"))
After you have mastered the strategies for two factors, you can adapt them to three or more factors as appropriate, based on how they interact and what you need.
emmeans()
and ref_grid()
for additional
arguments that may prove useful. Many of the most useful arguments are
passed to ref_grid()
.Most non-graphical functions in the emmeans package produce one of two classes of objects. The functions emmeans()
, emtrends()
, ref_grid()
, contrast()
, and pairs()
return emmGrid
objects (or lists thereof, class emm_list
). For example
EMM <- emmeans(mod, "Treatment")
The functions summary()
, confint()
, test()
, joint_tests()
, and others
return summary_emm
objects (or lists thereof, class summary_eml
):
SEMM <- summary(EMM)
If you display EMM
and SEMM
, they look identical; that's because emmGrid
objects are displayed using summary()
. But they are not identical. EMM
has all the ingredients needed to do further analysis, e.g. contrast(EMM, "consec")
will estimate comparisons between consecutive Treatment
means. But SEMM
is just an annotated data frame and we can do no further analysis with it. Similarly, we can change how EMM
is displayed via arguments to summary()
or relatives, whil;e in SEMM
, everything has been computed and those results are locked-in.
This is probably the most common issue, and it can happen when a treatment is coded as a numeric predictor rather than a factor. Instead of getting a mean for each treatment, you get a mean at the average of those numerical values.
treatment
with factor(treatment)
and re-fit the model.at = list(treatment = c(3,5,7))
to the emmeans()
call.emtrends()
pairwise ~ ...
recipe {#pairwise}
The basic object returned by emmeans()
and contrast()
is of class emmGrid
,
and additional emmeans()
and contrast()
calls can accept emmGrid
objects.
However, some options create lists of emmGrid
objects, and that makes things
a bit confusing. The most common case is using a call like
emmeans(model, pairwise ~ treat * dose)
, which computes the means and
all pairwise comparisons -- a list of two emmGrid
s. If you try to obtain additional
contrasts, say, of this result, contrast()
makes a guess that you want to
run it on just the first element.
This causes confusion (I know, because I get a lot of questions about it).
I recommend that you avoid using the pairwise ~
construct altogether:
Get your means in one step, and get your contrasts in separate step(s).
The pairwise ~
construct is generally useful if you have only one factor;
otherwise, it likely gives you results you don't want.
There are several of these vignettes that offser more details and more advanced topics. An index of all these vignette topics is available here.
The strings linked below are the names of the vignettes; i.e., they can
also be accessed via vignette("
name", "emmeans")
emmGrid
objects: "utilities"Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.