1
$\begingroup$

It's not uncommon to read a scientific article in which the authors perform an ANOVA followed by post-hoc t-tests, but don't report the results of the ANOVA at all, and instead skip straight into the results of the post-hoc t-tests. Assuming all post-hocs are reported, is there any way that a reader could back-calculate the results of the ANOVA? For simplicity, let's say it's a 2-way ANOVA.

Edit here is an example:

Load data

library(palmerpenguins)
library(emmeans)

data(package = 'palmerpenguins')
dat = penguins

Categorical x continuous interaction

model1 = lm(body_mass_g  ~ species*flipper_length_mm ,data=dat)
# anova
anova(model1)
Analysis of Variance Table

Response: body_mass_g
                           Df    Sum Sq  Mean Sq F value    Pr(>F)    
species                     2 146864214 73432107 534.660 < 2.2e-16 ***
flipper_length_mm           1  24776495 24776495 180.398 < 2.2e-16 ***
species:flipper_length_mm   2   1519564   759782   5.532  0.004327 ** 
Residuals                 336  46147424   137344                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# post-hoc tests
res1 = emtrends(model1, pairwise ~ species, var = "flipper_length_mm")
res1$contrasts
 contrast           estimate   SE  df t.ratio p.value
 Adelie - Chinstrap    -1.74 7.86 336  -0.222  0.9733
 Adelie - Gentoo      -21.79 6.94 336  -3.139  0.0052
 Chinstrap - Gentoo   -20.05 8.19 336  -2.448  0.0394

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

Categorical x categorical interaction

model2 = lm(body_mass_g  ~ species*sex,data=dat)
# anova
anova(model2)
Analysis of Variance Table

Response: body_mass_g
             Df    Sum Sq  Mean Sq F value    Pr(>F)    
species       2 145190219 72595110 758.358 < 2.2e-16 ***
sex           1  37090262 37090262 387.460 < 2.2e-16 ***
species:sex   2   1676557   838278   8.757 0.0001973 ***
Residuals   327  31302628    95727                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# post-hoc tests
res2 = emmeans(object = model2, pairwise ~ species|sex)
res2$contrasts
sex = female:
 contrast           estimate   SE  df t.ratio p.value
 Adelie - Chinstrap     -158 64.2 327  -2.465  0.0377
 Adelie - Gentoo       -1311 54.4 327 -24.088  <.0001
 Chinstrap - Gentoo    -1153 66.8 327 -17.246  <.0001

sex = male:
 contrast           estimate   SE  df t.ratio p.value
 Adelie - Chinstrap      105 64.2 327   1.627  0.2357
 Adelie - Gentoo       -1441 53.7 327 -26.855  <.0001
 Chinstrap - Gentoo    -1546 66.2 327 -23.345  <.0001

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

In either case, is it possible to compute any of the ANOVA results, using only the post-hoc test results?

$\endgroup$
2
  • 1
    $\begingroup$ Please post an example dataset.csv, some (simple!) ANOVA results, and t-test results. If you also include the code you used, that would be helpful to folks responding. $\endgroup$
    – J_H
    Jan 8 at 2:46
  • 1
    $\begingroup$ I've added some examples. $\endgroup$
    – David B
    Jan 8 at 16:28

1 Answer 1

2
$\begingroup$

In principle, this may be possible if you have a completely randomized design. When there are blocks, pairing, or mixed effects, it gets harder or impossible.

Simplest case is with a table of cell means and their standard errors:

> emmeans(model2, ~sex*species)
 sex    species   emmean   SE  df lower.CL upper.CL
 female Adelie      3369 36.2 327     3298     3440
 male   Adelie      4043 36.2 327     3972     4115
 female Chinstrap   3527 53.1 327     3423     3632
 male   Chinstrap   3939 53.1 327     3835     4043
 female Gentoo      4680 40.6 327     4600     4760
 male   Gentoo      5485 39.6 327     5407     5563

Confidence level used: 0.95 

Then you can use the emmobj function to construct an emmGrid object:

> eobj = emmobj(c(3369,4043,3527,3939,4680,5485), V=diag(c(36.2,36.2,53.1,53.1,40.6,39.6)^2),
+     levels = list(sex=c("f","m"), species=c("a","c","g")), df = 327)

> eobj # shows same table as before
 sex species estimate   SE  df lower.CL upper.CL
 f   a           3369 36.2 327     3298     3440
 m   a           4043 36.2 327     3972     4114
 f   c           3527 53.1 327     3423     3631
 m   c           3939 53.1 327     3835     4043
 f   g           4680 40.6 327     4600     4760
 m   g           5485 39.6 327     5407     5563

Confidence level used: 0.95 

> joint_tests(eobj)
 model term  df1 df2 F.ratio p.value
 sex           1 327 311.580  <.0001
 species       2 327 747.804  <.0001
 sex:species   2 327   8.735  0.0002

This is comparable to the Type III anova for the model:

> joint_tests(model2)
 model term  df1 df2 F.ratio p.value
 species       2 327 746.924  <.0001
 sex           1 327 311.838  <.0001
 species:sex   2 327   8.757  0.0002

But the anova you show is the Type I (sequential) anova, and you can't get that from eobj.

Now, with only the pairwise comparisons shown in the illustration, you do't have any information about the sex effects or sex:species effects. If you had also the results for pairwise comparisons of sex for each species, I think it may be possible to cobble a table of means and SEs, though with a different grand mean, by arbitrarily setting one of the cell means to zero and working from there. But it gets messy because each SE you see is $\sqrt{1/2}$ times the square root of the sum of two cell variances. In principle, it's possible to work through this but it won't be fun.

It is not possible to do the anova for model1 where we have comparisons of slopes. Again it may be possible if we know the slope estimates themselves, and some kind of information about the interaction. But it would be very ugly.

Again, I emphasize these ideas won't work for mixed models or others where cell means are in any way correlated.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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