I have data for three groups that I need to compare the means of. Currently, I am running a multiple linear regression with a single categorical feature to assess whether my test groups differ from my control group (and by how much).

import numpy as np
import pandas as pd
import statsmodels.formula.api as smf

n = 1000

control = pd.DataFrame({'group': ['control']*n,
           'target': np.random.normal(0.0, 1.0, size=(n,))

test1 = pd.DataFrame({'group': ['test1']*n,
           'target': np.random.normal(5.0, 1.0, size=(n,))

test2 = pd.DataFrame({'group': ['test2']*n,
           'target': np.random.normal(10.0, 1.0, size=(n,))

data = pd.concat([control, test1, test2], axis=0)

levels = ['control', 'test1', 'test2']

mod = smf.glm(formula="target~ C(group, levels=levels)", data=data).fit()

I can extend this to see how much the test groups differ with the following code:

levels = ['test1', 'test2']

mod = smf.glm(formula="target ~ C(group, levels=levels)", data=data[data['group'].isin(['test1', 'test2'])]).fit()

However, I'm not convinced that running a second regression is a good idea (I think I will have to correct for multiple testing).

My question is: can I run this all in a single regression? If the answer is yes, how would I do it?


1 Answer 1


Doing regressions only on 2-group subsets of the data throws away the information and the extra number of data points that you already have incorporated into the model based on all 3 groups.

The general case is a pretty standard situation: ANOVA as a first step, followed by some post-hoc test of particular comparisons, either pre-specified in the design or looking at all possible comparisons among groups. Those post-hoc tests are based on the results of the complete ANOVA, combining information from all groups.

In your particular case, if all you want to know is whether each of the two test groups differs from control, then do simple dummy coding (also called treatment coding) of your group categorical predictor, with the control group taken as reference. Then linear regression give you coefficients and standard errors for the differences of each of the 2 treatment groups from the control, to form the basis of statistical tests.

Technically one should correct for doing two comparisons but in practice that isn't always done. With a simple case like this a reader can always gauge for herself whether that correction would have made a substantive difference in the conclusions of your study.

In more general cases, there are many different approaches for testing particular comparisons among groups after ANOVA while taking multiple comparisons into account.

Dunnett's test was designed specifically for comparing multiple test groups against a control. Tukey's test examines all pairwise comparisons. The Newman-Keuls method is more powerful but less conservative than Tukey's test. Duncan's test goes even further in the direction of higher power but greater risk of false positives. Scheffé's method applies to all possible contrasts combining group means, not just pairwise comparisons.

When there are very many tests being performed as in gene-expression analysis, control of false-discovery rate is generally preferred to the above types of test that control family-wise error rate.

So stick with your 3-group model and get the comparisons and tests you need out of that.

  • $\begingroup$ Sorry, I wasn't as clear as I could be. 'In your particular case, if all you want to know is whether each of the two test groups differs from control' that is what I want, but in addition I want to know if the two test groups differ from each other, and can this be run in a single regression? $\endgroup$
    – Tom Kealy
    Commented Jun 8, 2020 at 18:28
  • $\begingroup$ @TomKealy for testing all pairwise comparisons Tukey's range test, linked in the answer, would be a well-accepted choice. When run on the full model it uses information from all groups to get estimates of standard errors, and inherently corrects for multiple comparisons. It does require assumptions that values within each group are normally distributed and that the within-group variance values are the same among groups, but those are already implicit in your use of linear regression. It's always good to check those assumptions and to consider alternate modeling strategies if they aren't met. $\endgroup$
    – EdM
    Commented Jun 9, 2020 at 14:37

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.