Which multiple comparison correction should I use?

Quick answer. Choose the correction from the comparisons you planned before seeing the p-values. All pairs of groups after a one-way ANOVA: Tukey’s test (Games–Howell if you will not assume equal variances — decide before seeing the data). Each treatment against one control: Dunnett’s test (if you will not assume equal variances, decided in advance: Welch t-tests against control, adjusted with Holm). A handful of separate hypotheses, any of which you would report as a finding: Holm — it is always at least as powerful as Bonferroni. Hundreds or thousands of tests in a screen (genes, proteins): Benjamini–Hochberg, which controls the false discovery rate. One comparison you named in advance as the primary question — in a protocol, registration or dated plan, not chosen after the results: no correction; label everything else exploratory.

Terms used below. A family is the set of tests you would report together in one conclusion. An adjusted p-value has already been corrected for the size of the family; compare it with 0.05 directly. Tukey’s, Dunnett’s and Games–Howell results come out adjusted. The familywise error rate is the chance of at least one false positive in the family. The false discovery rate is the expected share of false positives among the results you call significant. Tukey’s HSD (honestly significant difference) compares every pair of groups. Tests are positively associated when a high value on one tends to go with high values on the others; if you cannot tell, use Holm. Step-down and step-up describe the order in which a procedure compares the sorted p-values. Holm–Šidák is a close variant of Holm; state whichever your code ran.

Decide by the question you are asking

Your comparisons	Use	Controls	R function	Python function
One pre-specified primary comparison	No correction; label other tests exploratory	—	—	—
Every pair of groups after one-way ANOVA	Tukey’s HSD	Familywise error	TukeyHSD(aov(y ~ group))	scipy.stats.tukey_hsd(*groups)
Every pair, equal variances not assumed (decided in advance; see Welch’s ANOVA)	Games–Howell	Familywise error	rstatix::games_howell_test(d, y ~ group)	scipy.stats.tukey_hsd(*groups, equal_var=False)
Each treatment against one control, equal variances assumed	Dunnett’s test	Familywise error	multcomp::glht(fit, mcp(group = "Dunnett"))	scipy.stats.dunnett(*treatments, control=control)
Each treatment against one control, equal variances not assumed (decided in advance)	Welch t-tests against control, then Holm	Familywise error	t.test(...), then p.adjust(p, "holm")	ttest_ind(..., equal_var=False), then multipletests(p, method="holm")
A small set of separate tests (several outcomes, selected contrasts)	Holm	Familywise error	p.adjust(p, "holm")	multipletests(p, method="holm")
Many tests in a screen, where a list of candidates is the output	Benjamini–Hochberg	False discovery rate	p.adjust(p, "BH")	scipy.stats.false_discovery_control(p)

The table names the functions; complete code you can run is in the worked examples below. Two error rates sit behind the table. The familywise error rate is the chance of at least one false positive among all your tests; use it when a single false claim would matter. The false discovery rate is the expected share of false positives among the results you call significant; it is the better target when you will follow up a list of hits and can tolerate a few duds (Benjamini & Hochberg 1995).

Write the row you will use into the analysis plan before you see the data, and make one choice per family of tests.

If nothing was written down, treat every comparison you ran as the family: all pairs (Tukey or Games–Howell) for group comparisons, Holm for separate outcomes, Benjamini–Hochberg for a screen. Say in Methods that the comparisons were not pre-specified. If you already hold a p-value, check whether it is adjusted: output from TukeyHSD, tukey_hsd, dunnett or glht is; a p-value from t.test or ttest_ind is not.

Why uncorrected tests mislead

Each test at α = 0.05 has a 5% chance of a false positive when there is no real effect. Across m independent tests the chance of at least one is 1 − 0.95^m:

Tests (m)	Chance of at least one false positive
1	5%
3	14%
5	23%
10	40%
20	64%
50	92%

Put as frequencies: if 100 labs each compared 10 genes that truly do not differ, about 40 of them would report at least one “significant” gene. Corrections bring that back to about 5 labs (familywise control) or keep the share of false hits among reported hits near 5% (false discovery control).

Count every test you would have reported as a finding — every pair, outcome or time point — when you set m.

The same six p-values, four corrections

Six p-values from six separate outcome measures (constructed example), adjusted with R’s p.adjust and statsmodels’ multipletests. Each cell is the adjusted p-value; compare it with 0.05.

Raw p	Bonferroni	Holm	Hochberg (1988) step-up, independent or positively associated tests only	Benjamini–Hochberg (1995), false discovery rate
0.003	0.018	0.018	0.018	0.018
0.009	0.054	0.045	0.045	0.024
0.012	0.072	0.048	0.048	0.024
0.030	0.180	0.090	0.082	0.045
0.041	0.246	0.090	0.082	0.049
0.270	1.000	0.270	0.270	0.270
Below 0.05	1	3	3	5

p <- c(0.003, 0.009, 0.012, 0.030, 0.041, 0.27)
p.adjust(p, method = "holm")   # also "bonferroni", "hochberg", "BH"

from statsmodels.stats.multitest import multipletests
from scipy import stats

p = [0.003, 0.009, 0.012, 0.030, 0.041, 0.27]
print(multipletests(p, method="holm")[1])   # also "bonferroni", "simes-hochberg", "fdr_bh"
print(stats.false_discovery_control(p))     # Benjamini–Hochberg

What the example shows:

• Holm replaced Bonferroni and found two more effects with the same familywise protection. R’s documentation puts it plainly: “There seems no reason to use the unmodified Bonferroni correction because it is dominated by Holm’s method, which is also valid under arbitrary assumptions” (R p.adjust).
• Hochberg is at least as powerful as Holm (here 0.082 against 0.090) but is valid only when the tests are independent or positively associated (Hochberg 1988; R p.adjust). Holm stays the default: switch to Hochberg only if you planned it and are sure the tests are independent or positively associated; if unsure, use Holm.
• Benjamini–Hochberg flagged five, because it promises something different: that, on average, no more than 5% of the flagged results are false. It holds for independent tests and for positively dependent ones, including one-sided comparisons of many treatments with one control (Benjamini & Yekutieli 2001). For arbitrary dependence use the more conservative "BY" method.

Use Holm for a small family of tests (Hochberg only if planned and the tests are independent or positively associated) and Benjamini–Hochberg for a screen.

Treatments against one control only: Dunnett (for all pairs, Tukey)

Constructed example: a control and three doses, six animals each. If you planned to compare each dose only with control, use Dunnett’s test: Tukey’s test also spends power on dose-versus-dose comparisons you did not plan. If you also want dose-versus-dose comparisons, Dunnett does not cover them; use Tukey. Both pool the variances of all groups; if you decided in advance not to assume equal variances, use Welch t-tests against control with Holm’s adjustment (Welch’s ANOVA).

Dose vs control	Unadjusted p (Student’s t-test, dose vs control only)	Tukey (6 comparisons)	Dunnett (3 comparisons)
Low (+0.98)	0.027	0.064	0.037
Mid (+1.65)	0.0009	0.0012	< 0.001
High (+2.37)	0.00006	0.00002	< 0.0001

The low dose reaches significance with Dunnett’s test (0.037) but stays above 0.05 with Tukey’s (0.064), on the same data, with both controlling the familywise error rate for the comparisons they cover. R’s multcomp agrees: 0.0375 for the low dose with set.seed(2026), against SciPy’s 0.037. (Dunnett’s and Tukey’s tests use the pooled variance from all four groups, so their p-values can be smaller than a two-group t-test’s.) In the paper, give each dose’s difference from control with the adjusted p-value from the test you planned.

d <- data.frame(y = c(10.2, 11.5, 9.8, 10.9, 11.1, 10.4,   11.1, 12.2, 10.7, 11.9, 12.5, 11.4,
                      11.9, 12.8, 11.4, 12.6, 13.0, 12.1,   12.6, 13.5, 12.2, 13.1, 13.9, 12.8),
                dose = rep(c("control", "low", "mid", "high"), each = 6))
d$dose <- factor(d$dose, levels = c("control", "low", "mid", "high"))  # control must be the first level
fit <- aov(y ~ dose, data = d)
TukeyHSD(fit)                            # all pairs
library(multcomp)
set.seed(2026)                           # Dunnett p-values use random integration
summary(glht(fit, linfct = mcp(dose = "Dunnett")))   # each dose vs control (first level)

from scipy import stats
import numpy as np

control = [10.2, 11.5, 9.8, 10.9, 11.1, 10.4]
low = [11.1, 12.2, 10.7, 11.9, 12.5, 11.4]
mid = [11.9, 12.8, 11.4, 12.6, 13.0, 12.1]
high = [12.6, 13.5, 12.2, 13.1, 13.9, 12.8]

print(stats.tukey_hsd(control, low, mid, high))                # all pairs
print(stats.dunnett(low, mid, high, control=control,
                    rng=np.random.default_rng(2026)))          # each dose vs control

Dunnett’s p-values come from numerical integration with a random element in both SciPy and R’s multcomp (which prints the integration error with the result). Values near 0.05 are stable to about two decimals; very small values vary more: in three unseeded SciPy runs of this example the high-dose p-value ranged from 0.000005 to 0.00002. Report very small values as “< 0.001” and set a seed, as above, so the result is repeatable (SciPy dunnett; Dunnett 1955).

Three mistakes that change the answer

• Correcting twice. Tukey’s, Dunnett’s and Games–Howell p-values are already adjusted. Report them as they are; applying Bonferroni on top makes the analysis stricter than its stated error rate.
• Choosing the family after seeing the results. Decide in the plan which comparisons form one family (all pairs, or each dose against control). Dropping comparisons that “did not work” before correcting shrinks m after the fact.
• Reading “not significant after correction” as “no effect”. Report each difference with its confidence interval, adjusted or not, so readers see how large the effect could be.

See two-way ANOVA vs multiple t-tests for a two-factor design (genotype × treatment), where the comparisons you plan inside the model form the family.

Check what your code (or AI assistant) actually ran

Software defaults differ, and code from an AI assistant may keep them. Before you report a corrected result, find the method name in the code:

Function	Default	What to check
R p.adjust(p)	"holm"	Fine for a few tests; set "BH" explicitly for a screen.
statsmodels multipletests(p)	"hs" (Holm–Šidák)	Not Bonferroni or Holm; name the method you report.
SciPy false_discovery_control(p)	"bh"	Benjamini–Hochberg; use method="by" for arbitrary dependence.
SciPy tukey_hsd(*groups)	equal_var=True (Tukey HSD; Tukey–Kramer if group sizes differ)	Pass equal_var=False for Games–Howell.
SciPy dunnett(...)	random integration without a seed	Pass rng for repeatable p-values.

Defaults from the R p.adjust documentation, statsmodels multipletests and SciPy tukey_hsd, confirmed by running the installed versions listed above. In the six-value example, the statsmodels default gives 0.018, 0.044, 0.047, 0.087, 0.087 and 0.27 — close to Holm, yet different numbers, so a Methods section that says “Holm” while the code ran the default does not reproduce exactly.

Ask an assistant precisely: “Adjust these p-values with Holm’s method using p.adjust(p, method = "holm"), and print the method name with the result.” In Python: “…using multipletests(p, method="holm")[1] from statsmodels, and print the method name.” In R, also check that the control is the first factor level before glht(..., "Dunnett"): R sorts levels alphabetically, so “vehicle” can end up last.

When you do not need a correction

Adjustment is needed “in confirmatory studies whenever results from multiple tests have to be combined in one final conclusion and decision” (Bender & Lange 2001). Three common cases:

• One primary comparison stated in advance. Test it at 0.05. Report other tests as exploratory, with intervals.
• Exploratory screening that will be confirmed in a new experiment. Prefer Benjamini–Hochberg adjusted p-values. If you report uncorrected p-values instead, say so, label every hit a candidate and confirm it in a new experiment.
• Separate questions in one paper (a toxicity result and an efficacy result, say): each is its own family only if the questions lead to separate conclusions and the families were set in the plan before the data came in.

See how to answer a reviewer about multiple comparisons if a reviewer asks for a correction anyway.

Write it in Methods

Name the family, the method and the error rate:

Each of the three doses was compared with the vehicle control using Dunnett’s test, which controls the familywise error rate at 0.05 across the three comparisons (R 4.3.3, multcomp).

The six secondary outcomes were tested with Welch’s t-tests, and p-values were adjusted with Holm’s method (p.adjust, R 4.3.3) to control the familywise error rate at 0.05. Adjusted p-values are reported.

Differential expression p-values for 12,000 genes were adjusted with the Benjamini–Hochberg procedure; genes with an adjusted p-value below 0.05 (false discovery rate 5%) are reported as candidates.

References

• Bender R, Lange S (2001). Adjusting for multiple testing—when and how? Journal of Clinical Epidemiology 54(4):343–349. doi:10.1016/S0895-4356(00)00314-0
• Benjamini Y, Hochberg Y (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. Journal of the Royal Statistical Society: Series B 57(1):289–300. doi:10.1111/j.2517-6161.1995.tb02031.x
• Benjamini Y, Yekutieli D (2001). The control of the false discovery rate in multiple testing under dependency. Annals of Statistics 29(4):1165–1188. doi:10.1214/aos/1013699998
• Dunnett CW (1955). A multiple comparison procedure for comparing several treatments with a control. Journal of the American Statistical Association 50(272):1096–1121. doi:10.1080/01621459.1955.10501294
• Games PA, Howell JF (1976). Pairwise multiple comparison procedures with unequal n’s and/or variances: a Monte Carlo study. Journal of Educational Statistics 1(2):113–125. doi:10.3102/10769986001002113
• Hochberg Y (1988). A sharper Bonferroni procedure for multiple tests of significance. Biometrika 75(4):800–802. doi:10.1093/biomet/75.4.800
• Holm S (1979). A simple sequentially rejective multiple test procedure. Scandinavian Journal of Statistics 6(2):65–70. JSTOR 4615733
• R Core Team. p.adjust: Adjust P-values for Multiple Comparisons. R documentation, package stats. stat.ethz.ch
• SciPy developers. scipy.stats.tukey_hsd, scipy.stats.dunnett, scipy.stats.false_discovery_control. SciPy 1.18 documentation. docs.scipy.org
• statsmodels developers. statsmodels.stats.multitest.multipletests. statsmodels 0.15 documentation. statsmodels.org