Benjamini-Hochberg Procedure

The BH procedure operates through an ordered, step-up algorithm:

• Ranking: All individual p-values from the univariate tests are sorted from smallest to largest: p(1) ≤ p(2) ≤ ... ≤ p(m)
• Threshold Calculation: For each ranked p-value p(k), a critical value is computed as (k/m) × α, where m is the total number of tests and α is the desired FDR level
• Cutoff Identification: The procedure finds the largest rank k where p(k) ≤ (k/m) × α
• Rejection: All hypotheses with p-values smaller than or equal to this identified p(k) are rejected This sequential nature ensures the procedure adapts to the observed distribution of p-values rather than applying a rigid universal cutoff.

The BH procedure controls the False Discovery Rate (FDR), not the more stringent Family-Wise Error Rate (FWER).

• FWER (e.g., Bonferroni correction): Controls the probability of making any Type I error. In high-dimensional biology with 20,000 genes, this is often too conservative, leading to many missed true biomarkers.
• FDR: Controls the expected proportion of false positives among all rejected hypotheses. If you reject 100 genes at an FDR of 0.05, you expect at most 5 to be false discoveries. This tolerance for a small fraction of errors makes the BH procedure far more powerful for exploratory biomarker discovery, where follow-up validation is planned.

The original 1995 proof by Benjamini and Hochberg required that the test statistics be independent. This assumption is often violated in genomics, where gene expression levels are correlated.

• Positive Regression Dependency on a Subset (PRDS): Benjamini and Yekutieli later proved that the BH procedure also controls the FDR under the more relaxed PRDS condition, which includes multivariate normal data with all non-negative correlations.
• General Dependency: For arbitrary dependency structures, the Benjamini-Yekutieli procedure offers a modified, more conservative version that guarantees FDR control by replacing α with α / Σ(1/i). Understanding these assumptions is critical for valid biomarker selection in correlated omics datasets.

A direct extension of the BH procedure is the q-value, introduced by John Storey. While the BH procedure gives a binary reject/fail-to-reject decision, the q-value provides an FDR analogue to the p-value for each individual feature.

• Definition: The q-value of a feature is the minimum FDR at which that feature would be declared significant.
• Interpretation: A gene with a q-value of 0.013 means that if you select all genes with a q-value ≤ 0.013, the expected FDR among that set is 1.3%.
• Estimation: q-values incorporate an estimate of π₀, the proportion of true null hypotheses, making them more powerful than raw BH-adjusted p-values when many features are non-null.

The BH procedure is a standard filter step in high-dimensional biomarker workflows:

• Univariate Pre-Filtering: Before applying multivariate models like LASSO or Random Forest, analysts often perform a differential expression analysis (e.g., limma or DESeq2) and retain only features with a BH-adjusted p-value < 0.05.
• Dimensionality Reduction: This reduces the feature space from tens of thousands to a few hundred or thousand candidates, making downstream model training computationally feasible and reducing noise.
• Multi-Omics Integration: The procedure is applied independently to each omics layer (transcriptomics, proteomics, metabolomics) before data fusion, ensuring only statistically reliable features from each modality are combined.

The BH procedure can be visualized on a p-value rank plot:

• The x-axis represents the rank (k) from 1 to m.
• The y-axis represents the p-value.
• A straight line is drawn from the origin with a slope of α/m.
• The procedure rejects all hypotheses whose p-values fall below this line, up to the last point that crosses it.
• This visualization immediately reveals the procedure's adaptive nature: in datasets with a strong signal (many small p-values), the line is crossed at a higher rank, yielding more discoveries. In null datasets, the line may not be crossed at all, correctly yielding zero rejections.

The expected proportion of Type I errors (false positives) among all rejected null hypotheses. Unlike the family-wise error rate, FDR scales gracefully with the number of tests. The Benjamini-Hochberg procedure guarantees that FDR ≤ α under independence or positive regression dependency. In a genomics study testing 20,000 genes, an FDR of 0.05 means that, on average, 5% of the genes called significant are expected to be false discoveries.

The minimum FDR at which a test may be called significant, serving as an FDR-adjusted analogue to the p-value. Introduced by John Storey, q-values are computed by estimating the proportion of true null hypotheses (π₀) from the empirical p-value distribution. A feature with a q-value of 0.013 can be rejected at any FDR threshold ≥ 0.013. This provides a more powerful and interpretable alternative to raw p-values in high-dimensional settings.

A conservative multiple testing correction that controls the family-wise error rate (FWER)—the probability of making at least one Type I error. It divides the significance threshold α by the number of tests m: α_adjusted = α / m. For 10,000 tests at α = 0.05, the threshold becomes 5 × 10⁻⁶. While it provides strict control, it severely reduces statistical power in high-dimensional biomarker studies, making it overly stringent for exploratory discovery.

An adaptation of the Benjamini-Hochberg procedure that controls the FDR under arbitrary dependency structures among test statistics. It modifies the BH critical value by dividing by a harmonic sum constant: c(m) = Σ(1/i) for i = 1 to m. This makes it more conservative but mathematically valid when tests exhibit complex correlations—common in genomic data where genes operate in co-expression networks. Use when independence cannot be assumed.

A Bayesian refinement defining the FDR conditional on at least one rejection occurring: pFDR = E[V/R | R > 0]. Storey's method estimates π₀—the proportion of true null hypotheses—from the p-value distribution's right tail, leading to the q-value framework. This approach is more powerful than the original BH procedure when a substantial fraction of tested features are expected to be non-null, a common scenario in biomarker discovery.

A step-down procedure that uniformly improves power over the simple Bonferroni correction while still controlling the FWER. It orders p-values from smallest to largest and compares each to a sequentially adjusted threshold: α / (m - i + 1). The first non-rejection halts the procedure. Less powerful than BH for FDR control but provides stronger guarantees, making it suitable for confirmatory clinical trials where any false positive is unacceptable.