Multiple Comparison Corrections

The Multiple Comparisons Problem

When you run multiple hypothesis tests, the probability of at least one false positive increases with the number of tests. This is the multiple comparisons problem, and it is one of the most pervasive issues in data analysis.

The Mathematics of Inflated Error Rates

If you run m independent tests, each at significance level alpha, the probability of at least one false positive (the family-wise error rate, FWER) is:

FWER = 1 - (1 - alpha)^m

For alpha = 0.05:

With 100 tests, you are almost guaranteed to get at least one false positive. This is not an edge case -- it is the norm in modern data analysis, where dashboards track dozens of metrics and experiments test multiple variants.

Where Multiple Comparisons Arise

1. Multiple metrics. An A/B test tracks conversion rate, revenue per user, session duration, and bounce rate. Each metric is a separate test.

2. Multiple variants. Testing 5 variants against a control means 5 comparisons.

3. Multiple segments. Analyzing the effect across device type, geography, and user tenure creates many subgroup tests.

4. Multiple time points. Checking results weekly over a 3-month experiment creates repeated tests.

5. Feature selection. Testing hundreds of features for association with an outcome (common in genomics, ML feature selection).

Why This Matters in Practice

A product team runs an A/B test on 20 metrics. One metric shows p = 0.03. The team declares victory and ships the change. But with 20 tests, the probability of at least one false positive is 64%. The "significant" result is likely noise.

This is not a hypothetical problem. Published research across medicine, psychology, and economics has been undermined by uncorrected multiple comparisons. The replication crisis in psychology is partly attributed to this issue.

Two Frameworks for Control

There are two philosophies for controlling multiple comparisons:

1. Family-wise error rate (FWER). Control the probability of any false positive across all tests. This is the strictest standard. Methods: Bonferroni, Holm, Tukey's HSD.

2. False discovery rate (FDR). Control the proportion of false positives among all significant results. This is less strict but more powerful. Methods: Benjamini-Hochberg, Benjamini-Yekutieli.

The choice depends on the cost of false positives versus the cost of missing real effects. FWER is appropriate when any false positive is costly (e.g., medical decisions). FDR is appropriate when you can tolerate some false positives in exchange for more discoveries (e.g., exploratory analysis, feature screening).

A Simulation

import numpy as np
from scipy import stats

np.random.seed(42)
n_simulations = 10000
n_tests = 20
alpha = 0.05

false_positives = 0
for _ in range(n_simulations):
    # All null hypotheses are true (no real effects)
    p_values = np.random.uniform(0, 1, n_tests)
    if (p_values < alpha).any():
        false_positives += 1

print(f"FWER with {n_tests} tests at alpha = {alpha}: {false_positives/n_simulations:.3f}")
# Output: FWER ~ 0.64 (matches theoretical 1 - 0.95^20)

This simulation confirms the theory: with 20 tests and no real effects, you get at least one false positive 64% of the time.

Bonferroni and Holm Corrections

The Bonferroni and Holm corrections control the family-wise error rate (FWER). They ensure that the probability of any false positive across all tests does not exceed alpha.

Bonferroni Correction

The Bonferroni correction is the simplest and most conservative method. Divide alpha by the number of tests:

alpha_adjusted = alpha / m

A result is significant if its p-value is less than alpha_adjusted.

import numpy as np

p_values = np.array([0.001, 0.012, 0.034, 0.048, 0.067, 0.089, 0.15, 0.23])
m = len(p_values)
alpha = 0.05

# Bonferroni
bonferroni_alpha = alpha / m
significant_bonf = p_values < bonferroni_alpha

print(f"Bonferroni threshold: {bonferroni_alpha:.4f}")
print(f"Significant: {significant_bonf}")

Advantages:

• Simple to understand and implement
• Controls FWER strongly (even for dependent tests)
• No assumptions about test independence

Disadvantages:

• Very conservative; low power
• With many tests, almost nothing is significant
• The correction assumes the worst case (all tests independent)

Holm-Bonferroni Correction

The Holm correction (Holm's step-down procedure) is uniformly more powerful than Bonferroni while still controlling FWER. It works by sorting p-values and applying progressively less stringent thresholds.

Procedure:

1. Sort p-values from smallest to largest: p1 <= p2 <= ... <= pm
2. Compare p1 to alpha/m
3. If p1 < alpha/m, reject H1 and compare p2 to alpha/(m-1)
4. Continue until you fail to reject. Stop there.

import numpy as np

def holm_correction(p_values, alpha=0.05):
    m = len(p_values)
    sorted_indices = np.argsort(p_values)
    sorted_p = p_values[sorted_indices]
    
    adjusted = np.zeros(m)
    for i, p in enumerate(sorted_p):
        adjusted[i] = min(p * (m - i), 1.0)
    
    # Enforce monotonicity
    adjusted = np.maximum.accumulate(adjusted)
    
    # Restore original order
    result = np.zeros(m)
    result[sorted_indices] = adjusted
    return result

p_values = np.array([0.001, 0.012, 0.034, 0.048, 0.067, 0.089, 0.15, 0.23])
adjusted = holm_correction(p_values, alpha=0.05)

for p, adj in zip(p_values, adjusted):
    print(f"p = {p:.3f}, adjusted = {adj:.4f}, significant = {adj < 0.05}")

Why Holm Is Better Than Bonferroni

Holm's procedure is more powerful because it uses the data. If the smallest p-value is very small, it gets the full Bonferroni threshold (alpha/m). But the second smallest gets a less stringent threshold (alpha/(m-1)), and so on.

The key insight: once you reject the smallest p-value, you have fewer remaining tests to worry about, so you can afford a slightly less strict threshold for the next one.

Comparing Bonferroni and Holm

from statsmodels.stats.multitest import multipletests

p_values = [0.001, 0.012, 0.034, 0.048, 0.067, 0.089, 0.15, 0.23]

# Bonferroni
reject_bonf, pvals_bonf, _, _ = multipletests(p_values, alpha=0.05, method='bonferroni')

# Holm
reject_holm, pvals_holm, _, _ = multipletests(p_values, alpha=0.05, method='holm')

print("Bonferroni:", pvals_bonf)
print("Holm:      ", pvals_holm)

Holm will always reject at least as many hypotheses as Bonferroni, and usually more. Both control FWER, but Holm does so with greater power.

When to Use FWER Methods

Use Bonferroni or Holm when:

• The number of tests is small (less than 20)
• Any false positive is costly (medical trials, safety tests)
• You need a simple, defensible method
• Tests may be dependent (Bonferroni handles this; Holm assumes some independence structure)

For most product experiments with a handful of metrics, Holm is the recommended default. It controls FWER while being less conservative than Bonferroni.

False Discovery Rate: Benjamini-Hochberg

The Benjamini-Hochberg (BH) procedure controls the false discovery rate (FDR) -- the expected proportion of false positives among all significant results. It is less strict than FWER methods but substantially more powerful.

FDR vs. FWER

FWER controls P(at least one false positive). This is the probability of making any error.

FDR controls E[false positives / total positives]. This is the proportion of significant results that are false.

If you run 100 tests and 10 are significant, FDR = 0.10 means you expect 1 of those 10 to be a false positive. FWER = 0.05 means there is only a 5% chance that any of the 10 is false.

FDR is more permissive. It allows some false positives in exchange for detecting more real effects.

The Benjamini-Hochberg Procedure

Procedure:

1. Sort p-values from smallest to largest: p1 <= p2 <= ... <= pm
2. Find the largest k such that p_k <= (k/m) * alpha
3. Reject all hypotheses with p-values <= p_k

import numpy as np

def benjamini_hochberg(p_values, alpha=0.05):
    m = len(p_values)
    sorted_indices = np.argsort(p_values)
    sorted_p = p_values[sorted_indices]
    
    # Find the largest k where p_k <= (k/m) * alpha
    thresholds = np.arange(1, m+1) / m * alpha
    below = sorted_p <= thresholds
    
    if below.any():
        max_k = np.max(np.where(below)[0]) + 1
        reject = np.zeros(m, dtype=bool)
        reject[sorted_indices[:max_k]] = True
    else:
        reject = np.zeros(m, dtype=bool)
    
    return reject

p_values = np.array([0.001, 0.012, 0.034, 0.048, 0.067, 0.089, 0.15, 0.23])
reject = benjamini_hochberg(p_values, alpha=0.05)

for p, r in zip(p_values, reject):
    print(f"p = {p:.3f}, significant = {r}")

Adjusted P-Values

The BH procedure can also produce adjusted p-values (q-values) that are directly comparable to alpha:

from statsmodels.stats.multitest import multipletests

p_values = [0.001, 0.012, 0.034, 0.048, 0.067, 0.089, 0.15, 0.23]
reject, pvals_adj, _, _ = multipletests(p_values, alpha=0.05, method='fdr_bh')

for p, adj, r in zip(p_values, pvals_adj, reject):
    print(f"p = {p:.3f}, q = {adj:.4f}, significant = {r}")

Benjamini-Yekutieli Correction

The standard BH procedure assumes tests are independent or positively dependent. When tests may be negatively dependent or have arbitrary dependence structures, use the Benjamini-Yekutieli (BY) correction:

reject, pvals_adj, _, _ = multipletests(p_values, alpha=0.05, method='fdr_by')

BY is more conservative than BH but controls FDR under arbitrary dependence. Use it when you cannot assume independence or positive dependence.

When to Use FDR Methods

Use Benjamini-Hochberg when:

• You are running many tests (20+)
• You can tolerate some false positives
• You are screening for candidates to validate later
• You want more power to detect real effects
• The analysis is exploratory rather than confirmatory

Common applications:

• Genomics: testing thousands of genes for association with a disease
• Feature selection: screening hundreds of features for predictive power
• A/B testing: evaluating many metrics or variants
• Dashboard monitoring: tracking dozens of metrics for anomalies

A Practical Comparison

from statsmodels.stats.multitest import multipletests
import numpy as np

np.random.seed(42)
# Simulate 50 tests: 10 true effects, 40 null
p_values = np.concatenate([
    np.random.uniform(0, 0.01, 10),   # true effects (small p-values)
    np.random.uniform(0, 1, 40)       # null hypotheses
])

methods = ['bonferroni', 'holm', 'fdr_bh', 'fdr_by']
for method in methods:
    reject, _, _, _ = multipletests(p_values, alpha=0.05, method=method)
    n_rejected = reject.sum()
    n_false = reject[10:].sum()  # false positives (from the 40 nulls)
    print(f"{method:12s}: {n_rejected} rejected, {n_false} false positives")

Typical output:

bonferroni   : 8 rejected, 0 false positives
holm         : 9 rejected, 0 false positives
fdr_bh       : 12 rejected, 1 false positive
fdr_by       : 10 rejected, 0 false positives

FDR methods detect more true effects but allow a few false positives. FWER methods are stricter but miss some real effects.

Choosing and Applying Corrections

Selecting the right correction method requires understanding your goals, the number of tests, and the cost of each error type.

Decision Framework

How many tests?
|-- Few (< 10)
|   +-- Any method works; Holm is a safe default
|-- Moderate (10-50)
|   +-- Confirmatory analysis? --> Holm (FWER)
|   +-- Exploratory analysis? --> Benjamini-Hochberg (FDR)
+-- Many (50+)
    +-- Need strict control? --> Bonferroni or BY
    +-- Screening for candidates? --> Benjamini-Hochberg

What is the cost of a false positive?
|-- Very high (medical, safety) --> FWER (Bonferroni/Holm)
|-- Moderate (product decisions) --> Holm or BH
+-- Low (exploratory, follow-up planned) --> BH

Are tests independent?
|-- Yes --> BH is valid
|-- Positively dependent --> BH is valid
+-- Arbitrarily dependent --> Use BY (more conservative)

Practical Guidelines for Product Teams

1. Pre-specify your correction method. Decide before seeing the data. Choosing after seeing which method gives "significant" results is p-hacking.

2. For A/B tests with a handful of metrics: Use Holm. It controls FWER and is less conservative than Bonferroni.

3. For experiments with many variants or segments: Use Benjamini-Hochberg. The number of comparisons grows quickly, and FWER methods become too conservative.

4. For exploratory analysis: Use BH and label results as "hypothesis-generating." Follow up with confirmatory tests.

5. For monitoring dashboards: Use BH or no correction, but require sustained significance before acting. A single anomaly is likely noise.

Implementing Corrections in Python

from statsmodels.stats.multitest import multipletests
import pandas as pd

# Example: A/B test with 15 metrics
results = pd.DataFrame({
    'metric': ['conversion', 'revenue', 'session_time', 'bounce_rate',
               'pages_per_session', 'add_to_cart', 'checkout_start', 'signup',
               'retention_d1', 'retention_d7', 'nps', 'clicks', 'scroll_depth',
               'video_plays', 'shares'],
    'p_value': [0.001, 0.008, 0.012, 0.034, 0.041, 0.052, 0.067,
                0.089, 0.12, 0.15, 0.23, 0.31, 0.45, 0.56, 0.78]
})

# Apply Holm correction
reject_holm, pvals_holm, _, _ = multipletests(
    results['p_value'], alpha=0.05, method='holm'
)
results['p_holm'] = pvals_holm
results['sig_holm'] = reject_holm

# Apply Benjamini-Hochberg
reject_bh, pvals_bh, _, _ = multipletests(
    results['p_value'], alpha=0.05, method='fdr_bh'
)
results['p_bh'] = pvals_bh
results['sig_bh'] = reject_bh

print(results[['metric', 'p_value', 'p_holm', 'sig_holm', 'p_bh', 'sig_bh']])

Common Mistakes

1. Not correcting at all. Running 20 tests and using raw p-values. This is the most common mistake and leads to a 64% false positive rate.

2. Correcting only "significant" results. Apply the correction to all tests, not just the ones that look significant. Selective correction invalidates the error control guarantee.

3. Mixing correction methods. Using Bonferroni for some tests and BH for others without a principled reason. Pick one method and apply it consistently.

4. Ignoring dependence. Tests on the same data (e.g., conversion rate and revenue per user from the same experiment) are not independent. Bonferroni handles this; BH assumes positive dependence.

5. Over-relying on corrections. Corrections control error rates but do not fix fundamental design problems. A poorly designed experiment with corrected p-values is still poorly designed.

Reporting Corrected Results

Always report both raw and corrected p-values:

Metric          | Raw p | Holm p | BH p | Significant?
----------------|-------|--------|------|-------------
Conversion      | 0.001 | 0.015  | 0.003| Yes (both)
Revenue         | 0.008 | 0.112  | 0.020| Yes (BH only)
Session time    | 0.012 | 0.168  | 0.024| Yes (BH only)
Bounce rate     | 0.034 | 0.408  | 0.057| No

This transparency lets stakeholders evaluate the strength of evidence. A result significant under both methods is stronger than one significant only under FDR.