Movatterモバイル変換


[0]ホーム

URL:


Utilities and options for emmeans

emmeans package, Version 2.0.1

Contents

  1. Updating anemmGrid object
  2. Setting options
    1. Setting and viewing defaults
    2. Optimal digits to display
    3. Startup options
  3. Combining and subsettingemmGridobjects
  4. Accessing results to use elsewhere
  5. Adding grouping factors
  6. Re-labeling and re-leveling anemmGrid

Index of all vignette topics

Updating anemmGrid object

Several internal settings are saved when functions likeref_grid(),emmeans(),contrast(), etc. are run. Those settings can be manipulatedvia theupdate() method foremmGrids. Toillustrate, consider thepigs dataset and model yetagain:

pigs.lm <- lm(log(conc) ~ source + factor(percent), data = pigs)pigs.emm <- emmeans(pigs.lm, "source")pigs.emm
##  source emmean     SE df lower.CL upper.CL##  fish     3.39 0.0367 23     3.32     3.47##  soy      3.67 0.0374 23     3.59     3.74##  skim     3.80 0.0394 23     3.72     3.88## ## Results are averaged over the levels of: percent ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95

We see confidence intervals but not tests, by default. This happensas a result of internal settings inpigs.emm.s that arepassed tosummary() when the object is displayed. If we aregoing to work with this object a lot, we might want to change itsinternal settings rather than having to rely on explicitly callingsummary() with several arguments. If so, just update theinternal settings to what is desired; for example:

pigs.emm.s <- update(pigs.emm, infer = c(TRUE, TRUE), null = log(35),                     calc = c(n = ".wgt."))pigs.emm.s
##  source emmean     SE df  n lower.CL upper.CL null t.ratio p.value##  fish     3.39 0.0367 23 10     3.32     3.47 3.56  -4.385  0.0002##  soy      3.67 0.0374 23 10     3.59     3.74 3.56   2.988  0.0066##  skim     3.80 0.0394 23  9     3.72     3.88 3.56   6.130 <0.0001## ## Results are averaged over the levels of: percent ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95

Note that by adding ofcalc, we have set a default tocalculate and display the sample size when the object is summarized. Seehelp("update.emmGrid") for details on the keywords that canbe changed. Mostly, they are the same as the names of arguments in thefunctions that construct these objects.

Of course, we can always get what we want via calls totest(),confint() orsummary()with appropriate arguments. But theupdate() function ismore useful in sophisticated manipulations of objects, or calledimplicitly via the... oroptions argument inemmeans() and other functions. Those options are passed toupdate() just before the object is returned. For example,we could have done the above update within theemmeans()call as follows (results are not shown because they are the same asbefore):

emmeans(pigs.lm, "source", infer = c(TRUE, TRUE), null = log(35),        calc = c(n = ".wgt."))

Back to contents

Setting options

Speaking of theoptions argument, note that the defaultinemmeans() isoptions = get_emm_option("emmeans"). Let’s see what thatis:

get_emm_option("emmeans")
## $infer## [1]  TRUE FALSE

So, by default, confidence intervals, but not tests, are displayedwhen the result is summarized. The reverse is true for results ofcontrast() (and also the default forpairs()which callscontrast()):

get_emm_option("contrast")
## $infer## [1] FALSE  TRUE

There are also defaults for a newly constructed reference grid:

get_emm_option("ref_grid")
## $is.new.rg## [1] TRUE## ## $infer## [1] FALSE FALSE

The default is to display neither intervals nor tests whensummarizing. In addition, the flagis.new.rg is set toTRUE, and that is why one sees astr() listingrather than a summary as the default when the object is simply shown bytyping its name at the console.

Setting and viewing defaults

The user may have other preferences. She may want to see bothintervals and tests whenever contrasts are produced; and perhaps shealso wants to always default to the response scale when transformationsor links are present. We can change the defaults by setting thecorresponding options; and that is done via theemm_options() function:

emm_options(emmeans = list(type = "response"),            contrast = list(infer = c(TRUE, TRUE)))

Now, newemmeans() results and contrasts follow the newdefaults:

pigs.anal.p <- emmeans(pigs.lm, consec ~ percent)pigs.anal.p
## $emmeans##  percent response   SE df lower.CL upper.CL##        9     31.4 1.28 23     28.8     34.1##       12     37.5 1.44 23     34.7     40.6##       15     39.0 1.70 23     35.6     42.7##       18     42.3 2.24 23     37.9     47.2## ## Results are averaged over the levels of: source ## Confidence level used: 0.95 ## Intervals are back-transformed from the log scale ## ## $contrasts##  contrast              ratio     SE df lower.CL upper.CL null t.ratio p.value##  percent12 / percent9   1.20 0.0671 23    1.038     1.38    1   3.202  0.0112##  percent15 / percent12  1.04 0.0604 23    0.896     1.20    1   0.650  0.8613##  percent18 / percent15  1.09 0.0750 23    0.911     1.29    1   1.194  0.5202## ## Results are averaged over the levels of: source ## Confidence level used: 0.95 ## Conf-level adjustment: mvt method for 3 estimates ## Intervals are back-transformed from the log scale ## P value adjustment: mvt method for 3 tests ## Tests are performed on the log scale

Observe that the contrasts “inherited” thetype = "response" default from the EMMs.

NOTE: Setting the above options doesnot change how existingemmGrid objects are displayed; it only affects onesconstructed in the future.

There is one more option –summary – that overrides allother display defaults for both existing and future objects. Forexample, specifyingemm_options(summary = list(infer = c(TRUE, TRUE))) willresult in both intervals and tests being displayed, regardless of theirinternal defaults, unlessinfer is explicitly specified ina call tosummary().

To temporarily revert to factory defaults in a single call toemmeans() orcontrast() orpairs(), specifyoptions = NULL in the call.To reset everything to factory defaults (which we do presently),null-out all of theemmeans package options:

options(emmeans = NULL)

Optimal digits to display

When anemmGrid object is summarized and displayed, thefactory default is to display it with just enough digits as is justifiedby the standard errors or HPD intervals of the estimates displayed. Youmay use the"opt.digits" option to change this. If it isTRUE (the default), we display only enough digits as isjustified (but at least 3). If it is set toFALSE, thenumber of digits is set using the R system’s default,getOption("digits"); this is often much more precision thanis justified. To illustrate, here is the summary ofpigs.emm displayed without optimizing digits. Compare itwith the first summary in this vignette.

emm_options(opt.digits = FALSE)pigs.emm
##  source   emmean         SE df lower.CL upper.CL##  fish   3.394492 0.03668122 23 3.318612 3.470373##  soy    3.667260 0.03744798 23 3.589793 3.744727##  skim   3.796770 0.03938283 23 3.715300 3.878240## ## Results are averaged over the levels of: percent ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95
emm_options(opt.digits = TRUE)  # revert to optimal digits

By the way, setting this option doesnot round thecalculated values computed bysummary.emmGrid() or saved inasummary)emm object; it simply controls the precisiondisplayed byprint.summary_emm().

Startup options

The options accessed byemm_options() andget_emm_option() are stored in a list namedemmeans within R’s options environment. Therefore, if youdesire options other than the defaults provided on a regular basis, thiscan be easily arranged by specifying them in your startup script for R.For example, if you want to default to Satterthwaite degrees of freedomforlmer models, and display confidence intervals ratherthan tests for contrasts, your.Rprofile file could containthe line

options(emmeans = list(lmer.df = "satterthwaite",                        contrast = list(infer = c(TRUE, FALSE))))

Back to contents

Combining and subsettingemmGrid objects

Two or moreemmGrid objects may be combined using therbind() or+ methods. The most common reason(or perhaps the only good reason) to do this is to combine EMMs orcontrasts into one family for purposes of applying a multiplicityadjustment to tests or intervals. A user may want to combine the threepairwise comparisons of sources with the three comparisons above ofconsecutive percents into a single family of six tests with a suitablemultiplicity adjustment. This is done quite simply:

rbind(pairs(pigs.emm.s), pigs.anal.p[[2]])
##  contrast              estimate     SE df t.ratio p.value##  fish - soy             -0.2728 0.0529 23  -5.153  0.0002##  fish - skim            -0.4023 0.0542 23  -7.428 <0.0001##  soy - skim             -0.1295 0.0530 23  -2.442  0.1364##  percent12 - percent9    0.1796 0.0561 23   3.202  0.0238##  percent15 - percent12   0.0378 0.0582 23   0.650  1.0000##  percent18 - percent15   0.0825 0.0691 23   1.194  1.0000## ## Results are averaged over some or all of the levels of: percent, source ## Results are given on the log (not the response) scale. ## P value adjustment: bonferroni method for 6 tests

The default adjustment is"bonferroni"; we could havespecified something different via theadjust argument. Anequivalent way to combineemmGrids is via the additionoperator. Any options may be provided byupdate(). Below,we combine the same results into a family but ask for the “exact”multiplicity adjustment.

update(pigs.anal.p[[2]] + pairs(pigs.emm.s), adjust = "mvt")
##  contrast              ratio     SE df lower.CL upper.CL null t.ratio p.value##  percent12 / percent9  1.197 0.0671 23    1.022    1.402    1   3.202  0.0214##  percent15 / percent12 1.039 0.0604 23    0.881    1.224    1   0.650  0.9681##  percent18 / percent15 1.086 0.0750 23    0.894    1.320    1   1.194  0.7305##  fish / soy            0.761 0.0403 23    0.656    0.884    1  -5.153  0.0002##  fish / skim           0.669 0.0362 23    0.574    0.779    1  -7.428 <0.0001##  soy / skim            0.879 0.0466 23    0.756    1.020    1  -2.442  0.1109## ## Results are averaged over some or all of the levels of: source, percent ## Confidence level used: 0.95 ## Conf-level adjustment: mvt method for 6 estimates ## Intervals are back-transformed from the log scale ## P value adjustment: mvt method for 6 tests ## Tests are performed on the log scale

Also evident in comparing these results is that settings are obtainedfrom the first object combined. So in the second output, where they arecombined in reverse order, we get both confidence intervals and tests,and transformation to the response scale.

To subset anemmGrid object, just use the subscriptingoperator[]. For instance,

pigs.emm[2:3]
##  source emmean     SE df lower.CL upper.CL##  soy      3.67 0.0374 23     3.59     3.74##  skim     3.80 0.0394 23     3.72     3.88## ## Results are averaged over the levels of: percent ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95

Accessing results to use elsewhere

Sometimes, users want to use the results of an analysis (say, anemmeans() call) in other computations. Thesummary() method creates asummary_emm objectthat inherits from thedata.frame class; so one may use thevariables therein just as those in a data frame.

AnemmGrid object has its own internal structure and wecan’t directly access the values we see displayed. If follow-upcomputations are needed, usesummary() (orconfint() ortest()), creates asummary_emm object which inherits fromdata.frame – making it possible to access the values. Forillustration, let’s add the widths of the confidence intervals in ourexample.

CIs <- confint(pigs.emm)CIs$CI.width <- with(CIs, upper.CL - lower.CL)CIs
##  source emmean     SE df lower.CL upper.CL CI.width##  fish     3.39 0.0367 23     3.32     3.47    0.152##  soy      3.67 0.0374 23     3.59     3.74    0.155##  skim     3.80 0.0394 23     3.72     3.88    0.163## ## Results are averaged over the levels of: percent ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95

By the way, the values stored internally are kept to full precision,more than is typically displayed:

CIs$emmean
## [1] 3.394492 3.667260 3.796770

If you want to display more digits, specify so using theprint method:

print(CIs, digits = 5)
##  source emmean       SE df lower.CL upper.CL CI.width##  fish   3.3945 0.036681 23   3.3186   3.4704  0.15176##  soy    3.6673 0.037448 23   3.5898   3.7447  0.15493##  skim   3.7968 0.039383 23   3.7153   3.8782  0.16294## ## Results are averaged over the levels of: percent ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95

Back to contents

Adding grouping factors

Sometimes, users want to group levels of a factor into a smallernumber of groups. Those groups may then be, say, averaged separately andcompared, or used as aby factor. Theadd_grouping() function serves this purpose. The functiontakes four arguments: the object, the name of the grouping factor to becreated, the name of the reference factor that is being grouped, and avector of level names of the grouping factor corresponding to levels ofthe reference factor. Suppose for example that we want to distinguishanimal and non-animal sources of protein in thepigsexample:

pigs.emm.ss <- add_grouping(pigs.emm.s, "type", "source",                            c("animal", "vegetable", "animal"))str(pigs.emm.ss)
## 'emmGrid' object with variables:##     source = fish, soy, skim##     type = animal, vegetable## Nesting structure:  source %in% type## Transformation: "log"

Note that the new object has a nesting structure (see more about thisin the“messy-data” vignette),with the reference factor nested in the new grouping factor. Now we canobtain means and comparisons for each group

emmeans(pigs.emm.ss, pairwise ~ type)
## $emmeans##  type      emmean     SE df  n lower.CL upper.CL##  animal      3.60 0.0267 23 19     3.54     3.65##  vegetable   3.67 0.0374 23 10     3.59     3.74## ## Results are averaged over the levels of: percent, source ## Results are given on the log (not the response) scale. ## Confidence level used: 0.95 ## ## $contrasts##  contrast           estimate     SE df t.ratio p.value##  animal - vegetable  -0.0716 0.0455 23  -1.573  0.1295## ## Results are averaged over the levels of: percent, source ## Results are given on the log (not the response) scale.

Back to contents

Re-labeling or re-leveling anemmGrid

Sometimes it is desirable to re-label the rows of anemmGrid, or cast it in terms of other factor(s). This canbe done via thelevels argument inupdate().

As an example, sometimes a fitted model has a treatment factor thatcomprises combinations of other factors. In subsequent analysis, we maywell want to break it down into the individual factors’ contributions.Consider, for example, thewarpbreaks data provided with R.We will define a single factor and fit a non homogeneous-variancemodel:

warp <- transform(warpbreaks, treat = interaction(wool, tension))library(nlme)warp.gls <- gls(breaks ~ treat, weights = varIdent(form = ~ 1|treat), data = warp)( warp.emm <- emmeans(warp.gls, "treat") )
##  treat emmean   SE   df lower.CL upper.CL##  A.L     44.6 6.03 8.02     30.7     58.5##  B.L     28.2 3.29 8.00     20.6     35.8##  A.M     24.0 2.89 8.00     17.3     30.7##  B.M     28.8 3.14 8.00     21.5     36.0##  A.H     24.6 3.42 8.00     16.7     32.5##  B.H     18.8 1.63 8.00     15.0     22.5## ## Degrees-of-freedom method: satterthwaite ## Confidence level used: 0.95

But now we want to re-cast thisemmGrid into one thathas separate factors forwool andtension. Wecan do this as follows:

warp.fac <- update(warp.emm, levels = list(                wool = c("A", "B"), tension = c("L", "M", "H")))str(warp.fac)
## 'emmGrid' object with variables:##     wool = A, B##     tension = L, M, H

So now we can do various contrasts involving the separatefactors:

contrast(warp.fac, "consec", by = "wool")
## wool = A:##  contrast estimate   SE   df t.ratio p.value##  M - L     -20.556 6.69 11.5  -3.074  0.0203##  H - M       0.556 4.48 15.6   0.124  0.9899## ## wool = B:##  contrast estimate   SE   df t.ratio p.value##  M - L       0.556 4.55 16.0   0.122  0.9881##  H - M     -10.000 3.54 12.0  -2.824  0.0269## ## Degrees-of-freedom method: satterthwaite ## P value adjustment: mvt method for 2 tests

Note: When re-leveling to more than one factor, you have to becareful to anticipate that the levels will be expanded usingexpand.grid(): the first factor in the list varies thefastest and the last varies the slowest. That was the case in ourexample, but in others, it may not be. Had the levels oftreat been ordered asA.L, A.M, A.H, B.L, B.M, B.H, then we would have had tospecify the levels oftension first and the levels ofwool second.

Back to contents

Index of all vignette topics


[8]ページ先頭

©2009-2025 Movatter.jp