Benjamini-Hochberg Procedure

The BH procedure operates by ranking raw p-values in ascending order and comparing each to a linearly adjusted threshold. For the i-th ranked p-value, the critical value is (i/m) * α, where m is the total number of tests and α is the desired FDR level. The algorithm finds the largest index k where p(k) ≤ (k/m) * α, then rejects all hypotheses with rank ≤ k. This step-up approach is more powerful than family-wise error rate methods like Bonferroni.

Unlike the Bonferroni correction, which controls the probability of making any false positive (Family-Wise Error Rate), the BH procedure controls the False Discovery Rate—the expected proportion of false positives among discoveries. This distinction is critical in genomics:

• FWER: P(≥1 false positive) ≤ α — extremely conservative for 20,000+ genes
• FDR: E[false positives / total rejections] ≤ α — tolerates a small fraction of errors This tolerance dramatically increases statistical power when screening thousands of biomarkers.

The original 1995 proof by Benjamini and Hochberg assumed test statistics are independent. Subsequent work by Benjamini and Yekutieli (2001) extended the proof to cases of positive regression dependency—a condition commonly met in genomics where gene expression measurements exhibit block-correlation structures. For arbitrary dependency structures, a modified BH procedure multiplies the threshold by a harmonic sum constant, though this conservative adjustment is rarely needed for typical microarray and RNA-seq data.

John Storey extended the BH framework by introducing the q-value, the minimum FDR at which a given test would be called significant. While the BH procedure outputs a binary reject/accept decision at a fixed α, the q-value provides a continuous measure of significance analogous to the p-value but calibrated for FDR. For a gene with q-value = 0.05, up to 5% of genes with equal or lower q-values are expected to be false positives. The qvalue R/Bioconductor package implements this via bootstrap estimation of π₀, the proportion of true null hypotheses.

The BH procedure is the default multiple testing correction in major differential expression packages:

• DESeq2: results(dds, pAdjustMethod='BH') applies BH correction to Wald test p-values
• edgeR: topTags(et, adjust.method='BH') in the exact test workflow
• limma: topTable(fit, adjust='BH') after empirical Bayes moderation
• Seurat: FindMarkers(test.use='wilcox', p_val_adj='BH') for single-cell analysis Each implementation sorts raw p-values and applies the step-up threshold, typically with α = 0.05 or 0.10.

When interpreting BH-adjusted results:

• Adjusted p-value of 0.03: If all genes with adjusted p ≤ 0.03 are declared significant, the expected FDR among them is 3%
• FDR = 10% is common in exploratory biomarker discovery to prioritize candidates for validation
• Avoid double-correction: Do not apply BH after Bonferroni—they address different error rates
• Pre-filtering matters: Removing low-count genes before BH correction increases power by reducing the multiplicity burden
• The procedure controls the expected FDR, not the realized FDR in any single experiment

The expected proportion of Type I errors among all rejected null hypotheses. Unlike the family-wise error rate, FDR scales gracefully with the number of tests. In a typical RNA-seq experiment with 20,000 genes, an FDR threshold of 0.05 means that approximately 5% of genes called significant are expected to be false positives. This is the error metric directly controlled by the Benjamini-Hochberg procedure.

The probability of making at least one Type I error across all hypothesis tests. Controlled by conservative methods like the Bonferroni correction, which divides the significance threshold α by the number of tests. For genomic studies with tens of thousands of genes, FWER control is often too stringent, resulting in many true positives being missed. The Benjamini-Hochberg FDR approach was developed specifically to address this limitation.

The minimum FDR at which a test is declared significant, representing the FDR analogue of the p-value. Developed by John Storey as a Bayesian extension of the Benjamini-Hochberg framework. A q-value of 0.03 for a gene means that if you call all genes with q ≤ 0.03 significant, the expected FDR is 3%. Unlike BH-adjusted p-values, q-values incorporate an estimate of π₀, the proportion of true null hypotheses.

Under the global null hypothesis, p-values follow a Uniform(0,1) distribution. In genomic experiments, the observed p-value histogram is a mixture: a uniform component from true nulls and a spike near zero from true alternatives. Inspecting this distribution is a critical diagnostic step before applying Benjamini-Hochberg. A severely non-uniform distribution may indicate failed assumptions, poor normalization, or unmodeled batch effects.

The simplest and most conservative multiple testing adjustment. It sets the significance threshold to α / m, where m is the number of tests. For 20,000 genes at α = 0.05, the threshold becomes 2.5 × 10⁻⁶. While it guarantees strong control of the FWER, it dramatically reduces statistical power in high-dimensional settings. The Benjamini-Hochberg procedure was designed as a more powerful alternative for discovery-oriented genomics.

The class of statistical adjustments required when performing thousands of simultaneous hypothesis tests on microarray or RNA-seq datasets. Without correction, testing 20,000 genes at α = 0.05 yields approximately 1,000 false positives by chance alone. Methods span from conservative FWER approaches to adaptive FDR procedures. The choice of method represents a fundamental trade-off between sensitivity and specificity in biomarker discovery.