跳转至

Bayesian machine learning approach for recurrent events studies using Soft Bayesian Additive Regression Trees (SBART)

作者: MengXing Chen, Debajyoti Sinha, Antonio Linero
主题: 流行病学
相关性: 6/10
链接: https://arxiv.org/abs/2606.12701

一、领域脉络与小综述

这个方向是什么

本文所涉方向是复发事件数据的强度函数建模，特别关注如何在非齐次泊松过程（NHPP）的框架下，用非参数（贝叶斯集成学习）方法同时估计时间常数基线、个体间未观测异质性（frailty）以及协变量与时间的复杂非线性交互作用。该子方向服务于流行病学/临床随访研究中的核心问题：当研究对象可能经历多次相同类型事件（如多次住院）时，如何建模计数的条件强度，并给出个体水平的预测与不确定性量化。

当前成熟度：半参数和参数方法高度成熟（如 Sinha 1993 的贝叶斯比例强度模型、Andersen-Gill 计数过程方法），但完全非参数且能捕获时变效应的贝叶斯集成学习工具仍处于早期拓展阶段。

发展脉络（history）

  • 奠基工作 (1990s–2000s)：Oakes (1992) 提出 frailty 模型框架，将一个潜在随机效应（frailty）乘入强度函数，以解释同一受试者内事件间的相关性。Sinha (1993) 在贝叶斯框架下做半参数比例强度模型（即 λ_i(t|W_i) = λ_0(t) W_i exp(βᵀx_i)），奠定了贝叶斯复发事件建模的基础。Cook & Lawless (2021) 系统研究了间歇观测下的独立性条件，厘清了在何种假设下可以使用标准似然。

  • 主要进展——集成学习的引入 (2010–2018)：Chipman et al. (2010) 提出 BART（Bayesian Additive Regression Trees），将浅决策树的集成作为非参数回归工具，并在多种任务上展示出色预测能力。Hahn et al. (2020) 将其用于因果推断（BCF），Deshpande et al. (2020) 将其扩展到时变系数模型（VCBART）。然而，BART 的硬分割导致估计不光滑。

  • 当前 frontier——软分割与生存数据扩展 (2018–2022)：Linero & Yang (2018) 提出 SBART，用概率性软分割替换硬指示函数，使集成回归函数可自动适应未知光滑性，并将后验收缩率推到稀疏、可加函数上的近 minimax 最优。Basak et al. (2022; 本文引 [17]) 将 SBART 扩展到成簇且区间删失的生存数据。Linero et al. (2022; [19]) 进一步为删失生存时间开发了 MBART（将生存时间建模为 NHPP 首次事件）和 CoxBART，并证明简化版 MBART 的后验收缩率在高维稀疏渐近意义下接近 minimax 最优。

  • 本文的位置：往上述流水线上加一级——将 SBART 从“单次事件/首次事件”推向“同一受试者的多次复发事件”，允许强度函数中的非参数项 同时依赖时间与协变量（b(t,x)），而不仅仅是协变量。

子线索聚类

  1. 频率派复发事件建模（Andersen & Gill 1982; Lin et al. 2000; Murris et al. 2025 [RecForest]）——侧重边际均值和累积强度，不建模个体水平的 frailty。
  2. 贝叶斯复发事件建模（Sinha 1993; Manda & Meyer 2005; Pennell & Dunson 2006）——用 frailty 捕获异质性，但协变量效应通常被假设为 log-linear 且不随时间变化。
  3. BART / SBART 集成学习的发展（Chipman 2010 → Linero & Yang 2018 SBART → Hahn 2020 BCF → Deshpande 2020 VCBART → Basak 2022 区间删失 SBART → Linero 2022 MBART/CoxBART）——一条从硬分割到软分割、从简单回归到复杂生存数据的演进链。
  4. 复发事件的树集成方法（Murris et al. 2025 的 RecForest）——只做频率派边际估计，不做个体轨迹预测；本文与它构成直接竞争与互补。

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

  1. 如何同时建模基线风险、时变协变量效应与个体未观测异质性，而又不施加比例风险假设？
  2. 如何将 BART/SBART 的无参数学习能力从一个强度过程（每次观测一个事件）扩展到用同一个过程生成的多次事件（每名受试者含多条事件时间）？
  3. 如何在 MCMC 中高效处理似然中涉及的时间积分 ∫ Φ(b(t,x)) dt——这个积分对每个 MCMC 迭代、每个受试者都得算？
  4. 当模型假设被违反（如齐次泊松过程真、frailty 分布误设）时，方法是否仍稳健？

目前主流方法（比例强度贝叶斯模型）的瓶颈：① 协变量效应必须常数（不随时间变）；② 回归函数只能 log-linear；③ 无法自动检测交互。 RecForest 避开了①–③但只做边际、不做个体预测。

⚠️ 作者的 framing

作者把缺口 frame 成："Existing methods either lack subject-specific intensity (RecForest) or impose time-constant log-linear covariate effects (Sinha 1993, Pennell & Dunson)"。因此他们的 RecSBART 被表现为 "the obviously next step"——将 SBART 推入复发事件 NHPP 设定，自然继承软分割的光滑性、自适应稀疏性以及贝叶斯不确定性量化优势。

什么明显该被引 / 该存在、却没出现在 intro 里？ 未提及 Dirichlet 过程混合 (DPM) 模型做复发事件强度的完全非参数贝叶斯、也未提 Gaussian Process (GP) 泊松过程强度法（如 Adams et al. 2009 的 thinned PP 框架，本文引了 Adams 但只用了其扩增技巧，未提 GP 替代）。此外，因果推断中针对重复事件的处理效应估计（如 Lin et al. 2000 的半参数边际方法做率比/率差的估计）没有被拿来定位——本文只讨论了"fit"而非"causal effect"。

张力

未见明显对立引用。作者与其他模型的比较都在模拟中量化呈现，而非声称对方"wrong"。

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

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

符号表 | 记号 | 含义 | 分类 | |------|------|------| | \(N_i(t)\) | 受试者 i 在时间区间 (0, t] 内经历的事件次数（计数过程） | 随机变量 / 样本 | | \(y_{i1}<\dots<y_{in_i}\) | 受试者 i 的 \(n_i\) 次复发事件发生时刻（有序） | 可观测数据 | | \(a_i\) | 受试者 i 的终止时间（结束随访、死亡等，非信息性） | 可观测数据（已知常数） | | \(x_i \in \mathbb{R}^p\) | 受试者 i 的 p 维固定协变量向量 | 可观测数据（已知） | | \(W_i\) | 受试者 i 的潜在 frailty（未观测的随机效应，乘入强度） | 潜在变量 | | \(\lambda_i(t|W_i, x_i)\) | 条件强度函数，给定 \(W_i\) 和 \(x_i\) | 模型定义量 | | \(\lambda_0\) | 常数基线强度（不随时间变化的基线） | 参数（标量） | | \(b(t,x)\) | 非参数函数：捕获协变量及时间的非线性 / 交互效应 | 参数（函数，由树集表示） | | \(\Phi(\cdot)\) | 标准正态 CDF（连接函数，将 b 的值映射到 (0,1)） | 已知函数（固定） | | \(\eta\) | Frailty 分布的超参数：\(W_i \sim \text{Gamma}(\eta, \eta)\) | 参数（标量） | | \(M\) | 集成中的软决策树棵数 | 用户指定超参数 | | \(\Psi_m, \varpi_m\) | 第 m 棵树的树结构和叶节点参数 | 参数（树结构 + 叶均值向量） | | \(\alpha_b\) | 树分支节点 b 的软分割带宽参数 | 参数（控制光滑性） |

模型（数据生成机制） 对每个受试者 i，

\[N_i(t) \mid W_i, x_i \;\sim\; \text{NHPP}\big(\lambda_i(t|W_i, x_i)\big),\]
其中
\[\lambda_i(t|W_i, x_i) = \lambda_0\; W_i\; \Phi\big(b(t, x_i)\big).\]

  • \(\lambda_0\) = 常数基线（取对数后模型无异于普通 Poisson 拦截）；
  • \(W_i\) = 乘性 frailty，Gamma(η, η) 分布，E[W_i]=1 用以识别（因为 λ_0 已吸收整体尺度）；
  • \(\Phi(b(t,x))\) ∈ (0,1)：这个构造的用意是 把 0 到 1 之间的一个乘子与 frailty 和基线相乘，使得整个强度有上界 λ_0 W_i（因为 Φ ≤ 1），从而使用 thinned Poisson 扩增。

可观测数据 与 不可观测量 - 可观测：对每个受试者 i，我们能看到的只有 \(\{y_{ij}\}_{j=1}^{n_i}\)（事件时刻）、终止时间 \(a_i\)、协变量 \(x_i\)。 - 不可观测：\(W_i\)（frailty）、\(b(t,x)\)（函数）、\(\lambda_0, \eta\)（标量参数）、树结构及叶均值、带宽α。 - 关键识别手段：W_i 的分布（Gamma）和 E[W_i]=1 假设是用于区分“基线风险 λ_0”与“额外个体间变异性 W_i”的常见做法。模型假设（2）的乘积形式本身是一个不可检验的识别假设——没有额外的外部验证数据，无法区分真正由时间-协变量交互贡献的结构性变异与由遗漏随机效应贡献的方差。

第二步：最小内核

特例——最简例子 (p=1, b(t,x)=βx 的简化版)
如果拿走几乎所有灵活性与光滑性，只保留论文的核心创新思路，可以想象一个简化问题：假设 \(p=1\)，且 \(b(t,x)\) 仅仅是协变量 x 的一个未知非线性函数 \(f(x)\)（不依赖 t），但你并不知道 f 的形式。你要估计 λ_0, W_i, f(x)，并用贝叶斯集成（软决策树）来逼近 \(f(x)\)。

在这个特例下，作者要解决的问题退化为： - 给定 n 个受试者，每人有若干复发时间，每个时间点来自一个 NHPP，强度为 λ_0 W_i Φ(f(x_i))。 - 你想要自然地从数据中学习 f(x)，而不是先假设它是什么形式（多项式、样条等）。 - 但因为积分 \(\int_0^{a_i} \Phi(f(x_i)) dt = a_i \Phi(f(x_i))\) 在 t 不进入 f 时变得简单，thinned Poisson 扩增就变得直观：
你引入一个“虚假”的潜在事件流（强度 λ_0 W_i (1−Φ(f(x_i)))），这些事件的观测等价于把每个子 i 的事件分成两类（可观测到的真实事件 vs 被丢失的潜在事件），然后通过一个二项 probit 机率（Φ vs 1−Φ）来区分它们。这样，原本不便于 MCMC 的泊松似然就转化成一个等价的带有潜变量 Z的 probit 回归似然——这正是该文两层扩增的紧凑形式。

所以核心内核：
整篇文章在数学上干的事是：将一个带有未知函数 b(t,x) 的非齐次泊松过程似然，通过“加上一层受 thinning 控制的虚拟事件 + 加上一层 probit 潜变量”的两次扩增，转化为一个能用 BART/SBART 连续响应回归来处理的形式。原本难处是似然里面的 ∫Φ(b(t,x)) dt；扩增后，所有有关 b(·) 的信息都以“被观测数据点位置上的 Φ(b) 或 1−Φ(b)”的形式出现，而它们又与一个截断正态潜变量等价。这样，b(·)的 SBART 后验更新就变成了一个标准共轭的叶参数更新问题（在潜变量条件下，叶均值 µ_ms 的后验是高斯的）。

三、这篇论文做了什么

三句话

  1. 研究了什么：针对复发事件数据（如多次住院）的条件强度函数的非参数贝叶斯估计，允许协变量效应随时间非线性变化、并包含个体 frailty。
  2. 核心工具：Soft Bayesian Additive Regression Trees (SBART) + 两层数据扩增（thinned Poisson + probit latent variable），将强度函数的后验计算转化为潜变量下的高斯回归更新。
  3. 主要结论：在仿真中（即使模型假设被违反），RecSBART 对累积强度的 AMSE 始终低于 RecForest 和贝叶斯比例强度模型；在结直肠癌再入院数据上，MSMR 最小（8.727 vs 10.153 vs 14.514），交叉验证泛化差距最小（8.7% vs 28.4%）。通过新提出的 Bayesian Marginal Effect (BME) 揭示了协变量-时间交互证据。

关键设定与假设

完整设定（在第二节最小记号上补全）

  • NHPP 假设（模型核心，不可直接检验）：给定 \(W_i\) 和 \(b(t,x)\)，事件数按非齐次泊松过程生成。
  • 乘法结构：\(\lambda_i = \lambda_0 \cdot W_i \cdot \Phi(b(t,x))\)。其中 \(\Phi\) 作为连接函数将 b 局限在 (0,1)，好处是使强度有上界，利于 thinned PP；代价是对 b 的超大值，Φ(b)≈1 后失去对额外变化的敏感。
  • Frailty 分布：\(W_i \sim \text{Gamma}(\eta,\eta)\)，独立同分布。E[W_i] = 1 用于识别。这个假设可能脆弱——若真实 frailty 分布非 Gamma，模型仍可能通过 b(·) 吸收部分变异，但作者在仿真中测试了 Gamma 与 Uniform 的情形并发现稳健性。
  • SBART 的非参数结构：b(t,x) 建模为 M 棵软决策树的和，每棵树的叶均值先验为 N(0,σ_μ²)，σ_μ 被设定为 3/(2√M)（默认正则化——使 b(·) 先验集中在合理范围）。树的深度先验参数 (γ=0.95, β=2) 抑制深度分裂。带宽 α_b 有 Gam(1, r_α) 先验，r_α=10。
  • 终止时间非信息性：a_i 不携带关于复发过程的信息。

相比已有文献 - 放宽：相比 Sinha 1993 的比例强度模型，移除了 log-linear 协变量效应和时间常数效应，允许 b(t,x) 中 t 和 x 有任意交互。 - 保持：frailty Gamma 假设和乘积结构（这是与 Sinha 1993 共用的）。 - 不变：与 Basak et al. 2022 ([17]) 相比，它们处理的是成簇且单次事件区间删失，而本文处理的是同一受试者的复发事件（时间直至终止），且后者强度依赖时间。

主要结果

模拟结果（Table I, III） | 场景 | 方法 | AMSE(Λ̂) | 排序 | |------|------|-----------|------| | A: HPP, misspecified | RecSBART | 0.035 | 1st（最小） | | A | RecForest | 0.036 | 2nd | | A | Prop.Intensity | 0.778 | worst | | B: NHPP, misspecified frailty | RecSBART | 0.032 | 1st | | B | RecForest | 0.036 | 2nd | | B | Prop.Intensity | 1.035 | worst | | C: NHPP, correctly specified | RecSBART | 0.027 | 1st | | C | RecForest | 0.030 | 2nd | | C | Prop.Intensity | 0.484 | worst |

也报告了 Frailty 估计：RecSBART 的 AMSE(Ŵ) = 0.0447，比例强度模型 0.0476（略好）。

真实数据结果（Section 6） - 数据：结直肠癌患者再入院研究，n=403，共458次再入院（frailtypack R包）。 - 三法对比：MSMR(RecSBART) = 8.727, RecForest = 10.153, Prop.Intensity = 14.514。 - 交叉验证：RecSBART 泛化差距 8.7%，RecForest 28.4%。 - BME 分析：发现随时间递增的性别效应、化疗效应、Dukes’ stage效应，以及三者的交互效应（特别是在低 Dukes’ stage 下性别×化疗交互更强）；BME 在 log 累积强度上不随时间常数，提示比例强度假设可能不成立。

证明路线与技术技巧（MCMC算法，非经典理论证明——本文无收敛定理）

算法整体路线（Algorithm 1） 1. Thinned Poisson 扩增（外层扩增）：对每个受试者 i，从 λ₀ W_i a_i 强度的 Poisson 生成虚拟候选时间点 \(c_{im} \sim U(0,a_i)\)。以概率 \(1-\Phi(b(c_{im},x_i))\) 接收这些候选点为“虚拟丢失事件”（记作 G_{ik}）。这一步使得完整的似然（式 (4)）等价于“真实事件由 potion Φ(b) 鉴别，虚拟事件由 (1-Φ(b)) 鉴别”的二元响应结构。 2. Probit 潜变量扩增（内层扩增）：对每个真实事件 y_{ij} 和每个虚拟事件 G_{ik}，引入潜变量 Z_{i·} ~ N(b(·), 1)，对虚拟事件截断在 (−∞, 0)、真实事件截断在 (0, ∞)。这一步将式 (5) 的二元响应似然转化为标准高斯回归似然（式 (6)）。 3. 更新叶参数 µ：在上述潜变量条件下，b(·) 的更新变成一个高斯回归——对每棵树的每个叶节点，其后验均值是相应潜变量在落入该叶节点下的样本的均值（乘以适当的先验精度）。这直接tractable。 4. 更新 λ₀, W_i, η：λ₀ 和 W_i 在扩增后是共轭Gamma更新；η 用 slice sampling（没有调参负担）。 5. 更新 SBART 树结构：用标准的 BART MCMC 步骤（加树、删树、改变分裂变量和 cutpoint）来探索树的生长/死亡，提出被接受概率由数据的部分残差（当前潜变量−其他树集合）决定。

关键跳跃点/难点 - 积分 ∫Φ(b(t,x)) dt：如果不做扩增，每次 MCMC 迭代都要评估这个积分对每个 i，代价太大。用 thinned Poisson 扩增后，此积分不再出现——被等价于“从 λ₀W_i a_i 中抽取虚拟事件数 p_i”的 Poisson 抽样取代。这一技巧来自 Adams et al. (2009) 的 GP-Poisson 框架，此处被适配到 SBART。 - 第二个潜在交叉点：probit 潜变量需要把 Φ(b) 解释为 P(Z > 0 | b)，但 Z 截断在 ±∞；这一步是 Albert & Chib (1993) 的标准技巧，在贝叶斯 probit 回归中被广泛使用。 - 软分割与硬分割的 ↔：因为 SBART 用 ψ(x; c, α) 概率性分配，每个数据点不再唯一地分配给一个叶节点，而是分配到所有叶节点，仅权重不同。这使得“叶节点上的样本”在严格意义上不存在——于是作者不得不把 b(t,x) 对树结构与叶参数的后验更新写成一个加权回归的形式（每个数据点对每个叶节点贡献一个权重 ψ 的乘积）。这在计算上比 BART 更贵（每棵树的每次分裂要评估每个数据点），但文中未量化这个额外开销。

技术技巧点名 - Thinned Poisson（Adams et al. 2009）——用于绕过似然中积分。 - Probit 潜变量（Albert & Chib 1993）——将二元响应 (Φ vs 1−Φ) 转化为正态回归。 - Slice Sampling（Neal 2003）——更新非共轭参数 η（frailty 精度）。 - SBART 的软分割（Linero & Yang 2018）——使用 CDF（此处为 inv-logit）作为分裂决策的概率。 - Jackknife 置信区间——对 BME 的不确定性量化，非 Bootstrap，用删除一受试者重估。

真实例子与应用

已在上面展开（Section 6）。补充： - 使用的数据集是 frailtypack R 包中的 colorectal 数据集，含 403 名结直肠癌患者、458 次入院。四个协变量（性别、化疗、Dukes' stage、Charlson 指数）。 - 方法部署：50 棵 SBART 树，2500 burn-in + 2500 sampling，资源的 BME 用后验样本计算。 - 例子主要验证两点：① RecSBART 的拟合优于两种 baseline（MSMR + 残差密度）；② 通过 BME（新定义）可以直观展示时变协变量效应与交互，而传统比例强度模型无法揭示这些模式。 - 结论比证明窄的情况：本文为纯贝叶斯计算方法+实证，无任何关于后验收缩率、一致性或 minimax 性质的理论定理。唯一提到的理论结果是引用了 Linero & Yang 2018 中关于 SBART 后验收缩到近 minimax 速率的结论，但这个结果只适用于 “simple variant of MBART”（简化版），并非本文所实现的 RecSBART 模型（该模型含额外的 frailty 和 λ₀）。

🔎 结论是否比证明窄

  • 是的：作者在文中大量声称 “RecSBART 适应未知光滑性与稀疏性”、“后验收敛” 等。但这些都是在别处证明的（被引的 [1] 证明的是 SBART 在标准回归设定下的收缩率，而非在复发事件 NHPP + frailty 设定下）。本文自始至终没有给出任何理论支撑——所有的表现评估都来自于仿真和真实数据。所以读者应注意：尽管 RecSBART 在模拟时表现良好，其贝叶斯后验的一致性、速率、模型选择一致性均未建立。

四、开放问题（简短，扎根具体语句）

  1. 高维协变量下的可扩展性
    作者说：“Future work may extend RecSBART to accommodate ultrahigh-dimensional covariates by incorporating sparsity-inducing priors”（Section 7）。目前 SBART 自身的稀疏性（在协变量选择上）通过树分裂概率的性质获得，但该性质在高维 p 较大时的后验收缩率只在原 SBART 论文的 additively separable 设定下被证明（Linero & Yang 2018），并未在当前设定下验证。

  2. Tensor / 结构化数据的限制
    作者承认 “currently restricted to vector-type covariates and responses, which limits its applicability to more complex data structures such as tensor data”（Section 7）。这可以连接你的 Tensor-network / einsum 复杂度工作：如果你能将 SBART 的叶更新推广到张量响应（如稀疏收缩路径的张量核表示），则可提供一个计算复杂度的清晰表征。

  3. BME 的统计性质未被证明
    作者首次定义了 BME（式 10），并用 jackknife 给出了数值置信区间，但BME 的后验均值是否一致地估计一个明确定义的因果量 / 预测量？在非参数模型 b(t,x) 中，BME 依赖于积分交叉后的所有其他协变量的分布（被平均掉）。这个定义类似于“部分依赖”或“ALE”，但本文没有讨论它与参数模型中的回归系数是否可比较（见 Fig 4 对比例强度假设的检验——很实用，但只算适配诊断而非推断）。

  4. Frailty 分布的假设强度
    模型假设 W_i ~ Gamma(η, η)。如果实际数据是更复杂的共享 frailty 结构（如嵌套、加性、时空相关），RecSBART 的 b(t,x) 能否自动吸收这些结构？ 这一条可基于 Cook & Lawless (2021) 的依赖过程假设加以形式化表征，本文未涉及。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论