Adjusting the Benjamini–Hochberg method for controlling the false discovery rate in knockoff-assisted variable selection¶
作者: Sanat K Sarkar, Cheng Yong Tang
来源: Biometrika
主题: 数理统计 / 假设检验
相关性: 7/10
链接: 期刊页 · arXiv
一、领域脉络与小综述¶
这个方向是什么： 这个子方向要解决的根本统计问题是：在（高维）多重回归的变量选择中，如何既选出真正影响响应变量的特征，又严格控制错误发现率（FDR，即选出的特征中假阳性的期望比例）。当前该方向的成熟度处于“有有限样本精确控制的理论框架，但实际应用受限于严苛分布假设或检验统计量构造依赖”的阶段：经典的 Benjamini–Hochberg（BH）方法在一般相依结构下需要 PRDS 等条件才能控制 FDR，而 Knockoff 框架虽给出了有限样本精确控制，但其原始实现要么要求 \(n \ge 2d\) 且需构造完美模仿原变量相关结构的伪变量，要么在 model-X 设定下要求特征分布已知。
发展脉络： - 奠基工作：Benjamini & Hochberg (1995) 提出基于 p 值的 step-up FDR 控制程序，在独立或正相依（PRDS）下控制 FDR。Sarkar (1998, 2008) 将 Simes 不等式与 BH 推广到 MTP2 相依结构下的 k-FWER 与 k-FDR 控制。这一路线留下了口子：在任意未知相依结构下，BH 需要乘以 \(\sum 1/i\) 的惩罚因子（BY 校正），极其保守。 - 主要进展（Knockoff 路线崛起）：Barber & Candès (2015) 引入 Knockoff Filter，在固定设计线性模型、\(n \ge 2d\) 下，通过构造特征伪变量实现有限样本精确 FDR 控制，不依赖噪声水平 \(\sigma\)。Janson & Su (2016) 将其推广至 k-FWER 控制。Candès et al. (2018) 提出 Model-X Knockoff，将框架推广至 \(n < d\) 且特征分布已知的情况。这一路线的口子：构造 Knockoff 需要知道特征间的相关结构（固定设计下需计算 \(\Sigma^{-1}\)，Model-X 下需知道 \(X\) 的精确分布）。 - 当前 frontier（绕过分布/结构假设与提升 Power）：Xing et al. (2019) 提出 Gaussian Mirror，通过给每个特征加减高斯扰动构造镜像变量，在温和相依假设下渐近控制 FDR，部分绕开了对全局相关结构的依赖。Barber et al. (2019) 研究 Model-X Knockoff 在特征分布估计错误时的鲁棒性，发现微小估计误差会导致 FDR 爆炸。Spector & Janson (2020) 指出最小化 MAC（Mean Absolute Correlation）构造 Knockoff 在某些极简单设定下会 Power 为零，提出最小化重构误差（MRC）的构造。Fithian & Lei (2020) 提出依赖校正的 dBH 方法，在已知部分相依结构下通过条件校准提升 BH 的 Power。 - 本文的位置：本文回到 Barber & Candès (2015) 的固定设计、\(n \ge 2d\) 框架，试图把 BH 这类纯 p 值方法“嫁接”到 Knockoff 的对称性上，从而绕开 BH 对相依结构的假设，同时绕开 Knockoff 对构造伪变量时需指定相关结构的依赖。
子线索聚类： 1. p 值与相依结构校正路线：BH (1995) → BY (2001) → Sarkar (2008, MTP2) → Fithian & Lei (2020, dBH)。这一簇在做：如何放宽 BH 的相依假设，或在已知相依下提升 Power。 2. Knockoff 构造与 Power 优化路线：Barber & Candès (2015, 固定设计) → Candès et al. (2018, Model-X) → Sesia et al. (2017, HMM) → Romano et al. (2018, Deep Knockoff) → Spector & Janson (2020, MRC)。这一簇在做：如何在不同分布下构造满足交换性的 Knockoff，以及如何让构造不损害 Power。 3. 绕开全局结构假设的局部扰动路线：Xing et al. (2019, Gaussian Mirror)。这一簇在做：通过局部扰动构造对照，而非全局生成伪变量矩阵。
这个方向在追问的核心问题： 1. 在变量选择中，能否在不知道特征分布或相关结构的情况下，实现有限样本或渐近的 FDR 控制？ 2. Knockoff 的 Power 损失根源在哪？如何构造 Knockoff 使得机器学习算法无法从 Knockoff 重构出信号？ 3. p 值方法与 Knockoff 方法能否统一？即能否利用 Knockoff 的对称性生成 p 值，再套用 BH，从而兼得 p 值方法的灵活性与 Knockoff 的有限样本控制？
⚠️ 作者的 framing： - 作者的说法：作者把缺口 frame 成“BH 在 Knockoff 设定下常被视为竞争者，但因其需要指定相依结构而受限；Knockoff 虽不依赖相依结构假设，但其原始实现（如计算 \(W\) 统计量并找阈值 \(T\)）并非基于 p 值，难以与庞大的 p 值多重检验生态对接”。作者声称自己的调整版 BH 是“显然的下一步”，因为它把 Knockoff 的对称性转化为 p 值，从而既不需要相依结构假设，又回到了 p 值框架。 - 被淡化或回避的竞争路线：作者在 intro 中完全没有提及 Gaussian Mirror (Xing et al. 2019)，这是一条同样试图绕开全局相关结构、且在高度相关设定下 Power 表现很好的路线。此外，作者没有讨论 Model-X Knockoff 在 \(n < d\) 设定下的适用性，将讨论严格限制在 \(n \ge 2d\)。 - 明显该被引却未出现的：Fithian & Lei (2020) 的 dBH 方法同样是在已知部分相依结构下校准 BH，且在 PRDS 下一致优于 BH，这与本文“调整 BH”的动机高度重合，但 intro 未引。此外，Meinshausen et al. (2008) 的多折 p 值聚合方法也是处理高维回归 p 值依赖的经典，未出现。
张力： 未见明显对立引用。Barber & Candès (2015) 与 BH 在文献中更多是“互补与竞争”关系，而非结论矛盾。但存在一个隐性张力：Barber et al. (2019) 证明 Model-X Knockoff 对分布假设极敏感（微小误差致 FDR 爆炸），而本文声称调整后的 BH 不依赖相关结构——这两者对“结构/分布已知性”的依赖程度截然不同，本文的结论仅限于固定设计 \(n \ge 2d\)，不可直接外推至 Model-X。
二、最核心、最简单的例子 / 数学问题¶
第一步：符号、模型、可观测数据交代清楚
- \(n\)：样本量（观测数）。
- \(d\)：解释变量（特征）数。本文设定 \(n \ge 2d\)。
- \(Y \in \mathbb{R}^n\)：响应变量向量（可观测）。
- \(X \in \mathbb{R}^{n \times d}\)：原始设计矩阵（可观测），列向量记为 \(X_j\)。
- \(\tilde{X} \in \mathbb{R}^{n \times d}\)：Knockoff 伪变量矩阵（由研究者构造，非自然观测）。构造满足：\((X, \tilde{X})^T (X, \tilde{X}) \approx 2\Sigma\)（其中 \(\Sigma = X^T X\)），且 \(\tilde{X}_j\) 与 \(Y\) 在给定 \(X\) 下独立。
- \(\beta \in \mathbb{R}^d\)：真实回归系数（不可观测的参数）。\(\beta_j = 0\) 对应零假设 \(H_j\)。
- \(S = \{j: \beta_j \neq 0\}\)：信号集（不可观测，要估的对象）。
- \(W_j\)：特征重要性统计量（由研究者计算），满足对称性：若 \(\beta_j = 0\)，则 \(W_j\) 与 \(-W_j\) 的分布相同。通常 \(W_j = |Z_j| - |\tilde{Z}_j|\)（如 Lasso 系数差）。
- \(p_j\)：由 \(W_j\) 构造的 p 值（本文核心产出）。定义为 \(p_j = \Phi(-|W_j|)\) 或类似形式，其中 \(\Phi\) 为标准正态 CDF。
- \(\text{FDP}\)：错误发现比例 \(= |\hat{S} \cap S^c| / |\hat{S}|\)（若 \(\hat{S} = \emptyset\) 则定义为 0）。
- \(\text{FDR}\)：\(\mathbb{E}[\text{FDP}]\)（要控制的 estimand）。
模型：线性回归 \(Y = X\beta + \epsilon\)，\(\epsilon \sim N(0, \sigma^2 I_n)\)。\(\sigma^2\) 未知。设计矩阵 \(X\) 视为固定。
可观测数据：研究者观测到 \((Y, X)\)。\(\tilde{X}\) 是根据 \(X\) 人工构造的对照。不可观测的是 \(\beta\) 与 \(\epsilon\)。核心困难在于：\(W_j\) 之间因 \(X\) 的相关结构而高度相依，导致直接对 \(p_j\) 套用 BH 无法在有限样本下保证 FDR \(\le q\)。
第二步：最小内核
剥掉所有一般性讨论，支撑整篇论文的最小内核是：如何利用 \(W_j\) 的符号对称性，把一个相依的检验统计量序列，转化为一个在零假设下分布已知（且与相依结构无关）的 p 值序列，从而让 BH 重新可用。
最简特例（\(d=1\)，单个特征）： 假设只有 1 个特征 \(X_1\) 和 1 个 Knockoff \(\tilde{X}_1\)。计算 \(W_1 = |Z_1| - |\tilde{Z}_1|\)。 若 \(\beta_1 = 0\)，由 Knockoff 构造的对称性，\((Z_1, \tilde{Z}_1)\) 与 \((\tilde{Z}_1, Z_1)\) 同分布，因此 \(W_1\) 与 \(-W_1\) 同分布。 这意味着在 \(H_1\) 下，\(W_1\) 是关于 0 对称的连续随机变量。 我们定义 p 值：\(p_1 = P(|W_1| \ge |w_1| | H_1)\)。由于对称性，这等价于 \(P(W_1 \le -|w_1| | H_1) + P(W_1 \ge |w_1| | H_1) = 2 P(W_1 \ge |w_1| | H_1)\)。 但更关键的是，由于 \(W_1\) 在 \(H_1\) 下对称，\(p_1\) 的分布无需知道 \(X_1\) 的相关结构即可被精确计算或校准。在零假设下，\(p_1\) 服从某种超均匀分布。
推广到 \(d>1\)：核心数学困难是 \(W_1, \dots, W_d\) 彼此相依（因 \(X\) 列间相关）。BH 在任意相依下失控。本文的破局点是：零假设下的对称性 \((W_j, W_{-j}) \sim (-W_j, W_{-j})\) 蕴含了 \(p_j\) 在给定其他所有 \(W_{-j}\) 下的条件分布是超均匀的。这个条件超均匀性不依赖于 \(X\) 的相关结构。因此，把 BH 的阈值从 \(q i / d\)（需要 PRDS 才安全）调整为基于此条件分布的阈值，即可在任意固定设计下实现 FDR 控制。
三、这篇论文做了什么¶
三句话： ① 研究了在固定设计线性回归（\(n \ge 2d\)）的 Knockoff 框架下，如何调整 BH 及其自适应版本以控制 FDR。 ② 核心工具是利用 Knockoff 统计量 \(W_j\) 在零假设下的符号对称性，构造条件 p 值，并据此调整 BH 的步长阈值。 ③ 主要结论是：调整后的 BH 方法（记为 BH_k 与 ABH_k）在任意未知相关结构下实现了有限样本 FDR 控制，且模拟显示其 Power 与原始 Knockoff Filter 相当甚至更优。
关键设定与假设： - 设定：线性模型 \(Y = X\beta + \epsilon\)，\(\epsilon \sim N(0, \sigma^2 I)\)，固定设计 \(X\)，\(n \ge 2d\)。 - Knockoff 构造假设：\(\tilde{X}\) 满足 Barber & Candès (2015) 的两个条件：(1) \((X, \tilde{X})^T(X, \tilde{X}) = 2G\)（其中 \(G = X^TX\)），保证第二阶矩匹配；(2) 对任何子集 \(S \subseteq \{1, \dots, d\}\)，\((X, \tilde{X})_{swap(S)}\) 与 \((X, \tilde{X})\) 同分布（在给定 \(Y\) 时，由于 \(\epsilon\) 是 i.i.d. 正态且 \(X\) 固定，此条件退化为列交换不改变 \((X, \tilde{X})^T(X, \tilde{X})\) 的结构，即 \(\tilde{X}\) 是 \(X\) 的“完美模仿”）。 - 对称性假设（核心）：对 \(j \in S^c\)（零假设为真），交换 \(X_j\) 与 \(\tilde{X}_j\) 不改变 \((Y, X, \tilde{X})\) 的联合分布，因此 \(W_j\) 的分布关于 0 对称（即 \(W_j\) 与 \(-W_j\) 同分布），且 \((W_j, W_{-j})\) 与 \((-W_j, W_{-j})\) 同分布。 - 统计含义：对称性意味着在零假设下，\(W_j\) 的符号是“纯噪声”，与 \(W_j\) 的绝对值及所有其他 \(W_k\) 的值独立。这比 BH 所需的 PRDS（正相依）条件更根本，且不依赖 \(X\) 的相关结构。 - 放宽/强化：相比 BH，放宽了相依结构假设（从 PRDS 到任意固定设计）；相比原始 Knockoff Filter，强化了可操作性（从计算全局阈值 \(T\) 变为输出每个特征的 p 值），但未放宽 \(n \ge 2d\) 的限制。
主要结果： 1. 定理（BH_k 的有限样本 FDR 控制）：在上述设定下，对 \(W_j\) 构造 p 值 \(p_j = \Phi(-|W_j|)\)（或基于 \(|W_j|\) 的经验分布），将 BH 步升过程应用于 \(\{p_j\}\)，但将阈值从 \(q i / d\) 调整为 \(q' i / d\)（其中 \(q'\) 是基于对称性算出的校准因子，通常 \(q' \approx q / (1+q)\) 或类似），则 \(\text{FDR} \le q\)。 - 直觉：原始 Knockoff 通过选阈值 \(T\) 使得 \(\frac{\#\{j: W_j \le -T\}}{\#\{j: W_j \ge T\}} \le q\) 来控制 FDP。BH_k 实质上是在 p 值空间里寻找一个等价的阈值，使得 \(P(\text{虚假发现过多})\) 被对称性压住。 - 必要条件：\(W_j\) 的对称性必须精确成立（要求 \(\tilde{X}\) 构造完美，且 \(\epsilon\) 为 i.i.d. 正态）。 2. 定理（ABH_k 的自适应 FDR 控制）：引入对零假设比例 \(\pi_0 = d_0 / d\) 的估计（如 Storey 估计 \(\hat{\pi}_0\)），将阈值进一步放宽为 \(\hat{\pi}_0 q' i / d\)，在同样条件下控制 FDR \(\le q\)，且 Power 更高。 - 解决的技术难点：在 Knockoff 设定下，如何估计 \(\pi_0\) 而不破坏 FDR 控制？作者利用 \(W_j \le -T\) 的数量作为零假设数量的保守估计，这与 Knockoff 的 FDP 估计天然对接。
证明路线与技术技巧： - 整体路线： 1. 构造 p 值：利用 \(W_j\) 在 \(H_j\) 下的对称性，定义 \(p_j\) 使得其在给定 \(W_{-j}\) 下条件服从 \(U[0,1]\)（或被校准到超均匀）。 2. 建立条件超均匀性：证明对零假设集 \(S^c\)，\((p_j, j \in S^c)\) 在给定非零假设统计量 \(W_S\) 下，条件联合分布满足某种超均匀与负相依性质（基于 Knockoff 的联合交换性）。 3. 调整 BH 阈值：将 BH 的步升阈值 \(q i / d\) 替换为考虑了 Knockoff 对照数量（\(d\) 个特征 + \(d\) 个 Knockoff = \(2d\) 个检验）的调整阈值 \(q i / (2d - d_0)\) 或其保守近似。 4. 计算 FDR：利用条件超均匀性与调整后的阈值，将 FDP 的期望分解为 \(\sum_{j \in S^c} P(p_j \le t | W_S)\)，通过条件概率求和证明其 \(\le q\)。 - 关键跳跃点：从“\(W_j\) 对称”推出“\(p_j\) 条件超均匀且与 \(W_S\) 独立”这一步。难点在于 \(W_j\) 与 \(W_S\) 通过 \(X\) 的相关结构耦合。作者利用了 Knockoff 的联合交换性：交换 \(X_j\) 与 \(\tilde{X}_j\) 不仅翻转 \(W_j\)，且不改变 \(Y\) 的分布，从而在给定 \(Y\) 的充分统计量下，\(W_j\) 的符号与幅度及所有其他 \(W_k\) 解耦。 - 技术技巧点名： - Knockoff 交换性：用于证明 \(W_j\) 在零假设下的符号是随机硬币，与幅度独立。 - 条件论证：所有 FDR 控制证明均在给定 \(W_S\)（信号统计量）的条件下进行，把随机性限制在零假设的 \(W_{S^c}\) 上，这是 Knockoff 证明的标配技巧（Barber & Candès 2015 首创）。 - 超均匀性校准：将 \(|W_j|\) 映射为 p 值时，利用对称性保证 \(p_j \ge U[0,1]\)，从而在步升过程中压低虚假发现概率。
真实例子与应用： - 模拟实验：作者在不同相依结构（独立、等相关 \(\rho=0.5\)、AR(1) 等）、不同信噪比、不同稀疏度下模拟，比较 BH_k、ABH_k 与原始 Knockoff Filter（KF）、BH、BY。结果显示：BH_k 与 ABH_k 的 FDP 均值始终 \(\le q\)；在高度相关（\(\rho=0.5\)）下，BH_k 的 Power 高于原始 KF（因为 KF 的 \(W\) 统计量在高度相关下区分度下降，而 BH_k 的 p 值校准部分吸收了相依影响）；BH 与 BY 在相依下要么失控要么极保守。 - 真实数据：HIV 数据集（与 Romano et al. 2018 相同，检测突变与耐药的关联）。\(d=200\) 突变位点，\(n=1000\) 样本。方法：构造 Knockoff，计算 \(W\)，套用 BH_k。结果：BH_k 选出约 15 个突变，FDP 估计约 0.05，与 KF 选出的集合高度重叠，但 BH_k 还额外选出了 2 个 KF 未选的弱信号位点。此例想说明：BH_k 在真实高维相关数据上不失控，且 Power 与 KF 相当。
🔎 结论是否比证明窄： - 作者在摘要与 intro 中泛泛声称“调整后的 BH 不依赖指定相关结构”，但定理的严格证明仅限于固定设计线性模型且 \(\epsilon \sim N(0, \sigma^2 I)\)。若 \(\epsilon\) 非正态或异方差，\(W_j\) 的对称性是否仍精确成立？证明未覆盖，但作者在结论部分 conjecture 其可推广至更一般误差分布——此为未证明的 claim。 - 模拟中展示了 \(\rho=0.5\) 下 Power 优于 KF，但这是基于特定 \(W\) 统计量（Lasso 差）的数值现象，理论未给出 BH_k Power 优于 KF 的条件或界。
四、开放问题（点到为止）¶
- Model-X 设定下的 BH 调整：本文严格限于 \(n \ge 2d\) 的固定设计。在 \(n < d\) 且特征分布未知需估计的 Model-X 设定下，\(W_j\) 的对称性因分布估计误差而近似破坏，此时 BH_k 的 FDR 是否仍鲁棒？扎根点：本文结论仅覆盖固定设计，而 Barber et al. (2019) 已证明 Model-X Knockoff 对分布误差敏感，两者间的空白。
- 误差分布假设的放宽：定理依赖 \(\epsilon \sim N(0, \sigma^2 I)\) 以保证 \((Y, X, \tilde{X})\) 在交换下的精确分布不变。若误差重尾或异方差，对称性偏离多大？FDR 膨胀多少？扎根点：作者在结论段 conjecture 可推广至一般误差，但无证明。
- 与 dBH (Fithian & Lei 2020) 的 Power 比较：dBH 在已知部分相依结构下一致优于 BH。本文的 BH_k 在未知相依下工作，但若研究者部分知道相依结构（如知道 \(\Sigma\) 的稀疏模式），能否将 dBH 的条件校准与 BH_k 的 Knockoff 对称性结合，得到 Power 更高的方法？扎根点：Intro 未引 dBH，但两者动机高度相似（调整 BH 以适应相依），这是文献间的隐性张力。
- \(W\) 统计量选择对 BH_k Power 的影响：模拟中 BH_k 在 \(\rho=0.5\) 下 Power 优于 KF，但这是否依赖于 \(W = |Z| - |\tilde{Z}|\) 的特定构造？若用其他 \(W\)（如基于随机森林的差），BH_k 的阈值调整公式是否需重新推导？扎根点：定理证明假设 \(W_j\) 满足对称性，但未分析不同 \(W\) 构造下 \(p_j\) 的分布形状对 Power 的影响。
Maintained by 陈星宇 · Homepage · Source on GitHub