Consistency of emmeans tidiers (with other tidiers)

[crsh]

While working on a PR to add a new emmeans-tidier, I noticed that the emmeans-tidiers have some internal and external inconsistencies:

1. The lsmobj-method uses the common arguments conf.int and conf.level as defined in the param_confint template. The other methods (e.g., emmGrid) do not provide these arguments and instead rely on the argument names native to the emmeans summary()-methods (e.g., infer and level).

I haven't looked exhaustively at other methods, but I additionally noticed some inconsistencies compared to other contrast tidiers, specifically tidy.TukeyHSD():

2. tidy.TukeyHSD() reports the contrasted conditions in a column labelled comparison in the form of a-b. In contrast, the emmeans tidiers return the same information in two columns labelled level1 and level2 (containing a and b).

fm1 <- aov(breaks ~ wool + tension, data = warpbreaks)
thsd <- TukeyHSD(fm1, "tension", ordered = TRUE)
tidy(thsd)

emmp <- pairs(emmeans(fm1, ~ tension))
tidy(emmp)

3. In the tibble returned by tidy.TukeyHSD(), the column containing p-values is labelled adj.p.value. In contrast, the emmeans tidiers label this column p.value regardless of whether it has been adjusted for multiple comparisons or not (see code above). Unless I missed something, the use of adj.p.value is currently unique to tidy.TukeyHSD().

It seems desirable to try to keep things consistent across methods where possible but particularly within the set of tidiers for a given package. I would, therefore, suggest the following changes, that I'd be willing to implement in a PR:

1. Add the arguments conf.int and conf.level to all emmeans tidiers.
2. Change reporting of contrast pairs in either tidy.TukeyHSD() or the emmeans-methods. I'm not sure which of the two is preferable here.
3. Either use adj.p.value in emmeans tidiers whenever p-values are adjusted for multiple comparisons or use p.value in tidy.TukeyHSD(). Again I'm not sure which is preferable.

Any thoughts?

[alexpghayes]

1. Yes! Please do!
2. Your call. If you are uncertain why don't you post example output from both here and we can discuss more.
3. Use adj.p.value in the emmeans tidiers.

Thanks for taking this on! Let me know how I can help!

[crsh]

2. Your call. If you are uncertain why don't you post example output from both here and we can discuss more.

Sure. I've looked around a little and it seems tidy.TukeyHSD() used to have an argument to separate the comparison column into level1 and level2 but it was removed recently to be consistent with the multcomp tidiers. So might make sense to also look at those tidiers. Doing so, I've noticed a couple of other things:

fm1 <- aov(breaks ~ wool + tension, data = warpbreaks)
thsd <- TukeyHSD(fm1, "tension")
as.data.frame(tidy(thsd))

     term comparison  estimate   conf.low  conf.high adj.p.value
1 tension        M-L -10.000000 -19.35342 -0.6465793 0.447421021
2 tension        H-L -14.722222 -24.07564 -5.3688015 0.001121788
3 tension        H-M  -4.722222 -14.07564  4.6311985 0.033626219

library("emmeans")
emmp <- pairs(emmeans(fm1, ~ tension))
as.data.frame(tidy(emmp, infer = c(T, T), level = 0.95))

  level1 level2  estimate std.error df   conf.low  conf.high statistic     p.value
1      M      L -10.000000  3.872378 50 -19.35342 -0.6465793 -2.582393 0.033626219
2      H      L -14.722222  3.872378 50 -24.07564 -5.3688015 -3.801856 0.001121788
3      H      M  -4.722222  3.872378 50 -14.07564  4.6311985 -1.219463 0.447421021

4. Unless I'm mistaken, I think the column df should be renamed to parameter, as is the case in other methods, right?
5. In contrast to the tidy.TukeyHSD()-output, there is no term column. Should there be one?

library("multcomp")
glhs <- summary(glht(fm1, linfct = mcp(tension = "Tukey")))
as.data.frame(tidy(glhs))

    lhs rhs   estimate std.error statistic     p.value
1 M - L   0 -10.000000  3.872378 -2.582393 0.033699296
2 H - L   0 -14.722222  3.872378 -3.801856 0.001109948
3 H - M   0  -4.722222  3.872378 -1.219463 0.447461036

6. The multcomp tidiers require manually calling summary() prior to tidying, whereas emmeans tidiers do so for you. Consequently, it's not possible to obtain test statistics and confidence interval (tidiers are separate) and there are no conf.int or conf.level arguments. tidy.glht seems to be what quick = TRUE should do. It returns only lhs, rhs, and estimate.

7. As discussed with respect to emmeans, these tidiers should probably also use adj.p.value, if applicable, rather than p.value.

8. In contrast to the tidy.TukeyHSD()-output, there is no term column. Should there be one?

Regarding point 2.: I like the multcomp approach here. It specifies the contrast and the null value. I'm not sure about the column names. Are lhs and rhsused like this in other tidiers? Otherwise contrast and null.value might be more descriptive?

[crsh]

Sorry to bother, @alexpghayes, but do you have any thoughts on this? I'd like to tackle this soon; some developments in my own package depend on this.

[alexpghayes]

Apologies!

2. Your call. If you are uncertain why don't you post example output from both here and we can discuss more.

...it seems tidy.TukeyHSD() used to have an argument to separate the comparison column into level1 and level2 but it was removed recently to be consistent with the multcomp tidiers.

I am happy with or without this argument, although I think providing it would be nice. My main concern is a consistent interface.

4. Unless I'm mistaken, I think the column df should be renamed to parameter, as is the case in other methods, right?

I don't think so. We typically do this for htest objects but not regression / contrast tables.

5. In contrast to the tidy.TukeyHSD()-output, there is no term column. Should there be one?

I'm not sure. If there is a way to select between multiple contrasts, or in general get contrasts for multiple categories features, then yes.

6. The multcomp tidiers require manually calling summary() prior to tidying, whereas emmeans tidiers do so for you. Consequently, it's not possible to obtain test statistics and confidence interval (tidiers are separate) and there are no conf.int or conf.level arguments. tidy.glht seems to be what quick = TRUE should do. It returns only lhs, rhs, and estimate.

Let's move over to the emmeans approach then! Also, we're getting rid of quick, so no need to implement special behavior there.

7. As discussed with respect to emmeans, these tidiers should probably also use adj.p.value, if applicable, rather than p.value.

Agree.

Regarding point 2.: I like the multcomp approach here. It specifies the contrast and the null value. I'm not sure about the column names. Are lhs and rhsused like this in other tidiers? Otherwise contrast and null.value might be more descriptive?

I like contrast and null.value as well!

[crsh]

Thanks, so to summarize I'll do the following:

emmeans

• add the arguments conf.int and conf.level
• if applicable, change p.value to adj.p.value
• if applicable, add term column
• use contrast (with A - B) and null.value columns

multcomp

• include summary()-call
• if applicable, change p.value to adj.p.value
• if applicable, add term column
• use contrast (with A - B) and null.value columns

TukeyHSD

• use contrast (with A - B) and null.value columns

[alexpghayes]

Exactly!