跳转至

Robust Subgroup Analysis for Heterogeneous Censored Data

作者: Zhaohui Xu, Daoji Li, Zemin Zheng
主题: 因果推断
相关性: 7/10
链接: https://arxiv.org/abs/2607.11389

一、领域脉络与小综述

这个方向是什么

这个子方向解决的根本问题是：在结局变量存在删失（censoring） 且总体存在未知异质性（heterogeneity） 的情况下，如何在不预先知道个体子组归属的前提下，同时识别子组结构并估计子组特异性的协变量效应。其核心挑战在于：删失机制掩盖了真实的结局分布，而异质性使得传统的同质性模型产生有偏估计。当前该方向的成熟度处于方法快速扩展期——融合惩罚（fusion penalty）框架已从线性回归扩展到多种模型，但针对删失数据的稳健方法仍较少，且大多依赖较强的误差假设（如子高斯）或固定维数设定。

发展脉络（history）

  1. 奠基工作：融合惩罚框架的提出

    • Ma and Huang (2017)：首次提出用凹成对融合惩罚（concave pairwise fusion penalization）对线性回归模型进行子组分析。核心思想是用个体特异截距（subject-specific intercepts）表示异质性，并对截距的成对差异施加SCAD惩罚，从而自动将观测分成子组。这篇论文奠定了整个子领域的方法论基础——不需要预先知道子组数，且建立了Oracle性质。它留下的主要口子是：仅适用于完整数据和线性模型。

  2. 主要进展：向不同模型和数据类型扩展

    • 完整数据下的扩展：Ma and Huang (2017) 的框架被迅速推广到Poisson回归（Chen et al., 2019）、分位数回归（Zhang et al., 2019; Lu et al., 2021）、加性部分线性模型（Liu and Lin, 2019）、函数型线性回归（Li et al., 2021）、面板数据模型（Wang and Zhu, 2022, 2024）和多变量响应回归（Wu et al., 2026）。这些工作验证了融合惩罚框架的通用性，但都未处理删失数据。
    • 删失数据下的首次尝试：Yan et al. (2021) 将融合惩罚框架扩展到删失数据下的加速失效时间（AFT）模型。他们采用Buckley-James插补来处理删失，并假设误差服从子高斯分布。这是该方向的关键一步，但留下了几个口子：(a) 插补可能扭曲子组结构；(b) 子高斯假设排除了重尾误差；(c) 理论结果仅在固定维数下建立。

  3. 当前Frontier与本文位置

    • 当前Frontier：在删失数据下，如何同时实现稳健性（对重尾误差和异常值）和高维适应性（允许协变量维数随样本量增长）。此外，如何用逆概率加权（IPW） 替代插补来处理删失，以避免模型误设带来的偏差，也是一个活跃方向。
    • 本文位置：本文（Xu, Li, Zheng, 2026）直接站在Yan et al. (2021) 的肩膀上，针对其留下的三个口子提出解决方案：(a) 用IPW替代插补，仅需估计删失时间的生存函数，而非完整的误差分布；(b) 引入M-估计（Huber损失）以容纳重尾误差，放松子高斯假设；(c) 将理论结果扩展到发散维数（diverging dimension）设定，允许p和q随n增长。

子线索聚类

  1. 融合惩罚框架的模型扩展：这条线索的核心是“将Ma and Huang (2017) 的成对融合惩罚应用到不同的回归模型”。代表工作包括Chen et al. (2019, Poisson), Zhang et al. (2019, 分位数), Liu and Lin (2019, 加性部分线性), Li et al. (2021, 函数型线性), Wang and Zhu (2022, 2024, 面板数据), Wu et al. (2026, 多变量响应)。这些工作主要关注完整数据，验证了框架的通用性。
  2. 删失数据下的子组分析：这条线索专门处理删失结局。代表工作包括Yan et al. (2021, AFT模型+插补) 和 Hu et al. (2021, Cox模型+部分似然)。本文属于此线索，但引入了IPW和M-估计作为新的技术工具。
  3. 稳健估计方法：这条线索关注对异常值和重尾误差的稳健性。代表工作包括Sun et al. (2020, 自适应Huber回归) 和 Cai et al. (2021, 函数型回归的稳健M估计)。本文借鉴了M-估计的思想，但将其与融合惩罚和删失数据处理相结合。

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

  1. 如何在不依赖强分布假设（如子高斯）的情况下，对删失数据进行稳健的子组识别？
  2. 当协变量维数随样本量增长时，融合惩罚框架的Oracle性质是否仍然成立？
  3. IPW与插补两种处理删失的策略，在子组分析场景下各自的优劣是什么？
  4. 如何自适应地选择稳健化参数（如Huber损失的τ）以平衡偏差和方差？

⚠️ 作者的framing

  • 作者把缺口frame成什么：作者将Yan et al. (2021) 的插补方法、子高斯假设和固定维数设定定位为三个主要“缺口”，并声称自己的方法（IPW + M-估计 + 发散维数理论）是“显然的下一步”。作者在引言中明确列出了四点差异（(a)-(d)），直接对标Yan et al. (2021)。
  • 哪些竞争路线被淡化或回避：作者在Remark 1中提到了基于秩的方法（rank-based methods）和加权估计方程（weighted estimating equations）也是处理删失下稳健估计的常用工具，并指出Fu et al. (2021) 将其用于AFT模型的聚类，但需要预先知道子组数。作者以此淡化这些竞争路线，强调自己的方法能同时识别未知子组数和估计效应。然而，作者并未深入讨论这些方法在子组分析场景下的潜在优势或与融合惩罚结合的可能性。
  • 什么明显该被引/该存在、却没出现在intro里？：作者在讨论部分提到了Sun et al. (2020) 的自适应Huber回归，但仅在“未来工作”中提及。在引言中，作者没有引用任何关于数据自适应选择稳健化参数（如τ）的工作，而这对于M-估计的实际应用至关重要。此外，作者没有引用关于高维删失数据下变量选择的文献（如Song et al., 2014），尽管其理论框架允许p和q发散。

张力

未见明显对立引用。所有被引工作基本遵循“融合惩罚是有效的子组分析工具”这一共识，差异主要在于模型设定和数据类型。

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

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

  • 符号：

    • \(T_i\)：第\(i\)个个体的失效时间（潜在量，不可完全观测）。
    • \(C_i\)：第\(i\)个个体的删失时间（潜在量）。
    • \(T_i^* = \min(T_i, C_i)\)：观测到的生存时间（可观测）。
    • \(\Delta_i = I(T_i \le C_i)\)：删失指示符（1=事件发生，0=删失；可观测）。
    • \(x_i \in \mathbb{R}^p\)：协变量向量，其效应\(\beta\)是同质的（对所有个体相同；可观测）。
    • \(z_i \in \mathbb{R}^q\)：异质性变量向量，其效应\(\alpha_i\)是个体特异的（可观测）。
    • \(\alpha_i \in \mathbb{R}^q\)：第\(i\)个个体的异质性效应（待估参数，\(nq\)维）。
    • \(\beta \in \mathbb{R}^p\)：同质协变量效应（待估参数）。
    • \(K\)：未知的子组数量。
    • \(G_k\)：第\(k\)个子组的个体集合。
    • \(\rho_k \in \mathbb{R}^q\)：第\(k\)个子组的公共异质性效应，即对于所有\(i \in G_k\)，有\(\alpha_i = \rho_k\)。
    • \(\varepsilon_i\)：误差项，独立于\(x_i\)和\(z_i\)。
    • \(S(t) = P(C > t)\)：删失时间的生存函数（需估计）。
    • \(\hat{S}(t)\)：\(S(t)\)的Kaplan-Meier估计。
    • \(\rho(\cdot)\)：损失函数（如Huber损失）。
    • \(p(\cdot, \lambda)\)：凹惩罚函数（如SCAD），\(\lambda\)是调优参数。

  • 模型：

    • 异质性加速失效时间（AFT）模型：
      \[\log(T_i) = z_i^\top \alpha_i + x_i^\top \beta + \varepsilon_i, \quad i=1,\dots,n.\]
    • 子组结构假设：存在一个未知的划分\(G_1,\dots,G_K\)，使得\(\alpha_i = \rho_k\)对所有\(i \in G_k\)成立。\(K\)和\(G_k\)均未知。
    • 删失机制：\(C_i\)与\(T_i\)独立（条件于协变量），且其生存函数\(S(t)\)是连续的。

  • 可观测数据：

    • 研究者实际能观测到的是\(\{(T_i^*, x_i^\top, z_i^\top, \Delta_i)\}_{i=1}^n\)。
    • 想要但观测不到的量：
      1. 真实的失效时间\(T_i\)（被\(C_i\)删失）。
      2. 个体所属的子组标签（即\(i \in G_k\)）。
      3. 个体特异的异质性效应\(\alpha_i\)（待估参数）。
      4. 子组数\(K\)。

第二步：讲最小内核

本文的核心思路可以用一个最简特例来理解：单变量异质性效应（\(q=1\)）、两个子组（\(K=2\)）、无删失（\(\Delta_i=1\)）、无同质协变量（\(p=0\)）。

  • 最简特例设定：

    • 模型退化为：\(\log(T_i) = \alpha_i + \varepsilon_i\)。
    • 真实子组结构：存在两个子组\(G_1\)和\(G_2\)，使得\(\alpha_i = \rho_1\)（若\(i \in G_1\)）或\(\alpha_i = \rho_2\)（若\(i \in G_2\)），且\(\rho_1 \neq \rho_2\)。
    • 可观测数据：\(\{(\log(T_i), z_i)\}_{i=1}^n\)，其中\(z_i\)是异质性变量（此处为1，因为模型无\(z_i\)，但为了保留框架，可认为\(z_i=1\)）。实际上，在这个特例中，我们直接观测到\(\log(T_i)\)。
    • 目标：估计\(\rho_1, \rho_2\)，并识别每个个体属于哪个子组。

  • 核心思路（在这个特例下）：

    1. 目标函数：我们想最小化一个带惩罚的损失函数：

      \[Q_n(\alpha_1,\dots,\alpha_n) = \sum_{i=1}^n \rho(\log(T_i) - \alpha_i) + \lambda \sum_{1 \le i < j \le n} p(|\alpha_i - \alpha_j|, \lambda).\]
      其中\(\rho(\cdot)\)是Huber损失（提供稳健性），\(p(\cdot, \lambda)\)是SCAD惩罚（促进融合）。

    2. 惩罚的作用：SCAD惩罚对成对差异\(|\alpha_i - \alpha_j|\)施加“大棒加胡萝卜”策略：

      • 如果\(|\alpha_i - \alpha_j|\)很小（即两个个体可能属于同一子组），惩罚会强烈鼓励它们相等（融合），从而将这两个个体合并。
      • 如果\(|\alpha_i - \alpha_j|\)很大（即两个个体属于不同子组），惩罚会停止增长（SCAD的凹性），从而避免对真实差异进行有偏压缩。

    3. M-估计的作用：使用Huber损失\(\rho(\cdot)\)代替最小二乘。当误差\(\varepsilon_i\)是重尾时，最小二乘估计会被极端值严重拉偏。Huber损失在残差较小时表现为二次函数（高效），在残差较大时表现为线性函数（稳健），从而自动降权异常值。

    4. 求解过程（直觉）：

      • 算法（RISA-ADMM）通过迭代更新\(\alpha_i\)。在每一步，它都会“检查”所有成对差异\(|\alpha_i - \alpha_j|\)。
      • 如果某对差异被SCAD惩罚“融合”为0，则这两个个体被归入同一子组。
      • 最终，所有\(\alpha_i\)会收敛到少数几个不同的值（即\(\hat{\rho}_1, \hat{\rho}_2\)），从而自动完成子组识别和参数估计。

  • 为什么这个特例抓住了核心：即使在这个最简单的设定下，本文的核心挑战——同时进行稳健估计和自动子组发现——已经完整呈现。Yan et al. (2021) 的方法在这个特例下会使用Buckley-James插补（但此处无删失，所以退化为最小二乘）并假设误差子高斯。而本文的方法通过Huber损失直接处理重尾误差，通过SCAD惩罚自动发现子组。论文的一般情形（有删失、有\(x_i\)、\(q>1\)、\(p>0\)、\(K>2\)）只是在这个内核上“加壳”：用IPW处理删失，用\(z_i\)向量化异质性效应，用\(x_i\)加入同质协变量。

三、这篇论文做了什么

三句话

  1. 研究了什么问题：针对删失数据下的异质性加速失效时间（AFT）模型，提出一种稳健的子组分析方法，能在不预先知道子组归属的情况下，同时识别子组并估计协变量效应。
  2. 核心工具/方法：将逆概率加权（IPW）（处理删失）、M-估计（Huber损失，处理重尾误差）和凹成对融合惩罚（SCAD，自动发现子组）相结合，并开发了RISA-ADMM算法进行求解。
  3. 主要结论：在温和正则条件下，证明了所提估计量具有Oracle性质（即能以趋于1的概率与已知真实子组结构下的Oracle估计量相等），并在发散维数设定下建立了估计量的相合性和渐近正态性。模拟和真实数据例子展示了其对重尾误差的稳健性。

关键设定与假设

  • 模型：异质性AFT模型（公式2.2）：\(\log(T_i) = z_i^\top \alpha_i + x_i^\top \beta + \varepsilon_i\)。
  • 目标函数（公式2.1）：
    \[Q_n(\alpha, \beta; \lambda) = \sum_{i=1}^n \frac{\Delta_i}{\hat{S}(T_i^*)} \rho(Y_i^* - z_i^\top \alpha_i - x_i^\top \beta) + \sum_{1\le i < j \le n} p(\|\alpha_i - \alpha_j\|, \lambda).\]
    • IPW项：\(\Delta_i / \hat{S}(T_i^*)\)。通过Kaplan-Meier估计删失生存函数的倒数来加权，使得在删失下对完整数据损失函数的期望无偏。相比Yan et al. (2021)的插补，它仅需估计\(S(t)\)，而非整个误差分布，更稳健且更易实现。
    • M-估计项：\(\rho(\cdot)\)为Huber损失。相比Yan et al. (2021)的子高斯假设，它允许误差\(\varepsilon_i\)具有重尾分布（如t分布），通过调节参数\(\tau\)控制二次和线性区域的边界。
    • 融合惩罚项：\(p(\cdot, \lambda)\)为SCAD惩罚。相比L1惩罚，SCAD的凹性避免了有偏估计，能更准确地恢复子组结构。
  • 关键假设：
    • Condition 1（惩罚函数）：标准假设，SCAD满足。
    • Condition 2（删失时间）：存在\(\nu>0\)使得\(P(C=\nu)>0\)且\(P(C>\nu)=0\)。这保证了删失生存函数在支撑集的上界处有正质量，从而IPW权重不会趋于无穷。相比Yan et al. (2021)，这是一个更温和的假设，常用于行政删失场景。
    • Condition 3（损失函数）：Huber损失的一、二阶导数有界，且\(E[\dot{\rho}(\varepsilon_i)]=0\)。这确保了M-估计的相合性。
    • Condition 4（协变量）：\(\sup_i \|x_i\| \le c_1\sqrt{p}\)，\(\sup_i \|z_i\| \le c_2\sqrt{q}\)，且设计矩阵\(U\)的最小特征值有正下界。相比Yan et al. (2021)的固定上界假设，本文允许\(\|x_i\|\)和\(\|z_i\|\)随维数增长，这是发散维数理论的关键。

主要结果

  • 定理1（Oracle估计量的性质）：

    • 陈述：假设已知真实子组结构，Oracle估计量\((\hat{\beta}^{or}, \hat{\rho}^{or})\)的估计误差以\(\xi_n = \max\{n^{1-\gamma}/G_{\min}, n^{1/2+r}\sqrt{Kq+p}/G_{\min}\}\)为界（几乎必然收敛），且是渐近正态的。
    • 直觉：只要最小子组大小\(G_{\min}\)足够大（远大于\(\max\{n^{1-\gamma}, n^{1/2+r}\sqrt{Kq+p}\}\)），Oracle估计量就是相合的。渐近正态性为统计推断（如构造置信区间）提供了基础。
    • 必要条件：\(G_{\min} \gg \max\{n^{1-\gamma}, n^{1/2+r}\sqrt{Kq+p}\}\)。这比Yan et al. (2021)的固定维数条件更严格，但允许p和q发散（Remark 4指出，当K固定时，\(p=o(n^{1-2r})\)和\(q=o(n^{1-2r})\)）。
    • 解决的技术难点：在发散维数下，处理IPW带来的Kaplan-Meier估计误差，并建立M-估计方程的相合性。

  • 定理2（Oracle性质）：

    • 陈述：如果最小信号差异\(b_n = \min_{k\neq k'} \|\rho_{0k} - \rho_{0k'}\| > c\lambda\)且\(\lambda \gg \xi_n\)，则存在目标函数\(Q_n(\alpha, \beta; \lambda)\)的一个局部极小点，以趋于1的概率等于Oracle估计量\((\hat{\alpha}^{or}, \hat{\beta}^{or})\)。
    • 直觉：只要子组间的真实差异足够大（能被惩罚参数\(\lambda\)“检测”到），且\(\lambda\)的衰减速度慢于估计误差\(\xi_n\)，那么融合惩罚就能完美恢复真实子组结构，且参数估计与已知子组结构时一样好。
    • 必要条件：\(b_n > c\lambda \gg \xi_n\)。这要求信号足够强，且调优参数选择得当。

  • 推论1（本文估计量的渐近性质）：

    • 陈述：在定理2的条件下，本文的估计量\((\hat{\beta}, \hat{\rho})\)与Oracle估计量渐近等价，因此也是渐近正态的。
    • 意义：这为实际应用中的统计推断提供了理论依据。

证明路线与技术技巧

  • 整体路线：

    1. 建立Oracle估计量的性质（定理1）：首先假设已知真实子组结构，将问题简化为一个带IPW的M-估计问题。通过泰勒展开、经验过程理论和Lindeberg-Feller中心极限定理，证明其相合性和渐近正态性。关键步骤是处理Kaplan-Meier估计\(\hat{S}(t)\)带来的额外误差（公式B.1）。
    2. 证明Oracle性质（定理2）：证明在真实子组结构附近，目标函数\(Q_n\)的局部极小点就是Oracle估计量。这分为两步：
      • Step 1：证明在由Oracle估计量定义的“融合”子空间\(\Omega\)（即\(\alpha_i = \alpha_j\)若\(i,j\)同组）上，Oracle估计量是唯一极小点。
      • Step 2：证明在\(\Omega\)的邻域内，任何偏离\(\Omega\)的点（即试图将不同子组的个体融合）都会导致目标函数值增大。这通过泰勒展开和惩罚函数的性质（SCAD的凹性）来实现，关键不等式是公式B.5和B.6。
    3. ADMM算法收敛性（命题1）：证明RISA-ADMM算法的原始残差和对偶残差都收敛到0。

  • 关键跳跃点：

    • 处理IPW误差：在定理1的证明中，需要将基于\(\hat{S}(t)\)的估计函数\(\Phi_n(\phi)\)与基于真实\(S(t)\)的\(\Phi_n^S(\phi)\)联系起来。作者利用Kaplan-Meier估计的一致收敛速度（\(\sup_t |\hat{S}(t) - S(t)| = o(n^{-1/2+r})\) a.s.）来建立公式B.1，这是整个渐近理论的基础。
    • Oracle性质的局部性：定理2的证明并非全局性的，而是证明Oracle估计量是目标函数的一个局部极小点。这依赖于信号差异\(b_n\)足够大，使得在Oracle估计量的一个小邻域内，不同子组的个体不会被错误地融合。

  • 技术技巧点名：

    • IPRLS（迭代惩罚重加权最小二乘）：用于求解ADMM算法Step 1中的子问题（公式3.4）。通过将Huber损失的一阶导数转化为权重，将非二次优化转化为迭代重加权最小二乘问题，从而得到闭式解（公式3.5）。
    • ADMM（交替方向乘子法）：通过引入辅助变量\(\eta_{ij} = \alpha_i - \alpha_j\)，将原问题分解为可分别求解的子问题，特别是将非光滑的SCAD惩罚项分离出来，使其可以通过软阈值算子高效求解。
    • Kaplan-Meier估计：用于估计删失时间的生存函数，是IPW方法的核心。
    • 泰勒展开与经验过程：用于建立M-估计方程的渐近性质。

真实例子与应用

  • 数据/场景：德国信用数据集（German credit data）。目标变量是违约时间\(T\)，异质性变量\(z\)是支票账户状态（二值），同质协变量\(x\)包括信用金额、居住时间等9个变量。目的是识别不同信用风险的子组。
  • 方法应用：分别用本文的RISA-ADMM和Yan et al. (2021)的BJ-ADMM拟合异质性AFT模型。
  • 结果：
    • 子组识别：BJ-ADMM未能识别任何子组（\(\hat{K}=1\)），而RISA-ADMM识别出\(\hat{K}=2\)个子组。这一结果与Pei et al. (2024)的发现一致。
    • 模型诊断：同质性AFT模型的残差呈现双峰分布（左图），而RISA-ADMM拟合后的残差呈单峰分布（右图），表明异质性模型更好地捕捉了数据中的潜在结构。
    • 变量显著性：BJ-ADMM未能识别出“现居时间”、“其他债务人”等4个变量在5%水平下显著，而RISA-ADMM识别出所有这些变量均显著，且其符号与经济学直觉一致（如“现居时间长”意味着更稳定，与更长的违约时间相关）。
  • 这个例子想说明什么：该例子旨在展示本文方法在实际应用中的有效性和优越性。它说明：(a) 忽略异质性会导致模型误设（残差双峰）；(b) 本文的稳健方法（RISA-ADMM）能够成功发现BJ-ADMM未能发现的子组结构；(c) 正确识别子组有助于发现对结局有重要影响的协变量，从而提供更合理的解释。

🔎 结论是否比证明窄

  • 是。定理2的Oracle性质（公式4.2）要求最小信号差异\(b_n > c\lambda\)，且\(\lambda \gg \xi_n\)。这意味着该性质成立的前提是子组间的真实差异足够大，且调优参数\(\lambda\)的选择必须非常精确（既要足够大以区分信号，又要足够小以避免过度惩罚）。在实际应用中，这个条件是否满足很难验证。论文的结论（如“我们的估计量具有Oracle性质”）在表述上可能比证明所覆盖的条件更宽泛。
  • 此外，定理1和推论1中的渐近正态性依赖于发散维数下\(G_{\min}\)的增长率条件。在有限样本下，这个近似可能不准确，特别是当子组大小不平衡时。论文在模拟中验证了有限样本性能，但并未系统性地探讨这些理论条件被违反时的表现。

四、开放问题

  1. 高维扩展（\(p, q \gg n\)）：论文的理论允许p和q发散，但要求\(p+q < n\)。作者在讨论部分明确指出，将方法扩展到\(p, q \gg n\)的高维设定是一个有趣且具有挑战性的方向。扎根点：Section 7, "It is interesting and challenging to extend our robust subgroup analysis approach to the high-dimensional setting where both p and q can greatly exceed n."

  2. 数据自适应的稳健化参数选择：论文使用固定的Huber阈值\(\tau=1.345\)，并发现子组识别对\(\tau\)不敏感，但参数估计精度受影响。作者提到可以借鉴Sun et al. (2020)的自适应Huber回归来选择\(\tau\)。扎根点：Section 7, "one could investigate how to extend the idea of choosing a data-adaptive robustification threshold proposed in Sun et al. (2020) to our subgroup analysis framework."

  3. 更一般的损失函数与算法：论文的RISA-ADMM算法基于IPRLS框架，理论上可扩展到任何具有一阶导数的损失函数。但作者仅用Huber损失进行了验证。探索其他稳健损失（如Tukey's bisquare）的性能，或开发更高效的算法（如随机ADMM）是开放问题。扎根点：Remark 2, "any convex or nonconvex robust loss function with a well-defined score function... can be incorporated into the proposed framework."

  4. 更弱的删失假设：Condition 2假设删失时间在某个有限点\(\nu\)处有正质量，这适用于行政删失。对于更一般的删失机制（如随机删失且支撑集无界），IPW权重的方差可能爆炸，需要更复杂的截断或稳定化技术。扎根点：Condition 2及其后的讨论（"This condition is commonly adopted... and is often satisfied in clinical studies with administrative censoring."）。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论