[论文解读] Differentiated uniformization: A new method for inferring Markov chains on combinatorial state spaces including stochastic epidemic models

Kevin Rupp, Rudolf Schill+7 more|arXiv (Cornell University)|Jan 1, 2021

Complex Network Analysis Techniques被引用 2

一句话总结

本文提出了一种新型方法——微分统一化（differentiating uniformization），通过利用转移率矩阵 Q 中的张量积结构，高效计算连续时间马尔可夫链（CTMC）中矩阵指数及其导数，适用于具有组合状态空间的模型（如随机 SIR 流行病模型）。该方法实现了参数（如感染率与恢复率）的可扩展、数值稳定的推断，通过奥地利首波 COVID-19 大流行的完整贝叶斯分析得到验证，揭示了时变传播动力学及不确定性量化结果。

ABSTRACT

Motivation: We consider continuous-time Markov chains that describe the stochastic evolution of a dynamical system by a transition-rate matrix $Q$ which depends on a parameter $θ$. Computing the probability distribution over states at time $t$ requires the matrix exponential $\exp(tQ)$, and inferring $θ$ from data requires its derivative $\partial\exp\!(tQ)/\partialθ$. Both are challenging to compute when the state space and hence the size of $Q$ is huge. This can happen when the state space consists of all combinations of the values of several interacting discrete variables. Often it is even impossible to store $Q$. However, when $Q$ can be written as a sum of tensor products, computing $\exp(tQ)$ becomes feasible by the uniformization method, which does not require explicit storage of $Q$. Results: Here we provide an analogous algorithm for computing $\partial\exp\!(tQ)/\partialθ$, the differentiated uniformization method. We demonstrate our algorithm for the stochastic SIR model of epidemic spread, for which we show that $Q$ can be written as a sum of tensor products. We estimate monthly infection and recovery rates during the first wave of the COVID-19 pandemic in Austria and quantify their uncertainty in a full Bayesian analysis. Availability: Implementation and data are available at https://github.com/spang-lab/TenSIR.

研究动机与目标

• 解决高维 CTMC 中矩阵指数及其导数计算不可行的问题，尤其针对具有组合状态空间的模型。
• 克服在如随机 SIR 流行病模型中存储和计算巨大转移率矩阵 Q 所带来的挑战。
• 开发一种无需依赖状态演化单调性的数值稳定、可扩展的 CTMC 参数推断方法。
• 实现对流行病参数（如感染率与恢复率）的完整贝叶斯推断，包括不确定性量化。
• 将统一化方法扩展至参数导数计算，从而实现对出生过程以外模型的推断。

提出的方法

• 将转移率矩阵 Q 表示为张量积之和，从而在不显式存储 Q 的情况下实现矩阵-向量乘法。
• 应用统一化方法，通过带指数衰减的幂级数展开计算矩阵指数 exp(tQ)。
• 推导微分统一化算法，以计算矩阵指数关于模型参数 θ 的导数。
• 利用张量结构在级数展开的每一步高效计算 P^n p(0)，避免完整的矩阵运算。
• 采用哈密顿蒙特卡洛采样实现对时变参数的完整贝叶斯推断。
• 原则上利用低秩张量格式，将可扩展性扩展至包含数百个突变或隔室的模型。

实验结果

研究问题

• RQ1 在 Q 大到无法存储的组合状态空间 CTMC 中，能否高效计算矩阵指数及其导数？
• RQ2 将 Q 表示为张量积之和，是否能实现瞬态分布及其参数敏感性的可扩展计算？
• RQ3 微分统一化方法在通用性方面是否优于现有方法，尤其在 SIR 模型等非单调过程中的表现？
• RQ4 当状态空间呈指数级增长时，如何对时变流行病参数执行完整贝叶斯推断？
• RQ5 该方法能否扩展至其他具有类似组合复杂度的模型，如肿瘤进展或捕食者-猎物动力学？

主要发现

• 微分统一化方法的计算复杂度随种群规模呈立方增长，尽管状态空间呈指数增长，但该方法在大规模 SIR 模型中仍具实用性。
• 该方法实现了奥地利首波 COVID-19 大流行期间每月感染率与恢复率的完整贝叶斯推断。
• 后验估计显示，2020 年 4 月至 8 月期间，恢复率约为每日 0.07（平均恢复时间约 2 周）。
• 3 月的估计表明恢复率较低，约为每日 0.03（一个月均值），可能由于恢复病例报告延迟所致。
• 该方法揭示了数月间参数估计存在显著不确定性，提示模型与数据可能存在拟合不足或参数非恒定。
• 该方法可推广至其他具有张量结构 Q 的 CTMC，包括肿瘤进展和化学反应网络。