CORA: a COrrelation-Redundancy-Aware FDR adjustment with genomic applications

Technical Summary: CORA – A Correlation-Redundancy-Aware FDR Adjustment

Problem Statement
In high-throughput genomic studies, such as RNA-seq and microarray analyses, researchers simultaneously test thousands of genes for differential expression. The standard approach for controlling the False Discovery Rate (FDR) is the Benjamini–Hochberg (BH) procedure. However, the BH procedure treats all hypotheses as independent or equally informative, ignoring the inherent correlation structure among genes. In biological reality, genes within the same pathway often exhibit strong co-expression. Consequently, BH may report clusters of highly correlated genes as separate, independent discoveries. This inflates the list of Differentially Expressed Genes (DEGs) with redundant information, potentially overstating the effective dimensionality of the biological signal without adding proportional biological insight. Existing methods either ignore correlation (BH, adaptive BH), are overly conservative under arbitrary dependence (Benjamini–Yekutieli, BY), or rely on complex assumptions (multivariate Gaussian models, latent factors) that are difficult to apply to heterogeneous genomic test statistics.

Methodology
The authors propose CORA (COrrelation-Redundancy-Aware), a modification of the BH procedure that explicitly incorporates pairwise gene correlations into the rejection threshold.

• Core Mechanism: For a set of $m$ ordered p-values $p_{(1)} \leq \dots \leq p_{(m)}$, CORA calculates a correction coefficient $c_i$ for the gene at rank $i$. This coefficient represents the accumulated absolute correlation between the current gene and all previously ranked (more significant) genes:
  $$c_i = \sum_{k < i} |r_{\sigma(k), \sigma(i)}|$$
  where $r$ denotes the Pearson correlation coefficient between gene expression profiles.
• Adjusted Threshold: The rejection threshold for the $i$-th gene is adjusted from the standard BH value ($i\alpha/m$) to:
  $$\alpha_i^{CORA} = \frac{\alpha(i - c_i)}{m}$$
  (with the exception of $i=1$, where $c_1=0$).
• Interpretation: If a gene is highly correlated with genes already deemed significant (high $c_i$), its effective threshold becomes stricter ($i - c_i < i$), penalizing redundant discoveries. If a gene is uncorrelated with previous rejections ($c_i \approx 0$), the threshold remains identical to BH.
• Theoretical Properties:
  1. Reduction to BH: Under population-level independence, sample correlations converge to zero as sample size increases, causing CORA to converge to the standard BH procedure (Theorem 1).
  2. Conservativeness: The rejection set of CORA is always a subset of the BH rejection set ($R_{CORA} \subseteq R_{BH}$). CORA never rejects a hypothesis that BH does not reject (Theorem 2).
  3. FDR Control: While a rigorous theoretical proof of FDR control under Positive Regression Dependence on Subsets (PRDS) remains an open problem, extensive simulations confirm that $FDR_{CORA} \leq \alpha$ empirically across all tested scenarios.

Key Results
The paper evaluates CORA through simulations and benchmarks on eight real-world datasets (three microarray, three RNA-seq).

• Simulation Performance:
  • Under independence, CORA behaves identically to BH.
  • Under correlation structures (Block AR(1), Compound Symmetry, and biologically realistic Pathway-Clustered), CORA maintains FDR control at or below the nominal level $\alpha$.
  • In terms of statistical power, CORA sits between the conservative BY procedure and the standard BH procedure. It offers a modest power loss relative to BH (typically 5–15%) but substantially higher power than BY (3–5×).
• Real Data Benchmarks:
  • DEG Reduction: CORA systematically reduces the number of reported DEGs compared to BH.
    • Microarrays: 5–17% reduction (e.g., Ovarian cancer dataset reduced from 1,082 to 895 genes, a 17.3% reduction).
    • RNA-seq: 17–23% reduction. The effect is stronger in RNA-seq due to higher pairwise correlations and larger numbers of DEGs.
  • Nature of Removed Genes: The genes removed by CORA (accepted by BH but rejected by CORA) are concentrated near the significance boundary. They are not the most strongly significant genes but rather marginal discoveries that are highly correlated with stronger signals.
  • Downstream Impact: In classification benchmarks (using DLDA, Linear SVM, and Naive Bayes), the reduction in the DEG list did not alter misclassification rates. The top-N gene lists (e.g., top 50 or 100) showed 94–100% overlap between CORA and BH, indicating that the most significant biological signals are preserved.

Significance and Claims
The authors position CORA not as a method to outperform BH in raw statistical power or classification accuracy, but as a tool for parsimony and biological interpretability.

• Principle of Parsimony: The primary contribution is providing a mathematically grounded mechanism to penalize redundant discoveries. When a cluster of 50 correlated genes is detected, reporting all 50 as independent findings overstates the evidence. CORA operationalizes the intuition that a gene "explained" by its correlated peers carries less marginal information.
• Practical Utility: For experimentalists planning validation (e.g., qPCR) or biomarker panel selection, CORA offers a smaller, more focused candidate list by removing redundant genes near the significance threshold. This is particularly valuable in RNA-seq studies with thousands of correlated DEGs.
• Implementation: The method is implemented in a freely available R package (CORA), requiring only p-values and a correlation matrix, without the need for multivariate Gaussian assumptions, resampling, or complex tuning parameters.

In summary, CORA bridges the gap between formal FDR control and correlation-aware feature selection, offering a principled way to refine DEG lists in high-correlation genomic settings without sacrificing the control of false discoveries.