| 第600行： |
第600行： |
| | ==谱分解== | | ==谱分解== |
| | | | |
| − | 谱分解（Spectral Clustering）的思想，是利用矩阵的谱，也就是特征向量所展示的物理含义来进行聚类和划分。
| |
| − |
| |
| − | ===最大化有效信息的谱分解方法===
| |
| − |
| |
| − | 该方法先对转移矩阵做特征分解，用特征值和特征向量构造状态间的相似性，再在嵌入空间中聚类，最后据此构造宏观动力学。
| |
| − |
| |
| − | '''输入'''：原始包含 <math>N</math> 个状态的马尔科夫链 <math>P</math> 和距离超参 <math>\epsilon</math>；
| |
| − |
| |
| − | # 对 <math>P</math> 进行矩阵的[[特征值分解]]，得到特征值集合 <math>\Lambda=\{\lambda_i\}^N_{i=1}</math> 与特征向量集合 <math>E=\{e_i\}^N_{i=1}</math>；
| |
| − | # 构建新的有偏特征向量集合 <math>E’=\{\lambda_ie_i|\lambda_i≠0\}^N_{i=1}</math>，该集合中的每个向量都被其对应的特征值加权，同时去除了零特征值及其对应的特征向量（新的状态数量为 <math>N'=rank(P)</math>）。直观地说，忽略特征值为 0 的特征向量是有意义的，因为它对应了动力学中的简并性，而非零特征值及其对应特征向量则包含了丰富的网络拓扑结构信息；
| |
| − | # 依据 <math>E'</math> 计算状态间的相似矩阵 <math>D_{N'×N'}</math>：
| |
| − | ## 如果状态 <math>v_i</math> 和 <math>v_j</math> 分别在对方的马尔科夫毯中，则通过新的特征向量计算两个节点的 cosine 相似性，作为两者间的相似度 <math>d_{ij}</math> 和 <math>d_{ji}</math>；
| |
| − | ## 否则将两个节点间的相似度设为无穷大 <math>\infty</math>（可以设个比较大的值，如10000）；
| |
| − | # 基于相似度矩阵 <math>D_{N'×N'}</math> 和一个超参 <math>\epsilon</math>（需要线性搜索，可以选择EI最大的参数），使用[[OPTICS]]算法（是一种基于密度的划分算法，旨在识别数据集中不同密度的划分结构）进行划分，输出对应超参 <math>\epsilon</math> 下的划分方式。
| |
| − |
| |
| − | 该方案的时间复杂度为<math>O(N^3)</math>。
| |
| − |
| |
| − | '''缺点'''：
| |
| − | # 对大规模动力学进行粗粒化时，对转移矩阵进行特征值分解的计算复杂度也较高
| |
| − | # 聚类算法所需的距离超参 <math>\epsilon</math> 需要线性搜索，需要加入人为的先验知识
| |
| − |
| |
| − | <br/>
| |
| − |
| |
| − | ===寻找成块划分的谱分解方法===
| |
| − |
| |
| − | 当马尔科夫矩阵存在块状结构，或者状态明显可以被划分为几类时，根据这样的划分，该矩阵就是成块的。对于简单的 <math>P</math>，如图6中的 <math>\bar{P}</math> 所示，我们能够很快找出它的成块划分 <math>\Phi</math>。然而在实际应用中，我们经常遇到的矩阵并不完美。例如，成块矩阵的状态排序可能被打乱（如图（6）中的 <math>P_1</math>），或者矩阵中混入了如 <math>P_2</math> 这样的噪声（如图6中的 <math>P</math>，<math>P = P_1 + P_2</math>，<math>P_1^TP_2 = 0</math>）。在这种情况下，我们通常面临的是像 <math>P</math> 这样的矩阵：既无法确定它是否成块，也无法决定它的划分，甚至不知道它的马尔科夫秩。
| |
| − |
| |
| − | [[文件:Lump fig1.png|缩略图|398x398像素|图6：Zhang<ref name=":4" /> 文章中的示意图。图中左面四个矩阵都是成块的马尔科夫矩阵，而右面的P_2是一个噪声矩阵，(P_1)^T P_2 = 0|替代=]]
| |
| − |
| |
| − | 针对这一问题，Anru Zhang<ref name=":4">Zhang, Anru, and Mengdi Wang. "Spectral state compression of markov processes." ''IEEE transactions on information theory'' 66.5 (2019): 3202-3231.</ref><ref name=":4" />提出了一种寻找最优划分的方法，其核心是利用'''稳态流矩阵'''（frequency matrix）<math>F = \mathrm{diag}(\pi) \cdot P</math> 的奇异值结构来恢复成块划分。
| |
| − | 其中 <math>\pi</math> 是 <math>P</math> 的稳态分布，<math>F_{ij} = \pi_i P_{ij}</math> 可以理解为系统在稳态下从状态 <math>i</math> 流向状态 <math>j</math> 的长期概率流量。具体步骤如下：
| |
| − |
| |
| − | # 获取一个 <math>n\times n</math> 维的马尔科夫矩阵 <math>P</math>，或者马尔科夫链的频率采样，并由此得到频率矩阵 <math>F</math>；
| |
| − | # 获取马尔科夫秩 <math>r<n</math>，或指定一个 <math>r</math>；
| |
| − | # 对 <math>F</math> 进行 SVD 分解，<math>F = U \Sigma V^T</math>，其中 <math>U</math> 为左奇异向量，<math>V</math> 为右奇异向量，<math>\Sigma = diag(\sigma_{1},\sigma_{2},...,\sigma_{N})</math> 是一个对角矩阵，包含所有有序奇异值；
| |
| − | # 通过下列公式得到最优划分 <math>\Omega_1, \Omega_2, ... , \Omega_r</math>
| |
| − |
| |
| − | <math>
| |
| − | \Omega_1, \Omega_2, ... , \Omega_r = argmin_{\Omega_1, \Omega_2, ... , \Omega_r} min_{v_1, v_2, ... , v_r} \sum_{s=1}^r \sum_{i \in \Omega_s} || V_{[i,:r]} - v_s ||_2^2
| |
| − | </math>
| |
| − |
| |
| − | 其中，<math>v_s</math>是第<math>s</math>类<math>\Omega_s</math>的特征向量。这个算法的核心思想是：在最优的划分中，<math>\Omega_s</math>中的状态<math>i</math>的右奇异向量会和<math>v_s</math>距离最小，也就是说，让每个微观态的右特征向量和它对应的宏观态的右特征向量尽可能接近。
| |
| − |
| |
| − | 这个算法看起来非常熟悉，实际上它与 K-Means 聚类算法非常相似。让<math>v_s</math>和组里的<math>V_{[i,:r]}</math>距离最小的方法，就是让<math>v_s</math>成为这些点的中心。回顾一下 K-Means 算法把 n 个点聚成 r 类的目标函数，其中第<math>i</math>个点被归到第<math>c^i</math>类，而<math>\mu_j</math>为第<math>j</math>类点的中心点：
| |
| − |
| |
| − | <math>J(c^1, ... , c^n, \mu_1, ... , \mu_r) = \sum_{i=1}^n || x_i - \mu_{c^i} ||^2 </math>
| |
| − |
| |
| − | 我们看到这两个目标函数非常相似。因此，我们可以通过上述算法或 k-Means 聚类来获取最优的状态划分，进而判断马尔科夫矩阵的成块性。
| |
| − |
| |
| − |
| |
| − |
| |
| − | 为什么是对 <math>F</math> 做 SVD 而不是对 <math>P</math>？
| |
| − | 论文证明了，当马尔科夫链关于划分 <math>\Omega_1, \ldots, \Omega_r</math> 是成块的时，<math>F</math> 存在分解 <math>F = F_1 + F_2</math>，其中 <math>F_1</math> 是秩为 <math>r</math> 的成块部分，<math>F_2</math> 是组内局部动力学部分，两者行空间正交（<math>F_1 F_2^T = 0</math>）。<math>F_1</math> 的右奇异向量 <math>V_{F_1}</math> 具有'''分段常数'''结构：同一组 <math>\Omega_k</math> 内的所有状态，其行向量 <math>(V_{F_1})_{[i,:]}</math> 完全相同，不同组之间则有明确间距。因此，K-Means 可以精确恢复分组。
| |
| − |
| |
| − | 从物理上看，<math>F</math> 编码的是'''稳态下的概率流量结构'''，而 <math>P</math> 编码的是'''条件转移分布'''。对 <math>F</math> 做 SVD 相当于在稳态概率流网络中寻找主要的流量通道结构——低频状态对应的流量通道自然被压低权重，不会干扰大尺度结构的识别。
| |
| − |
| |
| − | 当系统存在明显的近成块结构时，<math>F</math> 的前 <math>r</math> 个奇异值主要反映的是组间的大尺度概率流通道，而后续的奇异值反映的是组内的局部流量波动。如果这两部分之间存在明显的谱间隙（singular value gap），即 <math>\sigma_r(F_1) \gg |F_2|</math>，那么 K-Means 聚类 <math>\hat{V}_F</math> 就能准确恢复成块划分。
| |
| − |
| |
| − | 这对应于一种特殊但常见的物理场景：系统各状态群之间的长期概率流量远弱于群内流量，形成类似“准吸引域”的结构，粒子长时间在群内游走，偶尔才跨越群边界。当 <math>P</math> 不是精确成块，而是近似成块（near-lumpable）时，该方法仍然可以工作，但需要满足一个关键条件：'''组间动力学 <math>F_1</math> 的奇异值必须主导 <math>F</math> 的整体奇异值谱'''，即 <math>\sigma_r(F_1) \gg \|F_2\|</math>。
| |
| − |
| |
| − | 需要注意的是，成块性的数学定义本身并'''不要求准吸引域结构（群内流量强于群间流量）'''。成块性只要求同组状态到各组的转移概率总量相同，并不约束转移的方向强弱。论文给出了分组错误率的上界：
| |
| − |
| |
| − | <math>
| |
| − | \text{Misclassification rate} \leq C \frac{\pi_{\max} r \tau^* \log^2(n)/n + (r\|F_2\|^2) \wedge \|F_2\|_F^2}{\sigma_r^2(F_1)}
| |
| − | </math>
| |
| − |
| |
| − | 分子中第一项是有限样本的统计误差，第二项是近似成块带来的结构偏差，分母是组间动力学的信号强度。因此，当系统虽然满足成块性、但组间转移概率反而很大时（即 <math>|F_2| \gtrsim \sigma_r(F_1)</math>），SVD 提取的前 <math>r</math> 个方向会反映组内的局部动力学模式而非组间的成块结构，此时该方法会给出错误的分组。
| |
| − |
| |
| − | <br/>
| |
| − |
| |
| − | ===基于动力学可逆性的谱分解方法===
| |
| − |
| |
| − | 张江等人关于动力学可逆性的文章<ref name="reversibility"/>提供了另一种基于 SVD 的粗粒化方法。该工作的出发点与上面寻找成块划分的方法不同：它不以恢复某种预设的结构性划分为目标，而是把马尔科夫动力学看作一个通信信道，把每个状态下的概率转移看作一条信息通路（information pathway）。当系统的 <math>N</math> 条通路中存在线性相关的行向量时，这些通路就是冗余的，它们要么传递极少的信息（对应接近零的奇异值），要么可以被其他通路的线性组合表示。
| |
| − |
| |
| − | 通过 SVD 可以识别出 <math>r_\epsilon</math> 条“核心通路”，其余 <math>N - r_\epsilon</math> 条则是可以消除的冗余通路。粗粒化就是把状态按照它们在核心通路上的表现进行归组，让宏观动力学只保留这些高效通道。其基本思想是将 <math>P</math> 中的行向量 <math>P_{i},\forall i \in [1,N]</math> 投影到 <math>P</math> 的几个较大奇异值对应的奇异向量所张成的子空间上，从而保留 <math>P</math> 的主要信息，并尽量保持 <math>\Gamma</math> 不变。
| |
| − |
| |
| − | 与 EI 相比，<math>\Gamma_\alpha</math> 能捕捉到 EI 无法感知的信息：EI 基于 KL 散度，只能度量行向量之间的'''相似性'''；而 <math>\Gamma_\alpha</math> 基于线性代数，还能度量行向量之间的'''线性依赖性'''。两个行向量可以彼此很不像（KL 散度大），但仍然线性相关（一个是另一个的线性组合）。后一种情况下，EI 看不出冗余，但 <math>\Gamma_\alpha</math> 能通过秩亏反映出来。
| |
| − |
| |
| − |
| |
| − |
| |
| − | 算法具体步骤如下：
| |
| − |
| |
| − | # 给定含有 <math>N</math> 个状态的马尔科夫矩阵 <math>P</math>，对其做 SVD 分解 <math>P = U \Sigma V^T</math>；
| |
| − | # 根据奇异值谱设定阈值 <math>\epsilon</math>，取前 <math>r_\epsilon</math> 个右奇异向量，构成 <math>V_{N \times r_\epsilon}</math>；
| |
| − | # 对 <math>P</math> 的行向量做 PCA 降维：<math>\tilde{P} = P \cdot V_{N \times r_\epsilon}</math>，得到每个状态在核心通路子空间上的投影；这一步的物理含义是：<math>P^T P</math> 是转移概率矩阵的自相关矩阵，其特征向量 <math>V</math> 给出了“哪些目的地模式是相关联的”这一信息的主轴。<math>\tilde{P}_i = P_i \cdot V_{N \times r_\epsilon}</math> 是状态 <math>i</math> 的出转移分布在这些主要模式上的投影，保留了该状态“在核心通路上的行为”，滤去了冗余通路上的波动。
| |
| − | # 对 <math>\tilde{P}</math> 的行向量用 K-Means 聚类为 <math>r_\epsilon</math> 组，得到划分 <math>\Phi</math>；
| |
| − | # 基于划分 <math>\Phi</math> 和稳态流矩阵 <math>F = \mathrm{diag}(\pi) P</math>，通过 <math>F' = \Phi^T F \Phi</math> 构建宏观稳态流矩阵，再对每行归一化得到宏观转移矩阵 <math>P'</math>。
| |
| − |
| |
| − | 最后用 <math>F' = \Phi^T F \Phi</math> 而非直接对 <math>P</math> 做行列聚合，可以保证粗粒化操作和时间演化算子具有可交换性，从而确保宏观动力学仍是合法的马尔科夫链。
| |
| − |
| |
| − |
| |
| − |
| |
| − | 虽然该工作与前面 Zhang & Wang 的方法都包含 SVD + K-Means 这一步，但两者之间仍存在差异：
| |
| − |
| |
| − | # Zhang & Wang的工作目标是在粗粒化过程中尽量保持马尔科夫性，并恢复成块结构；而动力学可逆性工作的目标是在粗粒化过程中减少 TPM 中的信道冗余性。
| |
| − | # Zhang & Wang方法的 SVD 对象是 <math>V_F</math>，稳态流矩阵的右奇异向量，而动力学可逆性的SVD对象是 <math>V_P</math>，马尔科夫转移矩阵的右奇异向量。当 <math>P</math> 恰好是低秩的时，<math>V_F</math> 和 <math>V_P</math> 张成完全相同的子空间（因为 <math>F = \mathrm{diag}(\pi) P</math>，两者秩相同时列空间一致）。当 <math>P</math> 是近似低秩（满秩但有奇异值谱间隔）时，稳态分布 <math>\pi</math> 的加权效应会使两者产生差异：对 <math>F</math> 做 SVD 天然压低了低频状态的影响，而直接对 <math>P</math> 做 SVD 则给所有状态平权。
| |
| − |
| |
| − | <br/>
| |
| | 谱分解（Spectral Clustering）的思想，是利用矩阵的谱，也就是特征向量所展示的物理含义来进行聚类和划分。 | | 谱分解（Spectral Clustering）的思想，是利用矩阵的谱，也就是特征向量所展示的物理含义来进行聚类和划分。 |
| | | | |