跳转至

Joint modeling of multistate and nonparametric multivariate longitudinal data

作者: Lu You, Falastin Salami, Carina Törn, Åke Lernmark, Roy Tamura
来源: Annals of Applied Statistics
主题: 非参数 / 半参数
相关性: 7/10
链接: https://doi.org/10.1214/24-aoas1889

一、领域脉络与小综述

这个方向是什么： 这个子方向要解决的根本统计问题是：当个体在时间维度上同时经历离散的状态转移（如疾病分期、死亡）与连续的纵向测量（如生物标志物轨迹）时，如何在一个统一框架下对这两类异质过程进行联合推断，以纠正纵向测量误差、借力纵向轨迹预测状态转移，并最终实现假设检验与未来状态预测。当前该方向已高度成熟，主流软件（如 JM / joineRML）已能处理参数纵向子模型与生存/多状态子模型的联合，但一旦纵向轨迹呈现非单调、多峰等复杂非线性形态，参数假设即成为瓶颈，非参数/半参数联合建模仍处于方法框架不断扩展、但理论性质（如渐近效率、收敛速度）尚不完善的阶段。

发展脉络（history）： 从 introduction 与文献线索看，该方向的发展可串成以下几步： - 奠基工作：联合模型的开端是单变量纵向与单终点生存的耦合（如 Wulfsohn & Tsiatis, 1997，引用句指出其用线性混合模型与 Cox 模型共享随机效应，奠定了 "shared parameter" 框架）；多状态模型的奠基则归于 Andersen & Keiding (2002) 与 Putter et al. (2007)，作者引用它们以确立多状态转移强度与占用概率的计数过程语言。 - 主要进展：从单终点向多终点/多状态扩展（Putter et al., 2007 将联合模型延拓至多状态转移）；从单变量纵向向多变量纵向扩展（Rizopoulos, 2012 与 Hickey et al., 2016 引入多变量联合模型，作者引用它们以说明 "已有软件可处理多变量参数纵向"）；从参数纵向向半参数/非参数纵向扩展（Rizopoulos et al., 2013 尝试在联合模型中放宽纵向参数假设，作者引用以承认 "半参数纵向已有人做，但多变量非参数纵向+多状态的联合仍空白"）。 - 当前 frontier：如何在多变量非参数纵向设定下，既避免参数误设，又维持共享潜在过程与多状态转移强度的可识别性与可计算性。作者在 intro 明确写道："Existing joint models ... assume parametric longitudinal submodels ... extension to multivariate nonparametric longitudinal data and multistate survival processes remains unexplored"。 - 本文的位置：填补 "多变量非参数纵向 + 多状态转移" 的联合建模空白，并在 TEDDY 队列上实证。

子线索聚类： 被引文献大致落在三条子线索上： 1. 共享参数联合模型（Shared parameter JM）：Wulfsohn & Tsiatis (1997) → Rizopoulos (2012) → Hickey et al. (2016)。这一簇用潜在随机效应 \(b_i\) 将纵向线性混合模型与 Cox 生存/多状态模型缝合，核心是 \(b_i\) 的分布假设与 EM/MCMC 估计。 2. 多状态计数过程模型：Andersen & Keiding (2002) → Putter et al. (2007)。这一簇用转移强度 \(\alpha_{rs}(t)\) 与占用概率 \(P_{rs}(s,t)\) 刻画状态流，核心是 Markov / semi-Markov 假设与 Aalen-Johansen 估计量。 3. 非参数/半参数纵向平滑：Rizopoulos et al. (2013) 与样条/核平滑文献。这一簇用基函数展开或核回归替代参数轨迹，核心是偏差-方差权衡与平滑参数选择。

这个方向在追问的核心问题： 1. 可识别性：当纵向子模型退化为非参数时，共享潜在过程 \(W_i(t)\) 的无穷维轨迹如何与有限维转移强度参数 \(\gamma\) 一起被识别？需要多少个纵向标志物、何种测量频率？ 2. 估计与计算：无穷维 \(W_i(t)\) 的非参数估计（样条/核）与有限维 \(\gamma\) 的似然估计如何联合优化？EM 算法在无穷维随机效应下的收敛性与计算可行性？ 3. 预测与校准：联合模型对未来状态占用概率的预测，其均方误差或校准度如何量化？非参数纵向部分引入的平滑偏差如何影响下游生存预测？

⚠️ 作者的 framing： - 作者把缺口 frame 成什么：作者将缺口框定为 "多变量非参数纵向 + 多状态转移" 的联合模型缺失，并强调 TEDDY 研究中生物标志物轨迹的非线性与多状态进展的复杂性，使得现有参数联合模型不适用，从而本文成为 "显然的下一步"。 - 竞争路线被淡化或回避：作者回避了半参数效率理论路线（如基于 influence function 的 debiased ML），也未提及因果中介/纵向因果推断路线（如 Robins 的 g-estimation 或 longitudinal marginal structural models），这些路线同样能处理纵向测量与状态转移的耦合，但作者将其完全排除在 intro 之外。 - 明显该被引却未出现的：半参数联合模型的理论文献（如基于 sieve / penalized spline 的联合模型渐近理论，如 Lu et al., 2009 或 Groll & Tutz, 2014 的渐近性质）；多状态因果推断（如 Young et al., 2020 对多状态 treatment effect 的识别与估计）。这些缺失是研究者值得去查的信号——作者可能刻意缩小 framing 以避开与理论更深的文献正面交锋。

张力： 未见明显对立引用。参数联合模型与非参数联合模型之间是互补而非矛盾；多状态与单终点之间是推广而非对立。但有一条隐性张力：非参数纵向的平滑偏差与多状态转移强度的参数假设之间，前者是无穷维、后者是有限维，二者在似然中的地位不对称，可能导致估计量的渐近偏差分布不均匀——这一张力作者未在 intro 中点明，但值得研究者深挖。

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

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

  • 符号：
  • \(i = 1, \dots, n\)：个体索引，\(n\) 为样本量。
  • \(s, r \in \{1, \dots, K\}\)：状态索引，\(K\) 为状态总数（如 1=无自身抗体，2=有单抗体，3=有多抗体，4=T1D）。
  • \(t\)：时间变量（如年龄，0-15 岁）。
  • \(Y_{ij}(t)\)：个体 \(i\) 的第 \(j\) 个纵向标志物在时间 \(t\) 的观测值，\(j = 1, \dots, J\)（\(J\) 为标志物数量，如 \(J=2\) 对应两种自身抗体）。
  • \(t_{ijm}\)：个体 \(i\) 的第 \(j\) 个标志物的第 \(m\) 次测量时间，\(m = 1, \dots, n_{ij}\)。
  • \(W_i(t) = (W_{i1}(t), \dots, W_{iJ}(t))\)：个体 \(i\) 的共享潜在过程，\(J\) 维无穷维随机函数，驱动纵向观测与状态转移。
  • \(\alpha_{rs}(t \mid W_i(t))\)：从状态 \(r\) 到状态 \(s\) 的转移强度（hazard），依赖于潜在过程 \(W_i(t)\)。
  • \(\gamma_{rs}\)：转移强度中的有限维参数（如 Cox 模型的回归系数）。
  • \(\beta_j\)：第 \(j\) 个标志物的非参数平滑函数的系数向量（样条基展开）。

  • \(\mathcal{H}_i(t)\)：个体 \(i\) 在时间 \(t\) 之前的观测历史，包含纵向测量与状态转移记录。

  • 模型：

  • 纵向子模型：\(Y_{ij}(t) = W_{ij}(t) + \epsilon_{ij}(t)\)，其中 \(W_{ij}(t)\) 是非参数平滑轨迹（用样条基 \(B(t)\) 展开：\(W_{ij}(t) = B(t)^T \beta_{ij} + b_{ij}\)，\(b_{ij}\) 为个体随机效应），\(\epsilon_{ij}(t) \sim N(0, \sigma_j^2)\) 为测量误差。
  • 多状态子模型：\(\alpha_{rs}(t \mid W_i(t)) = \alpha_{rs,0}(t) \exp(\gamma_{rs}^T W_i(t))\)，即 Cox 型比例风险模型，转移强度的 baseline \(\alpha_{rs,0}(t)\) 非参数，\(\gamma_{rs}\) 为有限维参数，链接到潜在过程 \(W_i(t)\)。

  • 共享机制：\(W_i(t)\) 同时出现在纵向观测的均值与多状态转移强度的指数项中，构成联合模型的耦合点。

  • 可观测数据：

  • 对每个个体 \(i\)，实际能观测到的是：
    1. 纵向测量 \(\{Y_{ij}(t_{ijm}) : j=1,\dots,J; m=1,\dots,n_{ij}\}\)，在离散时间点上的重复测量。
    2. 状态转移记录 \(\{(r_i(t), s_i(t)) : t \in \text{转移时间点}\}\)，即个体在哪些时间点从哪个状态转移到哪个状态，以及右删失时间 \(C_i\)。
  • 不可观测 / 潜在量：
    1. 连续时间潜在过程 \(W_i(t)\) 的完整轨迹（只在测量时间点有带误差的观测 \(Y_{ij}(t)\)）。
    2. Baseline 转移强度 \(\alpha_{rs,0}(t)\) 的完整函数。
    3. 个体随机效应 \(b_{ij}\)。
  • 识别的关键假设：条件独立性（纵向测量误差 \(\epsilon_{ij}(t)\) 与状态转移过程条件独立，给定 \(W_i(t)\)；且 \(W_i(t)\) 的当前值足以驱动转移强度，即无更长历史依赖）。

第二步：讲最小内核

剥掉多变量 (\(J>1\))、多状态 (\(K>4\))、非参数样条展开等一般性外壳，最小内核是： 单变量纵向 (\(J=1\)) + 两状态转移 (\(K=2\), 即 "健康" → "疾病") + 样条展开的潜在轨迹。

在这个最简特例下： - 纵向子模型：\(Y_i(t) = B(t)^T \beta_i + b_i + \epsilon_i(t)\)，其中 \(B(t)\) 是低维样条基（如 3 个基函数），\(\beta_i\) 为固定效应系数，\(b_i\) 为单维随机效应。 - 多状态子模型：\(\alpha_{12}(t \mid W_i(t)) = \alpha_{0}(t) \exp(\gamma W_i(t))\)，\(W_i(t) = B(t)^T \beta_i + b_i\)。 - 联合似然：\(L_i = \int \left[ \prod_m f(Y_i(t_m) \mid b_i) \right] \times \left[ \alpha_{12}(t \mid W_i(t)) \exp\left(-\int_0^{C_i} \alpha_{12}(u \mid W_i(u)) du\right) \right] \times f(b_i) db_i\)。

要证的命题退化成什么： 在这个特例下，核心数学困难不再是无穷维非参数轨迹的估计（因为样条基已将无穷维投影到有限维 \(\beta_i\)），而是随机效应 \(b_i\) 的积分与非参数 baseline \(\alpha_0(t)\) 的联合估计。具体地，EM 算法的 E-step 需要在给定观测 \((Y_i, \text{转移记录})\) 下计算 \(b_i\) 的后验期望，而 \(b_i\) 的后验受转移强度的非参数部分 \(\alpha_0(t)\) 影响，导致 E-step 无法闭式计算，必须依赖数值积分或 MCMC。

证明怎么走、为什么成立： 本文的估计路线（在一般情形下）是： 1. 用 penalized spline 拟合纵向子模型，得到 \(W_i(t)\) 的初步估计 \(\hat{W}_i(t)\)（第一阶段）。 2. 将 \(\hat{W}_i(t)\) 代入多状态子模型，用 Cox partial likelihood 估计 \(\gamma\)，用 Nelson-Aalen 估计 \(\alpha_{rs,0}(t)\)（第二阶段）。 3. 若采用全似然方法，则联合优化纵向样条系数、随机效应分布参数与多状态参数，通过 EM 算法迭代，E-step 用 Laplace 近似或 MCMC 计算 \(b_i\) 后验，M-step 更新参数。

在最简特例下，两阶段方法的合理性来自测量误差的渐近可忽略性（当纵向测量频率 \(n_{i1} \to \infty\) 时，\(\hat{W}_i(t) \to W_i(t)\)，第二阶段估计量渐近等价于真实 \(W_i(t)\) 代入的结果）。但这一渐近可忽略性在有限测量频率下不成立，且非参数平滑偏差会渗入转移强度估计——这是最小内核暴露的核心数学困难。

三、这篇论文做了什么

三句话： ① 研究了多变量非参数纵向数据与多状态转移过程的联合建模问题，以 TEDDY 队列中自身抗体轨迹与 T1D 进展为实证驱动。 ② 核心方法是用共享潜在过程 \(W_i(t)\) 将非参数样条纵向子模型与 Cox 型多状态子模型耦合，采用两阶段或全似然 EM 算法估计。 ③ 主要结论是：提出的联合模型在模拟中能恢复真实转移强度与纵向轨迹，在 TEDDY 数据上能检验组间转移率差异并预测未来状态占用概率，优于分离建模。

关键设定与假设： - 设定：\(n\) 个独立个体，每个个体有 \(J\) 个纵向标志物（\(J \geq 2\)）与 \(K\) 个状态（\(K \geq 4\)），纵向测量在离散时间点，状态转移在连续时间观测（可能右删失）。 - 假设 1：条件独立性：\(Y_{ij}(t)\) 的测量误差 \(\epsilon_{ij}(t)\) 与状态转移过程条件独立，给定 \(W_i(t)\) 与随机效应 \(b_i\)。这是联合模型可识别的关键，作者在 Section 2 明确写出。 - 假设 2：比例风险链接：\(\alpha_{rs}(t \mid W_i(t)) = \alpha_{rs,0}(t) \exp(\gamma_{rs}^T W_i(t))\)，即潜在过程 \(W_i(t)\) 对转移强度的效应是乘法的且与时间无关。相比已有文献（如 Putter et al., 2007 的 Markov 假设），这里放宽了 Markov 性（允许 \(W_i(t)\) 依赖历史），但强化了比例风险假设。 - 假设 3：样条展开的潜在轨迹：\(W_{ij}(t) = B_j(t)^T \beta_j + b_{ij}\)，其中 \(B_j(t)\) 为样条基，\(\beta_j\) 为固定效应，\(b_{ij}\) 为随机效应。这相比 Rizopoulos (2012) 的线性混合模型假设，放宽了轨迹的参数形式，但引入了平滑参数选择与偏差问题。 - 假设 4：随机效应分布：\(b_i \sim N(0, D)\)，\(D\) 为未结构化协方差矩阵。这是经典假设，与 Hickey et al. (2016) 一致。

主要结果： - 定理/命题 1（联合似然构造）：在上述假设下，个体 \(i\) 的联合似然可分解为纵向部分、多状态部分与随机效应部分的乘积积分，具体为 \(L_i(\theta) = \int f(Y_i \mid b_i; \theta_Y) f(S_i \mid b_i; \theta_S) f(b_i; \theta_b) db_i\)，其中 \(\theta = (\theta_Y, \theta_S, \theta_b)\) 分别为纵向、多状态与随机效应参数。直觉：共享 \(b_i\) 使得纵向与多状态不再独立，联合似然通过积分消除 \(b_i\) 的潜在性。 - 定理/命题 2（两阶段估计的渐近性质）：作者未给出严格的渐近定理（如收敛速度或效率界），但在 Section 3 声称两阶段估计量在纵向测量频率足够高时渐近等价于全似然估计量。必要条件：纵向测量次数 \(n_{ij} \to \infty\)（测量误差渐近可忽略）与样条基维数随样本量适当增长（控制平滑偏差）。技术难点：非参数平滑偏差与多状态参数估计的交互——偏差如何渗入 \(\gamma\) 的估计，作者未严格量化。 - 定理/命题 3（预测公式）：给定个体 \(i\) 在时间 \(s\) 的纵向历史 \(\mathcal{H}_i(s)\)，未来时间 \(t\) 的状态占用概率为 \(P(S_i(t) = r \mid \mathcal{H}_i(s)) = \int P(S_i(t) = r \mid W_i) f(W_i \mid \mathcal{H}_i(s)) dW_i\)，通过 Monte Carlo 积分计算。直觉：预测需要将 \(W_i\) 的后验分布（从纵向历史推断）前推到未来时间，再代入多状态模型计算转移概率。

证明路线与技术技巧： - 整体路线： 1. 构造联合似然（纵向 × 多状态 × 随机效应的积分）。 2. 两阶段估计：先用 penalized spline 拟合纵向子模型，提取 \(\hat{W}_i(t)\)；再将 \(\hat{W}_i(t)\) 代入多状态 Cox 模型，用 partial likelihood 与 Nelson-Aalen 估计 \(\gamma\) 与 \(\alpha_{rs,0}(t)\)。 3. 全似然估计：EM 算法，E-step 用 Laplace 近似或 MCMC 计算 \(E[b_i \mid Y_i, S_i; \theta^{(k)}]\)，M-step 分别更新纵向样条系数、多状态参数与随机效应协方差 \(D\)。 4. 预测：从 \(\hat{W}_i(t)\) 的后验抽取样本，代入多状态模型模拟未来状态轨迹，计算占用概率。 - 关键跳跃点： - E-step 中 \(b_i\) 的后验分布受多状态似然的影响，无法闭式计算（因为多状态似然包含 \(\exp(-\int \alpha_{rs,0}(u) \exp(\gamma^T W_i(u)) du)\) 的积分项）。作者用 Laplace 近似 绕过这一困难，将后验近似为高斯分布，从而闭式计算条件期望与方差。这一近似在 \(\gamma\) 较小或测量频率较高时较精确，但在 \(\gamma\) 较大或测量稀疏时可能失效——作者未给出 Laplace 近似误差的严格界。 - Penalized spline 的平滑参数选择：作者用 restricted maximum likelihood (REML) 或 cross-validation 选择惩罚参数，但未讨论平滑偏差对下游多状态估计的影响——这是两阶段方法的已知弱点（偏差传递）。 - 技术技巧点名： - Penalized spline smoothing：用于纵向子模型的非参数轨迹估计，惩罚项控制样条系数的粗糙度，起作用在于避免过拟合同时允许非线性轨迹。 - Laplace approximation：用于 EM 算法的 E-step，近似 \(b_i\) 的后验分布，起作用在于将不可积的后验转化为可计算的高斯近似。 - Cox partial likelihood：用于多状态子模型的 \(\gamma\) 估计，起作用在于避免 baseline \(\alpha_{rs,0}(t)\) 的非参数估计干扰有限维参数的推断。 - Nelson-Aalen estimator：用于 baseline 转移强度的非参数估计，起作用在于提供 \(\alpha_{rs,0}(t)\) 的步函数估计，无需参数假设。 - Monte Carlo integration：用于预测步骤，从 \(b_i\) 的后验抽取样本并模拟未来状态轨迹，起作用在于计算不可闭式的状态占用概率。

真实例子与应用： - 数据 / 场景：TEDDY（The Environmental Determinants of Diabetes in the Young）出生队列，追踪从出生到 15 岁的儿童，观测两种自身抗体（GAD65, IA-2）的纵向出现与消失，以及 T1D 的发病。状态定义为：1=无抗体，2=有单抗体（GAD65 或 IA-2），3=有多抗体，4=T1D。样本量约 8,676 个儿童，纵向测量约 3-6 个月一次。 - 怎么把本文方法用上去： 1. 纵向子模型：对每个儿童的 GAD65 与 IA-2 测量值（连续量，如抗体水平），用 penalized spline 拟合非线性轨迹（抗体水平随年龄的上升与下降）。 2. 多状态子模型：对状态转移（1→2, 2→3, 3→4 等），用 Cox 模型参数化转移强度，链接到样条拟合的潜在抗体轨迹 \(\hat{W}_i(t)\)。 3. 预测：给定儿童在年龄 5 岁的抗体历史，预测其在 10 岁时处于各状态的概率。 - 得到什么结果： - 假设检验：检验 HLA 基因型组（高风险 vs 低风险）的转移强度差异，得到 \(\gamma\) 的显著估计，表明高风险组的 1→2 与 2→3 转移率更高。 - 预测：联合模型的 AUC（预测 T1D 发病）高于仅用基线信息的模型，表明纵向抗体轨迹的动态信息提升了预测精度。 - 这个例子想说明什么：验证联合模型在真实复杂队列中的可行性，展示非参数纵向轨迹捕捉非线性动态（如抗体短暂出现后消失）的能力，以及相对于分离建模（先拟合纵向再独立拟合多状态）的预测优势。

🔎 结论是否比证明窄： - 作者在 Section 3 声称两阶段估计量 "asymptotically equivalent to full likelihood estimator when measurement frequency is high"，但未给出严格定理或收敛速度证明。这一声称比实际证明宽——实际证明可能需要额外的平滑偏差界与测量误差渐近可忽略的严格条件（如 \(n_{ij} \to \infty\) 的速度相对于 \(n \to \infty\) 的速度）。 - 作者在 Section 4 的预测步骤中声称 "prediction error can be calibrated by simulation"，但未给出预测误差的渐近分布或均方误差界，这一声称比实际证明宽。 - 全似然 EM 算法的收敛性：作者未证明 EM 在此无穷维设定下的收敛性（如是否收敛到局部最大、收敛速度如何），仅声称 "EM algorithm converges in practice"——这是泛泛 claim，缺乏严格保证。

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

  1. 非参数平滑偏差对多状态参数估计的渗入与修正：两阶段方法中，纵向样条的平滑偏差如何定量影响 \(\gamma\) 的估计偏差与渐近分布？是否可以构造 debiased 估计量（如基于 influence function 的偏差修正）恢复 \(\sqrt{n}\)-收敛？扎根在 Section 3 的 "asymptotically equivalent" 声称与未给出的严格渐近定理。
  2. 联合模型的半参数效率界：在非参数纵向轨迹与多状态转移强度的联合设定下，\(\gamma\) 的半参数效率界是什么？全似然估计量是否达到此界？扎根在 Section 2 的联合似然构造与未讨论的效率理论。
  3. Laplace 近似在稀疏测量下的误差控制：当纵向测量频率较低（如 TEDDY 中某些儿童测量稀疏）时，EM 算法 E-step 的 Laplace 近似误差如何量化？是否有替代方案（如 adaptive Gaussian quadrature 或 MCMC 的收敛诊断）？扎根在 Section 3.2 的 Laplace 近似描述与未给出的误差界。
  4. 因果中介视角下的多状态转移效应：本文的联合模型是纯关联性的，若将状态转移视为中介（如抗体出现 → 多抗体 → T1D），如何在此非参数纵向框架下定义与识别自然直接/间接效应？扎根在 intro 的 "simultaneous inference" framing 与完全未提及的因果推断文献。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论