跳转至

Isotonic subgroup selection

作者: Manuel M Müller, Henry W J Reeve, Timothy I Cannings, Richard J Samworth
来源: Journal of the Royal Statistical Society Series B
主题: 数理统计 / 假设检验
相关性: 7/10
链接: 期刊页 · arXiv

一、领域脉络与小综述

这个方向是什么
子群选择（subgroup selection）问题：给定协变量‑响应样本 \((X_i, Y_i) \in \mathcal{X} \times \mathbb{R}\)，目标是识别回归函数 \(f(x) = \mathbb{E}(Y \mid X = x)\) 超过预设阈值 \(\tau\) 的协变量子集 \(R = \{x \in \mathcal{X} : f(x) > \tau\}\)。这是后选择推断的一个特例，在临床试验（例如找出受益于新疗法的患者亚群）和实证科学中广泛出现。核心统计挑战有二：一是必须对最终报告的子集提供可验证的 Type I error 控制（避免报告实际上不存在的亚群），二是希望检验的 power 尽可能高（即真子集被报告出来的比例）。该方向当前成熟度不高：已有大量应用层面的警示，但在非参数框架下同时保证有限样本误差控制和最优 minimax power 的方法仍属空白。

发展脉络（history）
按时间与逻辑顺序梳理：

  • 临床/应用层面的警示（约 1997 – 2020）: Senn & Harrell (1997); Feinstein (1998); Rothwell (2005); Wang et al. (2007); Kaufman & MacLehose (2013); Altman (2015); Gabler et al. (2016); Lipkovich et al. (2017); Watson & Holmes (2020)。这些工作反复论证不加调整的亚组检索会严重夸大假阳性率，但未给出可控的方法论。

  • 形状约束回归的估计理论（1955 起，近十年尤其活跃）：

  • 奠基：Ayer et al. (1955), Brunk (1955), van Eeden (1956) 提出并发展了 isotonic 回归。

  • 风险界与 oracle 不等式：Meyer & Woodroofe (2000), Zhang (2002), Chatterjee (2014), Chatterjee et al. (2015), Bellec (2018), Han et al. (2019), Deng & Zhang (2020), Fokianos et al. (2020), Pananjady & Samworth (2022)。这些工作确定了多元 isotonic 回归的 minimax 速率为 \(n^{-1/d}\)（在固定格点设计下）以及分段常数信号的适应速率。但它们关注的是回归函数的整体 MSE，而非其中哪部分大于阈值这个子集推断问题。

  • 多重检验与逻辑约束（1969 – 2019）：

  • Gabriel (1969) 提出“coherent”多重检验原则：拒绝集应构成一个上集（upper set）。

  • Meijer & Goeman (2015) 和 Ramdas et al. (2019) 将假设组织成有向无环图（DAG），并按顺序检验以控制 FWER / FDR。Ramdas et al. (2019) 进一步处理了在线设定。但这些方法假设 p‑value 已经事先计算好，没有讨论如何从数据构造这些 p‑value 并同时利用非参数结构。

  • Anytime‑valid / martingale 检验（1939 – 2021）：

  • Ville (1939), Wald (1947), Robbins (1970) 开创了基于 martingale 的序贯检验。

  • Duan et al. (2020) 基于 “masking” p‑value 提出交互式全域零检验；Howard et al. (2021) 给出了非参数下的 time‑uniform 置信序列。Wasserman et al. (2020) 提出通用推断方法（split LRT）。这些工作提供了对单个复合假设的时时有效 p‑value，但并未将其放入结构化多重检验框架下解决子群选择。

  • 异质性治疗效应的子集识别（近五年）：Ballarini et al. (2018), Watson & Holmes (2020) 等尝试用预测个体治疗效应（PITE）或 ML 检测交互，但多依赖离散化或缺乏有限样本误差控制。

本文位置：首次将上述四条线索缝合——在多元 isotonic 回归下，用单调性自然诱导 DAG 结构假设，用 martingale confidence sequence 构造每个点的 anytime‑valid p‑value，再用 Gabriel 的 coherent 原则设计多重检验程序 RISS，同时实现了非渐近均匀 Type I error 控制以及匹配 minimax 下界（至多对数因子）的 power。

子线索聚类：

线索 代表工作 核心关注
形状约束回归的理论估计 Han et al. (2017); Deng & Zhang (2020); Bellec (2018); Chatterjee (2014) 对单调/凸函数的 MSE / oracle 界
多重检验与 DAG 逻辑约束 Gabriel (1969); Meijer & Goeman (2015); Ramdas et al. (2019); Goeman & Solari (2010) 给定 p‑value 后如何控制 FWER/FDR
Anytime‑valid martingale 检验 Robbins (1970); Duan et al. (2020); Howard et al. (2021); Wasserman et al. (2020) 构造有限样本有效的、可随时停止的 p‑value
子群选择的临床/应用层面 Kaufman & MacLehose (2013); Gabler et al. (2016); Ballarini et al. (2018) 强调问题严重但缺少理论工具

核心问题：① 如何在连续协变量空间上对“\(f(x) > \tau\)”这一族假设提供 uniform Type I error 控制？② 在非参数形状约束下，能达到的最优 power（被拒绝区域的 Lebesgue 测度）是多少？③ 如何将 anytime‑valid 理念与结构化多重检验结合实现可在线计算的算法？④ 能否推广到分类响应、分位数回归、异质性治疗效应？

⚠️ 作者的 framing
作者在引言中把缺口框定为：“现有子群选择方法要么牺牲 Type I error（例如不加调整的逐个检验），要么只针对有限个预先指定的子组；在 isotonic 假设下，自然的偏序结构使我们能同时解决 error 控制和 power 最优。”

  • 被本文淡化/回避的竞争路线：
  • 基于树/分治的方法（SIDES, GUIDE, PRIM）——理由可能是这些方法需要离散化且缺乏整体误差控制。
  • 直接通过函数估计 + bootstrap 校准阈值的方法——本文未讨论。
  • 一个明显的缺失引用：Yang & Barber (2019) 关于 isotonic 回归一致置信带的工作；以及 Deng, Han, Zhang (2020) 关于多元 isotonic 点态置信区间的 pivot 极限理论。后者的推断思路（利用块最大‑最小估计量）与本文构造 p‑value 的视角可能互补，值得研究者去查是否可嫁接。

张力：未发现被引文献之间存在显式矛盾。不同线索间多属互补而非竞争。

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

第一步：符号、模型与可观测数据

  • 协变量：\(X \in [0,1]^d\)，设计为固定 \(r\)-regular 格点（等距网格点），记 \(n\) 为总样本量，各坐标方向格点数 \(m \approx n^{1/d}\)。
  • 响应：\(Y = f(X) + \varepsilon\)，\(\varepsilon\) 为 sub‑Gaussian 噪声，方差参数 \(\sigma\) 已知（或给定上界）。
  • 未知回归函数：\(f\) 在 \([0,1]^d\) 上对每个坐标单调不减（non‑decreasing in each coordinate）。
  • 阈值：\(\tau\) 是预先指定的常数（例如主观设定的“有效”水平）。
  • 目标子集：\(R = \{x \in [0,1]^d : f(x) > \tau\}\)。由单调性，\(R\) 必为 上集（upper set）：若 \(x \in R\)，则对所有坐标分量都大于或等于 \(x\) 的点也属于 \(R\)。
  • 可观测数据：设计点 \(x_i\)（格点位置）及其对应的 \(y_i\)。不可观测的是 \(f\) 在非设计点上的值以及噪声的每一实现。本文假设设计点位置已知，且可重复采样（条件于设计）。
  • 参数：样本量 \(n\)，协变量维数 \(d\)，显著性水平 \(\alpha \in (0,1)\)。

第二步：最小内核——一元分段常数特例

取 \(d=1\)，设计点为 \(x_1 < x_2 < \cdots < x_n\) 均匀分布在 \([0,1]\) 上（格点间距 \(1/n\)）。真回归函数设为

\[f(x) = \begin{cases} 0, & x \leq \theta, \\ 1, & x > \theta, \end{cases}\]
其中 \(\theta \in (0,1)\) 未知。阈值 \(\tau = 0.5\)，则真上集 \(R = (\theta, 1]\)。
噪声 \(\varepsilon_i \sim \text{subG}(\sigma^2)\) 独立。

任务：报告一个上集 \(\widehat{R} = (\widehat{\theta}, 1]\)，使得
- Type I error 控制：\(\Pr(\widehat{R} \subseteq R) \geq 1 - \alpha\)，即拒绝集不含“假阳性”点（没有报告 \(x \notin R\) 的点）。
- Power 最大化：在满足该误差约束下，最小化期望缺失测度 \(\mathbb{E}[|R \setminus \widehat{R}|] = \mathbb{E}[|\min(\theta, \widehat{\theta}) - \theta|]\)（当 \(\widehat{\theta} > \theta\) 时缺失为 0）。

最小内核做法：

  1. 对每个点 \(x_i\) 构造 p‑value：固定 \(i\)，检验 \(H_i : f(x_i) \leq 0.5\)。由单调性，\(H_i\) 蕴含所有 \(j \leq i\) 的 \(f(x_j) \leq 0.5\)。构造检验统计量 \(S_i = \sum_{j=i}^{n} (Y_j - 0.5)\)。在 \(H_i\) 下，\(S_i\) 的期望 \(\leq 0\)。利用 Howard et al. (2021) 的 sub‑Gaussian 时间均匀置信序列得到 p‑value

    \[\widetilde{p}_i = \exp\!\left(-\frac{(S_i)^2}{2 \sum_{j=i}^n (Y_j \vee \tau) + \cdots}\right)\]
    （细节略，核心是这个 p‑value 对任何停止时间皆有效：若 \(H_i\) 真，\(\Pr(\widetilde{p}_i \leq \alpha) \leq \alpha\)。）

  2. 多重检验：由于单调性，\(R\) 是上集，因此若拒绝 \(x_i\)（报告 \(f(x_i) > 0.5\)），则对所有 \(j > i\) 也必须拒绝。反之，若不拒绝 \(x_j\)（认为 \(f(x_j) \leq 0.5\)），则对所有 \(i < j\) 也不能拒绝。这就是 Gabriel (1969) 的 coherent 拒绝原则。

  3. RISS 在该特例中的实现：将所有 p‑value \(\{ \widetilde{p}_1, \ldots, \widetilde{p}_n \}\) 从大到小排序（实际上根据 DAG 的偏序），然后找到一个最小的阈值 \(c\) 使得所有满足 \(\widetilde{p}_i \leq c\) 的点组成的集合恰为某个上集。等价地，可定义 \(\widehat{\theta} = \min\{ x_i : \widetilde{p}_i \leq \alpha \text{ 且对所有 } j < i, \widetilde{p}_j \leq \alpha \}\)，再拒绝 \((\widehat{\theta}, 1]\)。

要证明的事：
- 对任意 \(\alpha\)，\(\Pr(\widehat{\theta} > \theta) \leq \alpha\)（即拒绝集不含假阳性）。
- 若 \(1-\theta\) 远离零（真集较大），则 \(\mathbb{E}[|\widehat{\theta} - \theta|] \leq C n^{-1} (\log n)^{\gamma}\)。其中 \(n^{-1}\) 即为 minimax 最优速率（因为一维下通过估计一个跳点，参数速率可达 \(O(n^{-1})\)）。

为什么这是一个“最小内核”：
- 一般 \(d \geq 2\) 的情形不过是把一维区间推广为多维上集，p‑value 构造变成对每个点取“所有坐标分量≥它的”子格点的均值，minimax 速率变为 \(n^{-1/d}\)。
- DAG 从一条链变成更复杂的偏序格，但 coherent 拒绝算法（Meijer‑Goeman 的 DAG 树）本质上是将一维的阈值规则推广到多维（通过“局部最小值”逻辑）。

因此，读者抓住一元例子的思想，就掌握了整篇论文的核心数学操作。

三、这篇论文做了什么

三句话：
1. 研究在多元 isotonic 回归中，给定样本 \((X_i,Y_i)\) 和阈值 \(\tau\)，如何识别回归函数超过 \(\tau\) 的协变量子集 \(R\)，并同时控制 Type I error 和达到最优 power。
2. 提出的 RISS (Reject Isotonic Subgroup Selection) 程序将每个设计点处的假设 \(H_x: f(x) \leq \tau\) 构建为 DAG，对每个点用 martingale 构造 anytime‑valid p‑value，再按 Gabriel coherent 原则进行多重检验。
3. 主要理论贡献：(i) 非渐近均匀 Type I error 控制——\(\Pr(\widehat{R} \subseteq R) \geq 1-\alpha\)；(ii) power 边界 \(\mathbb{E}[|R \setminus \widehat{R}|] \leq C (\log n)^{\gamma} n^{-1/d}\)；(iii) minimax 下界 \(c n^{-1/d}\)，从而证明 RISS 在该意义下最优（至多对数因子）。

关键设定与假设（在第二节记号基础上补全）： - 设计：主要是固定 \(r\)-regular 格点（各方向格点数 \(m = \lfloor n^{1/d} \rfloor\)），也扩展到随机设计（假设协变量密度有界且 Lipschitz）。 - 噪声：\(\varepsilon\) 是 sub‑Gaussian，参数 \(\sigma\) 已知；对分类扩展假设二项噪声。 - 单调性：\(f\) 对每个坐标单调不减（弱增）。 - 信号强度条件：Power bound 要求存在 \(\delta > 0\) 使得在真上集 \(R\) 中，\(f(x) \geq \tau + \delta\)，且下界证明在较弱的条件下仍成立。
- 与已有文献的差异：相比现有 isotonic 估计工作（如 Han et al. 2017）只关心 \(\ell_2\) 风险，本文首次关注子集推断。相比多重检验工作（如 Ramdas et al. 2019），本文提供了如何从数据构造 p‑value 并保证同时有效性的具体操作。

主要结果（理论型，选最关键三个定理）：

  • Theorem 1 (Type I error control)：对任意 \(\alpha \in (0,1)\) 和任意单调 \(f\)，

    \[\Pr\bigl( \widehat{R} \subseteq R \bigr) \geq 1 - \alpha,\]
    其中概率基于数据生成。证明关键：每个点的 p‑value 在各自的零假设下是 super-uniform（\(\Pr(\widetilde{p}_x \leq u) \leq u\)），而 RISS 的拒绝集定义保证了只有当所有被拒绝点的 p‑value 都 \(\leq \alpha\) 时才可能产生假阳性；通过 union bound over 所有可能的上集（其实是在 DAG 结构上做 union bound），得到整体控制。

  • Theorem 2 (Power bound)：设 \(\delta > 0\) 且 \(f(x) \geq \tau + \delta\) 对几乎所有 \(x \in R\)。则存在常数 \(C, \gamma > 0\) 使得

    \[\mathbb{E}\bigl[|R \setminus \widehat{R}|\bigr] \leq C\, (\log n)^{\gamma}\, \frac{1}{n^{1/d}}.\]
    这里 \(|\cdot|\) 是 Lebesgue 测度。证明思路：对离 \(R\) 边界有一定距离（至少 \(\ell\)）的点，构造 p‑value 以指数速率衰减；然后通过测度论得到整体 bound。需要用到 Bullwinkle 不等式（Howard et al. 2021）来控制 p‑value 的尾部。

  • Theorem 4 (Minimax lower bound)：对任何满足 Type I error 条件（\(\Pr(\widehat{R} \subseteq R) \geq 1-\alpha\)）的估计量 \(\widehat{R}\)，存在单调函数 \(f\) 使得

    \[\mathbb{E}\bigl[|R \setminus \widehat{R}|\bigr] \geq \frac{c}{n^{1/d}}.\]
    证明采用构造性硬例子：取 \(f\) 为分段常数，\(R\) 是一个体积为 \(1/2\) 的上集，并利用 Le Cam’s method 或 Fano’s inequality。这一步表明 RISS 的 power 速率是最优的（至多对数因子差）。

证明路线与技术技巧（理论型必写，具体）：

  • 整体路线：
  • 对每个格点构造 p‑value：利用 Howard et al. (2021) 的 empirical Bernstein martingale。具体地，对给定 \(x\)，定义其“上锥” \(C(x) = \{ x' \in [0,1]^d : x' \geq x \text{ componentwise} \}\)。取所有位于 \(C(x)\) 内的设计点，计算累积和 \(S_n(x) = \sum_{i: X_i \in C(x)} (Y_i - \tau)\)。在零假设 \(H_x\) 下该和的上偏差通过一个 time‑uniform bound 控制，得到 p‑value \(\widetilde{p}_x\)。
  • 组织 DAG：将所有格点按坐标偏序构成 DAG（\(x \preceq y\) 当且仅当每一坐标分量前者不大于后者）。由单调性，拒绝某一点蕴含拒绝其所有后继（DAG 中的“后代”）。这对应 Gabriel 的 coherent 性质。
  • 多重检验算法（RISS 实现）：采用 Meijer & Goeman (2015) 的 DAG 序贯 reject 算法——从局部最小节点（没有前驱）开始，依次测试其 p‑value 是否小于临界值，若小于则拒绝该节点及其整个后代，否则继续。临界值通过 Bonferroni 类型调整保证 FWER。关键验证：由于 martingale p‑value 的 super-uniform 性质，该算法满足 Type I error 控制。

  • Power 分析：对每个在真集 \(R\) 内且离边界 \(\epsilon\) 以上的点，根据 Hoeffding 型不等式证明其对偶 p‑value 几乎肯定小于任意多项式阈值；然后对边界层进行积分，得到期望缺失测度 bound。

  • 关键跳跃点：

  • Lemma 5：证明了 \(\widetilde{p}_x\) 的确是 p‑value（即使对依赖数据构造的上锥随机停时）。这依赖于 Howard et al. (2021) 的 empirical Bernstein 停时引理。难点在于噪声方差未知时需用样本方差替代，但仍保持 valid。
  • Proposition 1（多重检验的 FWER 控制）：论证 RISS 序列拒绝算法保持 coherence 且 FWER \(\leq \alpha\)。这里用到了 Gabriel (1969) 的“coherent test imply FWER control” 这一经典结论。

  • Lemma 8（幂尾部的指数衰减）：对于真集内的点，构造了下鞅并证明 \(\widetilde{p}_x\) 以指数速率趋于零。

  • 技术技巧点名：

  • martingale / time‑uniform confidence sequence（Howard et al. 2021）——用于构造每个点的 p‑value。
  • empirical Bernstein 不等式——在未知方差下替代 Hoeffding。
  • DAG 闭包与拓扑排序——用于多重检验算法。
  • Gabriel (1969) coherent 原则——证明拒绝集自然形成上集。
  • Le Cam’s method——用于 minimax 下界。

真实例子：
论文第四节的“Fuel consumption dataset”——来自 UCI 的 Auto MPG 数据。使用两个协变量（重量 weight、加速度 acceleration），响应为 MPG（每加仑英里数）。设置阈值 \(\tau = 20\)（即希望找出 MPG > 20 的车对应的特征区域）。RISS 程序在格点设计下（20×20 网格）运行，输出一个被拒绝的上集。结果表明：该上集主要由重量轻且加速度大的车构成，符合直觉。对比方法：如果忽略多重检验，直接对每个点进行单次检验（无校正），会得到更大但含许多假阳性的区域。RISS 则更保守但保证了错误拒绝概率低于 5%。作者还展示了在协变量空间中画出 RISS 的拒绝边界，并与简单阈值法对比，凸显其形状约束带来的平滑性和可解释性。

🔎 结论是否比证明窄：
- Theorem 2 的 power bound 要求 \(f(x) \geq \tau + \delta\) 对几乎所有 \(x \in R\)，若 \(\delta\) 很小（信号弱），bound 中的常数 \(C\) 可能随 \(\delta\) 变大。本文未给出紧依赖，只说了“存在某个常数”。实践中若信号弱到噪声水平附近，power 会显著下降。
- 所有结果针对的是固定格点设计。虽然 Section 3.3 讨论了随机设计的推广，但条件更强（要求协变量概率密度有界且 Lipschitz），且未提供相同程度的 power 最小最大化证明。文中明确说“这些是初步结果，有待进一步研究”。
- 对分类响应和二项噪声的扩展（Section 3.4）只给出了 FWER 控制，没有给出 power 的minimax 界。作者在讨论中承认了这一点（“Extending the power analysis to these settings is an interesting direction for future work”）。

四、开放问题（点到为止，扎根具体语句）

  1. 弱信号下的 power 紧界：Theorem 2 的常数依赖于 \(\delta\)，但当 \(\delta \to 0\) 时，bound 退化很快。论文 Theorem 2 注释中写道“the constant \(C\) may depend on \(\delta\) in a way that we have not made explicit”。寻找在信号刚好等于 \(\tau\) 时的精确 minimax 率（可能涉及 \(n^{-1/(d+1)}\)？）尚属开放。

  2. 随机设计的更紧结果：Section 3.3 的随机设计扩展仅给出 FWER 控制，但 power bound 未与固定设计同样紧。论文明确说“we leave a full power analysis under random design for future work”。这需要一个关于协变量密度光滑性的更精细处理。

  3. 异质性治疗效应的理论保证：虽然 Section 5 展示了如何将 RISS 用于 HTE 设定（替代响应为估计的 CATE），但该部分本身不含任何有限样本理论。论文末段写道“the theoretical development for the HTE setting is beyond the scope of this work”。补上完整的 FWER 和 power 分析是一个自然方向。

  4. 与其他形状约束的结合：论文只利用了单调性。若假设函数为凸函数或 Lipschitz，是否可以类似构造 DAG？凸函数不诱导偏序，但可诱导其他形式的结构化假设。本文在结论中提到“we believe similar ideas could be applied to other shape‑restricted classes such as convex functions”，但未给出任何具体操作。

此外，读者可去核实：Deng, Han, Zhang (2020) 的工作（关于多元 isotonic 点态置信区间）是否可用于直接构造子群的置信区域？若可，其 power 如何与 RISS 比较？这没有出现在本论文的讨论中。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论