Joint mirror procedure: controlling false discovery rate for identifying simultaneous signals¶

作者: Linsui Deng, Kejun He, Xianyang Zhang
来源: Biometrics
主题: 数理统计 / 假设检验
相关性: 7/10
机构绿灯: Chinese University of Hong Kong（US News 前 50，免分进入精读）
链接: https://doi.org/10.1093/biomtc/ujae142

一、领域脉络与小综述¶

这个方向是什么¶

这个子方向解决的根本问题是：在同时检验多个“联合假设”（joint hypothesis）时，如何控制错误发现率（FDR）。这里的“联合假设”是指，一个特征（如一个基因、一个SNP）是否“有意义”，取决于它是否同时满足多个子假设（如中介分析中的暴露-中介效应与中介-结局效应同时非零，或可重复性分析中同一个信号在多个独立研究中同时显著）。传统的多重检验方法（如Benjamini-Hochberg过程）是为单个假设设计的，直接应用于联合假设会面临两个问题：一是检验统计量的分布复杂（因为联合假设的p值通常由多个子检验的p值组合而成），二是错误定义不清晰（一个“假发现”可能部分子假设为真、部分为假）。这个方向当前处于方法快速发展期，已有若干专门方法，但有限样本FDR控制的理论保证仍是一个核心挑战。

发展脉络（history）¶

奠基工作：从单假设到联合假设的初步尝试
- Benjamini & Hochberg (1995)：提出了经典的BH过程，为独立或正相依的单个假设提供了FDR控制。这是所有后续工作的基石。
- Benjamini & Yekutieli (2001)：将BH过程推广到任意相依结构（通过BY过程），代价是更保守的界。这为处理特征间相依的联合假设提供了理论起点。
- 早期中介分析/可重复性分析：早期工作（如Sobel检验、Baron & Kenny方法）通常分别检验两个子假设，然后通过某种方式（如乘积p值、最大值p值）组合，但缺乏严格的FDR控制。
主要进展：专门为联合假设设计的FDR控制方法
- Huang et al. (2015)：提出了“联合FDR”（joint FDR）的概念，并开发了基于BH过程的“联合BH”方法。该方法将每个联合假设的p值定义为两个子p值的最大值（或乘积），然后直接应用BH过程。作者在引言中评价：“However, the joint BH procedure may be conservative when the two p-values are correlated, and its FDR control is only guaranteed under independence or positive dependence (PRDS) conditions.” 这指出了第一个主要缺口：对相依性的处理有限。
- Barber & Candès (2015)：提出了“镜像统计量”（mirror statistics）方法，用于在有限样本下控制FDR，无需对p值分布做任何假设。该方法通过构造一个在零假设下对称的统计量，并利用其对称性来估计假发现数量。这是本文的核心技术灵感来源。作者在引言中明确说：“Our work is inspired by the mirror statistics approach of Barber and Candès (2015) for FDR control in the single-hypothesis setting.”
- Dai et al. (2022)：提出了“联合镜像”（joint mirror）过程的早期版本，但作者在引言中指出：“Dai et al. (2022) proposed a joint mirror procedure for mediation analysis, but their method does not provide finite-sample FDR control and relies on asymptotic approximations.” 这指出了第二个主要缺口：缺乏有限样本保证。
当前Frontier：有限样本FDR控制与更精细的错误度量
- 本文（Deng, He & Zhang, 2024）：站在Barber & Candès (2015)和Dai et al. (2022)的肩膀上，提出了Joint Mirror (JM) 过程。其核心贡献是：
  - 有限样本FDR控制：通过迭代收缩拒绝域和留一法（leave-one-out）技术，首次在有限样本下（而非渐近）证明了FDR控制。
  - 引入复合FDR（cFDR）：提出了一种更严格的错误度量，为每个假发现按其空分量数量赋予权重，这比传统的“全有或全无”的FDR更精细。
  - 融入偏序信息：算法设计上可以整合偏序信息（如某些联合假设的“重要性”排序），以提高统计功效。

子线索聚类¶

这些被引文献大致落在两条子线索上：

基于p值的方法：如Huang et al. (2015)的联合BH过程。这类方法简单易行，但FDR控制依赖于p值的分布假设（独立性或PRDS），且在相依性较强时可能保守或失效。
基于镜像统计量的方法：如Barber & Candès (2015)的单假设镜像过程，以及本文的JM过程。这类方法不依赖p值的精确分布，通过构造对称统计量来估计假发现，能在有限样本下控制FDR，但统计量的构造和算法设计更复杂。

这个方向在追问的核心问题¶

如何实现有限样本FDR控制？ 大多数现有方法（如联合BH、Dai et al. 2022）的FDR控制是渐近的，或依赖于特定分布假设。能否在有限样本下、对任意相依结构都提供严格的FDR控制？
如何定义“假发现”？ 对于一个联合假设，如果它被拒绝，但只有部分子假设为真，这算不算假发现？传统的“全有或全无”定义可能过于粗糙。复合FDR（cFDR） 提供了一种更精细的加权方案，但如何设计能控制cFDR的过程？
如何利用结构信息提高功效？ 在联合假设中，特征之间可能存在偏序关系（如中介分析中，暴露-中介效应和中介-结局效应在生物学上可能有先后顺序）。如何利用这种信息来更有效地分配“检验预算”，从而提高统计功效？
如何处理特征间的相依性？ 在基因组学等应用中，成千上万个特征（如基因）的检验统计量是高度相依的。如何设计一个对任意相依结构都稳健的FDR控制过程？

⚠️ 作者的Framing¶

作者把缺口Frame成什么？ 作者将缺口明确地frame为“现有联合镜像方法（Dai et al. 2022）缺乏有限样本FDR控制”，以及“现有方法对假发现的定义不够精细”。因此，本文的JM过程被呈现为“显然的下一步”：它既提供了有限样本保证，又引入了更精细的cFDR度量。
哪些竞争路线被淡化或回避了？ 作者淡化了基于p值的联合BH方法，主要批评其对相依性的限制。对于更复杂的、基于贝叶斯的方法（如局部FDR），作者完全没有提及。这可能是因为贝叶斯方法通常需要先验假设，与本文的“无分布假设”框架不兼容。
什么明显该被引/该存在、却没出现在intro里？ 作者没有引用任何关于高维中介分析或高维可重复性分析中FDR控制的最新工作。例如，在基因组学中，有大量工作研究如何在高维（p >> n）情况下进行中介分析，并控制FDR。这些工作通常使用惩罚回归或去偏Lasso，其FDR控制是渐近的。本文的JM过程能否推广到高维设定？这是一个明显的、未被讨论的缺口。

张力¶

未见明显对立引用。所有被引工作都指向一个共识：联合假设的FDR控制是一个重要且困难的问题，现有方法各有优缺点，需要更好的有限样本保证。

二、最核心、最简单的例子 / 数学问题¶

第一步：把符号、模型、可观测数据交代清楚¶

符号：
- 特征：i = 1, ..., m，表示m个待检验的特征（如m个基因）。
- 联合假设：对于每个特征i，我们检验一个联合零假设H_{0,i}，它由K个子零假设组成。在本文中，K=2（中介分析：H_{0,i}^{(1)}：暴露-中介效应为零；H_{0,i}^{(2)}：中介-结局效应为零。可重复性分析：H_{0,i}^{(1)}：研究1中效应为零；H_{0,i}^{(2)}：研究2中效应为零）。
- 检验统计量：对于每个子假设k，我们有一个检验统计量T_i^{(k)}（如z-score、t-statistic）。T_i = (T_i^{(1)}, T_i^{(2)})是特征i的二维检验统计量向量。
- 镜像统计量：M_i^{(k)}，是子假设k的镜像统计量。其构造方式为：M_i^{(k)} = sign(T_i^{(k)}) * (|T_i^{(k)}| - c)，其中c是一个常数（通常取0或某个分位数）。关键性质：在子零假设H_{0,i}^{(k)}下，M_i^{(k)}的分布关于0对称。
- 联合镜像统计量：W_i = min(M_i^{(1)}, M_i^{(2)})。这是本文的核心统计量。直觉：如果一个特征在两个子假设下都是信号，那么M_i^{(1)}和M_i^{(2)}都应该很大且为正，因此W_i也会很大且为正。如果至少一个子假设为零，那么对应的镜像统计量会以等概率为正或负，导致W_i可能为负。
- 拒绝域：R(t) = {i: W_i > t}，即所有联合镜像统计量大于阈值t的特征集合。
- 假发现数量：V(t) = |{i in R(t): H_{0,i} is true}|，即被错误拒绝的联合零假设的数量。
- 假发现比例：FDP(t) = V(t) / max(|R(t)|, 1)。
- 目标FDR水平：q（如0.05或0.1）。
- 复合假发现数量（cV）：cV(t) = sum_{i in R(t)} w_i，其中w_i是特征i的权重，等于其被拒绝的联合假设中，真正为零的子假设的数量。例如，如果H_{0,i}^{(1)}为真但H_{0,i}^{(2)}为假，且特征i被拒绝，则w_i = 1。如果两个子假设都为真，则w_i = 2。
- 复合FDP（cFDP）：cFDP(t) = cV(t) / max(|R(t)|, 1)。
- 目标cFDR水平：q（与FDR相同，但控制的是加权后的错误）。
模型：
- 数据生成机制：对于每个特征i，我们观测到T_i = (T_i^{(1)}, T_i^{(2)})。这些统计量可以来自任何模型（如线性回归、逻辑回归），只要我们能构造出在零假设下对称的镜像统计量。关键假设：在子零假设H_{0,i}^{(k)}下，M_i^{(k)}的分布关于0对称。这个假设通常可以通过构造M_i^{(k)} = sign(T_i^{(k)}) * (|T_i^{(k)}| - c)来满足，其中c是某个常数，且T_i^{(k)}在零假设下是对称分布（如均值为0的正态分布）。
- 要估的对象：我们想要找到一组被拒绝的特征R，使得E[cFDP] <= q（控制复合FDR），同时最大化被拒绝的特征数量（即统计功效）。
可观测数据：
- 研究者实际能观测到的是：对于每个特征i，两个检验统计量T_i^{(1)}和T_i^{(2)}。这些是样本数据。
- 想要但观测不到的是：每个子零假设H_{0,i}^{(k)}的真假状态。我们只能通过统计推断来估计。

第二步：讲最小内核¶

最简特例：假设K=2，且所有特征的检验统计量T_i^{(1)}和T_i^{(2)}是独立同分布的，且都服从标准正态分布N(0,1)。我们想要控制传统的FDR（即w_i = 1如果H_{0,i}为真，否则w_i = 0）。

在这个特例下，JM过程的核心思路可以简化为：

构造联合镜像统计量：对于每个特征i，计算W_i = min(|T_i^{(1)}|, |T_i^{(2)}|)。注意，这里我们忽略了符号，因为T_i^{(k)}在零假设下对称，所以|T_i^{(k)}|的分布是半正态分布。关键：如果H_{0,i}为真（两个子假设都为零），那么W_i的分布是已知的（两个独立半正态分布的最小值）。如果H_{0,i}为假（至少一个子假设非零），那么W_i倾向于取更大的值。
估计假发现数量：对于任意阈值t，我们想要估计V(t)。在单假设镜像过程中，Barber & Candès (2015) 使用V_hat(t) = |{i: W_i < -t}|来估计V(t)。这是因为，在零假设下，W_i的分布关于0对称，所以W_i < -t和W_i > t的概率相等。但在联合假设下，这个性质不再成立，因为W_i是两个镜像统计量的最小值，其分布不一定对称。
JM过程的迭代收缩：JM过程通过迭代来解决这个问题。
- 初始化：设置一个初始阈值t_0（如0），并计算初始拒绝集R(t_0)。
- 迭代：在第s步，我们有一个候选拒绝集R_s。我们使用一个保守的估计量cV_hat(R_s)来估计cV(R_s)。这个估计量是通过“留一法”或“数据分裂”技术构造的，它保证E[cV_hat(R_s)] >= E[cV(R_s)]。
- 停止规则：如果cV_hat(R_s) / max(|R_s|, 1) <= q，则停止，并输出R_s作为最终拒绝集。否则，我们收缩拒绝域：选择一个更严格的阈值t_{s+1} > t_s，得到更小的拒绝集R_{s+1}，然后重复。

为什么这个迭代能工作？ 关键在于，随着拒绝域收缩，假发现的数量cV会减少。而保守的估计量cV_hat会随着cV的减少而减少（或至少不增加）。因此，最终会达到一个点，使得cV_hat / |R| <= q。由于cV_hat是保守的，我们有E[cV / |R|] <= q，即FDR被控制。

这个特例的证明：在这个特例下，证明的核心是构造一个有效的cV_hat。一个简单的方法是使用数据分裂：将数据随机分成两部分，一部分用于构造镜像统计量，另一部分用于估计假发现。但数据分裂会损失功效。JM过程使用更精巧的留一法，在不损失太多功效的情况下构造保守估计量。

三、这篇论文做了什么¶

三句话¶

研究了什么问题：本文研究了在同时检验多个联合假设（如中介分析、可重复性分析）时，如何在有限样本下控制错误发现率（FDR）和一种更精细的加权错误度量——复合FDR（cFDR）。
核心工具/方法：提出了Joint Mirror (JM) 过程，这是一种迭代算法，通过逐步收缩拒绝域，并利用留一法（leave-one-out）技术构造一个保守的假发现数量估计量，直到该估计量低于目标FDR水平。
主要结论：JM过程在有限样本下控制cFDR（从而也控制FDR），且能融入偏序信息以提高统计功效。模拟和真实数据应用验证了其有效性。

关键设定与假设¶

设定：
- 有m个特征，每个特征对应一个联合零假设H_{0,i}，由K个子零假设H_{0,i}^{(k)}组成。本文主要关注K=2的情况，但方法可推广到K>2。
- 对于每个子假设k，我们有一个检验统计量T_i^{(k)}，并构造镜像统计量M_i^{(k)}。
- 联合镜像统计量定义为W_i = min_{k in [K]} M_i^{(k)}。
- 目标：控制复合FDR（cFDR）在水平q以下。
假设：
- 假设1（镜像对称性）：对于每个特征i和每个子假设k，在子零假设H_{0,i}^{(k)}下，M_i^{(k)}的分布关于0对称。这是镜像方法的核心假设，通常可以通过构造M_i^{(k)} = sign(T_i^{(k)}) * (|T_i^{(k)}| - c)来满足，其中c是某个常数，且T_i^{(k)}在零假设下是对称分布。
- 假设2（相依性结构）：JM过程对特征间的相依性没有限制。它适用于任意相依结构，包括强相依性。这是相比基于BH过程的方法（需要PRDS假设）的一个关键优势。
- 假设3（偏序信息）：JM过程可以融入偏序信息，但这不是必须的。偏序信息用于在迭代过程中更有效地选择要拒绝的特征。

主要结果¶

定理1（有限样本cFDR控制）：在假设1下，JM过程输出的拒绝集R满足E[cFDP(R)] <= q，其中cFDP(R) = cV(R) / max(|R|, 1)。这是本文的核心理论结果。它保证了在有限样本下，无论特征间的相依性如何，复合FDR都能被控制在目标水平以下。
- 直觉：证明的关键在于构造了一个保守的cFDP估计量cFDP_hat(R)，使得E[cFDP_hat(R)] >= E[cFDP(R)]。然后，JM过程的停止规则保证了cFDP_hat(R) <= q，从而E[cFDP(R)] <= q。
- 必要条件：假设1（镜像对称性）是必须的。如果这个假设不成立，那么镜像统计量的对称性被破坏，保守估计量可能不再保守。
- 解决的技术难点：如何构造一个在有限样本下、对任意相依结构都保守的cFDP估计量？作者使用了留一法：对于每个特征i，在计算其镜像统计量时，暂时将其从数据中移除，然后用剩余数据来估计其假发现的可能性。这种方法巧妙地避免了数据分裂带来的功效损失，同时保证了保守性。
定理2（融入偏序信息的cFDR控制）：如果特征之间存在一个已知的偏序关系（如某些特征“更重要”），JM过程可以优先拒绝这些特征，同时仍然控制cFDR。这个定理保证了在利用先验信息时，FDR控制不会被破坏。
算法细节：
- 初始化：计算所有特征的联合镜像统计量W_i，并按降序排列。
- 迭代：从最大的W_i开始，逐步将特征加入拒绝集。在每一步，计算保守的cFDP估计量。如果估计量超过q，则停止，并输出当前拒绝集（不包括最后加入的特征）。
- 留一法估计量：对于每个候选拒绝集R，其保守的cFDP估计量cFDP_hat(R)的计算方式如下：
  1. 对于每个特征i in R，计算其“留一”镜像统计量W_i^{(-i)}，即用除i外的所有数据重新计算W_i。
  2. 对于每个特征i in R，估计其假发现的可能性p_i。p_i的构造利用了镜像统计量的对称性：如果H_{0,i}为真，那么W_i和W_i^{(-i)}应该具有相似的分布，且W_i为正的概率与W_i^{(-i)}为负的概率相关。
  3. 保守估计量cV_hat(R) = sum_{i in R} p_i，然后cFDP_hat(R) = cV_hat(R) / max(|R|, 1)。

证明路线与技术技巧¶

整体路线：
1. 构造保守估计量：证明的核心是构造一个cFDP_hat(R)，使得对于任何拒绝集R，都有E[cFDP_hat(R)] >= E[cFDP(R)]。作者使用留一法来构造这个估计量。
2. 建立停时性质：JM过程是一个停时规则：它在第一个满足cFDP_hat(R_s) <= q的步骤s停止。证明的关键是，对于这个停时tau，有E[cFDP(R_tau)] <= q。
3. 利用鞅或超鞅：证明通常利用cFDP_hat(R_s)和cFDP(R_s)之间的某种鞅或超鞅关系。通过证明cFDP_hat(R_s) - cFDP(R_s)是一个非负的鞅（或上鞅），可以推出停时下的期望不等式。
4. 处理相依性：留一法估计量的一个关键性质是，它不依赖于特征间的独立性。即使特征高度相依，留一法估计量仍然保持保守性。这是因为留一法只依赖于单个特征的条件分布，而不是联合分布。
关键跳跃点：
- 从单假设镜像到联合镜像：Barber & Candès (2015)的镜像统计量是标量，其对称性直接用于估计假发现。在联合假设下，W_i = min(M_i^{(1)}, M_i^{(2)})的分布不再对称。作者的关键想法是：不再直接利用W_i的对称性，而是利用留一法来构造一个条件期望，这个条件期望在零假设下是保守的。
- 留一法估计量的构造：如何定义p_i（特征i是假发现的概率）？作者使用了一个巧妙的技巧：p_i = P(W_i^{(-i)} < -W_i | data)。这个概率的直觉是：如果H_{0,i}为真，那么W_i和W_i^{(-i)}应该“差不多”，所以W_i^{(-i)}小于-W_i的概率应该接近1/2。如果H_{0,i}为假，那么W_i会很大，所以这个概率会很小。通过证明E[p_i] >= P(H_{0,i} is true | i in R)，作者证明了cV_hat的保守性。
技术技巧点名：
- 留一法（Leave-one-out）：用于构造保守的假发现数量估计量，避免数据分裂带来的功效损失。
- 条件期望与概率不等式：用于证明留一法估计量的保守性。
- 停时理论：用于证明迭代过程的FDR控制。
- 镜像统计量（Mirror Statistics）：核心工具，用于在无分布假设下构造对称统计量。

真实例子与应用¶

中介分析：
- 数据：使用一个公开的基因表达数据集，研究DNA甲基化（中介变量）在基因型（暴露）与基因表达（结局）之间的中介作用。有m个CpG位点（特征），每个位点对应一个联合假设（基因型→甲基化效应，甲基化→表达效应）。
- 方法应用：对每个CpG位点，计算两个效应的检验统计量（如z-score），构造联合镜像统计量，然后应用JM过程。
- 结果：JM过程在cFDR水平q=0.1下识别出R个显著的CpG位点。与联合BH过程相比，JM过程识别出更多的显著位点（更高的统计功效），同时cFDP的估计值低于q。
- 例子想说明什么：验证JM过程在实际中介分析中的有效性，展示其相比现有方法的功效优势。
可重复性分析：
- 数据：使用两个独立的GWAS（全基因组关联研究）数据集，研究同一个性状（如身高）的遗传关联。有m个SNP（特征），每个SNP对应一个联合假设（在第一个研究中显著，在第二个研究中显著）。
- 方法应用：对每个SNP，计算两个研究中的效应估计和标准误，构造镜像统计量，然后应用JM过程。
- 结果：JM过程识别出在两个研究中都显著相关的SNP集合。与联合BH过程相比，JM过程在控制cFDR的同时，识别出更多可重复的SNP。
- 例子想说明什么：展示JM过程在跨研究可重复性分析中的实用性，验证其在不同数据类型（GWAS vs. 基因表达）上的通用性。

🔎 结论是否比证明窄¶

证明的严格性：定理1的证明是严格的，它依赖于假设1（镜像对称性）和留一法估计量的构造。结论“有限样本cFDR控制”是精确的。
潜在的泛化问题：
- 高维情况：证明中假设m是固定的。当m远大于样本量n时（高维情况），检验统计量的构造（如通过Lasso估计的系数）可能不再满足假设1（镜像对称性），或者留一法估计量的计算变得不可行。作者在结论中提到了这一点：“Extending the JM procedure to high-dimensional settings is an important direction for future research.” 这是一个明确的、比证明更窄的结论。
- K>2的情况：虽然方法可以推广到K>2，但证明中主要关注K=2。对于K>2，联合镜像统计量的构造（如min）可能不是最优的，且留一法估计量的计算复杂度会急剧增加。作者没有给出K>2时的理论保证。
- 偏序信息的利用：定理2保证了融入偏序信息后cFDR控制仍然成立，但没有给出如何最优地利用偏序信息（如如何定义“重要性”权重）的指导。这是一个开放问题。

四、开放问题¶

高维推广：如何将JM过程推广到高维设定（p >> n）？在高维情况下，检验统计量（如通过惩罚回归或去偏Lasso得到的）的分布可能不再对称，且留一法估计量的计算成本过高。扎根点：作者在结论中明确提到“Extending the JM procedure to high-dimensional settings is an important direction for future research.”
K>2的最优联合镜像统计量：对于K>2，使用min作为联合镜像统计量是否最优？是否存在其他组合方式（如max、乘积、加权和）能提供更高的统计功效？扎根点：作者在引言中提到了“The JM procedure can be extended to handle more than two hypotheses”，但没有给出K>2时的理论分析和最优性讨论。
偏序信息的最优利用：如何系统地、数据自适应地利用偏序信息？例如，能否将偏序信息编码为一个先验分布，然后通过贝叶斯方法整合到JM过程中？扎根点：作者在算法部分提到了“can incorporate partial ordering information”，但没有给出如何选择“重要性”权重的通用准则。
与计算-统计权衡的联系：JM过程的计算复杂度是多少？对于大规模m（如m=10^6），留一法估计量的计算是否可行？是否存在计算上更高效的近似方法，同时保持FDR控制？扎根点：这是一个隐含的开放问题，源于方法本身的计算需求。对于关注计算-统计权衡的研究者，这是一个自然的切入点。

Maintained by 陈星宇 · Homepage · Source on GitHub