The emmeans package provides a variety of post hoc analyses such as obtaining estimated marginal means (EMMs) and comparisons thereof, displaying these results in a graph, and a number of related tasks.
This vignette illustrates basic uses of emmeans with
lm_robust objects. For more details, refer to the emmeans package itself and its vignettes.
warpbreaks dataset provided in base R has the results of a two-factor experiment. We start by fitting a model
Typical use of
emmeans() is to obtain predictions, or marginal means thereof, via a formula of the form
~ primary.variables | by.variables:
These results may be plotted as side-by-side intervals or as an interaction-style plot:
This particular example has a response transformation. That transformation is detected and we may back-transform to the original scale:
confint(emm, type = "response")
We may do comparisons and other contrasts:
Note that with a log transformations, it is possible to back-transform comparisons, and they become ratios. With other transformations, back-transforming is not possible.
Let’s create a variation on this example where one cell is omitted:
Note that the empty cell is detected and flagged as non-estimable.
Some additional explanation here. EMMs are based on a reference grid, defined as the grid created by all possible combinations of factor levels, together with the mean of each numerical predictor. The reference grid here (
rgi) is also an
"emmGrid" object just like the previous
emm. The grid itself is available as a data frame via the
grid member, and you can verify that the above results match those of the
predict function for the model:
predict(warpi.rlm, newdata = rgi@grid, se.fit = TRUE)
There is one exception for the empty cell. I will leave it as a user exercise to demonstrate that if we were to use different contrasts when fitting
warpi.rlm, the predictions will be the same except for the empty cell.
If there is a multivariate response, it is treated as another factor that is crossed with the other factors in the model. To illustrate, consider the dataset
MOats, provided in emmeans:
By default, the pseudo-factor is named
rep.meas, but we can change it if we like:
emmeans(MOats.rlm, pairwise ~ nitro, mult.name = "nitro")
This illustrates an additional feature of
emmeans that we can put a contrast method in the left side of a formula.
There are numerous capabilities of emmeans not illustrated here. See that package’s help files and vignettes. Using
vignette("basics", "emmeans") is a good starting point.