2
$\begingroup$

What test for pairwise comparisons does emmeans uses by default, when executing

    emmeans(model, pairwise~predictor)?

As far as I can understand the Tukey method (Tukey HSD) is used by default just for p-values adjustment, not for pairwise comparisons by themselves. My R knowledge is too poor to deconstruct the raw code of emmeans on Github, so hope someone will shed light on the issue.

$\endgroup$
2
  • $\begingroup$ I don't understand your question. The Tukey method IS implementable by simply adjusting the P values for each t ratio. The P values are obtained using the Studentized range distribution instead of the t distribution. $\endgroup$
    – Russ Lenth
    Sep 18, 2020 at 21:16
  • $\begingroup$ Thank you for a clarification. What i meant is that the Tukey test is used to adjust the P values when 'method = "tukey"' flag is noted in emmeans command, what is the default option. However, you can mention 'method="none"' what will lead to anadjusted P values. So, my question was about the default test performed, that emits the initial P values, that later will be adjusted by Tukey, Bonferroni or any other way. $\endgroup$ Sep 21, 2020 at 8:08

1 Answer 1

4
$\begingroup$

OK, I think I understand it now. The answer is that the term "adjusted P value" doesn't necessarily mean that a set of P values is obtained, and then we adjust those values somehow. There are some cases where this is true, but others where we compute them differently to begin with.

To illustrate, let's look at one of the built-in datasets in emmeans

> fiber.lm <- lm(strength ~ diameter + machine, data=fiber)
> library(emmeans)
> emm <- emmeans(fiber.lm, "machine")
> pairs(emm, adjust = "none")
 contrast estimate   SE df t.ratio p.value
 A - B       -1.04 1.01 11 -1.024  0.3280 
 A - C        1.58 1.11 11  1.431  0.1803 
 B - C        2.62 1.15 11  2.283  0.0433 

Let's do these manually. First, get the t ratios:

> t.rat <- c(-1.024, 1.431, 2.283)

Calculate the unadjusted P values; these are twice the right-hand tail areas:

> (upv <- 2 * (1 - pt(abs(t.rat), 11)))
[1] 0.32782764 0.18022061 0.04330727 

These match the results from pairs().

Now, if we want a Bonferroni adjustment, we adjust these by multiplying by the number of tests:

> 3 * upv
[1] 0.9834829 0.5406618 0.1299218

You can verify this using pairs(emm, adjust = "bonf") (results not shown).

If we want a Tukey adjustment, we use the Studentized Range distribution and compute the right-hand tail areas; we also have to multiply the t ratios by the square root of 2 because of how the distribution is scaled:

> 1 - ptukey(sqrt(2) * abs(t.rat), 3, 11)
[1] 0.5778266 0.3595113 0.1006032

You can verify this with the software:

> pairs(emm, adjust = "tukey")
 contrast estimate   SE df t.ratio p.value
 A - B       -1.04 1.01 11 -1.024  0.5781 
 A - C        1.58 1.11 11  1.431  0.3596 
 B - C        2.62 1.15 11  2.283  0.1005 

P value adjustment: tukey method for comparing a family of 3 estimates

Note that the Bonferroni adjustment is indeed an adjustment to the regular P values upv. But the Tukey-adjusted P values are an entirely separate calculation based on the t ratios but not on upv. So there is no "default" set of P values underlying the Tukey-adjusted P values.

I hope this adds clarity.

$\endgroup$
3
  • $\begingroup$ Thank for an exellent explanation! $\endgroup$ Sep 29, 2020 at 11:57
  • $\begingroup$ @Russ Lenth, I have a question regarding your answer. Since the results above using Tukey adjustment and none adjustment give very different p values for B-C comparison (significant p value without adjustment and nonsignificant p value with tukey adjustment), which one should I report to describe the difference between B and C? $\endgroup$
    – Chloe
    Jan 10, 2021 at 10:05
  • $\begingroup$ People have different opinions about that. But I suggest that in ordinary situations, one should use the Tukey adjustment for pairwise comparisons. People have different rules based on what was planned ahead and orthogonality, but too often those are abused and it's best to go with an accepted standard. $\endgroup$
    – Russ Lenth
    Jan 10, 2021 at 14:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.