跳转至

Joint Group-Based Trajectory Modeling for Paired Repeated Measures: An Application to Audiometric Phenotypes and Dietary Associations

作者: Ying Chen, Sharon Curhan, Kenneth I. Vaden Jr, Judy R. Dubno, Molin Wang
主题: 因果推断
相关性: 6/10
链接: https://arxiv.org/abs/2607.23858

一、领域脉络与小综述

这个方向是什么

这个子方向是组基轨迹建模（Group-Based Trajectory Modeling, GBTM），一种基于有限混合模型的半参数方法，用于从重复测量数据中识别具有不同轨迹模式的潜在亚组。其根本科学问题是：在存在个体异质性的纵向/有序测量数据中，如何将总体分解为若干个具有不同发展轨迹的、有临床或生物学意义的潜在类别，并研究这些类别与协变量（如暴露、干预）的关联。当前成熟度较高，在临床、犯罪学、发展心理学等领域有广泛应用，但其核心假设——给定潜在组后轨迹条件独立——在处理配对器官（如双耳、双肾）数据时经常被违反。

发展脉络（history）

  • 奠基工作：Nagin (1999) 提出半参数组基方法（GBTM），将轨迹建模为有限个多项式函数，并用最大似然估计组参数和成员概率。Jones, Nagin and Roeder (2001) 开发了SAS过程 PROC TRAJ，使方法可被广泛使用。Nagin (2005) 出版了专著，系统化该框架。这些工作确立了GBTM作为“严格分类”工具的地位——它假设组内无随机变异，所有变异由组均值轨迹和独立误差解释。
  • 主要进展：Nagin and Tremblay (2001) 提出双轨迹建模（dual-trajectory modeling），用于分析两个相关行为的发展轨迹，但依赖关系仅通过连接概率（linking probabilities） 建模，即两个轨迹的组分配是独立的，仅通过一个转移矩阵关联。Nagin et al. (2018) 提出多轨迹建模（multi-trajectory modeling），将多个相关轨迹联合建模，但强制它们共享同一个潜在组分配，且仍假设条件独立。这两条路线都未能同时满足：(i) 处理配对轨迹间的层次相关性，(ii) 允许配对轨迹属于不同潜在组。
  • 当前frontier：处理层次依赖性的主流方法是增长混合模型（Growth Mixture Model, GMM）（Muthén and Shedden, 1999; Proust-Lima, Philipps and Liquet, 2017），它通过引入个体水平随机效应来捕捉组内变异和相关性。但作者明确指出，GMM“inflate within-group variability and blur distinct phenotypic shapes”（第2页），这与GBTM“严格分类”的临床目标相悖。Yang et al. (2024) 尝试处理耳间依赖，但依赖预定义表型，且仅在分类后边际地考虑相关性。
  • 本文的位置：作者声称“no GBTM framework currently exists that can simultaneously (i) handle the hierarchical correlation inherent in paired trajectories, and (ii) jointly model the paired trajectories while permitting them to be assigned to different latent groups”（第3页）。本文填补了这一空白，提出联合GBTM，显式建模配对轨迹的层次协方差结构，同时允许它们属于不同潜在组。

子线索聚类

  1. 传统GBTM及其扩展：Nagin (1999, 2005), Jones, Nagin and Roeder (2001), Nagin and Odgers (2010a,b)。核心：条件独立假设，组内无随机效应。适用于单条轨迹。
  2. 处理相关轨迹的GBTM扩展：Nagin and Tremblay (2001)（双轨迹，独立分配+连接概率），Nagin et al. (2018)（多轨迹，共享组分配+条件独立）。这些方法要么依赖关系建模不充分，要么强制共享组。
  3. 处理层次依赖性的混合模型：GMM (Muthén and Shedden, 1999; Proust-Lima, Philipps and Liquet, 2017)，通过随机效应捕捉相关性，但模糊了组间差异。Ng and McLachlan (2014) 和 Sheng et al. (2022) 处理多水平增长轨迹的聚类，但属于GMM框架。
  4. 听力表型分类的应用：Dubno et al. (2013), Vaden et al. (2017) 使用监督学习；Parthasarathy et al. (2020) 使用高斯混合模型；Yang et al. (2024) 使用软分类+回归。这些方法要么分别处理左右耳，要么依赖预定义表型，未充分利用联合层次结构。

核心问题与瓶颈

  • 核心问题1：如何处理配对轨迹间的层次依赖性（如耳间、耳内相关性），同时保持GBTM的“严格分类”哲学？
  • 核心问题2：如何允许配对轨迹属于不同潜在组（如左耳和右耳可能呈现不同听力表型）？
  • 核心问题3：当存在稀有潜在组时，如何稳定估计组参数和协变量-组关联？
  • 已知瓶颈：传统GBTM的条件独立假设在配对数据中不成立，导致有偏估计和错误推断（Goldstein, 2011; Davies, Glonek and Giles, 2017）。GMM虽能处理依赖性，但会模糊组间差异，不适合需要严格分类的临床场景。

⚠️ 作者的framing

  • 作者把缺口frame成：现有GBTM扩展要么不能处理层次相关性（双轨迹模型），要么强制共享组分配（多轨迹模型），因此本文是“显然的下一步”——同时解决这两个问题。
  • 被淡化或回避的竞争路线：GMM被明确批评为“inflate within-group variability and blur distinct phenotypic shapes”（第2页），但作者未讨论GMM中是否可以通过约束随机效应方差来近似GBTM的严格分类。此外，作者未讨论贝叶斯非参数方法（如狄利克雷过程混合模型）作为替代方案。
  • 什么明显该被引/该存在、却没出现在intro里：作者未引用任何关于混合模型可识别性的经典文献（如Frühwirth-Schnatter, 2006），尽管本文的联合似然函数涉及大量参数，可识别性是一个关键问题。此外，未引用Davies, Glonek and Giles (2017) 关于协方差误设对GBTM影响的论文，尽管该论文直接支持了本文的动机。（值得研究者去查的问题：本文的联合模型在什么条件下是可识别的？）

张力

未见明显对立引用。所有被引工作基本一致地认为：GBTM的条件独立假设在配对数据中不成立，需要处理层次依赖性。GMM和GBTM之间的分歧是哲学性的（连续分布vs.严格分类），而非技术性矛盾。

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

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

  • 符号：
  • \( i = 1, \dots, n \)：参与者索引。
  • \( j = 1, 2 \)：耳索引（左、右）。
  • \( q = 1, \dots, Q \)：频率索引（如500, 1000, ..., 8000 Hz），\( Q = 7 \)。
  • \( Y_{ijq} \)：参与者 \( i \) 的耳 \( j \) 在频率 \( q \) 处的听力阈值（可观测）。
  • \( \mathbf{Y}_{ij} = (Y_{ij1}, \dots, Y_{ijQ})^\top \)：耳 \( j \) 的 \( Q \)-维听力阈值向量（可观测）。
  • \( \mathbf{Y}_i = (\mathbf{Y}_{i1}^\top, \mathbf{Y}_{i2}^\top)^\top \)：参与者 \( i \) 的 \( 2Q \)-维配对听力阈值向量（可观测）。
  • \( K_{ij} \in \{1, \dots, K\} \)：耳 \( j \) 的潜在组标签（不可观测，要估计）。
  • \( \phi_k = (\phi_{k0}, \dots, \phi_{kp})^\top \)：组 \( k \) 的轨迹参数（\( p+1 \) 维，如二次多项式系数）（要估计）。
  • \( \mathbf{W} \)：\( Q \times (p+1) \) 设计矩阵，由频率 \( h_1, \dots, h_Q \) 的多项式基构成（已知）。
  • \( \varepsilon_i = (\varepsilon_{i1}^\top, \varepsilon_{i2}^\top)^\top \)：\( 2Q \)-维误差向量（不可观测）。
  • \( \Sigma_{k_1 k_2} \)：给定联合组 \( (K_{i1}=k_1, K_{i2}=k_2) \) 下的 \( 2Q \times 2Q \) 协方差矩阵（要估计）。
  • \( \pi_k \)：边际组分配概率，\( P(K_{ij}=k) = \pi_k \)。
  • \( \mathbf{X}_{ij} = [1, \mathbf{x}_{ij}^\top]^\top \)：耳 \( j \) 的协变量向量（可观测）。
  • \( \beta_k \)：组 \( k \) 的 logistic 回归系数（要估计）。
  • \( b_i = (b_{i2}, \dots, b_{iK})^\top \)：参与者 \( i \) 的随机效应向量（不可观测），\( b_i \sim \text{MVN}(0, \Gamma) \)，\( \Gamma = \text{diag}\{\gamma_2, \dots, \gamma_K\} \)（要估计）。

  • \( \pi_{k_1 k_2}(\mathbf{X}_i | b_i) = P(K_{i1}=k_1, K_{i2}=k_2 | b_i) \)：给定随机效应下的联合组分配概率。

  • 模型：

  • 数据生成机制：给定联合组 \( (K_{i1}=k_1, K_{i2}=k_2) \)，配对听力阈值向量 \( \mathbf{Y}_i \) 服从多元正态分布：
    \[\mathbf{Y}_i | (k_1, k_2) \sim \text{MVN}\left( (\mathbf{I}_2 \otimes \mathbf{W}) \phi_{k_1, k_2}, \Sigma_{k_1 k_2} \right),\]
    其中 \( \phi_{k_1, k_2} = (\phi_{k_1}^\top, \phi_{k_2}^\top)^\top \)，\( \mathbf{I}_2 \) 是 \( 2 \times 2 \) 单位矩阵，\( \otimes \) 是Kronecker积。
  • 协方差结构：\( \Sigma_{k_1 k_2} = (\mathbf{I}_2 \otimes \mathbf{V}) \Lambda_{k_1 k_2} (\mathbf{I}_2 \otimes \mathbf{V}) \)，其中 \( \mathbf{V} = \text{diag}(\sigma_0 + \sigma_1, \dots, \sigma_0 + Q\sigma_1) \) 是频率特异标准差的对角矩阵，\( \Lambda_{k_1 k_2} \) 是 \( 2Q \times 2Q \) 相关矩阵，具有交换性结构（见公式(4)）。

  • 组分配机制：给定随机效应 \( b_i \)，两耳的组分配条件独立：

    \[\pi_{k_1 k_2}(\mathbf{X}_i | b_i) = \pi_{k_1}(\mathbf{X}_{i1} | b_i) \pi_{k_2}(\mathbf{X}_{i2} | b_i),\]
    其中 \( \pi_k(\mathbf{X}_{ij} | b_i) \) 由混合效应多项logistic回归定义（公式(18)）。

  • 可观测数据：\( \{ (\mathbf{X}_i, \mathbf{Y}_i) : i=1,\dots,n \} \)，即每个参与者的协变量和配对听力阈值向量。

  • 不可观测/潜在量：组标签 \( K_{ij} \)、随机效应 \( b_i \)、误差项 \( \varepsilon_i \)。

第二步：最小内核

最简特例：假设只有 \( Q=2 \) 个频率点（如500 Hz和4000 Hz），\( K=2 \) 个潜在组（如“正常”和“受损”），且协方差结构为最简单的形式：耳内相关 \( \rho_1 \) 相同，耳间相关 \( \rho_2 \) 相同，无频率间交叉相关（即 \( \rho_3=0 \)）。此时，模型退化为：

  • 每个耳 \( j \) 的 \( 2 \)-维观测 \( \mathbf{Y}_{ij} \) 在给定组 \( k \) 下服从二元正态分布，均值 \( \mathbf{W} \phi_k \)，协方差 \( \mathbf{V} \Lambda_k \mathbf{V} \)，其中 \( \Lambda_k = \begin{pmatrix} 1 & \rho_1 \\ \rho_1 & 1 \end{pmatrix} \)。
  • 配对数据 \( \mathbf{Y}_i = (\mathbf{Y}_{i1}^\top, \mathbf{Y}_{i2}^\top)^\top \) 在给定联合组 \( (k_1, k_2) \) 下服从四元正态分布，均值 \( (\mathbf{W} \phi_{k_1}^\top, \mathbf{W} \phi_{k_2}^\top)^\top \)，协方差：
    \[\Sigma_{k_1 k_2} = \begin{pmatrix} \mathbf{V} \Lambda_{k_1} \mathbf{V} & \rho_2 \mathbf{V} \mathbf{1}_{2 \times 2} \mathbf{V} \\ \rho_2 \mathbf{V} \mathbf{1}_{2 \times 2} \mathbf{V} & \mathbf{V} \Lambda_{k_2} \mathbf{V} \end{pmatrix},\]
    其中 \( \mathbf{1}_{2 \times 2} \) 是全1矩阵。

核心思路：传统GBTM假设给定组后两耳独立，即 \( \rho_2 = 0 \)。本文的核心是显式估计 \( \rho_2 \)（以及更一般的层次相关结构），而不是假设它为0。这通过将协方差矩阵 \( \Sigma_{k_1 k_2} \) 参数化并纳入似然函数来实现。在EM算法中，E步计算后验概率 \( \tau_{ik_1 k_2} \) 时，使用完整的 \( \Sigma_{k_1 k_2} \)（包含 \( \rho_2 \)），从而校正了因忽略耳间相关导致的偏差。M步中，\( \phi_k \) 和 \( \Sigma_{k_1 k_2} \) 的更新也依赖于这些后验概率，形成一个自洽的迭代过程。

为什么成立：当 \( \rho_2 > 0 \) 时，忽略它会导致：(i) 似然函数误设，使得参数估计有偏；(ii) 后验概率 \( \tau_{ik_1 k_2} \) 计算错误，导致组分配偏差。本文通过显式建模 \( \rho_2 \)，恢复了似然的正确形式，从而得到一致估计。模拟结果（表1）证实：当 \( \rho_2 = 0.8 \) 时，传统GBTM的覆盖概率（CP）降至87-92%，而本文方法（TwoStg）的CP回升至93-95%。

三、这篇论文做了什么

三句话

  1. 研究了什么问题：针对配对重复测量数据（如双耳听力阈值）中传统GBTM条件独立假设被违反的问题，提出联合GBTM框架，显式建模配对轨迹间的层次依赖关系，同时允许两条轨迹遵循不同的潜在模式。
  2. 核心工具/方法：基于有限混合模型和EM算法，开发了两种估计程序——单阶段（OneStg） 全最大似然估计（理论基线）和两阶段（TwoStg） 稳健估计（先拟合无协变量混合模型，再估计协变量-组关联），并采用交换性或非结构化协方差结构。
  3. 主要结论：模拟表明，忽略层次依赖会导致显著偏差（相对偏差增大、覆盖概率下降），而所提方法能有效校正；在CHEARS-AAA数据应用中，两阶段方法识别出4种听力表型，并发现更高DASH饮食依从性与更低概率的较差表型（表型3和4）相关。

关键设定与假设

  • 配对数据：每个参与者贡献两条相关轨迹（如左右耳），每条轨迹在有序域（频率）上测量。
  • 层次协方差结构：误差协方差 \( \Sigma_{k_1 k_2} \) 分解为耳内相关（\( \rho_1 \)）、耳间同频相关（\( \rho_2 \)）和耳间异频相关（\( \rho_2 \rho_3 \)），见公式(4)。假设所有组共享相同的 \( \rho_1, \rho_2, \rho_3 \)（“common correlation structure”），或允许组特异（“group-specific”）。
  • 组分配机制：给定随机效应 \( b_i \)，两耳组分配条件独立（公式(7)下方）。随机效应 \( b_i \) 服从对角协方差 \( \Gamma \) 的多元正态分布。
  • 与已有文献的对比：相比传统GBTM（假设条件独立），本文放松了该假设；相比GMM（引入随机效应到均值），本文将相关性建模在残差协方差中，保持均值轨迹的“严格分类”性质。
  • 可识别性假设：未明确讨论，但两阶段方法通过先拟合无协变量混合模型来缓解稀有组导致的不可识别问题。

主要结果

  • 模拟结果（表1和表2）：
  • 当耳间相关 \( \rho_2 = 0 \) 时，三种方法（IndGBTM, OneStg, TwoStg）表现相似，偏差和覆盖概率接近名义水平。
  • 当 \( \rho_2 > 0 \) 时，OneStg和TwoStg的相对偏差（RB） 显著低于IndGBTM，覆盖概率（CP） 更接近95%。例如，\( n=500, \rho_2=0.8 \) 时，IndGBTM对 \( \beta_{12} \) 的CP为87.3%，而TwoStg为92.9%（表1）。
  • TwoStg通常具有最小的经验标准误（ESE），表明其估计更稳定。
  • 当测量变异性较小时（低 \( \sigma_0, \sigma_1 \)），所有方法表现相似，因为此时耳间相关的影响被噪声淹没。
  • 真实数据结果（表3和图1）：
  • 使用TwoStg方法（交换性协方差结构，\( K=4 \)），识别出4种听力表型：表型1（73.6%的耳，接近正常）、表型2（6.8%）、表型3（14.3%）、表型4（5.3%，最差）。
  • 更高DASH饮食依从性（Q5 vs Q1）与更低概率的表型3（OR=0.77, 95% CI: [0.48, 1.06]）和表型4（OR=0.89, [0.05, 1.73]）相关，但置信区间较宽（可能因稀有组）。
  • 图1显示，识别的表型与Dubno et al. (2013) 的4种表型（older-normal, metabolic, sensory, metabolic plus sensory）相似。

证明路线与技术技巧

  • 整体路线：
  • 模型设定：写出联合似然函数（公式(2)），包含对随机效应 \( b_i \) 的积分。
  • EM算法框架：引入潜在变量 \( Z_{ik_1 k_2} = I(K_{i1}=k_1, K_{i2}=k_2) \)，写出完全数据对数似然（公式(3)）。
  • E步：计算后验概率 \( \tau_{ik_1 k_2}^{(ν)} \)（公式(7)），使用当前参数估计和随机效应的经验贝叶斯估计 \( \hat{b}_i^{(ν)} \)。
  • M步：将参数分为 \( \Theta_1 = \{\phi, \Sigma\} \) 和 \( \Theta_2 = \{\beta, \Gamma\} \)，分别最大化 \( U_1 \) 和 \( U_2 \)。
    • 更新 \( \phi_k \)：通过最大化 \( U_1 \) 得到加权最小二乘解（细节在补充材料S2）。
    • 更新 \( \Sigma_{k_1 k_2} \)：使用矩估计（公式(14)-(15)）或数值优化（L-BFGS-B，公式(13)）。
    • 更新 \( \beta \) 和 \( b_i \)：使用Zhang (2021) 的IRLS算法，避免数值积分。
  • 两阶段变体：第一阶段拟合无协变量混合模型（使用因子化边际概率 \( r_{k_1} r_{k_2} \) 代替 \( \pi_{k_1 k_2} \)），得到稳定的后验概率 \( \tau_{ik_1 k_2}^* \)；第二阶段固定这些概率，仅估计 \( \Theta_2 \)。
  • 关键跳跃点：
  • 稀有组问题：当某些联合组 \( (k_1, k_2) \) 稀有时，直接估计 \( \pi_{k_1 k_2} \) 不稳定。两阶段方法通过因子化 \( r_{k_1} r_{k_2} \) 来稳定估计，因为边际概率 \( r_k \) 估计更稳健。
  • 随机效应积分：直接最大化 \( U_2 \) 需要数值积分，计算昂贵。作者采用Zhang (2021) 的IRLS方法，将积分近似为条件期望，避免重复数值积分。
  • 协方差估计：当每组数据有限时，组特异非结构化协方差矩阵估计不稳定。作者提出“common correlation structure”（公式(14)-(15)），共享所有组的耳内和耳间协方差，减少参数数量。
  • 技术技巧点名：
  • EM算法：处理潜在变量和随机效应的标准框架。
  • 经验贝叶斯：在E步中用 \( \hat{b}_i^{(ν)} \) 近似随机效应，避免积分。
  • IRLS（迭代重加权最小二乘）：用于GLMM中的 \( \beta \) 更新（Zhang, 2021），避免数值积分。
  • L-BFGS-B：用于数值优化协方差参数（公式(13)）。
  • Bootstrap：用于两阶段方法的标准误估计，仅重抽样最高层次（参与者），以考虑第一阶段估计的不确定性。

真实例子与应用

  • 数据：CHEARS-AAA子队列（NHS II），3047名女性参与者，每人有左右耳在7个频率的听力阈值，以及DASH饮食依从性评分（分5个五分位数）和其他协变量。
  • 方法应用：
  • 使用二次多项式设计矩阵 \( \mathbf{W} \)（公式(17)）建模频率-阈值轨迹。
  • 使用混合效应多项logistic回归（公式(18)）建模DASH评分与潜在组关联，调整年龄、种族、吸烟、BMI等。
  • 比较交换性和非结构化协方差结构，通过BIC选择 \( K=4 \) 和交换性结构。
  • 由于稀有组（表型2和4仅占6.8%和5.3%），单阶段方法遇到可识别性问题，因此主要报告两阶段结果。
  • 结果：见表3和图1。更高DASH依从性与更低概率的较差表型相关，但置信区间较宽（可能因稀有组和样本量有限）。
  • 这个例子想说明：验证方法在真实数据中的可行性，展示其能识别有临床意义的表型（与先前文献一致），并揭示饮食与听力表型的关联趋势。

🔎 结论是否比证明窄

  • 窄结论：模拟中，当测量变异性较小时（低 \( \sigma_0, \sigma_1 \)），所有方法表现相似，本文方法未显示显著优势（见补充材料表S1-S2）。作者在讨论中承认：“when the correlation between trajectories within a cluster is negligible, the conventional GBTM serves an appropriate and robust analytic choice”（第18页）。这意味着本文方法的优势仅在高信噪比和高耳间相关时体现。
  • 泛化claim：作者声称框架可扩展到“more complex hierarchical designs, such as further longitudinal clustering”（第18页），但未提供任何理论或模拟证据。这是一个conjecture，而非严格证明。
  • 可识别性：作者未证明联合模型（公式(2)）的参数可识别性，仅通过模拟和实证表明算法收敛。对于稀有组，两阶段方法通过因子化边际概率来缓解，但未给出理论保证。

四、开放问题

  1. 可识别性理论：本文的联合GBTM模型（公式(2)）在什么条件下参数可识别？特别是当某些联合组 \( (k_1, k_2) \) 稀有时，是否存在非可识别性？扎根于：作者未讨论可识别性，仅通过两阶段方法“stabilize”估计（第8页）。
  2. 扩展到纵向数据：作者提出未来工作可扩展到“longitudinal data involving multiple time points”（第18页），但未给出具体模型或算法。如何将时间维度（多个时间点）的随机效应与频率维度的层次相关结合？扎根于：Discussion最后一段。
  3. 理论性质：本文仅通过模拟评估估计性能，未提供任何渐近理论（如一致性、渐近正态性、收敛速度）。对于两阶段方法，第一阶段估计的不确定性如何影响第二阶段推断？扎根于：作者承认两阶段方差估计“do not take into account the extra variation introduced by fitting the mixture model in the first stage”（第11页），并建议用bootstrap校正。
  4. 稀有组的更优处理：两阶段方法通过因子化边际概率来稳定稀有组估计，但这是否是最优的？能否引入惩罚似然或贝叶斯先验来进一步改善？扎根于：作者指出“when group sizes were small or highly imbalanced, two-stage approach is highly recommended”（第18页），但未探索其他方法。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论