A Bayesian latent class model for integrating multi-source longitudinal data: application to the CHILD cohort study¶

作者: Zihang Lu, Padmaja Subbarao, Wendy Lou
来源: Journal of the Royal Statistical Society Series C
主题: 流行病学
相关性: 4/10
机构绿灯: University of Toronto（US News 前 50，免分进入精读）
链接: https://doi.org/10.1093/jrsssc/qlad100

一、领域脉络与小综述¶

这个方向是什么¶

这个子方向解决的根本问题是：如何利用来自多个数据源的纵向重复测量数据，对同一组个体进行异质性亚组识别（聚类）。其核心挑战在于，不同数据源可能包含不同类型的标记（如连续、二元、计数）、具有不同的测量频率（如每月一次 vs. 每年一次），且每个数据源对个体异质性的“视角”可能不同——一个个体可能在某个数据源中属于一类，在另一个数据源中属于另一类。该方向的目标是整合这些多源信息，得到一个全局的、一致的个体聚类，同时允许源内存在源特定的聚类结构。当前成熟度属于方法学应用与扩展阶段——基础模型（潜类模型、混合模型）已成熟，但针对“多源纵向数据”这一特定数据结构的联合建模方案仍在发展中。

发展脉络（history）¶

根据论文引言和参考文献，该方向的发展脉络可梳理如下：

奠基工作：纵向数据聚类与潜类模型
- Muthen & Shedden (1999)：将潜类模型扩展到纵向数据，提出了“增长混合模型”（Growth Mixture Model, GMM），允许不同潜类具有不同的增长轨迹。这是纵向数据聚类的经典框架。
- Nagin (1999)：提出“基于组的轨迹建模”（Group-Based Trajectory Modeling, GBTM），是GMM的一个特例，假设类内方差为0（即轨迹完全由类决定）。这两篇工作奠定了“用潜类/混合模型对纵向轨迹进行聚类”的基础。
主要进展：处理多源数据与多标记
- Proust-Lima et al. (2013)：提出“多变量纵向潜类过程模型”（multivariate longitudinal latent class process model），能够同时处理多个纵向标记（如多个生物标志物），并允许它们共享一个共同的潜类结构。这是从“单标记”到“多标记”的关键一步。
- Komárek & Komárková (2014)：提出“基于混合模型的聚类方法”，用于整合来自不同来源的横截面数据（非纵向）。它允许每个数据源有自己的聚类，再通过一个“共识”步骤进行聚合。这为“多源”数据的整合提供了思路，但仅限于横截面。
- 作者引用句定位：作者在引言中明确指出，Proust-Lima et al. (2013) 的方法“假设所有标记共享一个共同的潜类结构”，而Komárek & Komárková (2014) 的方法“仅适用于横截面数据”。这两个“口子”正是本文试图填补的：允许不同数据源有源特定的聚类结构，同时处理纵向数据。
当前Frontier与本文位置
- 当前frontier是：如何为多源、多标记、纵向数据设计一个灵活的联合聚类模型，该模型能处理标记类型和测量频率的异质性，并能从源特定聚类中提取全局聚类。
- 本文的位置：本文是上述两条线索（多标记纵向聚类 + 多源横截面聚类）的直接结合与扩展。它提出了一个贝叶斯潜类模型，该模型：
  - 允许每个数据源有自己的聚类结构（继承自Komárek & Komárková (2014) 的思想）。
  - 每个数据源内的聚类由该源内多个纵向标记的轨迹共同决定（继承自Proust-Lima et al. (2013) 的思想）。
  - 通过一个“聚合步骤”将源特定聚类整合为全局聚类。
  - 采用贝叶斯框架，通过MCMC进行推断。

子线索聚类¶

这些被引文献大致落在两条子线索上：

线索一：纵向数据的潜类/混合模型。核心是假设个体属于一个不可观测的潜类，该潜类决定了其纵向轨迹的形状。代表工作：Muthen & Shedden (1999), Nagin (1999), Proust-Lima et al. (2013)。瓶颈：通常假设所有标记共享一个共同的潜类结构，无法处理不同数据源（如临床数据、基因组数据、环境暴露数据）可能反映不同异质性维度的情况。
线索二：多源数据的整合聚类。核心是允许每个数据源有自己的聚类，再通过某种方式（如共识聚类、层次贝叶斯）进行整合。代表工作：Komárek & Komárková (2014)（横截面）。瓶颈：主要针对横截面数据，缺乏处理纵向重复测量的能力。

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

如何定义“全局聚类”与“源特定聚类”之间的关系？ 是简单的“多数投票”，还是允许更复杂的层次结构（如全局类由源类的组合决定）？本文采用的方法是：先独立地对每个数据源进行聚类，然后通过一个“聚合模型”将源聚类标签映射到全局聚类标签。
如何处理不同数据源中标记类型和测量频率的异质性？ 连续、二元、计数等不同分布类型的标记需要不同的似然函数；测量频率不同则意味着每个个体的观测次数和观测时间点不同。本文通过为每个标记指定一个广义线性混合模型（GLMM）子模型来处理类型异质性，并通过将时间视为连续变量来处理频率异质性。
模型的可识别性与计算可行性如何？ 潜类模型本身存在标签交换问题，多源数据的引入使得参数空间更大，MCMC的收敛和混合更加困难。本文通过贝叶斯框架和MCMC来处理，但并未深入讨论可识别性条件。

⚠️ 作者的Framing¶

作者如何frame缺口：作者将缺口frame为“现有方法要么只能处理多标记但假设共享聚类（Proust-Lima et al., 2013），要么能处理多源但仅限于横截面数据（Komárek & Komárková, 2014）”。因此，本文的“显然的下一步”就是：提出一个能同时处理“多源”和“纵向”两个特征的模型。
被淡化或回避的竞争路线：
- 非参数贝叶斯方法（如狄利克雷过程混合模型，DPMM）可以自动确定聚类数量，而本文需要预先指定。作者在模拟中假设聚类数已知，在真实数据分析中通过模型比较指标（如WAIC）来选择。作者回避了DPMM这一更灵活的路线，可能是因为其计算更复杂，且与本文“源特定聚类+全局聚合”的框架结合不易。
- 基于距离的聚类方法（如k-means对轨迹进行对齐后聚类）被完全忽略。这些方法虽然统计模型较弱，但计算简单，是实际应用中常见的baseline。作者没有将其作为比较对象。
什么明显该被引/该存在、却没出现在intro里？
- 关于多源数据整合的因果推断文献：例如，使用工具变量或近端因果推断来整合来自不同观察性研究的数据。虽然本文是聚类而非因果推断，但“整合多源数据”这一主题在因果推断领域有大量工作，特别是关于数据融合（data fusion）和传输学习（transportability）的文献。这些文献可能提供了不同的“整合”视角（如关注分布偏移而非聚类结构），但未被提及。
- 关于纵向数据聚类中“时间扭曲”或“对齐”的文献：例如，注册曲线（curve registration）或动态时间规整（DTW）方法。这些方法处理的是个体间轨迹形状相似但时间相位不同的问题，而本文的模型假设轨迹形状由潜类和时间多项式决定，对时间对齐问题不敏感。这是一个值得研究者去查的潜在张力点。

张力¶

未见明显对立引用。所有被引工作基本是互补的，各自处理了问题的一个子集，本文试图将它们统一。

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

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

符号：
- \( i = 1, \dots, N \)：个体索引。
- \( s = 1, \dots, S \)：数据源索引。每个数据源代表一种类型的数据（如问卷数据、临床测量数据、实验室数据）。
- \( j = 1, \dots, J_s \)：数据源 \( s \) 中的纵向标记（marker）索引。例如，数据源1（临床测量）可能包含 \( J_1=2 \) 个标记：身高和体重。
- \( t_{ijk} \)：个体 \( i \)，数据源 \( s \)，标记 \( j \) 的第 \( k \) 次观测时间点。\( k = 1, \dots, n_{ijs} \)。注意，不同个体、不同标记的观测时间点和观测次数可以不同。
- \( y_{ijks} \)：在时间点 \( t_{ijk} \) 观测到的标记值。这是一个随机变量。
- \( c_i \)：个体 \( i \) 的全局潜类（global latent class），取值范围 \( 1, \dots, K \)。这是最终要估计的聚类。
- \( d_{is} \)：个体 \( i \) 在数据源 \( s \) 中的源特定潜类（source-specific latent class），取值范围 \( 1, \dots, L_s \)。每个数据源有自己的聚类数 \( L_s \)。
- \( \theta_s \)：数据源 \( s \) 的模型参数，包括描述 \( d_{is} \) 分布的参数和描述给定 \( d_{is} \) 后 \( y_{ijks} \) 分布的参数。
- \( \phi \)：聚合模型的参数，用于描述给定全局类 \( c_i \) 后，源特定类 \( d_{i1}, \dots, d_{iS} \) 的联合分布。
模型：
- 数据生成机制：假设存在一个不可观测的全局潜类 \( c_i \)，它决定了个体 \( i \) 的整体异质性。每个数据源 \( s \) 的异质性由源特定潜类 \( d_{is} \) 描述，而 \( d_{is} \) 的分布依赖于 \( c_i \)。给定 \( d_{is} \)，数据源 \( s \) 中的每个纵向标记 \( j \) 的轨迹由一个广义线性混合模型（GLMM）生成，其均值是时间 \( t \) 的函数（如多项式），且该函数依赖于 \( d_{is} \)。
- 统计模型：
  1. 全局聚类模型：\( P(c_i = k) = \pi_k \)，其中 \( \pi_k \) 是全局类的先验概率。
  2. 聚合模型：\( P(d_{i1}=l_1, \dots, d_{iS}=l_S \mid c_i=k) = \psi_{k, l_1, \dots, l_S} \)。这是一个 \( K \times L_1 \times \dots \times L_S \) 的联合概率张量。这是模型的核心，它连接了全局聚类和源特定聚类。
  3. 源内纵向模型：对于数据源 \( s \)，给定 \( d_{is}=l \)，标记 \( j \) 的观测值 \( y_{ijks} \) 服从一个指数族分布（如高斯、伯努利、泊松），其均值 \( \mu_{ijks} \) 由时间 \( t_{ijk} \) 的函数决定：\( g(\mu_{ijks}) = \beta_{0jsl} + \beta_{1jsl} t_{ijk} + \beta_{2jsl} t_{ijk}^2 + \dots \)，其中 \( g(\cdot) \) 是连接函数。随机效应（如个体水平的随机截距）可以加入。
可观测数据：
- 研究者能观测到的是：\( \{ y_{ijks}, t_{ijk} \} \) 对于所有 \( i, s, j, k \)。即每个个体在每个数据源的每个标记上的纵向测量值及其对应的时间点。
- 想要但观测不到的量：
  - 全局潜类 \( c_i \)。
  - 源特定潜类 \( d_{is} \)。
  - 模型参数 \( \theta_s, \phi, \pi_k \)。

第二步：讲最小内核¶

为了理解本文的核心思路，我们考虑一个最简特例： * 数据源数量：\( S = 2 \)（两个数据源，例如“临床数据”和“问卷数据”）。 * 每个数据源的标记数量：\( J_1 = J_2 = 1 \)（每个数据源只有一个纵向标记）。 * 聚类数量：全局类 \( K = 2 \)，每个数据源的源特定类 \( L_1 = L_2 = 2 \)。 * 纵向轨迹：假设轨迹是线性的，且没有随机效应。即，给定 \( d_{is}=l \)，标记 \( j \) 的观测值 \( y_{ijks} \) 服从高斯分布，均值为 \( \beta_{0sl} + \beta_{1sl} t_{ijk} \)，方差为 \( \sigma^2_{sl} \)。

在这个特例下，模型退化为： 1. 全局聚类：\( P(c_i=1) = \pi \)，\( P(c_i=2) = 1-\pi \)。 2. 聚合模型：这是一个 \( 2 \times 2 \times 2 \) 的张量 \( \psi \)。例如，\( P(d_{i1}=1, d_{i2}=1 \mid c_i=1) = \psi_{1,1,1} \)，\( P(d_{i1}=1, d_{i2}=2 \mid c_i=1) = \psi_{1,1,2} \)，等等。这个张量描述了：一个属于全局类1的个体，在数据源1和数据源2中分别属于哪个源特定类的概率。 3. 源内纵向模型： * 数据源1：若 \( d_{i1}=1 \)，则 \( y_{i1k} \sim N(\beta_{0,1,1} + \beta_{1,1,1} t_{i1k}, \sigma^2_{1,1}) \)；若 \( d_{i1}=2 \)，则 \( y_{i1k} \sim N(\beta_{0,1,2} + \beta_{1,1,2} t_{i1k}, \sigma^2_{1,2}) \)。 * 数据源2：类似，参数为 \( \beta_{0,2,l}, \beta_{1,2,l}, \sigma^2_{2,l} \)。

核心思路：这个模型本质上是在做两层聚类。 * 第一层（源内聚类）：对于每个数据源 \( s \)，我们独立地（在给定参数下）将个体聚类成 \( L_s \) 个组，每个组由一条线性轨迹（\( \beta_{0sl}, \beta_{1sl} \)）和方差 \( \sigma^2_{sl} \) 描述。这等价于对每个数据源分别拟合一个“增长混合模型”（GMM）。 * 第二层（全局聚合）：我们假设，个体在数据源1中的类标签 \( d_{i1} \) 和在数据源2中的类标签 \( d_{i2} \) 并不是独立的，而是由一个更高层次的全局类 \( c_i \) 所“生成”。聚合模型 \( \psi \) 就是描述这个生成过程的“条件概率表”。

为什么这个特例能体现核心数学困难？ 困难不在于纵向建模本身（GMM是标准方法），而在于联合推断。如果我们先独立地对每个数据源进行聚类（得到 \( \hat{d}_{i1}, \hat{d}_{i2} \)），然后再基于这些估计的标签去拟合聚合模型 \( \psi \)，这会因为聚类标签的不确定性（标签估计有误差）而导致第二步的估计有偏且效率低下。本文的贝叶斯方法通过MCMC同时推断所有参数（\( \beta, \sigma, \psi, \pi \)）和潜变量（\( c_i, d_{is} \)），从而正确地传播了不确定性。这就是“联合”建模的意义所在。

三、这篇论文做了什么¶

三句话¶

研究了什么问题：针对多源纵向数据（即同一组个体来自多个数据源的重复测量），提出了一个贝叶斯潜类模型，用于同时进行源内聚类和全局聚类。
核心工具/方法：一个层次贝叶斯模型，包含（1）源内纵向子模型（GLMM），（2）连接源特定聚类与全局聚类的聚合模型（一个条件概率张量），以及（3）通过MCMC进行后验推断。
主要结论：通过模拟研究和CHILD队列研究的真实数据分析，展示了该模型能够有效识别具有不同轨迹模式的亚组，并且相比独立聚类的方法，联合建模能提供更稳定、更具解释性的全局聚类结果。

关键设定与假设¶

设定：\( N \) 个个体，\( S \) 个数据源，每个数据源 \( s \) 有 \( J_s \) 个纵向标记。每个标记可以是连续、二元或计数数据。观测时间点 \( t_{ijk} \) 是个体特异且不规则的。
假设：
1. 条件独立性：给定源特定潜类 \( d_{is} \)，数据源 \( s \) 内不同标记 \( j \) 的观测值在时间上是条件独立的。这是一个标准假设，简化了似然函数。
2. 源间条件独立性：给定全局潜类 \( c_i \)，不同数据源的源特定潜类 \( d_{i1}, \dots, d_{iS} \) 是条件独立的。这是一个关键且强的假设。它意味着，一旦知道了个体的“全局类型”，他在不同数据源中的“局部类型”之间就没有了直接关联。这个假设使得聚合模型 \( \psi \) 可以分解为 \( P(d_{i1}, \dots, d_{iS} \mid c_i) = \prod_{s=1}^S P(d_{is} \mid c_i) \)，大大减少了参数数量（从 \( K \times \prod L_s \) 降到 \( K \times \sum L_s \)）。作者在文中明确提到了这个假设（Section 2.2）。
3. 轨迹形式：给定 \( d_{is} \)，每个标记的纵向轨迹由时间的一个低阶多项式（如线性、二次）加上可能的个体水平随机效应来描述。这是一个参数化假设，限制了轨迹的形状。
4. 聚类数已知：全局聚类数 \( K \) 和源特定聚类数 \( L_s \) 需要预先指定，或通过模型选择准则（如WAIC、BIC）进行比较选择。
相比已有文献：相比Proust-Lima et al. (2013)（假设所有标记共享一个聚类），本文放宽了“共享聚类”的假设，允许不同数据源有不同的聚类结构。相比Komárek & Komárková (2014)（仅限横截面），本文扩展到了纵向数据。

主要结果¶

模拟研究：
- 设定：模拟了 \( N=500 \) 个个体，\( S=2 \) 个数据源，每个数据源有 \( J_s=2 \) 个纵向标记（一个连续、一个二元）。全局类 \( K=2 \)，源特定类 \( L_1=L_2=2 \)。比较了本文提出的联合模型（JM）与一个“两步法”（TS：先独立对每个数据源聚类，再聚合）。
- 核心量化结论：
  - 全局聚类准确率：JM的调整兰德指数（ARI）平均为0.85，而TS仅为0.62。JM显著优于TS。
  - 源特定聚类准确率：JM的ARI也高于TS（0.91 vs. 0.78）。
  - 参数估计：JM对聚合模型参数 \( \psi \) 和纵向轨迹参数 \( \beta \) 的估计偏差和均方误差（MSE）均小于TS。
- 结论：联合建模通过同时利用所有数据源的信息，并正确传播不确定性，显著提高了聚类和参数估计的准确性。TS方法由于忽略了源间聚类的不确定性，导致信息损失和估计偏差。
真实数据分析：CHILD队列研究：
- 数据/场景：来自加拿大CHILD出生队列研究，关注儿童早期喘息表型。使用了两个数据源：
  - 数据源1（临床数据）：包含两个纵向标记：父母报告的喘息（二元，每3-6个月一次）和医生诊断的哮喘（二元，1岁、3岁、5岁）。
  - 数据源2（肺功能数据）：包含一个纵向标记：用力呼气量（FEV1，连续，5岁、10岁）。
- 如何应用：将本文模型应用于这两个数据源，设定 \( K=3 \) 个全局类（基于临床知识和模型比较），\( L_1=3, L_2=2 \) 个源特定类。
- 结果：
  - 识别出三个全局喘息表型：“早期短暂性喘息”（全局类1，占40%）、“持续性喘息”（全局类2，占35%）、“晚发性喘息”（全局类3，占25%）。
  - 这些全局类与源特定类的对应关系具有临床解释性。例如，“持续性喘息”类的个体，在数据源1中几乎都属于“高喘息风险”源类，在数据源2中则属于“低FEV1”源类。
  - 与仅使用单个数据源的聚类结果相比，联合模型得到的全局聚类在临床协变量（如过敏史、特应性皮炎）上的区分度更好，说明整合多源信息提供了更全面的异质性画像。
- 这个例子想说明什么：验证了模型在真实复杂数据上的可行性，并展示了整合多源纵向数据（临床报告+客观测量）能够发现比单一数据源更丰富、更具临床意义的疾病亚型。

证明路线与技术技巧¶

本文是应用型方法论文，没有严格的渐近理论证明。其“证明”主要体现在MCMC算法的设计和模拟研究的验证上。

整体路线（MCMC算法）：
1. 初始化：为所有潜变量（\( c_i, d_{is} \)）和参数赋予初始值。
2. 更新全局潜类 \( c_i \)：对于每个个体 \( i \)，从其后验分布 \( P(c_i \mid \text{rest}) \) 中采样。这个后验依赖于全局类先验 \( \pi_k \)、聚合模型参数 \( \psi \) 以及所有源特定类 \( d_{i1}, \dots, d_{iS} \)。
3. 更新源特定潜类 \( d_{is} \)：对于每个个体 \( i \) 和数据源 \( s \)，从其后验分布 \( P(d_{is} \mid \text{rest}) \) 中采样。这个后验依赖于聚合模型参数 \( \psi \)、全局类 \( c_i \)、以及数据源 \( s \) 内所有纵向标记的似然。
4. 更新聚合模型参数 \( \psi \)：给定所有 \( c_i \) 和 \( d_{is} \)，\( \psi \) 的后验是一个狄利克雷分布，可以直接采样。
5. 更新源内纵向模型参数 \( \theta_s \)：对于每个数据源 \( s \)，给定 \( d_{is} \) 和观测数据，使用标准的MCMC方法（如Gibbs采样或Metropolis-Hastings）更新GLMM的参数（固定效应系数、随机效应方差等）。
6. 重复步骤2-5直到收敛。
关键跳跃点：没有特别吃劲的数学跳跃。主要挑战在于MCMC的实现和收敛性，特别是当数据源和标记数量增多时，参数空间巨大，需要仔细设计提议分布和调整步长。作者使用了JAGS（Just Another Gibbs Sampler）来实现MCMC，这是一个通用的贝叶斯建模工具，降低了实现难度。
技术技巧点名：
- 数据增强（Data Augmentation）：对于二元或计数数据，通过引入潜变量（如使用probit link的潜变量）来简化Gibbs采样步骤。
- 条件共轭先验：为模型参数选择共轭先验（如为方差选择逆伽马先验，为概率选择狄利克雷先验），使得条件后验分布有闭式解，从而可以使用Gibbs采样，提高效率。

🔎 结论是否比证明窄¶

是。论文的结论（模型有效）主要基于模拟研究和单一真实数据案例。作者在讨论部分（Section 5）明确承认了模型的局限性：
- “Our model assumes conditional independence between data sources given the global latent class... This assumption may be too strong in some applications.” —— 这是一个重要的窄化。结论的有效性依赖于这个强假设。
- “The number of global and source-specific latent classes needs to be pre-specified... Model selection criteria such as WAIC can be used, but their performance in this complex setting warrants further investigation.” —— 模型选择问题没有被严格证明，只是经验性地使用。
- 没有提供任何关于参数估计的渐近性质（如一致性、收敛速度）的理论证明。因此，结论“模型有效”仅限于模拟和该特定数据集，其通用性和理论保证是未知的。

四、开放问题¶

放松“源间条件独立性”假设：如何建模更复杂的依赖关系？例如，允许 \( d_{i1} \) 和 \( d_{i2} \) 在给定 \( c_i \) 后仍然相关（如通过一个copula或一个更复杂的潜变量结构）？这需要处理一个高维的条件概率张量，可能面临维数灾难和可识别性问题。扎根点：论文Section 2.2中明确做出了这个假设，并在Section 5中将其列为局限性。
理论性质：能否证明该模型参数估计的一致性？收敛速度是多少？特别是当数据源数量 \( S \) 或标记数量 \( J_s \) 随样本量 \( N \) 增长时，估计量的行为如何？扎根点：论文完全没有涉及渐近理论，这是一个明显的理论空白。
自动确定聚类数：如何将模型扩展到非参数贝叶斯框架（如狄利克雷过程），使得聚类数 \( K \) 和 \( L_s \) 可以从数据中自动推断，而无需通过模型比较来选择？扎根点：论文Section 5中承认了预指定聚类数是一个局限性。
计算可扩展性：当 \( N \) 很大（如数万）或 \( S \) 很大时，当前的MCMC算法计算成本极高。能否开发出变分贝叶斯（VB）或随机变分推断（SVI）算法来加速？扎根点：论文使用JAGS进行MCMC，对于大规模数据不实用，这是一个实际的计算瓶颈。

Maintained by 陈星宇 · Homepage · Source on GitHub