0
$\begingroup$

I want to assess whether the differences I find for different groups (of a factor variable) are significantly different from each other. Let's keep it hypothetical and assume a data structure like that:

Value     Group     A_binary     B_binary     C_binary
9         A         1            0            0
5         A         1            0            0
7         A         1            0            0
4         B         0            1            0
3         B         0            1            0
1         C         0            0            1

Now the underlying question is: "Does group affiliation affect the value in a significant way; that is, are the means found for A, B, C significant from each other?".

Now I probably run a t-test to see whether the mean "value" for one group differs from the one of another and so forth. I'd have to do this for every pair: A and B, B and C, A and C.

I'm wondering whether there's a neat, easy and fast way to do this in R?

$\endgroup$
3
$\begingroup$

Analysis of variance

If you want to obtain just a single global p-value then this would simply be a one-way analysis of variance (ANOVA). This can be carried out in a number of ways in R including:

oneway.test(Sepal.Width ~ Species, data = iris, var.equal = TRUE)
##  One-way analysis of means
## 
## data:  Sepal.Width and Species
## F = 49.16, num df = 2, denom df = 147, p-value < 2.2e-16

or

m <- lm(Sepal.Width ~ Species, data = iris)
summary(m)
## [...]
## F-statistic: 49.16 on 2 and 147 DF,  p-value: < 2.2e-16

But there are many other routes to obtain this test or variations of it, e.g., adjusting for heteroscedasticity etc.

Pairwise t-tests

If you want to conduct all pairwise t-tests then this is a multiple testing problem and the corresponding p-values should be adjusted for that. The pairwise.t.test() function in base R has a number of methods for this, defaulting to Holm's.

pairwise.t.test(iris$Sepal.Width, iris$Species)
##  Pairwise comparisons using t tests with pooled SD 
## 
## data:  iris$Sepal.Width and iris$Species 
## 
##            setosa  versicolor
## versicolor < 2e-16 -         
## virginica  9.1e-10 0.0031    
## 
## P value adjustment method: holm 

A more general approach would be to consider this as a set of linear hypothesis tests in a linear model, using Tukey contrasts. This can be done carried out in package multcomp supporting both the global test (ANOVA), the adjusted pairwise t-tests, and a number of further adjustment methods:

library("multcomp")
g <- glht(m, linfct = mcp(Species = "Tukey"))
summary(g, test = adjusted("holm"))
##   Simultaneous Tests for General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## Fit: lm(formula = Sepal.Width ~ Species, data = iris)
## 
## Linear Hypotheses:
##                             Estimate Std. Error t value Pr(>|t|)    
## versicolor - setosa == 0    -0.65800    0.06794  -9.685  < 2e-16 ***
## virginica - setosa == 0     -0.45400    0.06794  -6.683 9.08e-10 ***
## virginica - versicolor == 0  0.20400    0.06794   3.003  0.00315 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- holm method)
## 
summary(g, test = Ftest())
##   General Linear Hypotheses
## 
## Multiple Comparisons of Means: Tukey Contrasts
## 
## Linear Hypotheses:
##                             Estimate
## versicolor - setosa == 0      -0.658
## virginica - setosa == 0       -0.454
## virginica - versicolor == 0    0.204
## 
## Global Test:
##       F DF1 DF2    Pr(>F)
## 1 49.16   2 147 4.492e-17
| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ (+1) I think the post-hoc tests section would be better placed under the ANOVA heading. $\endgroup$ – mkt - Reinstate Monica Apr 19 '18 at 8:54
  • $\begingroup$ Okay that makes a lot of sense. Quick question abou the pairwise.t.test() output though: what does the "-" sign tell me in this case? Highest significance of p<0.000? Or the opposite? Also, I was wondering if I could simply run a linear model like that on my (hypothetical) data: summary(lm(testdat$Value ~ testdat$Group)) This would make Group_A my reference category to which I can then compare the other groups. In this case we would see that group C is significantly different from A (alpha = 5%), but this is not the case for group B. Does that make sense? $\endgroup$ – Fabian Habersack Apr 19 '18 at 9:07
  • $\begingroup$ @FabianHabersack (1) The - signals that the corresponding pair is not tested (for versicolor vs. versicolor there is nothing to test). (2) Never use lm(testdat$Value ~ testdat$Group), use lm(Value ~ Group, data = testdat) instead. (3) Instead of all pairwise tests (= Tukey contrasts) you can also use all against one (= Dunnett contrasts). The uncorrected tests are reported by summary(m). For multiple testing adjustment summary() for glht(m, linfct = mcp(Species = "Dunnett")) can be used. $\endgroup$ – Achim Zeileis Apr 19 '18 at 19:57
  • $\begingroup$ @mkt Most of these tests can be understood as either post-hoc tests or ANOVAs or "normal" significance tests. There are sometimes subtle differences wrt variance estimation, finite-sample or degrees-of-freedom correction, conditional vs. unconditional justifications etc. But for the examples chosen in my answer, the test statistics and p-values from the different views coincide exactly. $\endgroup$ – Achim Zeileis Apr 19 '18 at 20:02
  • 1
    $\begingroup$ @AchimZeileis I don't know - probably doesn't matter too much either way. $\endgroup$ – mkt - Reinstate Monica Apr 20 '18 at 4:48

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.