跳转至

Built-in Selection Bias in Proportional Hazards Models with Omitted Covariates: Simulation Evidence and Alternative Approaches

作者: Ayoub Bifenzi, Helene Jacqmin-Gadda
主题: 流行病学
相关性: 7/10
链接: https://arxiv.org/abs/2606.19982

一、领域脉络与小综述

  • 这个方向是什么：本方向关注在时间-事件分析（survival analysis）中，Cox比例风险模型估计的危险比（Hazard Ratio, HR） 因 “非可压缩性（non-collapsibility）” 而固有的选择偏倚问题。核心科学问题是：即使在随机对照试验（RCT）中，若存在影响结局但独立于处理的遗漏协变量，基于Cox模型估计的HR是否会系统性地偏离真实条件因果效应？如果是，偏离的方向和程度如何，以及有哪些替代方法可以避免或纠正这种偏倚。

  • 发展脉络（history）：

    1. 奠基工作：Cox (1972) 提出Cox比例风险模型，定义了危险比作为时间-事件分析中最常用的效应度量。Vaupel, Manton & Stallard (1979) 引入“脆弱性（frailty）”概念，指出未观测的个体间异质性会影响群体水平的风险率，这是理解该偏倚的理论雏形。
    2. 偏倚的识别与形式化：Hernán (2010) 在《Epidemiology》发表“The hazards of hazard ratios”，明确指出即使来自RCT的危险比也可能对因果推断产生误导。这是本文引用的关键催化剂之一。Greenland (2003) 在《Epidemiology》上区分了经典混杂与“碰撞分层偏倚”（collider-stratification bias），为理解选择偏倚提供了因果图语言。Stensrud & Hernán (2025,2020) 进一步呼吁研究者考虑使用不需要比例风险假设的方法，并提出使用生存差异作为替代。Post, Van den Heuvel & Putter (2024) 在《Lifetime Data Anal》中正式使用结构因果模型形式化了Cox模型的内建选择偏倚——这是最新的理论形式化工作。
    3. 替代方法与模拟评估：Balan & Putter (2020) 在《Stat Meth Med Res》上提供了frailty模型的教程，但作者明确提到其在这一场景下的鲁棒性未曾通过模拟评估（原文："their robustness in this context has not been evaluated through simulations"）。Aalen, Cook & Roysland (2015) 在《Lifetime Data Anal》上讨论了Cox分析是否给出因果效应。Struthers & Kalbfleisch (1986) 在《Biometrika》上建立了误指定比例风险模型中偏最大似然估计的渐近性质，证明了偏倚方向及其对删失分布的依赖性。Daniel, Zhang & Farewell (2021) 在《Biometr J》上给出了非可压缩效应估计量比较的经验结果。Vansteelandt et al. (2024) 在《J Am Stat Assoc》上提出了“assumption-lean” Cox回归，明确指出系数解释的复杂性。
    4. 本文的位置：作者声称，尽管已有大量理论文献指出危险比的非可压缩性和选择偏倚，但据他们所知，尚无研究通过模拟量化Cox模型与Weibull PH模型的偏倚程度，系统比较frailty模型、AFT模型和生存差异这些替代方法的校正能力，并在真实RCT数据上验证。 本文旨在填补这一“系统的模拟评估+真实数据验证”的缺口。

  • 子线索聚类：被引文献可大致分为以下子线索：

    1. 偏倚的理论形式化（因果图与数理统计）：Greenland (2003)，Aalen et al. (2015)，Post et al. (2024)。这簇工作用DAG或结构方程说明选择偏倚为什么发生，以及它是碰撞偏倚的一种形式。
    2. 偏倚的数学性质（点估计、渐近、假设检验）：Struthers & Kalbfleisch (1986)，Vansteelandt et al. (2024)，Daniel et al. (2021)。这簇工作推导了误指定模型的偏估计极限、符号、方向及其对删失分布的依赖，以及假设检验的有效性。
    3. 替代方法：Balan & Putter (2020)（frailty模型教程），Crowther, Royston & Clements (2023)（灵活参数AFT），Stensrud & Hernán (2025,2020)（生存差异）。这簇介绍或测试了可避免/校正偏倚的统计模型。
    4. 应用/模拟示范（较少）：该方向上系统模拟工作较少，本文本身就是填补这一缺口的尝试。唯一具体被引的真实数据来源是Lawton et al. (2017)（RTOG 9202试验）。

  • 这个方向在追问的核心问题（2-4个），以及当前主流方法与已知瓶颈：

    1. 偏倚的大小与方向：在给定RCT设定中，遗漏一个与处理独立但影响结局的协变量，会使处理效应（HR）朝零方向衰减（attenuation）？衰减程度取决于什么（遗漏协变量的效应β_U、基线风险、删失分布）？
    2. 何时假设检验有效？：当真实处理效应β_C=0时，Cox模型对β_unadj的检验是否控制第一类错误？当β_C≠0时，检验功效损失多大？
    3. 替代方法的鲁棒性：Frailty模型是否能正确恢复条件HR exp(β_C)？哪种AFT误差分布设定对遗漏协变量最不敏感？生存差异是否真正collapsible且估计偏差可忽略？
    4. 实际应用中的权衡：在真实数据中，当有部分协变量已被观测（如Gleason评分、PSA）、但仍可能存在未观测异质性时，调整与不调整的估计差异多大？哪种替代方法给出最稳定、最可信的因果结论？

  • ⚠️ 作者的 framing（必须明确标注成“这是作者的说法”）： 作者把缺口frame成：“虽然理论早已指出HR的非可压缩性和选择偏倚，但从未有系统的模拟工作量化其大小，也未曾比较过几种替代方法（frailty, AFT, 生存差异）在真实设定中的表现。” 作者将此缺口定性为一个“亟需的模拟与验证工作”。竞争路线被淡化的包括：使用边际结构模型（Marginal Structural Models）或逆概率加权——这在因果推断中是处理时依混杂的常用方法，但本文完全未提及；同样狭义的工具变量方法或比值比（Odds Ratio）的可压缩性讨论也未出现。明显该被引但未出现的是：关于风险比因果解释的广泛争论（如针对“危险比悖论”的多篇论述），以及对AFT模型受遗漏协变量影响的具体理论结果。 本文的所有替代方法都依赖参数假设（AFT误差分布、frailty分布、spline自由度），这些假设的轻微违背可能导致多少偏倚？作者并未充分探讨。这一未被讨论的鲁棒性边界值得研究者注意。

  • 张力：未见明显对立引用。已引用工作对“危险比非可压缩性”的基本判断高度一致，分歧主要在于替代方法的选择与解释。

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

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

    • 符号：

      • T：真实的（潜在）事件时间（随机变量）。
      • C：右删失时间（随机变量），假定与T独立。
      • T_i* = min(T_i, C_i)，观测到的追踪时间。
      • δ_i = I(T_i ≤ C_i)，事件发生指示变量（1=事件，0=删失）。
      • X：处理变量（可观测），在RCT中为二值变量，0=对照组，1=治疗组。
      • U：未观测的协变量（潜在的、未被测量或未知的影响结局的变量）。在RCT设定中，X ⊥ U（基线独立）。
      • β_C：条件对数危险比（log-hazard ratio for X after adjusting for U，即真正的因果效应）。
      • β_U：U的对数危险比（log-hazard ratio for U，影响事件风险的强度）。
      • λ_C(t|X, U) = λ_0,C(t) exp(β_C X + β_U U)：真实的条件危险函数，假设为比例危险结构。
      • β_unadj：未调整模型（仅含X的Cox模型）中的对数危险比估计量。
      • λ_M(t|X)：真实边际危险函数（由真实条件模型边缘化U得到）；HR_M(t) = λ_M(t|X=1)/λ_M(t|X=0)，时变边际危险比。
      • θ：frailty模型的方差参数（对Gamma frailty，Var(exp(U)) = θ）。

    • 模型： 真实的数据生成机制（DGM）为方程(1)的比例危险模型： λ_C(t_i | X_i, U_i) = λ_0,C(t_i) exp(β_C X_i + β_U U_i) 基线风险λ_0,C(·)通常是未知的（半参数）或参数化的（如Weibull）。X和U是独立生成的。C独立于T。我们观测到的是(T_i, δ_i, X_i)，但U_i没有被观测到*。

    • 可观测数据： 对于每个个体i，研究者观测到：

      • 处理分配X_i（已知，二值）。
      • 追踪时间T_i* (min(T_i, C_i))。
      • 事件指示符δ_i（0/1）。 不可观测的：U_i（存在于真实模型但缺失）；T_i本身（若δ_i=0则被删失）。

  • 第二步：讲最小内核——把原文的许多假设和一般性都剥掉，找出最简特例。

    最简特例：d=1（一个遗漏的连续协变量U），二值处理X，只有两个时间点？不，原文本身就是演示了最简单的连续情形。我们直接聚焦于最简特例：X为二值，U为连续N(0,1)，β_C＝-0.6，β_U＝1，Weibull基线风险，Cox PH模型。在此特例下，核心命题是“即使X⊥U，未调整Cox模型估计的HR也会偏离真实β_C，并且偏离为衰减（绝对值变小）”。

    在这个特例下，真实的条件危险函数是： λ_C(t|X,U) = λ_0,Weibull(t) exp(-0.6X + 1·U)

    可观测数据：对每个个体i，观测到X_i, T_i*, δ_i；U_i未观测到。研究者用仅含X的Cox模型（公式2）去拟合，得到一个估计值β̂_unadj（一个标量常数）。

    核心思路：为什么β̂_unadj ≠ β_C？因为危险函数λ_C(t|X,U)在每一个时刻t都条件于存活到t（即T>t）。存活压力是一个“碰撞变量”（collider），受X和U共同影响。在基线，X⊥U；但在存活者子集（T>t）中，由于高危者在早期更易失败，不同处理组的U分布会变得不平衡（参见原文Figure 2）。这种动态不平衡等同在模型中遗漏了一个随时间的混杂因子U，因此Cox模型的估计被推向零——即β_unadj均值-0.44而不是-0.6（表S1 β_U=1）。因此，本文的核心数学困难在于：在比例危险结构下解析地刻画这个被推向零的偏倚大小，以及寻找能恢复β_C的方法（frailty、AFT、生存差异），并检验这些方法在模拟中的表现以确认哪个最靠谱。

三、这篇论文做了什么

  • 三句话：①研究了Cox比例风险模型中因遗漏与处理独立的协变量U而导致的“内建选择偏倚”（危险比的非可压缩性），量化了该偏倚的大小，并展示了它对删失分布的依赖性；②比较了三种替代方法：frailty模型（直接对异质性建模）、加速失效时间模型（AFT，效应可达可压缩）以及生存差异（利用Kaplan-Meier或时变系数的Cox模型）；③通过广泛的模拟（跨越三种U分布、四种β_U强度）和RTOG 9202 RCT的真实数据，验证了：未调整的Cox/Weibull PH偏倚显著且随β_U增大而增大（标准化偏倚可达178%） ，而参数Weibull frailty模型（当基线正确指定时）和具有对数正态/对数逻辑误差的AFT模型最为鲁棒，生存差异估计（KM或时变Cox）也几乎无偏。

  • 关键设定与假设：

    • PH真模型：整个偏倚定义和推导依赖于真实数据生成机制遵循比例危险结构（方程1）。这是一个关键假设，但文中并未检验其合理性——此假设本身若违反，则会同时影响真模型和替代方法的表现。
    • X⊥U：处理与遗漏协变量在基线时独立（RCT设定）。这一假设衍生了“collider偏倚”的核心故事。
    • U的代表性：在模拟中，U分别取N(0,1)（连续）、log-gamma(1,1)（右偏）、Ber(0.5)（离散）。这覆盖了常见的连续和离散异质性类型。
    • 基线风险参数化：在参数PH模型中，假设Weibull基线分布。在非参数frailty估计中不做此假设（但作者发现半参数frailty因估计问题表现不佳）。
    • AFT模型假设：误差与(X,U)独立（ε ⟂ (X,U)），这一假设保证了AFT系数的可压缩性（方程16成立）。若ε与U相关（如U是ε的一部分），AFT的鲁棒性将退化——这正是作者在讨论中提到需要谨慎检查拟合优度的原因。
    • 删失独立于事件时间：C ⊥ T。

  • 主要结果（理论型的核心定理或量化结果）：

    1. 偏倚的数学形式（Section 2.2-2.3）：证明了真实边际危险比HR_M(t)是时变的，且与真实条件HR exp(β_C)不相等。给出了在Gamma frailty特殊情形下的闭式解（公式10），显示lim_{t→∞} HR_M(t)=1，即偏倚完全将效应衰减为零。这一点与模拟中观察到的β_unadj向0衰减一致（图3, 4-6）。
    2. 偏倚的方向与对称性（公式13及参考Struthers & Kalbfleisch 1986）：|β_unadj| ≤ |β_C|，即未调整Cox估计的效应总是被衰减（绝对值缩小）。这一结果对RCT中任何遗漏的相关U都成立。但方向性（是向零还是反向） 取决于U的类型和系数符号，原文给定β_U>0且β_C<0，偏倚是使效应向零（绝对值变小）。
    3. 假设检验的稳健性（Section 3.2）：当β_C=0（零假设），未调整Cox模型中的Wald检验对错误指定仍保留名义第一类错误率。当β_C≠0，检验功效受损（因为估计β_unadj被衰减，t统计量变小）。这是文中为数不多的理论结果，直接从Struthers & Kalbfleisch (1986) 引用得出。
    4. 核心模拟数值（Figure 4, Table S1）：以U~N(0,1)场景为例，β_U从0.2增至1时：
      • 半参数Cox：标准化偏倚从13.6%升至178.5%；覆盖率从94.1%降至59.3%。
      • Weibull frailty：标准化偏倚从5.9%升至26.3%；覆盖率始终在95%左右（除β_U=1时为95.0%）。参数frailty模型显著优于未调整模型和半参数frailty模型。
      • AFT（对数正态、对数逻辑）：标准化偏倚不超过26%（但对数正态在β_U=0.2时稍大，为33.9%），覆盖率皆在94%以上。
      • 生存差异（KM或时变Cox）：所有β_U水平下，标准化偏倚绝对值≤18%，覆盖率在92-96.5%之间。
    5. 真实数据（RTOG 9202）结果（Table 1）：模型估计的效应差异明显。未调整Cox（-0.173） vs. Weibull frailty（-0.233），后者更稳健且核实调整协变量后变化很小。AFT对数正态/对数逻辑估计在两个设定间几乎不变（例如对数正态：0.167 vs 0.163）。这些结果与模拟结论高度一致，验证了参数frailty和AFT的鲁棒性。

  • 证明路线与技术技巧：

    • 整体路线：本文不是纯理论证明，而是模拟实证+少数推导。证明方面集中在：
      1. 在2.2节通过对真实条件危险函数的边缘化推导，给出了HR_M(t)的表达式并证明其不等于exp(β_C)。关键步骤：交换积分与极限（引用Post et al. 2024），并利用存活作为collider的分析（附录A）。
      2. 在2.3节通过Gamma frailty模型推导HR_M(t)的闭式解（附录B）。关键技巧：变量变换u=log(v)，将积分化为Gamma分布密度积分（完整共轭结构）。这一推导展示了当U是Gamma分布时，生存函数的积分有闭式。
    • 关键跳跃点：最难的部分是证明未调整Cox的β_unadj不等于β_C并且总是衰减（|β_unadj| ≤ |β_C|）。文中直接引用Struthers & Kalbfleisch (1986) 的定理，未给出自己的推导。另一难点是在Gamma frailty情形下推导HR_M(t)的解析表达式（附录B），其中涉及复杂Gamma函数积分。
    • 技术技巧点名：
      • 边缘化积分：在式(4)-(5)中，将条件危险函数关于U的条件分布（U|T>t,X=x）积分，得到真实边际危险函数。关键点是这个条件分布是时变的。
      • 变量代换（u=log(v)）：在附录B中将Gamma分布的积分转化为标准Gamma密度积分。
      • Gamma密度恒等式：利用了Γ(1/θ+1) / (θ^{1/θ} Γ(1/θ)) = 1/(1+θA_{0,C})等恒等式简化闭式解。
      • M-estimation / Partial Likelihood理论：从引用Struthers & Kalbfleisch (1986) 的结果出发，依赖于部分似然估计的渐近理论。
      • Bootstrap / Greenwood公式：用于生存差异的标准误估计（Section 6末尾）。
      • Cubic Spline：对时变系数的Cox模型，用自然三次样条参数化β(t)（式(11)）。

  • 真实例子与应用：

    • 用的什么数据/场景：RTOG 9202随机对照试验（Lawton et al., 2017）。n=1116名局部晚期前列腺癌患者。处理X：长期间质性雄激素剥夺（24个月 vs 不额外治疗，在原2个月AD+R的基础上）。结局：从放疗结束到局部/远处复发或死亡的复合事件（51.35%事件）。协变量：Gleason评分（3水平）、肿瘤分期（T2 vs T3-T4）、log(PSA+0.1)、年龄（中心化/10）。
    • 怎么把本文方法用上去：作者估计了8个模型（两种不调整的PH模型、两种frailty、四种AFT），每个模型分“仅调整X”和“调整X+4个基线协变量”两种情况。对生存差异，用KM和时变Cox模型估计S(t|X=1)-S(t|X=0)在三个分位时间点。
    • 得到什么结果：未调整Cox（-0.173）比Weibull frailty（-0.233）更接近0，差异明显。当加入4个基线协变量（Gleason, 分期, PSA, 年龄）后，Cox的估计移动至-0.216，frailty变化很小（-0.224），说明这4个变量捕捉了U中的相当一部分。AFT对数正态/对数逻辑在两个场景中极其稳定（如0.167 vs 0.163）。时变Cox生存差异估值为正且显著（p<0.02）。这个例子想说明：在真实数据中，未调整Cox模型确实存在偏倚（估计偏向零），且其解释受调整与否显著影响；而frailty（尤其Weibull frailty）和AFT（尤其对数正态/对数逻辑）表现出跨设定稳健性，更适于因果推断。

  • 🔎 结论是否比证明窄：是的。作者在讨论中（Section 9）将结果推广到了“在RCT中，若遗漏了影响结局的U，则Cox PH估计总是有偏且衰减”。但这一结论的严格证明只适用于U是连续且独立于X、且真实数据生成比例危险的设定（Struthers&Kalbfleisch,1986）。然而，文中对AFT鲁棒性的模拟结果是在特定误差分布假设下得出的，作者也提醒“需仔细评估其基本假设”。另外，对半参数frailty模型，模拟显示其标准误的估计不稳定，但作者注明结论仅针对本文的个体frailty设定，而此类模型最初是为集群数据设计的——这意味着不能用本文结果一概否定所有半参数frailty模型。引用原句："However, it requires to well characterize the baseline risk function because the semi-parametric frailty model exhibits poor results, due to the under-estimation of the standard error of the estimates and also a bias for the targeted parameter estimate when the impact of the missing covariate and thus the variance of the frailty increase."

四、开放问题

  • 问题1（半参数frailty模型的估计改进）：本文发现半参数frailty模型（个体frailty）的渐近标准误被低估，导致覆盖率低（例如在β_U=1时覆盖率仅69.5%）。作者在Discussion中提到“it would be of interest to improve the inference procedure for the semi-parametric frailty model to better accommodate subject-specific frailty”。 这是一个明确的算法/推断问题：如何利用bootstrap或其他方法获得准确的标准误估计？这对偏倚校正的实际应用很重要。

  • 问题2（AFT模型的假设检验与稳健性边界）：AFT模型（特别是对数正态/对数逻辑误差）在模拟中极其鲁棒，但这是基于ε与(X,U)独立的假设。若U为某些特定类型（如离散且与ε相关），AFT的偏倚会如何？作者在Section 4.5只给出了鲁棒性的经验证据，未提供解析结果。 一个开放问题是：在给定真实PH模型下，当AFT误差分布被误指定时，γ̂_1的渐近偏倚的公式是什么？

  • 问题3（多协变量下的生存差异估计）：生存差异估计在模拟中表现很好。但作者仅在二值X且无协变量的场景中验证了KM方法。原文提到：对于连续X，需使用时变Cox模型，并选择“a finite and clinically relevant set of X values”。 这是一个开放的计算/实用问题：如何在低维核估计或机器学习生态中，自动估计多协变量条件下时变生存差异？

  • 问题4（理论证明更精确的偏倚界限）：现有引用（Struthers & Kalbfleisch 1986）只给出|β_unadj| ≤ |β_C|，但未给出|β_unadj|作为β_U和基线效应的显式函数。本文的模拟给出了可视化趋势，但尚未发展一个通用的解析偏差近似公式。 开放问题：在Gamma frailty之外的一般U分布下，β_unadj的极限公式是什么？能否用一个简单的校正因子恢复β_C？（类似“回归稀释偏倚”的校正公式。）

Maintained by 陈星宇 · Homepage · Source on GitHub

评论