跳转至

Model-based clustering of time-dependent observations with common structural changes

作者: Riccardo Corradin, Luca Danese, Wasiur R. KhudaBukhsh, Andrea Ongaro
来源: Statistics and Computing
主题: 流行病学
相关性: 6/10
链接: 期刊页 · arXiv

一、领域脉络与小综述

这个方向是什么

本文研究的子方向是基于模型的时间序列聚类，其核心问题是：给定一组时间序列（如不同国家的疫情曲线），如何根据它们内部结构变化（structural changes / change points）的同步性进行分组。与传统的基于形状相似性（如欧氏距离、DTW）或基于模型参数相似性（如混合回归模型）的聚类不同，本文的聚类标准是结构变化发生的时间点是否一致——两条序列若在同一时刻发生均值、方差或趋势的突变，则归为一组。这个方向处于时间序列聚类与变化点检测的交叉地带，成熟度中等：已有若干方法（如ClustSeg、funFEM）处理“共享变化点”的聚类，但大多假设变化点数量已知或变化模式简单（如多项式分段），且缺乏一个统一、灵活的贝叶斯框架来同时处理变化点位置的不确定性和聚类的不确定性。

发展脉络（history）

从introduction和参考文献中，可以梳理出以下发展脉络：

  1. 奠基工作：单序列变化点检测（1960s-2010s）

    • Chernoff & Zacks (1964)：早期工作，基于假设检验检测单个均值变化。这是变化点检测的经典起点。
    • Niu, Hao & Zhang (2016)：对多变化点检测（特别是高斯模型中的均值变化）进行了选择性综述，总结了从假设检验到一致性估计的进展。本文引用它作为“一般性回顾”。

  2. 主要进展：多序列聚类与共享变化点（2010s-2020s）

    • Same et al. (2011) - ClustSeg：提出一个混合模型，每个簇内的序列共享一个分段多项式回归模型（变化点即多项式系数的变化点）。本文引用它作为“与我们方法相关”的工作，但指出其聚类标准是“局部行为的相似性”，而非“变化点同步性”。
    • Bouveyron et al. (2015) - funFEM：提出一个基于函数型数据的判别式混合模型（FunFEM），用于聚类时间序列。本文引用它，但同样指出其聚类标准不是变化点同步。
    • Martínez & Mena (2014)：提出一个基于Pitman抽样公式的非参数变化点检测模型，用于马尔可夫过程。本文引用它作为贝叶斯非参数变化点检测的一个例子。
    • Pedroso, Loschi & Quintana (2023) - Multipartition Model：提出一个“多分区模型”，允许不同参数（如均值、方差）在不同时间点发生变化，从而识别哪些参数变了、何时变的。本文引用它作为“不同参数有不同变化点”的模型。
    • Page, Quintana & Dahl (2022)：提出一个依赖随机分区的模型，用于处理时间序列中分区的时序依赖性。本文引用它作为“变化点检测假设随机分区之间存在时间依赖”的模型。

  3. 当前Frontier：贝叶斯非参数与流行病学应用（2020s）

    • Bhatt et al. (2020) - Semi-mechanistic Bayesian model：提出一个半机制贝叶斯模型，用于估计COVID-19的时变再生数\(R_t\)，并分析非药物干预的效果。本文将其作为主要动机来源，但指出其模型是“半机制”的，而非纯数据驱动的聚类。
    • KhudaBukhsh & Rempała (2024)：证明了SEIR模型可以用一个具有时变感染率和恢复率的SIR模型来近似，并提供了大种群极限下的泛函大数定律。本文引用它来支撑其流行病学模型（SIR）的数学基础。
    • Hong & Li (2020), Gleeson et al. (2021)：讨论了具有时变感染率的流行病模型。本文引用它们作为“时变参数”建模的例子。
    • 本文的位置：作者将本文定位为填补一个空白：现有方法要么聚类标准不是变化点同步（如ClustSeg, funFEM），要么变化点检测模型不够灵活（如假设变化点数量已知），要么缺乏一个统一的贝叶斯框架来处理聚类和变化点的不确定性。本文提出一个通用的、基于随机序（random orders）的贝叶斯建模策略，可以灵活地与各种时间依赖模型（如SIR、ARIMA）结合，专门用于“共享变化点”的聚类。

子线索聚类

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

  1. 变化点检测与分割（Change Point Detection & Segmentation）：核心是检测单条或多条时间序列中的突变点。

    • 代表工作：Chernoff & Zacks (1964), Niu et al. (2016), Martínez & Mena (2014), Pedroso et al. (2023), Page et al. (2022)。
    • 核心问题：如何准确、一致地估计变化点的位置和数量？如何处理不同参数的不同变化点？如何建模变化点之间的依赖关系？

  2. 时间序列聚类（Time Series Clustering）：核心是根据某种相似性度量将多条时间序列分组。

    • 代表工作：Same et al. (2011) - ClustSeg, Bouveyron et al. (2015) - funFEM, Pamminger & Frühwirth-Schnatter (2010)。
    • 核心问题：如何定义时间序列之间的“相似性”？如何同时进行聚类和模型估计？如何处理高维或函数型数据？

  3. 流行病学建模与推断（Epidemiological Modeling & Inference）：核心是使用数学模型（如SIR、SEIR）来理解和预测疾病传播。

    • 代表工作：Bhatt et al. (2020), KhudaBukhsh & Rempała (2024), Hong & Li (2020), Gleeson et al. (2021), Seymour et al. (2022), Chitwood et al. (2022)。
    • 核心问题：如何从观测数据（如病例数、死亡数）中估计关键参数（如\(R_t\)）？如何处理时变参数和报告延迟？如何量化不确定性？

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

  1. 如何定义“结构变化”的同步性？ 是变化点位置完全相同，还是允许一定的时间偏移？本文采用“完全同步”的假设（同一组内的序列共享相同的变化点时间集）。
  2. 如何同时处理变化点数量和聚类数量的不确定性？ 贝叶斯非参数方法（如狄利克雷过程混合模型）是自然的选择，但如何将其与变化点检测模型结合是一个挑战。本文通过随机序和联结（ties）机制来解决。
  3. 如何将聚类模型与复杂的时间依赖结构（如流行病学模型）结合？ 大多数聚类方法假设序列独立或使用简单的自回归模型。本文声称其策略是“通用的”，可以结合SIR等非线性模型。
  4. 聚类结果的统计性质是什么？ 聚类是否一致（即随着样本量增加，聚类结果收敛到真实分组）？变化点估计是否一致？本文没有提供任何理论保证。

⚠️ 作者的 framing（必须明确标注成“这是作者的说法”）

  • 作者把缺口 frame 成什么？ 作者声称现有方法要么聚类标准不是“变化点同步”（如ClustSeg, funFEM），要么变化点检测模型不够灵活（如假设变化点数量已知），要么缺乏一个统一的贝叶斯框架。因此，本文的“显然的下一步”是：提出一个通用的、基于随机序的贝叶斯建模策略，专门用于“共享变化点”的聚类，并能与各种时间依赖模型结合。
  • 哪些竞争路线被他淡化或回避了？
    • 非贝叶斯方法：作者完全回避了频率学派或基于优化的聚类方法（如k-means with DTW, 或基于距离的层次聚类）。这些方法虽然不直接建模变化点，但计算上可能更高效，且在某些场景下也能达到类似目的。作者没有讨论为什么贝叶斯方法是唯一或更好的选择。
    • 理论性质：作者完全回避了聚类一致性和变化点估计一致性的理论问题。对于一个方法学论文，这通常是一个重要的考量。作者将其定位为“应用导向的建模创新”，从而回避了理论证明的责任。
    • 计算复杂度：作者没有讨论其MCMC算法的计算复杂度，特别是当时间序列数量（N）和长度（T）都很大时，算法的可扩展性如何。这对于一个实际应用（如聚类所有欧盟国家）是至关重要的。
  • 什么明显该被引 / 该存在、却没出现在 intro 里？
    • 基于距离的层次聚类方法：例如，使用DTW距离或基于特征（如变化点位置）的距离进行层次聚类。这是最直接的竞争方法，但作者完全没有提及。
    • 频率学派的多序列变化点检测：例如，基于CUSUM或似然比的方法，用于检测多条序列是否在相同时间点发生变化。这些方法更简单、计算更快，且可能有理论保证。作者没有讨论为什么需要复杂的贝叶斯模型。
    • 关于聚类一致性的理论文献：例如，关于贝叶斯聚类后验一致性（posterior consistency）的文献。作者没有引用任何相关理论工作，这暗示了本文的方法缺乏理论支撑。

张力

未见明显对立引用。所有被引工作都在各自的子线索内发展，没有出现彼此矛盾或在略不同条件下得相反结论的情况。作者通过精心选择引用，构建了一个“现有方法都不完美，本文是完美解决方案”的叙事。

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

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

  • 符号：

    • \(i = 1, \dots, N\)：时间序列的索引（如不同的国家）。
    • \(t = 1, \dots, T\)：时间点的索引（如天数）。
    • \(y_{it}\)：第\(i\)条时间序列在时间\(t\)的观测值（如每日新增病例数）。这是可观测的。
    • \(\lambda_i\)：第\(i\)条时间序列的潜在分区（latent partition）。它是一个将时间点\(\{1, \dots, T\}\)划分成若干连续区块的划分。每个区块对应一个“结构稳定”的阶段，区块的边界就是变化点。这是潜在/不可观测的。
    • \(\lambda\)：所有\(N\)条序列的潜在分区的集合，即\(\lambda = \{\lambda_1, \dots, \lambda_N\}\)。这是要估计的。
    • \(\hat{\lambda}\)：对\(\lambda\)的点估计，通过最小化Binder损失函数得到。
    • \(\theta_{it}\)：第\(i\)条序列在时间\(t\)的模型参数（如SIR模型中的感染率\(\beta\)、恢复率\(\gamma\)，或AR模型中的均值\(\mu\)）。它依赖于分区\(\lambda_i\)：在同一区块内，\(\theta_{it}\)是常数（或遵循一个简单的函数形式）。
    • \(\Theta\)：所有模型参数的集合。
    • \(\pi(\lambda)\)：潜在分区的先验分布。本文的核心创新在于定义这个先验。
    • \(\pi(\Theta | \lambda)\)：给定分区下模型参数的先验分布。
    • \(L(y | \lambda, \Theta)\)：给定分区和参数下的似然函数。

  • 模型：

    • 数据生成机制：每条时间序列\(y_i\)由一个时间依赖模型生成，该模型的参数\(\theta_{it}\)在时间上分段常数。分段常数意味着存在变化点，将时间轴分成若干段。不同序列\(i\)和\(j\)可能共享相同的变化点位置（即属于同一组）。
    • 统计模型：这是一个贝叶斯层次模型：
      1. 先验：\(\lambda \sim \pi(\lambda)\)（由随机序诱导的联结先验）。
      2. 先验：\(\Theta | \lambda \sim \pi(\Theta | \lambda)\)。
      3. 似然：\(y | \lambda, \Theta \sim L(y | \lambda, \Theta)\)。
    • 目标：推断后验分布\(p(\lambda, \Theta | y)\)，并从中得到对\(\lambda\)的点估计\(\hat{\lambda}\)。

  • 可观测数据：

    • 研究者能观测到的是：\(N\)条时间序列，每条长度为\(T\)，即\(\{y_{it}\}_{i=1, t=1}^{N, T}\)。
    • 研究者想要但观测不到的是：
      1. 每条序列的潜在分区\(\lambda_i\)（即变化点的位置和数量）。
      2. 哪些序列共享相同的变化点（即聚类结构）。
      3. 模型参数\(\Theta\)。

第二步：讲最小内核

本文的核心思路可以浓缩为一个最简特例：假设只有\(N=2\)条时间序列，每条序列的长度为\(T=3\)（时间点1, 2, 3）。我们想知道这两条序列是否在同一个时间点发生了结构变化。

  • 最简特例下的记号：

    • \(y_1 = (y_{11}, y_{12}, y_{13})\)，\(y_2 = (y_{21}, y_{22}, y_{23})\)。
    • 潜在分区\(\lambda_1\)和\(\lambda_2\)。由于\(T=3\)，可能的连续分区只有几种：
      • 无变化点：\(\{ \{1,2,3\} \}\)（记为分区A）
      • 一个变化点在时间1之后：\(\{ \{1\}, \{2,3\} \}\)（记为分区B）
      • 一个变化点在时间2之后：\(\{ \{1,2\}, \{3\} \}\)（记为分区C）
      • 两个变化点：\(\{ \{1\}, \{2\}, \{3\} \}\)（记为分区D）
    • 聚类结构：我们关心的是\(\lambda_1\)和\(\lambda_2\)是否“联结”（tied）在一起。如果它们联结，意味着它们必须是同一个分区。例如，如果它们联结，那么\(\lambda_1 = \lambda_2\)，要么都是A，要么都是B，等等。如果不联结，它们可以独立地取任何分区。

  • 核心思路（随机序诱导联结）：

    1. 随机序：作者引入一个“随机序”（random order）的概念。简单来说，就是给所有可能的变化点位置（这里是时间点1和2之间，以及时间点2和3之间）一个随机的顺序。这个顺序决定了哪些变化点“存活”下来，从而形成分区。
    2. 联结（Ties）：关键创新在于，作者允许不同序列的随机序共享一部分。如果两条序列的随机序完全共享，那么它们的变化点位置就完全一样，从而分区也完全一样。这就是“联结”。
    3. 最简例子：
      • 假设我们只考虑一个变化点的情况（分区B或C）。
      • 对于序列1，我们生成一个随机序，决定变化点是在时间1之后（B）还是时间2之后（C）。
      • 对于序列2，我们生成另一个随机序。
      • 联结机制：我们引入一个概率\(p\)，表示两条序列的随机序是“联结”的（即共享同一个随机序）。以概率\(1-p\)，它们是独立的。
      • 如果联结，那么\(\lambda_1 = \lambda_2\)：要么都是B，要么都是C。
      • 如果不联结，那么\(\lambda_1\)和\(\lambda_2\)可以独立地是B或C。
    4. 这个例子说明了什么：
      • 聚类标准：聚类（即联结）直接对应于“共享变化点”。如果两条序列联结，它们的变化点位置必然相同。
      • 不确定性：\(p\)是一个先验参数，它控制了我们对“两条序列是否属于同一组”的先验信念。后验推断会基于数据\(y_1, y_2\)来更新这个信念。
      • 一般化：这个“随机序 + 联结”的机制可以很容易地推广到\(N\)条序列和\(T\)个时间点的情况。通过定义更复杂的联结结构（如狄利克雷过程混合模型），可以自动推断出有多少个组，以及每组内有哪些序列。

这个最小内核清晰地展示了本文的核心数学思想：通过随机序的共享（联结）来诱导分区之间的共享，从而将“共享变化点”的聚类问题转化为一个贝叶斯非参数模型中的联结结构推断问题。

三、这篇论文做了什么

三句话

  1. 研究了什么问题：提出一种新的基于模型的时间序列聚类方法，其核心假设是：若两条时间序列的结构变化（change points）发生在同一时刻，则它们属于同一组。
  2. 核心工具/方法：利用随机序（random orders） 对每条序列的结构变化点进行潜在表示，并通过诱导不同观测之间的联结（ties） 实现聚类。该策略是通用的，可与多种已知的时间依赖模型（如SIR、AR）结合。
  3. 主要结论：通过模拟研究和COVID-19真实数据应用，展示了该方法能够有效地根据结构变化的同步性对时间序列进行聚类，并提供了有意义的流行病学解释。

关键设定与假设

  • 设定：有\(N\)条时间序列，每条长度为\(T\)。每条序列由一个时间依赖模型生成，该模型的参数在时间上是分段常数（即存在变化点）。目标是基于变化点位置的同步性将序列聚类。
  • 假设：
    1. 变化点同步性假设：同一组内的所有序列共享完全相同的变化点时间集。这是本文的核心假设，也是聚类的标准。这是一个很强的假设，在现实中可能不完全成立（例如，不同国家的疫情政策变化可能相差几天）。
    2. 连续分区假设：每个时间序列的分区\(\lambda_i\)是一个连续的分区（contiguous partition），即每个区块由连续的时间点组成。这是时间序列变化点问题的标准假设。
    3. 模型通用性：本文的方法框架不依赖于特定的时间依赖模型。作者在模拟和实例中使用了SIR模型和AR模型，但声称可以扩展到其他模型。
    4. 贝叶斯框架：所有推断都在贝叶斯框架下进行，需要为分区、模型参数等指定先验分布。
  • 相比已有文献的强化/放宽：
    • 强化：相比ClustSeg（Same et al., 2011）和funFEM（Bouveyron et al., 2015），本文强化了聚类标准，明确要求组内序列共享完全相同的变化点位置，而不是仅仅共享相似的局部行为。
    • 放宽：相比传统的单序列变化点检测方法，本文放宽了“序列独立”的假设，允许不同序列之间通过共享变化点而相互关联。

主要结果

本文是应用导向的，没有提供新的理论收敛性或效率结果。主要结果是方法本身和其在模拟与真实数据上的表现。

  • 方法贡献：提出了一个基于随机序的贝叶斯建模框架，用于“共享变化点”的时间序列聚类。该框架的核心是：

    • 随机序表示：每个时间序列的分区\(\lambda_i\)由一个随机序\(\sigma_i\)和一个截断机制生成。随机序决定了变化点出现的顺序，截断机制决定了有多少个变化点“存活”。
    • 联结机制：不同序列的随机序\(\sigma_i\)通过一个狄利克雷过程（DP） 或Pitman-Yor过程（PYP） 混合模型联结在一起。这意味着，属于同一组的序列共享同一个随机序，从而共享相同的变化点位置。这自动实现了聚类，并允许组数未知。
    • 通用性：该框架可以“包裹”在任何时间依赖模型（如SIR、AR、高斯过程）之上，只需将模型参数\(\theta_{it}\)定义为分区\(\lambda_i\)的函数即可。

  • 模拟研究：

    • 设定：生成\(N=10\)条时间序列，每条长度\(T=50\)。数据由分段常数均值的高斯模型生成。真实分组为3组，每组内序列共享相同的变化点位置。
    • 结果：该方法能够以高概率正确识别出真实的组数和分组。与不进行聚类的单序列变化点检测方法相比，聚类方法显著提高了变化点估计的准确性（特别是对于噪声较大的序列），因为组内序列可以“借用”彼此的信息。
    • 对比：与ClustSeg（Same et al., 2011）进行了比较。ClustSeg倾向于将序列聚成更多、更小的组，而本文的方法更准确地恢复了真实的分组结构。

  • 真实数据应用：COVID-19疫情聚类

    • 数据：欧盟27个国家在2020年3月1日至2021年12月31日期间的每日新增COVID-19病例数。
    • 模型：使用一个随机SIR模型作为时间依赖模型。该模型假设每个国家内部遵循SIR动力学，但感染率\(\beta\)和恢复率\(\gamma\)是分段常数，其变化点由本文的聚类模型决定。
    • 结果：该方法将27个国家聚类成若干组。作者发现，这些组与各国实施非药物干预（NPIs）的时间点高度相关。例如，较早实施严格封锁的国家（如意大利、西班牙）被聚为一组，而较晚实施或实施较宽松的国家（如瑞典）被聚为另一组。这验证了该方法能够捕捉到由政策变化驱动的结构变化同步性。
    • 这个例子想说明什么：展示了该方法在真实、复杂场景下的应用价值，并验证了其聚类结果具有流行病学上的可解释性。

证明路线与技术技巧

本文没有提供任何理论证明（如聚类一致性、变化点估计的收敛速度）。因此，没有“证明路线”可言。技术技巧主要体现在建模和计算上。

  • 建模技巧：

    • 随机序的构造：作者使用一个截断的泊松过程或Beta-Bernoulli过程来生成随机序。具体来说，每个时间点\(t\)（除了最后一个）都有一个“存活”的指示变量\(z_{it}\)，表示该时间点是否是一个变化点。这些\(z_{it}\)的联合分布由一个随机序决定。这种构造方式使得变化点的数量是随机的，且变化点位置是连续的。
    • 联结的诱导：通过将不同序列的随机序\(\sigma_i\)视为来自一个狄利克雷过程混合模型（DPMM） 的样本，实现了聚类。DPMM的基分布是随机序的先验，浓度参数控制了聚类的倾向。这样，属于同一组的序列自然共享同一个随机序，从而共享相同的变化点位置。

  • 计算技巧：

    • MCMC算法：由于模型是贝叶斯的，后验推断通过MCMC进行。作者设计了一个Gibbs采样器，交替采样：
      1. 每个序列的随机序\(\sigma_i\)（或等价地，其分区\(\lambda_i\)）。
      2. 序列之间的联结结构（即分组分配）。
      3. 给定分区下的模型参数\(\Theta\)。
    • 分裂-合并（Split-Merge）MCMC：为了更有效地探索分区空间，作者采用了分裂-合并MCMC（Jain & Neal, 2007）。这种算法允许在一次迭代中同时改变多个序列的分组分配，避免了逐个序列调整的低效性。
    • 点估计：从MCMC样本中，通过最小化Binder损失函数（Wade & Ghahramani, 2018）得到对潜在分区\(\hat{\lambda}\)的点估计。Binder损失函数惩罚了将两个时间点错误地分到同一组或不同组的行为。

🔎 结论是否比证明窄

是的，结论比证明窄。作者在introduction和结论中声称该方法具有“通用性”和“灵活性”，但全文没有任何理论证明来支持这些声称。具体来说： * 没有聚类一致性证明：没有证明随着\(N\)或\(T\)的增长，聚类结果会收敛到真实分组。 * 没有变化点估计一致性证明：没有证明变化点位置的估计是相合的。 * 没有MCMC收敛性证明：虽然使用了MCMC，但没有证明算法收敛到目标后验分布。 * 没有计算复杂度分析：没有分析算法的时间复杂度，特别是对于大规模数据（大\(N\)、大\(T\)）的可扩展性。

因此，本文的结论（“该方法有效”）主要基于模拟和真实数据上的表现，而非严格的数学证明。作者在文中也承认了这一点，将其定位为“应用导向的建模创新”。

四、开放问题

  1. 理论性质：本文方法是否具有聚类一致性和变化点估计一致性？在什么条件下成立？这是最直接的开放问题，扎根于本文完全没有提供任何理论保证这一事实。研究者可以尝试证明，在一定的正则性条件下，后验分布会集中在真实的分组和变化点位置附近。

  2. 计算可扩展性：本文的MCMC算法对于大规模数据（例如，\(N=1000\)条序列，每条\(T=1000\)个时间点）的计算成本如何？能否设计更高效的变分推断（VI）算法或在线算法？扎根于本文没有讨论计算复杂度，且MCMC通常难以扩展到大数据集。

  3. 放松“完全同步”假设：本文假设同一组内的序列共享完全相同的变化点时间集。在现实中，不同序列的变化点可能只相差几天（例如，不同国家实施封锁的时间相差一周）。如何将模型扩展到允许组内序列的变化点存在时间偏移（time shift）？扎根于本文的核心假设（变化点同步性假设） 过于严格。

  4. 与其他聚类标准的比较：本文的聚类标准是“共享变化点”。在什么情况下，这个标准比基于形状相似性（如DTW）或基于模型参数相似性的标准更优？是否存在一个统一的框架，可以同时考虑多种聚类标准？扎根于本文回避了与其他聚类标准的比较，只与ClustSeg进行了有限的对比。

Maintained by 陈星宇 · Homepage · Source on GitHub

评论