马尔科夫链的粗粒化

来自集智百科 - 复杂系统|人工智能|复杂科学|复杂网络|自组织
跳到导航 跳到搜索

马尔科夫链是一类特殊的随机过程，在给定当前状态下，对于未来状态的预测便不再依赖于过去状态。它可以看作一个n维的状态的序列[math]\displaystyle{ s^{(t)}\ = \{1, ..., n\}^{(t)} }[/math]，每一步的状态转换都有马尔科夫矩阵（满足每一行和为1的条件的方阵）[math]\displaystyle{ P }[/math] 决定，即[math]\displaystyle{ p(s^{(t+1)}) = p(s^{(t)}) P }[/math]。[math]\displaystyle{ P }[/math] 的每一行对应的每个状态转移到其他状态的概率，即 [math]\displaystyle{ p_{ij} = p(s^{(t+1)} = j | s^{(t)} = i) }[/math]。比如，当 [math]\displaystyle{ s^{(t)} }[/math] 等于第一个状态的时候，[math]\displaystyle{ P }[/math] 的第一行展示了 [math]\displaystyle{ s^{(t+1)} }[/math] 状态的概率。

对马尔科夫链做粗粒化的意义是什么呢？我们看到文献中着重强调这三点：

  1. 我们在研究一个超大的系统的时候，比如复杂城市系统时，并不会关注每一个微观的状态，在粗粒化中我们希望能过滤掉一些我们不感兴趣的噪声和异质性，而从微观尺度中总结一些中尺度或宏观尺度的规律；
  2. 有些状态的转移概率非常相似，所以可以被看成同一类状态，对这种马尔科夫链做partitioning可以减少系统表示的冗余性；
  3. 在用到马尔科夫决策过程的强化学习里，对马尔科夫链做粗粒化可以减少状态空间的大小，提高训练过程的样本效率和任务表现。[1]

在粗粒化过程中，我们有很多目标。

  1. 首先，我们想尽量降低压缩后的动力学维度；
  2. 其次，我们想在压缩过程中尽量保留动力学中的信息；
  3. 最后，我们想要压缩后的过程和原本的过程尽量保持一致。

在这篇文章里，我们会围绕着这三个维度展开讨论，并介绍一系列针对马尔科夫链的粗粒化方法。


粗粒化

马尔科夫链不只是一个马尔科夫矩阵那么简单，还涉及到了不同的状态，也会有状态之间的转移概率，所以我们可以从不同的角度来看待粗粒化：

  1. 对马尔科夫矩阵做粗粒化：在看待一个方阵时，我们可以把它想象成一个网络的邻接矩阵。节点为状态，连边为状态转移概率，状态转移的过程为一个粒子在节点间跳转，而粗粒化则是对这个网络做降维。
  2. 对状态空间做粗粒化：我们这里讨论的是离散状态的马尔科夫链，状态空间的大小等于马尔科夫链的大小。粗粒化可以看作对状态空间做分组并投影到新的状态空间。
  3. 对概率空间做粗粒化：从概率论出发，把状态的发生看作事件，通过计算事件发生的频率来构建概率空间。在降维表达事件的同时也构建了新的概率空间。

这三种出发点看似互有联系，但是具体关心的东西不同。第一种关注状态之间的连边拓扑结构，第二种关注不同状态的历史分布和相似性，第三种关注如何降维表达整个系统中的事件概率分布。


对马尔科夫链做粗粒化并没有一个统一的定义，目前该题目也没有一个完整的文献综述。在许多文献中，粗粒化coarse-graining，聚类/聚合partitioning/lumping/clustering/aggregation和降维dimension reduction是重叠等价的。

本词条只讨论离散时间马尔科夫链，且动力学[math]\displaystyle{ P }[/math]是时齐（不随时间变化）的，不讨论连续马尔科夫过程。我们可以定义微观状态空间为[math]\displaystyle{ S=\{s_1, s_2, ... ,s_n\} }[/math]，微观转移函数[math]\displaystyle{ P \in \mathcal{R}_{n \times n} }[/math]是一个大小为n的马尔科夫矩阵。同样的，粗粒化后的宏观状态空间为[math]\displaystyle{ S'=\{s'_1, s'_2, ... ,s'_r\} }[/math]，宏观转移函数[math]\displaystyle{ P' \in \mathcal{R}_{r \times r} }[/math]。 多数文献中会出现如图1的定义[2]：

对状态空间[math]\displaystyle{ S }[/math]和转移函数[math]\displaystyle{ P(t):S \rightarrow S }[/math]，我们需要一个投影函数[math]\displaystyle{ \Phi:S \rightarrow S' }[/math]，和一个新的转移函数[math]\displaystyle{ P'(t):S' \rightarrow S' }[/math]，而且[math]\displaystyle{ S' }[/math]和[math]\displaystyle{ P' }[/math]要符合马尔科夫链定义。

我们能看出[math]\displaystyle{ \Phi }[/math]是一个对状态空间[math]\displaystyle{ S }[/math]的分组函数（得到[math]\displaystyle{ S' }[/math]）；而且，我们还需要把[math]\displaystyle{ P }[/math]降维成[math]\displaystyle{ P' }[/math]。如果我们对这个降维/粗粒化过程进行分解的话，我们可以看到这个过程可以被分为三步。

  1. 划分（对状态空间进行粗粒化）：对状态进行分组，决定[math]\displaystyle{ S }[/math]空间中的状态要如何被投影到[math]\displaystyle{ S' }[/math]空间上，也就是得出[math]\displaystyle{ \Phi }[/math]。
  2. 聚合（对转移概率进行粗粒化）：基于分组函数[math]\displaystyle{ \Phi }[/math]，我们要决定如何把原来的状态转移矩阵[math]\displaystyle{ P }[/math]降维成新的[math]\displaystyle{ P' }[/math]。
  3. 判断粗粒化方案[math]\displaystyle{ \Phi }[/math]和[math]\displaystyle{ P' }[/math]符合什么样的性质。


划分：对状态空间进行粗粒化

我们假设投影函数[math]\displaystyle{ \Phi }[/math]是一个[math]\displaystyle{ n \times r }[/math]的线性投影矩阵，且符合下列两个条件：(1)所有元素非负；(2)每行加和等于1。

马尔科夫链的粗粒化可以分成硬划分(Hard partitioning)和软划分(Soft partitioning)。硬划分是指把若干个微观状态分成一个宏观状态类，且一个微观状态不能同时属于多个宏观状态类，而软划分则会有可能出现一个微观状态被分在属于多个宏观类的情况。

在硬划分的情况下，[math]\displaystyle{ \Phi_{i, j} = \begin{cases} 1\text{, if}\ s_i \in s_j', \\ 0\text{, otherwise} \end{cases} }[/math][3]，即一个表示[math]\displaystyle{ i }[/math]属于[math]\displaystyle{ j }[/math]列的零一分组矩阵。如果是软划分的话，[math]\displaystyle{ \Phi }[/math]则不一定会是零一矩阵。

在这个词条里，我们主要讨论硬划分[4]。对软划分感兴趣的朋友可以自行参阅这两篇论文[5] [6]。


聚合：对转移概率进行粗粒化

得出[math]\displaystyle{ \Phi }[/math]后，对微观马尔科夫矩阵[math]\displaystyle{ P }[/math]的降维可以表示为

[math]\displaystyle{ \begin{aligned} P' = W P \Phi, \end{aligned} }[/math]

 

 

 

 

(1)

其中[math]\displaystyle{ W \in \mathbb{R}^{r \times n} }[/math]。Buchholz[3]的文章中把[math]\displaystyle{ W }[/math]称作distributor matrix，而[math]\displaystyle{ \Phi \in \mathbb{R}^{n \times r} }[/math]称作collector matrix。

我们假设了[math]\displaystyle{ W, \Phi, P, P' }[/math]都是线性的，所以我们可以看到这个公式描述了图1中从左上角的t时刻宏观状态出发，经过逆向粗粒化[math]\displaystyle{ W }[/math]得到了t时刻微观状态，走微观动力学[math]\displaystyle{ P }[/math]得到了t+1时刻的微观状态，再经粗粒化[math]\displaystyle{ \Phi }[/math]得到t+1时刻的宏观状态。

而

[math]\displaystyle{ \begin{aligned} W = D^{-1}(\Phi)^T, D = diag(\alpha \Phi), \end{aligned} }[/math]

 

 

 

 

(2)

其中[math]\displaystyle{ \alpha }[/math]是一个非负对角矩阵。

所以，我们能看到，对马尔科夫链做粗粒化的时候，不仅要对状态空间做降维([math]\displaystyle{ \Phi }[/math])，还要对马尔科夫矩阵都要做降维，当中也包含了对概率空间的降维([math]\displaystyle{ W }[/math])。 [math]\displaystyle{ W }[/math] 式子1中我们发现，[math]\displaystyle{ \Phi }[/math]和[math]\displaystyle{ P }[/math]都是给定的，所以我们在这一步要做的是决定[math]\displaystyle{ W }[/math]。[math]\displaystyle{ W }[/math]的基础限制是要使得[math]\displaystyle{ P' }[/math]符合马尔科夫矩阵定义，所以有无限多种[math]\displaystyle{ W }[/math]的选项。从这些选项中，我们也能挑出某些符合我们想要的特征的[math]\displaystyle{ W }[/math]。


判断粗粒化方案符合的性质

最基本的粗粒化方法的约束是要确保粗粒化后的新的[math]\displaystyle{ S' }[/math]和[math]\displaystyle{ T' }[/math]要符合马尔科夫链定义。除此之外也有一些其他的规则，比如：

宏观最大预测性

宏观最大预测性是指，在将微观状态粗粒化为宏观状态后，我们只用宏观状态就能预测后续的宏观状态。

具体来说，只需要知道宏观状态[math]\displaystyle{ A^{(t)} }[/math]，就可以最大化[math]\displaystyle{ A^{(t+1)} }[/math]的预测效果，而微观状态[math]\displaystyle{ s^{(t+1)} }[/math]虽然能提供一些额外信息，但是对于预测[math]\displaystyle{ A^{(t+1)} }[/math]并无帮助。这意味着条件互信息[math]\displaystyle{ I(A^{(t+1)}; s^{(t+1)} | A^{(t)}) = 0 }[/math]。宏观动力学不需要任何来自微观状态的修正，就能独立且正确地反映系统在宏观层面上的演化。


微观最大预测性

微观最大预测性是指，在将微观状态粗粒化为宏观状态后，我们只用宏观状态就能预测后续的微观状态。

具体来说，只需要知道宏观状态[math]\displaystyle{ A^{(t)} }[/math]，就可以最大化[math]\displaystyle{ s^{(t+1)} }[/math]的预测效果。


粗粒化的交换律commutative[2]

图（1）展示了从某个微观状态[math]\displaystyle{ s^{(t)} }[/math]出发，到达下一刻的宏观状态[math]\displaystyle{ A^{(t+1)} }[/math]的两条路径[7]：

  1. 先经过微观动力学达到[math]\displaystyle{ s^{(t+1)} }[/math]，再进行粗粒化得到[math]\displaystyle{ A^{(t+1)} }[/math]
  2. 先进行粗粒化得到[math]\displaystyle{ A^{(t)} }[/math]，再经过宏观动力学到达[math]\displaystyle{ A^{(t+1)} }[/math]

如果我们把粗粒化和（宏观/微观）动力学看作是两个算子，那这两个算子的可交换性指的是：无论从哪个微观状态分布出发，走这两条路径得到的[math]\displaystyle{ A^{(t+1)} }[/math]的分布是一样的，即，[math]\displaystyle{ s^{(t)} \Phi P' = s^{(t)} P \Phi }[/math]。


动力学一致性[2]

好的粗粒化方案是能够尽可能保证原始的网络（或转移矩阵）和粗粒化后的网络（或转移矩阵）尽可能地相似，动力学一致性衡量了在原始网络和粗粒化后网络上演化的状态的轨迹的一致性。详情请见复杂网络中的因果涌现词条。

以下关于动力学的一致性检验的描述节选自复杂网络中的因果涌现：

定义微观网络[math]\displaystyle{ G }[/math]与宏观网络[math]\displaystyle{ G_M }[/math]，分别在两个网络上进行随机游走。在某个时间[math]\displaystyle{ t }[/math]，[math]\displaystyle{ G }[/math]上的预期分布为[math]\displaystyle{ P_m(t) }[/math]，而[math]\displaystyle{ G_M }[/math]上的预期分布为 [math]\displaystyle{ P_M(t) }[/math]。

将[math]\displaystyle{ P_m(t) }[/math]分布叠加到宏观宏观[math]\displaystyle{ G_M }[/math]的相同节点上，得到[math]\displaystyle{ P_{M|m}(t) }[/math]分布。通过计算[math]\displaystyle{ P_M(t) }[/math]和[math]\displaystyle{ P_{M|m}(t) }[/math]之间的KL散度来衡量其不一致性（inconsistency）。如果结果为零，则表示动力学一致。公式如下：

[math]\displaystyle{ inconsistency=\sum_{t=0}^T D_{KL}[P_M(t)||P_{M|m}(t)] }[/math]

也就是说，不一致性 表示了 宏观和微观网络节点在不同时刻的概率分布的KL散度之和。


注意，这是我们从不同来源[2]中总结出来的规则。目前我们还没找到一个严格的完备的马尔科夫粗粒化应该符合的定义。

我们下面会介绍不同的概念，它们分别满足上述的不同指标。


划分标准

基础划分

马尔科夫状态中有一些基础的定义[8]，可以帮助我们对状态进行划分。

滑过Transient和常返Recurrent

滑过Transient状态是指从该状态出发后，经过任意n步也不存在回到该状态的可能路径。

常返Recurrent状态是指从该状态出发后，经过n步后能够回到该状态。

我们可以定义一个常返状态组为：组中任意两个状态都互相可达。

所以，我们就能把一条马尔科夫链分解成多个常返状态组和若干个滑过态。


周期性

如果一个常返状态组存在步数[math]\displaystyle{ n\gt 1 }[/math]，且每个常返状态[math]\displaystyle{ i }[/math]只会经过[math]\displaystyle{ n }[/math]步才会回到[math]\displaystyle{ i }[/math]，那么就可以把这个常返状态组分解成[math]\displaystyle{ n }[/math]个子集。这里我们借用了知乎上的例子[8]进行说明：

假设一个马尔科夫链的状态按照这样的周期（步数[math]\displaystyle{ n=3 }[/math]）循环，

[math]\displaystyle{ \{s_1, s_2\} \rightarrow \{s_3, s_4, s_5\} \rightarrow \{s_6\} \rightarrow \{s_1, s_2\} }[/math]

我们就可以把状态划分成[math]\displaystyle{ \{s_1, s_2\} }[/math]，[math]\displaystyle{ \{s_3, s_4, s_5\} }[/math]，[math]\displaystyle{ \{s_6\} }[/math]三个子集。


成块性

马尔科夫链的成块性（Lumpability）是一种描述马尔科夫链是否可聚类或可粗粒化性质的概念，最早由Kemeny和Snell在其1969年著作Finite Markov Chains[9]中阐述。Lump一词有‘块’的意思。判断一个马尔科夫链是否可成块（lumpable），就是要判定是否存在一种状态划分方式，能够将原始的马尔科夫链压缩为一个保留其动力学特征的简化链。这种划分需确保：虽然抹除了聚类内部状态的微观差异，但能精确保持块间转移机制。这一特性为状态空间的合理划分提供了严格的数学判据。

符合成块性的粗粒化方案有很多很好的性质，如保证了对宏观层面上的最大预测性、保证了微观与宏观过程的一致性等。这使得成块性成为马尔科夫链维度约简研究中最常用的标准之一。

成块性的初始定义比较复杂，详见马尔科夫链的成块性。在这里，我们介绍其定义的充分必要条件，一个更直观的解析。在某些文献中，作者们会直接使用这个通俗版的解释，而且这两者之间是存在充分必要条件关系的[9]。

成块性的充分必要条件

设[math]\displaystyle{ p_{s_k \rightarrow s_m} = p(s^{(t)} = s_m | s^{(t-1)} = s_k) }[/math]，[math]\displaystyle{ p_{s_k \rightarrow A_i} = p(s^{(t)} \in A_i | s^{(t-1)} = s_k) = \sum_{s_m \in A_i} p_{s_k \rightarrow s_m} }[/math]，[math]\displaystyle{ p_{A_i \rightarrow A_j} = p(A^{(t)} = A_j | A^{(t-1)} = A_i) }[/math]，作者[9]提出了判断一个马尔科夫链对给定划分 [math]\displaystyle{ A }[/math] 是否成块的充分必要条件为：

定理1：

给定划分 [math]\displaystyle{ A }[/math]是成块的充分必要条件为：

对于任意一对[math]\displaystyle{ A_i, A_j }[/math]，每一个属于[math]\displaystyle{ A_i }[/math]的状态[math]\displaystyle{ s_k }[/math]的[math]\displaystyle{ p_{s_k \rightarrow A_j} }[/math]都是一样的。也就是说，[math]\displaystyle{ A_i }[/math]里的每一个状态[math]\displaystyle{ s_k }[/math]转移到[math]\displaystyle{ A_j }[/math]中的所有状态的转移概率之和是相同的，

[math]\displaystyle{ \forall i, j, \forall s_k, s_m \in A_i, p_{s_k \rightarrow A_j} = p_{s_m \rightarrow A_j} }[/math]

 

 

 

 

(3)

满足这个条件后，我们就可以将[math]\displaystyle{ A_j }[/math]中的所有微观状态[math]\displaystyle{ s_l }[/math]合并为一个宏观状态[math]\displaystyle{ A_j }[/math]，因为它们从其他宏观状态的转入概率是一致的。同样的，我们也可以将[math]\displaystyle{ A_i }[/math]中的所有微观状态[math]\displaystyle{ s_k }[/math]合并为一个宏观状态[math]\displaystyle{ A_i }[/math]，因为这些微观状态[math]\displaystyle{ s_k }[/math]（和[math]\displaystyle{ A_i }[/math]）到其他宏观状态的转出概率是一致的。

例如，在下图中，[math]\displaystyle{ s_1 }[/math]和[math]\displaystyle{ s_2 }[/math]都属于[math]\displaystyle{ A_1 }[/math]，它们到[math]\displaystyle{ A_1 }[/math]和[math]\displaystyle{ A_1 }[/math]的转移概率是一样的；同样，[math]\displaystyle{ s_3 }[/math]和[math]\displaystyle{ s_4 }[/math]都属于[math]\displaystyle{ A_2 }[/math]，它们到宏观状态的转移概率也是一样的。因此，[math]\displaystyle{ [[1,2], [3, 4]] }[/math]是[math]\displaystyle{ P }[/math]的一个成块划分。

成块性充分必要条件实例


基于这个条件，我们可以定义成块划分的宏观动力学[math]\displaystyle{ A^{(t)} = A^{(t-1)} P' }[/math]：

定义1：成块性的宏观动力学

对于给定的微观状态[math]\displaystyle{ S }[/math]，微观动力学[math]\displaystyle{ p_{s_k \rightarrow s_m} }[/math]和成块的划分 [math]\displaystyle{ A }[/math]，宏观动力学[math]\displaystyle{ p_{A_i \rightarrow A_j} }[/math]等于[math]\displaystyle{ A_i }[/math]中任意状态[math]\displaystyle{ s_k }[/math]到[math]\displaystyle{ A_j }[/math]的概率：

[math]\displaystyle{ \begin{aligned} p_{A_i \rightarrow A_j} = p_{s_k \rightarrow A_j} = \sum_{s_m \in A_j} p_{s_k \rightarrow s_m}, k \in A_i \end{aligned} }[/math]

 

 

 

 

(4)

这个公式表明，群组[math]\displaystyle{ A_i }[/math]到群组[math]\displaystyle{ A_j }[/math]的转移概率 = 群组[math]\displaystyle{ A_i }[/math]中任意状态[math]\displaystyle{ s_k }[/math]到群组[math]\displaystyle{ A_j }[/math]的转移概率 = 群组[math]\displaystyle{ A_i }[/math]中任意状态[math]\displaystyle{ s_k }[/math]到群组[math]\displaystyle{ A_j }[/math]中的所有状态的转移概率之和。

我们也可以直接定义[math]\displaystyle{ P }[/math]和[math]\displaystyle{ P' }[/math]之间的关系：

[math]\displaystyle{ P' = \Phi^T P (\Phi^T)^{-1} = \Phi^T P \Phi (\Phi^T \Phi)^{-1} }[/math]，

其中，[math]\displaystyle{ \Phi^T \Phi }[/math]是一个对角阵，对角线上的元素为每个[math]\displaystyle{ A_i }[/math]包括了多少个微观状态的计数，其逆矩阵就是把对角线上的元素全部取倒数。

对于一个马尔科夫链，有些分组是成块的，有些分组是不成块的。这里用一个简单的例子提供正例和反例。

我们用到的是这样一个转移矩阵：

[math]\displaystyle{ \begin{aligned} P = P_{s_i \rightarrow s_j} = \left [ \begin{array}{ccc:c} 0.3 & 0.3 & 0.3 & 0.1 \\ 0.3 & 0.3 & 0.3 & 0.1 \\ 0.3 & 0.3 & 0.3 & 0.1 \\ \hdashline0.1 & 0.1 & 0.1 & 0.7 \end{array} \right ] \end{aligned} }[/math]


我们可以看到，123行都是一样的，所以把他们[math]\displaystyle{ \{s_1, s_2, s_3\} }[/math]分在一组[math]\displaystyle{ A_1 }[/math]，并把状态[math]\displaystyle{ s_4 }[/math]分在另外一组[math]\displaystyle{ A_2 }[/math]，看起来像是一个合理的正例。接下来我们来计算一下。

第1组[math]\displaystyle{ A_1 }[/math]的状态[math]\displaystyle{ s_1 }[/math]到第1组[math]\displaystyle{ A_1 }[/math]，[math]\displaystyle{ p_{s_1 \rightarrow A_1} = p_{A_1 \rightarrow A_1} = \sum_{s_m \in A_1} p_{s_1 \rightarrow s_m} = p_{s_1 \rightarrow s_1} + p_{s_1 \rightarrow s_2} + p_{s_1 \rightarrow s_3} = 0.3 \times 3 = 0.9 }[/math]

第1组[math]\displaystyle{ A_1 }[/math]的状态[math]\displaystyle{ s_1 }[/math]到第2组[math]\displaystyle{ A_2 }[/math]，[math]\displaystyle{ p_{s_1 \rightarrow A_2} = p_{A_1 \rightarrow A_2} = \sum_{s_m \in A_2} p_{s_1 \rightarrow s_m} = p_{s_1 \rightarrow s_4} = 0.1 }[/math]

第2组[math]\displaystyle{ A_2 }[/math]的状态[math]\displaystyle{ s_4 }[/math]到第1组[math]\displaystyle{ A_1 }[/math]，[math]\displaystyle{ p_{s_4 \rightarrow A_1} = p_{A_2 \rightarrow A_1} = \sum_{s_m \in A_1} p_{s_4 \rightarrow s_m} = p_{s_4 \rightarrow s_1} + p_{s_4 \rightarrow s_2} + p_{s_4 \rightarrow s_3} = 0.1 \times 3 = 0.3 }[/math]

第2组[math]\displaystyle{ A_2 }[/math]的状态[math]\displaystyle{ s_4 }[/math]到第2组[math]\displaystyle{ A_2 }[/math]，[math]\displaystyle{ p_{s_4 \rightarrow A_2} = p_{A_2 \rightarrow A_2} = \sum_{s_m \in A_2} p_{s_4 \rightarrow s_m} = p_{s_4 \rightarrow s_4} = 0.7 }[/math]

[math]\displaystyle{ \begin{aligned} P = P_{A_k \rightarrow A_l} = \left [ \begin{array}{c:c} 0.9 & 0.1 \\ \hdashline0.3 & 0.7 \end{array} \right ] \end{aligned} }[/math]

而反例就是[math]\displaystyle{ \{\{s_1,s_2\}, \{s_3,s_4\}\} }[/math]。


请注意：成块性有一个非常重要的前提：我们需要给定一个对马尔科夫链的状态划分，然后根据马尔科夫状态转移矩阵来判断这个划分对这个马尔科夫链是否成块。

所以，成块性并不是一个马尔科夫链本身的属性，而是对马尔科夫链的某个状态划分的属性。直接说某个马尔科夫链成块或是不成块是不严谨的。相反，我们应该说，对于某个马尔科夫链的某个状态划分是成块或不成块的。然而，为了简略表达，当我们说某个马尔科夫链成块时，我们通常是指它存在某种成块的状态划分。

这一前提也提醒我们，在讨论成块性时，必须明确状态划分的方式。不同的划分方式可能导致不同的结论。一个马尔科夫链可能在某种划分下是成块的，而在另一种划分下则不成块。因此，成块性更多地反映了状态划分与马尔科夫链的适配性，而不是马尔科夫链本身的特性。


成块性的主要性质

成块性作为一种较为严格的划分方式，虽然在压缩维度上可能不如其他方法显著，但它能够保留在宏观层面上所有有用的信息，并且确保粗粒化前后过程的一致性。

宏观最大预测性

宏观最大预测性是指，在将微观状态粗粒化为宏观状态后，我们就不再需要微观状态，只用宏观状态就能预测后续的宏观状态。

具体来说，只需要知道宏观状态[math]\displaystyle{ A^{(t)} }[/math]，就可以最大化[math]\displaystyle{ A^{(t+1)} }[/math]的预测效果，而微观状态[math]\displaystyle{ s^{(t+1)} }[/math]虽然能提供一些额外信息，但是对于预测[math]\displaystyle{ A^{(t+1)} }[/math]并无帮助。这意味着条件互信息[math]\displaystyle{ I(A^{(t+1)}; s^{(t+1)} | A^{(t)}) = 0 }[/math]。宏观动力学不需要任何来自微观状态的修正，就能独立且正确地反映系统在宏观层面上的演化。

尽管成块性划分能够准确预测宏观动力学，但它无法确保像计算力学中的因果态那样的微观最大预测性。因为在聚类过程中会损失微观状态的信息，使得我们无法辨别某个[math]\displaystyle{ A^{(t)} }[/math]具体是哪个微观态[math]\displaystyle{ s^{(t)} }[/math]。这个信息的损失是不可逆的，所以从微观粗粒化到宏观之后，就无法再还原到微观了。

成块划分和粗粒化的可交换性

如果我们把粗粒化和（宏观/微观）动力学看作是两个算子，那这两个算子的可交换性指的是：无论从哪个微观状态分布出发，走这两条路径得到的[math]\displaystyle{ A^{(t+1)} }[/math]的分布是一样的，即，[math]\displaystyle{ s^{(t)} \Phi P' = s^{(t)} P \Phi }[/math]。

我们在马尔科夫链的成块性中推导了结论：基于成块性划分的粗粒化方法满足粗粒化和动力学这两个算子的可交换性。

一致性Consistency

我们在马尔科夫链的成块性中推导了结论：基于成块性划分的粗粒化方法构建的宏观网络与微观网络是一致的。


因果态

定义

因果态相关详情请参照计算力学词条。在计算力学中，状态可以被理解为智能体所接收到的全部历史信息，即历史信息即状态，所以它将系统从过去到未来的变化用一个离散的稳定随机过程描述，状态的取值空间则为双无限序列可数集合[math]\displaystyle{ \overleftrightarrow{S}=\cdots s_{-2} s_{-1} s_0 s_1 s_2\cdots }[/math]，也就是说，一个状态是一个无限长的时间序列。基于当前的时刻[math]\displaystyle{ t }[/math]，我们可以将[math]\displaystyle{ \overleftrightarrow{S} }[/math]分为单侧前向序列[math]\displaystyle{ \overleftarrow{s_t}=\cdots s_{t-3} s_{t-2} s_{t-1} }[/math]和单侧后向序列[math]\displaystyle{ \overrightarrow{s_t}=s_t s_{t+1} s_{t+2} s_{t+3}\cdots }[/math]两个部分，所有可能的未来序列[math]\displaystyle{ \overrightarrow{s_t} }[/math]形成的集合记作[math]\displaystyle{ \overrightarrow{S} }[/math]，所有可能的历史序列[math]\displaystyle{ \overleftarrow{s_t} }[/math]形成的集合记作[math]\displaystyle{ \overleftarrow{S} }[/math]。某一个时刻的状态指的是截止到当前时刻的历史序列。

而因果态相当于是对状态的一种粗粒化描述。这样的状态定义没有假设系统的马尔科夫性，而我们为了把因果态和马尔科夫链的粗粒化理论进行比较，可以定义马尔科夫链的因果态。首先，我们假设系统是一个[math]\displaystyle{ K }[/math]阶的马尔科夫链，即系统在[math]\displaystyle{ t }[/math]时刻的状态分布是由前K个时刻的状态决定的。所以，我们能够重新定义[math]\displaystyle{ \overleftarrow{S} }[/math]为[math]\displaystyle{ S^K }[/math]，而[math]\displaystyle{ \overleftarrow{S} }[/math]为[math]\displaystyle{ S }[/math]。

[math]\displaystyle{ Pr(s^{(t)} | s^{(t-1)}, ... , s^{(0)} ) = Pr(s^{(t)} | s^{(t-1)}, ... , s^{(t-K)} ) }[/math]

对于集合[math]\displaystyle{ \overset{\leftarrow}{S} }[/math]的划分可以有很多种，若某一种划分能够在预测能力最强的同时又非常简洁，那么它肯定是最优的划分，我们把这种用最优的划分方法得到的宏观态称为因果态。下面，我们给出因果态的形式化定义：

对于任意时刻[math]\displaystyle{ t }[/math] 和[math]\displaystyle{ t^{'} }[/math]，在给定过去状态[math]\displaystyle{ \overleftarrow{s_t} }[/math]的条件下，未来状态[math]\displaystyle{ \overrightarrow{s} }[/math] 的分布与给定过去状态[math]\displaystyle{ \overleftarrow{s_{t^{'}}} }[/math]的条件下，未来状态[math]\displaystyle{ \overrightarrow{s} }[/math] 的分布完全相同。那么我们就称[math]\displaystyle{ \overleftarrow{s_t} }[/math]与[math]\displaystyle{ \overleftarrow{s_{t^{'}}} }[/math]等价，记为：[math]\displaystyle{ \overleftarrow{s_t}\sim \overleftarrow{s_{t^{'}}} }[/math]，或简记为[math]\displaystyle{ t\sim t^{'} }[/math]，其中“[math]\displaystyle{ ∼ }[/math] ” 表示由等效未来状态所引起的等价关系，也叫预测等价性（predictive equivalence），可以用公式表示为：

[math]\displaystyle{ t\sim t^{'}\equiv \overleftarrow{s_{t^{'}}}\sim \overleftarrow{s_t} \triangleq P(\overrightarrow{s}|\overleftarrow{s_t} )=P(\overrightarrow{s} |\overleftarrow{s_{t^{'}}} ) }[/math]，

也就是说，若[math]\displaystyle{ t }[/math] 和[math]\displaystyle{ t^{'} }[/math]对未来状态预测的分布相同，则它们对应的状态就都具有一种等价性，这就是因果态定义的来源。

我们可以直接将因果态定义为一个从状态集到状态子集的映射，即[math]\displaystyle{ \epsilon{:}\overleftarrow{S}\mapsto2^{\overset{\leftarrow}{S}} }[/math]，其中[math]\displaystyle{ 2^{\overset{\leftarrow}{S}} }[/math]是[math]\displaystyle{ \overleftarrow{S} }[/math]的幂集，且[math]\epsilon[/math]需要满足如下条件：

[math]\displaystyle{ \epsilon(\stackrel{\leftarrow}{s})\equiv\{\stackrel{\leftarrow}{s}^{\prime}|\mathrm{P}(\overrightarrow{S}=\overrightarrow{s}\mid\stackrel{\leftarrow}{S}=\stackrel{\leftarrow}{s})=\mathrm{P}(\overrightarrow{S}=\overrightarrow{s}\mid\stackrel{\leftarrow}{S}=\stackrel{\leftarrow}{s}^{\prime})，\mathrm{for~all~}\overrightarrow{s}\in\overrightarrow{S},\stackrel{\leftarrow}{s}^{\prime}\in\stackrel{\leftarrow}{S}\} }[/math]

其中，[math]\displaystyle{ \stackrel{\leftarrow}{s} }[/math]为历史序列的随机变量。

所有的因果态的集合，可以由预测等价关系"[math]\sim[/math]"诱导的对状态空间的一种划分来定义，即状态集与等价关系的商集：[math]\mathcal{S}=\stackrel{\leftarrow}{S}/\sim[/math]，其中[math]\displaystyle{ \mathcal{S} }[/math]为因果态划分。

为了说明因果态的概念，我们绘制下面的示意图以举例说明：

如上图所示[10]，左侧的数字代表[math]\displaystyle{ t }[/math]时刻的状态序列，右侧的箭头形状代表对未来状态预测的分布，可以观察到[math]\displaystyle{ t_9 }[/math]和[math]\displaystyle{ t_{13} }[/math]时刻的箭头形状完全相同，说明它们对未来状态预测的分布相同，则处于相同的因果态；同样的道理，在[math]\displaystyle{ t_{11} }[/math]时刻，它的箭头形状与[math]\displaystyle{ t_9 }[/math]和[math]\displaystyle{ t_{13} }[/math]时刻不同，则处于不同的因果态。


因果态的主要性质

为什么我们要关注因果态这样一种特定的宏观态呢？下面是因果态的两个主要性质，能够体现出因果态是一种最理想的宏观态：

最大预测性

性质1（因果态具有最大预测性）：对于宏观态[math]\displaystyle{ \mathcal{∀R} }[/math]和正整数[math]\displaystyle{ L }[/math]，都有[math]\displaystyle{ H[\stackrel{\rightarrow}{S}^L|\mathcal{R}]\geq H[\stackrel{\rightarrow}{S}^L|\mathcal{S}] }[/math]，[math]\displaystyle{ \stackrel{\rightarrow}{S}^L }[/math]为长度为[math]\displaystyle{ L }[/math]的未来序列集合，[math]\displaystyle{ \mathcal{S} }[/math]为因果态划分，所有可能的历史序列记作[math]\displaystyle{ \overleftarrow{S} }[/math]。[math]\displaystyle{ H[\stackrel{\rightarrow}{S}^L|\mathcal{R}] }[/math]和[math]\displaystyle{ H[\stackrel{\rightarrow}{S}^L|\mathcal{S}] }[/math]是[math]\displaystyle{ \stackrel{\rightarrow}{S}^L }[/math]的条件熵。可以理解为因果态集合[math]\displaystyle{ \mathcal{S} }[/math]是在所有的划分宏观态集合[math]\displaystyle{ \mathcal{R} }[/math]中，预测能力最强的一种宏观态，证明过程如下：

根据因果态的定义可知：[math]\displaystyle{ \epsilon(\stackrel{\leftarrow}{s})\equiv\{\stackrel{\leftarrow}{s}^{\prime}|\mathrm{P}(\stackrel{\rightarrow}{S}=\stackrel{\rightarrow}{s}\mid\stackrel{\leftarrow}{S}=\stackrel{\leftarrow}{s})=\mathrm{P}(\stackrel{\rightarrow}{S}=\stackrel{\rightarrow}{s}\mid\stackrel{\leftarrow}{S}=\stackrel{\leftarrow}{s}^{\prime})\} }[/math]

上式等价于[math]\displaystyle{ \mathrm{P}(\stackrel{\rightarrow}{S}=\stackrel{\rightarrow}{s} |\mathcal{S}=\epsilon(\stackrel{\leftarrow}{s}))=\mathrm{P}(\stackrel{\rightarrow}{S}=\stackrel{\rightarrow}{s}\mid\stackrel{\leftarrow}{S}=\stackrel{\leftarrow}{s}) }[/math]

因为[math]\displaystyle{ H[\stackrel{\rightarrow}{S}^L|\mathcal{S}]~=~H[\stackrel{\rightarrow}{S}^L|~\stackrel{\leftarrow}{S}] }[/math]，且[math]\displaystyle{ H[\stackrel{\to}{S}^L|\mathcal{R}] \geq H[\stackrel{\to}{S}^L|\stackrel{\leftarrow}{S}] }[/math]

所以[math]\displaystyle{ H[\stackrel{\rightarrow}{S}^L|\mathcal{R}]\geq H[\stackrel{\rightarrow}{S}^L|\mathcal{S}] }[/math]


最小随机性

性质2（因果态具有最小随机性）：设[math]\displaystyle{ \hat{\mathcal{R}} }[/math]和[math]\displaystyle{ \hat{\mathcal{R}}^{\prime} }[/math]为满足性质1中不等式等号成立的宏观态，则对于所有的[math]\displaystyle{ \hat{\mathcal{R}} }[/math]和[math]\displaystyle{ \hat{\mathcal{R}}^{\prime} }[/math]，都有[math]\displaystyle{ H[\hat{\mathcal{R}}^{\prime}|\hat{\mathcal{R}}]\geq H[\mathcal{S}^{\prime}|\mathcal{S}] }[/math]，其中[math]\displaystyle{ \hat{\mathcal{R}}^{\prime} }[/math]和[math]\displaystyle{ \mathcal{S}^{\prime} }[/math]分别是该过程的下一时刻状态和下一时刻因果态。

该性质可以被理解为在相同预测能力的前提下，因果态集合[math]\displaystyle{ \mathcal{S} }[/math]在所有的宏观态[math]\displaystyle{ \mathcal{R} }[/math]的集合中，随机性最小。

用互信息的角度去理解的话，上式等价于[math]\displaystyle{ I(\mathcal{S}^{\prime};\mathcal{S})\geq I(\hat{\mathcal{R}}^{\prime};\hat{\mathcal{R}}) }[/math]，可以理解为因果态为所有的与它自己的下一时刻的互信息中的最大的一个宏观态。

若想更深入的理解因果态的性质可以阅读Cosma Rohilla Shalizi 和James Crutchfield合写的一篇论文[11]，里面有因果态更多的性质和对应的形式化证明过程。


尺度闭合性

当研究跨尺度的粗粒化时，判断粗粒化方案的一个指标是，得出的宏观动力学是否能够独立反映微观动力学在宏观层面的演化行为，而不需要使用微观动力学的状态信息对其进行修正。也就是说，判断微观和宏观动力学之间是不是独立计算和演化的。

这个问题由来已久，比如对于一个人类群体，既可以从社会学经济学角度解释系统的现象，也可以从生物学生理过程的层次去描述，还可以从分子尺度的物理学去描述。这些尺度上的各自规律可以独立成立，由此诞生了各个学科。这种分层次的现象一直被哲学家们用因果闭合或因果闭包的概念进行讨论。Anil Seth在他的动力学解耦[12]文章中认为，在一个复杂系统中，比如一个动物的身体，当宏观动力学（肢体的活动）能够自主演化，不需要回到微观动力学（细胞的活动）去提供预测信息的时候，我们就认为存在涌现。

而Fernando E. Rosas等人的工作中[13]研究了三种闭合方式：因果闭合、信息闭合和计算闭合，作为尺度间是否独立演化的指标。这里的闭合是指宏观尺度对微观尺度的闭合。当达到闭合时，即代表宏观尺度达到了某种程度上的‘独立演化’。


基础概念

软件和硬件

作者们使用‘软件’和‘硬件’来类比系统。当我们使用一台电脑计算“1+2”的时候，我们知道这台电脑会得出“3”的结果。这是因为我们知道软件的输入是“1+2”，而不会这样想：把“1+2”翻译成一条二进制的机器码，交给软件的底层代码进行操作，得到“3”的机器码再翻译回来；更不会去想，计算机的硬件，比如CPU里各个组成中的晶片经过了什么样的变化来处理这条指令的。

所以，这代表了，这类计算机的软件和硬件是可以分开讨论的。虽然‘软件’需要依托在‘硬件’之上运行，所以‘软件’代表了某个宏观尺度，‘硬件’代表了某个微观尺度。软件中每一个计算都离不开硬件上的特定变化，但是当我们去理解‘软件’的运行规律的时候，却不需要回到‘硬件’层面上来解释。这就是‘软件’和‘硬件’之间的计算独立性：‘软件’（对‘硬件’）闭合了。

计算机中的软件和硬件是人类定义的，那自然界中的其他非人造系统中会不会也存在这样的软件硬件的区分呢？这就是这篇文章主要探讨的事情。


单尺度因果态

词条的上一部分已经介绍过了计算力学中因果态的概念，为了与原文符号对应，我们重新提供因果态的定义。在一个时序系统中，我们想要有一个ϵ-machine（ϵ-机器），能够根据当前时间步的状态和历史记忆，来预测未来状态。

无论系统是确定的或概率生成的，当我们已知 [math]\displaystyle{ t }[/math] 时刻的状态和历史记忆 [math]\displaystyle{ \overleftarrow{X_t}=\{... , X_{t-2}, X_{t-1}, X_t\} }[/math] 时，我们能够预测未来状态的轨迹 [math]\displaystyle{ p(\overrightarrow{X}^L_t| \overleftarrow{X_t}) }[/math]，其中[math]\displaystyle{ \overrightarrow{X}^L_t=\{X_{t+1}, X_{t+2}, ... , X_{t+L}\} }[/math]，简写成[math]\displaystyle{ \overrightarrow{X_t} }[/math]。

不同时刻 [math]\displaystyle{ t }[/math] 和 [math]\displaystyle{ t' }[/math] 的状态和历史记忆虽然不完全一样，但是他们可能会导致同样的未来状态，[math]\displaystyle{ p(\overrightarrow{X_t}| \overleftarrow{X_t}) = p(\overrightarrow{X_t'}| \overleftarrow{X_t'}) }[/math]。所以，计算力学中把这些具有相同未来的 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 划为等价euqivalent的状态，即[math]\displaystyle{ \overleftarrow{x_t} \equiv \overleftarrow{x_{t'}}\ \ \text{iff}\ \ p(\overrightarrow{x_t}| \overleftarrow{x_t} ) = p(\overrightarrow{x_{t'}}|\overleftarrow{x_{t'}}) }[/math]，并定义因果态[math]\displaystyle{ E = \{E_t\}_{t \in \mathbb{N} }, E_t = \epsilon( \overleftarrow{X_t} ) }[/math]

所以，ϵ-机器的工作是根据 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 而判断 [math]\displaystyle{ t }[/math] 时刻系统处于哪个因果态 [math]\displaystyle{ E_t }[/math] ，并根据该因果态预测未来轨迹 [math]\displaystyle{ \overrightarrow{X_t} }[/math] 。

计算力学认为，包括人类在内的智能体，每时每刻都在对周围环境观察，得到大量的观测数据，并在自己的内部（大脑中）构建简洁的机制来预测观测数据的变化（一种隐马尔科夫模型），而因果态和ϵ-机器因其简洁性和有效性，被认为是一种最理想的预测机制模型。

我们可以从最小解释模型的角度出发来理解因果态。最小解释模型是指，我们要如何去用最简洁的宏观状态来解释和描述系统发生的现象，而且这种解释模型要在微观层面扰动时保持‘稳定’。所以它追求的是两方面：首先需要稳定，其次尺度越高越好。而这个宏观状态就是因果态。

因果态示例，出自Rosas于因果涌现读书会第五季第一期的分享


跨尺度因果态

设我们对微观动力学 [math]\displaystyle{ X_t }[/math] 粗粒化成宏观动力学 [math]\displaystyle{ Z_t = f(X_t) }[/math]，[math]\displaystyle{ f }[/math] 为粗粒化函数。上述的ϵ-机器只能处理在 [math]\displaystyle{ X }[/math] 空间上单一尺度的因果态，而作者们提出了另一种机器，υ-machine（υ-机器）来尝试定义跨尺度的计算机器和因果态。他们把因果态定义中的未来状态从微观的 [math]\displaystyle{ X }[/math] 改为宏观的 [math]\displaystyle{ Z }[/math]，并定义微观状态 [math]\displaystyle{ X }[/math] 中的等价关系为

[math]\displaystyle{ \overleftarrow{x_t} \equiv \overleftarrow{x_{t'}} \ \ \text{iff}\ \ p( \overrightarrow{z_t} | \overleftarrow{x_t}) = p(\overrightarrow{z_{t'}}| \overleftarrow{x_{t'}}) }[/math]

也就是说，υ-机器的目标不再是预测 [math]\displaystyle{ X }[/math] 的未来轨迹，而是其粗粒化后的 [math]\displaystyle{ Z }[/math] 的未来轨迹，并把拥有同样 [math]\displaystyle{ p(\overrightarrow{z_t}| \overleftarrow{x_t}) }[/math] 的 [math]\displaystyle{ \overleftarrow{x_t} }[/math] 和 [math]\displaystyle{ \overleftarrow{x_{t'}} }[/math] 定义为υ-机器的因果态，[math]\displaystyle{ U = \{U_t\}_{t \in \mathbb{N} },\ U_t = υ( \overleftarrow{X_t} ) }[/math]

它的意义是什么呢？作者给出的解释是，υ-机器的因果态是一个最优的信息瓶颈：υ-机器是一种对 [math]\displaystyle{ \overleftarrow{x_t} }[/math] 最粗的粗粒化粒度，并满足对 [math]\displaystyle{ \overrightarrow{z_t} }[/math]的最大预测性的粗粒化方案。

所以，我们现在有三种机器：

  1. 微观动力学 [math]\displaystyle{ X_t }[/math] 本身的ϵ-机器，因果态为 [math]\displaystyle{ E }[/math]；
  1. 宏观动力学 [math]\displaystyle{ Z_t }[/math] 本身的ϵ-机器，因果态为 [math]\displaystyle{ E' }[/math]；
  1. 微观与宏观动力学 [math]\displaystyle{ Z_t = f(X_t) }[/math] 之间的υ-机器，因果态为 [math]\displaystyle{ U }[/math]；
三种因果态的示意图，摘录自文章原文

这三种机器的因果态有什么区别呢？首先，[math]\displaystyle{ U_t }[/math] 要比 [math]\displaystyle{ E'_t }[/math] 更强（powerful），因为后者只包括了宏观的信息，而前者同时包括了宏观和微观的信息。其次，[math]\displaystyle{ E_t }[/math] 要比 [math]\displaystyle{ U_t }[/math] 更强，因为 [math]\displaystyle{ U_t }[/math] 只保留了宏观的最大预测性，而 [math]\displaystyle{ E_t }[/math] 能保留微观的最大预测性（同时因为我们知道 [math]\displaystyle{ f }[/math] ，所以 [math]\displaystyle{ E_t }[/math] 也能保留宏观的最大预测性）。

我们能够了解一些关于它们的推断。

推断1：微观 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 对宏观 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 提供的信息量 等于 υ-机器的因果态 [math]\displaystyle{ U_t }[/math] 提供的信息量，因为

[math]\displaystyle{ I(\overleftarrow{X_t}; \overrightarrow{Z_t}) = \sum_{\overleftarrow{x_t}} p(\overleftarrow{x_t}) D_{KL} \left ( p(\overrightarrow{z_t}|\overleftarrow{x_t}) || p(\overrightarrow{z_t}) \right ) = \sum_{u_t} p(u_t) D_{KL} \left ( p(\overrightarrow{z_t}|u_t) || p(\overrightarrow{z_t}) \right ) = I(U_t; \overrightarrow{Z_t}) }[/math]

推断2：υ-机器的因果态 [math]\displaystyle{ U_t }[/math] 对宏观 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 提供的信息量 大于等于 宏观 ϵ-机器的因果态 [math]\displaystyle{ E'_t }[/math]。证明如下：

首先，[math]\displaystyle{ I(E'_t ; \overrightarrow{Z_t}) = I( \overleftarrow{Z_t} ; \overrightarrow{Z_t}) }[/math]；其次，[math]\displaystyle{ I(U_t ; \overrightarrow{Z_t}) = I( \overleftarrow{X_t} ; \overrightarrow{Z_t}) }[/math]；

因为 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 是 [math]\displaystyle{ \overrightarrow{X_t} }[/math] 的粗粒化，所以 [math]\displaystyle{ I( \overleftarrow{X_t} ; \overrightarrow{Z_t}) }[/math] 肯定比 [math]\displaystyle{ I( \overleftarrow{Z_t} ; \overrightarrow{Z_t}) }[/math] 大。因此，[math]\displaystyle{ I(U_t ; \overrightarrow{Z_t}) = I( \overleftarrow{X_t}; \overrightarrow{Z_t}) \geq I( \overleftarrow{Z_t} ; \overrightarrow{Z_t}) = I(E'_t ; \overrightarrow{Z_t}) }[/math]。证毕。


了解了跨尺度因果态的定义后，我们就能逐个了解这篇文章中的三种闭合定义。


因果闭合 Causal Closure

因果闭合的定义为：[math]\displaystyle{ Z_t }[/math] 的ϵ-机器 与 [math]\displaystyle{ Z_t = f(X_t) }[/math] 之间的υ-机器 之间存在等价关系，即 [math]\displaystyle{ E' }[/math] 和 [math]\displaystyle{ U }[/math] 之间存在双射（bijective）投影关系（两个空间内的因果态分组相互一一对应）。

在这种情况下，不同 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 对 [math]\displaystyle{ Z }[/math] 轨迹预测的差异也会在 [math]\displaystyle{ \overleftarrow{Z_t} }[/math] 上反映出来。也就是说，这些差异在粗粒化 [math]\displaystyle{ f }[/math] 过程中会被完全保留，不会被更改或抹除。

当我们对 [math]\displaystyle{ Z }[/math] 进行干预的时候，微观 [math]\displaystyle{ X }[/math] 的状态也会跟着改变，就像当我们加热了一锅水的时候，其中的水分子的能量也会改变。因果闭合的物理意义为：在宏观上的干预，和微观上与之对应的干预，作用在宏观层面上产生的结果完全相同。对 [math]\displaystyle{ Z }[/math] 空间进行的干预影响到 [math]\displaystyle{ X }[/math] 空间时，按照 [math]\displaystyle{ X }[/math] 的动力学进行演化后再粗粒化得到的 [math]\displaystyle{ Z }[/math]，与直接在 [math]\displaystyle{ Z }[/math] 空间上进行演化没有区别。

这说明了，当我们对 [math]\displaystyle{ Z }[/math] 进行干预的时候，仅凭 [math]\displaystyle{ Z }[/math] 自身的动力学就能正确的演化，不用再回到微观 [math]\displaystyle{ X }[/math] 层溯源。所以这就是宏观和微观之间因果闭合了，即宏观动力学在因果性上能做到独立的演化。因此，我们可以把宏观的结果 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 视为只来自于宏观的干预 [math]\displaystyle{ \overleftarrow{Z_t} }[/math] ，与微观 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 无关。


信息闭合 Information Closure

上面所述的因果闭合的建模基于动力学的因果态，强调因果的概念。我们也可以关注计算力学的另一方面，即计算与预测。

因果态对 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 预测能力可以写成信息论的形式，即变量对 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 提供的互信息。[math]\displaystyle{ Z_t }[/math] 本身的ϵ-机器和 [math]\displaystyle{ Z_t = f(X_t) }[/math] 之间的υ-机器的等价关系可以写成：

[math]\displaystyle{ I( \overleftarrow{Z_t} ; \overrightarrow{Z_t}) = I( \overleftarrow{X_t} ; \overrightarrow{Z_t}) }[/math]

这说明了，宏观的 [math]\displaystyle{ \overleftarrow{Z_t} }[/math] 和微观的 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 给 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 提供的信息量是一样的。所以，信息闭合是指，对于一个粗粒化后的宏观过程，已知其自身过去的信息 [math]\displaystyle{ \overleftarrow{Z_t} }[/math] 就足以进行最优预测，进一步引入底层（微观）[math]\displaystyle{ \overleftarrow{X_t} }[/math] 信息也不能提升预测能力——也就是说，微观的信息 [math]\displaystyle{ \overleftarrow{X_t} }[/math] 对宏观未来 [math]\displaystyle{ \overrightarrow{Z_t} }[/math] 的预测是冗余的。所谓一个尺度在某某方面闭合，便是在说这个尺度上系统可以“自给自足”，不需要再考虑其他尺度。

虽然定义方式不一样，但是我们不难看出，根据跨尺度因果态介绍中的推断2，因果闭合和信息闭合这两个对 [math]\displaystyle{ Z_t = f(X_t) }[/math] 的条件是等价的：由于[math]\displaystyle{ I(U_t ; \overrightarrow{Z_t}) = I(E'_t ; \overrightarrow{Z_t}) }[/math]，所以[math]\displaystyle{ I(U_t ; \overrightarrow{Z_t}) = I( \overleftarrow{X_t}; \overrightarrow{Z_t}) = I( \overleftarrow{Z_t} ; \overrightarrow{Z_t}) = I(E'_t ; \overrightarrow{Z_t}) }[/math]。


另外，当信息闭合时也满足： [math]\displaystyle{ I(\overrightarrow{Z_t}; \overleftarrow{X_t} | \overleftarrow{Z_t}) = I(\overleftarrow{Z_t}, \overleftarrow{X_t} ; \overrightarrow{Z_t}) - I(\overrightarrow{Z_t} ; \overleftarrow{Z_t}) = I(\overleftarrow{X_t} ; \overrightarrow{Z_t}) - I(\overrightarrow{Z_t} ; \overleftarrow{Z_t}) = 0 }[/math]


计算闭合 Computational Closure

定义 [math]\displaystyle{ Z }[/math] 的ϵ-机器为计算闭合，当满足：存在粗粒化函数 [math]\displaystyle{ f^* }[/math]使得[math]\displaystyle{ E'_t = f^*(E_t) }[/math]，其中[math]\displaystyle{ E_t }[/math] 和 [math]\displaystyle{ E'_t }[/math] 分别为 [math]\displaystyle{ X }[/math] 和 [math]\displaystyle{ Z }[/math] 的ϵ-机器的因果态。

也就是说， [math]\displaystyle{ E_t }[/math] 和 [math]\displaystyle{ E'_t }[/math] 之间存在一个粗粒化的（一对一或多对一）映射关系。同样的，[math]\displaystyle{ X }[/math] 和 [math]\displaystyle{ Z }[/math] 之间也存在粗粒化的映射关系，所以我们可以看到这四个变量之间的关系可以归结于如下两条从 [math]\displaystyle{ X }[/math] 到 [math]\displaystyle{ E' }[/math] 的路径：

  1. [math]\displaystyle{ X \stackrel{\epsilon}{\rightarrow} E \stackrel{f^*}{\rightarrow} E' }[/math]
  2. [math]\displaystyle{ X \stackrel{f}{\rightarrow} Z \stackrel{\epsilon'}{\rightarrow} E' }[/math]

我们把 [math]\displaystyle{ f }[/math] 和 [math]\displaystyle{ f^* }[/math] 视作粗粒化算子，[math]\displaystyle{ \epsilon }[/math] 和 [math]\displaystyle{ \epsilon' }[/math] 视作因果态算子，那么我们看到两条路径中的粗粒化算子和因果态算子是可交换的。

计算闭合可交换性示意图，摘录自文章原文


上面提到过ϵ-机器可以看作最小解释模型，描述了动力学的机制。所以，因果态空间 [math]\displaystyle{ E }[/math] 和 [math]\displaystyle{ E' }[/math] 可以看作是动力学的‘机制空间’。而计算闭合的一个意义就是指出，这种动力学不仅在本来的变量/状态空间中能够被粗粒化 [math]\displaystyle{ f }[/math] ，而且动力学的机制也能够被粗粒化 [math]\displaystyle{ f^* }[/math]。所以，当计算闭合时，宏观 [math]\displaystyle{ Z }[/math] 空间的机制能通过粗粒化微观 [math]\displaystyle{ X }[/math] 空间的机制而获得。

对于一台机器，我们肯定希望它的宏观的计算过程要更加的简洁，也就是说宏观的ϵ-机器最好是微观的ϵ-机器的一种粗粒化。难道宏观机器还可以不是微观机器的粗粒化吗？粗粒化意味着宏观变量是微观变量的单射。这个映射过程可能引入更多的噪音，导致宏观过程的因果态反而比微观的还要多。我们看下面这样一个例子：

考虑一个马尔可夫链 [math]\displaystyle{ X }[/math]，有三个状态 [math]\displaystyle{ a }[/math]、[math]\displaystyle{ b }[/math] 和 [math]\displaystyle{ c }[/math]，[math]\displaystyle{ p(a,i)=p(b,i) \gt 0, \forall i }[/math]，[math]\displaystyle{ p(c,b)=0 }[/math]，且 [math]\displaystyle{ p(b,a) \neq p(c,a) }[/math]。这样，微观ϵ-机器只有两个因果态，[math]\displaystyle{ \{a,b\} }[/math] 和 [math]\displaystyle{ c }[/math]。

假设我们的粗粒化分组方案为[math]\displaystyle{ \{a, \{b,c\}\} }[/math]，即[math]\displaystyle{ Z_t = A\ \text{if}\ X_t = a; Z_t = B\ \text{if}\ X_t = b\ or\ c }[/math]。 可以证明当 [math]\displaystyle{ Z }[/math] 连续处于 [math]\displaystyle{ B }[/math] 状态的时间越长，[math]\displaystyle{ X_t = c }[/math] 的概率就越高（因为 的转移被禁止），所以 [math]\displaystyle{ Z }[/math] 不是马尔可夫的。由于 [math]\displaystyle{ p(b,a) \neq p(c,a) }[/math]，因此 [math]\displaystyle{ Z }[/math] 的ϵ-机器需要区分以不同长度的 [math]\displaystyle{ B }[/math] 结尾的历史，例如 [math]\displaystyle{ \{A, (BA), (BBA), (BBBA)\} }[/math] 。所以， [math]\displaystyle{ E' }[/math] 的因果态数量多于 [math]\displaystyle{ E }[/math] 的因果态数量，因此 [math]\displaystyle{ E' }[/math] 不可能是 [math]\displaystyle{ E }[/math] 的粗粒化。

我们本以为，把底层系统用机器描述后，宏观层的ϵ-机可能只需把底层机器的状态合并即可。但事实远非如此。在这个例子中，原始系统有3个状态，但因为和看起来完全一样，所以微观机器只需要区分 和 两种状态。但宏观的粗粒化选择却是，把和映射到一个宏观状态上。因为底层动力学的特殊性，想预测下一个何时出现，宏观ϵ-机需要记住最近连续出现了几个。也就是说，宏观行为实际上增加了记忆长度，需要无穷多个不同的因果态，而不仅仅是简单的两种状态。所以，这并不是一个好的、满足计算闭合要求的粗粒化方案。

如果满足计算闭合的要求，我们就可以从计算模型的角度对一个系统进行多尺度建模。在考虑多尺度建模的时候，我们只需要关注各个尺度机器上的粗粒化过程即可。


计算闭合与前面的因果闭合/信息闭合之间的关系是：如果一个粗粒化函数 [math]\displaystyle{ Z_t = f(X_t) }[/math] 是因果闭合/信息闭合的，那它一定是计算闭合的。证明如下：

首先，υ-机器的因果态 [math]\displaystyle{ U_t }[/math] 是 微观 ϵ-机器的因果态 [math]\displaystyle{ E_t }[/math] 的一种粗粒化方案（证明请参照论文[13] Proposition 2）。而且，因果闭合的条件说明了，[math]\displaystyle{ U_t }[/math] 等于宏观 ϵ-机器的因果态 [math]\displaystyle{ E'_t }[/math]。所以，[math]\displaystyle{ E'_t }[/math] 和 [math]\displaystyle{ E_t }[/math] 之间也应该存在粗粒化方案，满足计算闭合的条件。


尺度闭合 与 成块性

在论文[14]中，作者给出了以下定理：

  1. 一个不可约irreducible，和非周期性aperiodic 的马尔科夫链 [math]\displaystyle{ X }[/math] ，[math]\displaystyle{ Z_t=f(X_t) }[/math] 是对 [math]\displaystyle{ X }[/math] 的状态空间上定义的粗粒化函数，那 [math]\displaystyle{ X }[/math] 的强成块性 在其达到稳态stationary分布时 等价于 [math]\displaystyle{ H(Z_{t+1}|Z_t) = H(Z_{t+1}|X_t) }[/math]。
  1. [math]\displaystyle{ X }[/math] 的弱成块性 在其达到稳态分布时 等价于 [math]\displaystyle{ H(Z_{t+1}|Z_t) = H(Z_{t+1}|Z_{\lt t}) }[/math]。


作者们[13]认为，当使用 [math]\displaystyle{ Z_t }[/math] 来预测 [math]\displaystyle{ Z_{t+1} }[/math] 时，可以拆分为两步：先利用宏观 [math]\displaystyle{ Z_t }[/math] 推断微观 [math]\displaystyle{ X_t }[/math]，再利用微观 [math]\displaystyle{ X_t }[/math] 估算宏观 [math]\displaystyle{ Z_{t+1} }[/math]。所以，

推断1：当一个马尔科夫链 [math]\displaystyle{ X_t }[/math] 存在粗粒化函数 [math]\displaystyle{ Z_t=f(X_t) }[/math]，而 [math]\displaystyle{ X_t }[/math] 对 [math]\displaystyle{ f }[/math] 强成块时，[math]\displaystyle{ Z_t }[/math] 满足因果闭合/计算闭合。（但是满足因果闭合/计算闭合的 [math]\displaystyle{ Z_t=f(X_t) }[/math] 不一定都是强成块的，详情请参考论文[13]反例6）。证明如下：

[math]\displaystyle{ \begin{aligned} I(\overrightarrow{Z_t}; \overleftarrow{X_t} | \overleftarrow{Z_t}) & = H(\overrightarrow{Z_t} | \overleftarrow{Z_t}) - H(\overrightarrow{Z_t} | \overleftarrow{X_t} , \overleftarrow{Z_t}) \\ & = H(\overrightarrow{Z_t} | \overleftarrow{Z_t}) - H(\overrightarrow{Z_t} | \overleftarrow{X_t}) \\ & = \sum_{l = t+1}^{t+L} \left ( H(Z_l | Z_{l-1}) - H(Z_l | \overrightarrow{Z_{t+1}}, X_t) \right ) \\ &\leq = \sum_{l = t+1}^{t+L} \left ( H(Z_l | Z_{l-1}) - H(Z_l | \overrightarrow{Z_{t+1}}, X_t) \right ) \\ &= 0 \end{aligned} }[/math]


推断2：假设存在 [math]\displaystyle{ Z=f(E)=f(\epsilon(X)) }[/math]，[math]\displaystyle{ E }[/math] 为 [math]\displaystyle{ X_t }[/math] 的ϵ-机器的因果态，那么 [math]\displaystyle{ E }[/math] 是弱成块的条件 等价于 [math]\displaystyle{ Z_t }[/math] 计算闭合的条件。证明如下：

首先，计算闭合的条件说明，[math]\displaystyle{ Z }[/math] 的ϵ-机器的因果态 [math]\displaystyle{ E' }[/math]与 [math]\displaystyle{ E }[/math]之间存在粗粒化函数 [math]\displaystyle{ f' }[/math]，所以 [math]\displaystyle{ E }[/math] 可以被粗粒化成具有马尔科夫性的 [math]\displaystyle{ E' }[/math]，满足了弱成块的定义。

而且， [math]\displaystyle{ E }[/math] 对 [math]\displaystyle{ f }[/math] 弱成块时，我们就能推断 [math]\displaystyle{ X }[/math] 也具有马尔科夫性，而 [math]\displaystyle{ Z }[/math] 的ϵ-机器的因果态 [math]\displaystyle{ E' }[/math]就是一个对 [math]\displaystyle{ Z }[/math] 马尔科夫链的粗粒化划分，写作 [math]\displaystyle{ E'_t = \epsilon' (Z_t) = \epsilon' ( f (E_t)) }[/math]，所以因果态 [math]\displaystyle{ E' }[/math]与 [math]\displaystyle{ E }[/math]之间存在粗粒化函数，满足了计算闭合的条件。


下图为原文中提供的图例[13]，总结了上述概念之间的关系：

闭合性的概念.png


  1. 粗粒化函数 [math]\displaystyle{ f }[/math] 的因果闭合性 等价于 其信息闭合性。
  2. 当 [math]\displaystyle{ f }[/math] 因果/信息闭合时，它一定是计算闭合的。
  3. 当一个稳态的马尔科夫链 [math]\displaystyle{ X_t }[/math] 对 [math]\displaystyle{ Z_t = f(X_t) }[/math] 强成块时，它一定是因果/信息闭合的。
  4. 当稳态的 [math]\displaystyle{ X_t }[/math]的因果态 [math]\displaystyle{ E_t }[/math] 对 [math]\displaystyle{ Z_t = f(E_t) }[/math] 弱成块时，[math]\displaystyle{ Z_t }[/math] 是计算闭合的。


因果涌现：有效信息

因果涌现（causal emergence）是指动力系统中的一类特殊的涌现现象，即系统在宏观尺度会展现出更强的因果特性。特别的，对于一类马尔科夫动力系统来说，在对其状态空间进行适当的粗粒化以后，所形成的宏观动力学会展现出比微观更强的因果特性——即更大的有效信息，那么我们称该系统发生了因果涌现[15][16]。

定义

在马尔科夫链中，任意时刻的状态变量[math]X_t[/math]都可以看作是原因，而下一时刻的状态变量[math]X_{t+1}[/math]就可以看作是结果，这样马尔科夫链的状态转移矩阵就是它的因果机制。

有效信息（Effective Information，简称EI）是因果涌现 Causal Emergence理论中的一个核心概念，它可以用来衡量一个马尔科夫动力学的因果效应的强度。这里，因果效应是指把动力学看做一个黑箱，那么不同的输入分布就会导致不同的输出分布，二者之间联系的紧密程度就是因果效应。有效信息通常可以分解为两个部分：确定性（Determinism）和简并性（Degeneracy）。确定性是指，在动力学的作用下，我们根据系统前一时刻的状态会以多大程度预测它下一时刻状态；简并性是指：我们能够以多大程度从下一时刻的状态预测上一时刻的状态。如果确定性越大，或简并性越小，则系统的有效信息就会越大。其他相关详情请参照有效信息词条。

[math]\displaystyle{ \begin{aligned} EI &= I(X_t,X_{t+1}|do(X_t\sim U(\mathcal{X})))=I(\tilde{X}_t,\tilde{X}_{t+1}) \\ &= \sum^N_{i=1}\sum^N_{j=1}Pr(\tilde{X}_t=i,\tilde{X}_{t+1}=j)\log \frac{Pr(\tilde{X}_t=i,\tilde{X}_{t+1}=j)}{Pr(\tilde{X}_t=i)Pr(\tilde{X}_{t+1}=j)}\\ &= \sum^N_{i=1}Pr(\tilde{X}_t=i)\sum^N_{j=1}Pr(\tilde{X}_{t+1}=j|\tilde{X}_t=i)\log \frac{Pr(\tilde{X}_{t+1}=j|\tilde{X}_t=i)}{Pr(\tilde{X}_{t+1}=j)}\\ &= \frac{1}{N}\sum^N_{i=1}\sum^N_{j=1}p_{ij}\log\frac{N\cdot p_{ij}}{\sum_{k=1}^N p_{kj}} \end{aligned} }[/math]

其中[math]\displaystyle{ \tilde{X}_t,\tilde{X}_{t+1} }[/math]分别为把t时刻的[math]X_t[/math]干预为均匀分布后，前后两个时刻的状态。[math]\displaystyle{ p_{ij} }[/math]为第i个状态转移到第j个状态的转移概率。从这个式子，不难看出，EI仅仅是概率转移矩阵[math]P[/math]的函数。

我们也可以将转移概率矩阵[math]P[/math]写成[math]N[/math]个行向量拼接而成的形式，即：[math]\displaystyle{ P=(P_1^T,P_2^T,\cdots,P_N^T)^T }[/math]。其中，[math]P_i[/math]矩阵[math]P[/math]的第[math]i[/math]个行向量，且满足条件概率的归一化条件：[math]||P_i||_1=1[/math]，这里的[math]||\cdot||_1[/math]表示向量的1范数。那么EI可以写成如下的形式：

[math]\displaystyle{ \begin{aligned} EI &= \frac{1}{N}\sum^N_{i=1}\sum^N_{j=1}p_{ij}\log\frac{N\cdot p_{ij}}{\sum_{k=1}^N p_{kj}}\\ &=\frac{1}{N}\cdot \sum_{i=1}^N\left(P_i\cdot \log P_i - P_i\cdot\log \bar{P}\right)\\ &=\frac{1}{N}\sum_{i=1}^N D_{KL}(P_i||\bar{P}) \end{aligned} }[/math]

 

 

 

 

(1)

将矩阵每列求均值，可得到平均转移向量[math]\displaystyle{ \overline{P}=\sum_{k=1}^N P_k/N }[/math]。[math]D_{KL}[/math]便是两个分布的KL散度。因此，EI是转移矩阵每个行转移向量[math]P_i[/math]与平均转移向量[math]\bar{P}[/math]的KL散度的均值。

下表展示的是三个不同的转移概率矩阵，以及它们的EI数值：

马尔科夫链示例
[math]\displaystyle{ P_1=\begin{pmatrix} &0 &0 &1 &0& \\ &1 &0 &0 &0& \\ &0 &0 &0 &1& \\ &0 &1 &0 &0& \\ \end{pmatrix} }[/math],

[math]\displaystyle{ P_2=\begin{pmatrix} &1/3 &1/3 &1/3 &0& \\ &1/3 &1/3 &1/3 &0& \\ &0 &0 &0 &1& \\ &0 &0 &0 &1& \\ \end{pmatrix} }[/math],

[math]\displaystyle{ P_3=\begin{pmatrix} &1/4 &1/4 &1/4 &1/4& \\ &1/4 &1/4 &1/4 &1/4& \\ &1/4 &1/4 &1/4 &1/4& \\ &1/4 &1/4 &1/4 &1/4& \\ \end{pmatrix} }[/math].

[math]EI(P_1)=2[/math] bits [math]EI(P_2)=1[/math] bits [math]EI(P_3)=0[/math] bits

上面所列的三个状态转移矩阵的EI为：2比特、1比特和0比特。由此可见，如果转移概率矩阵中出现更多的0或1，也就是行向量多是独热向量（也叫做one-hot向量，即某一个位置为1，其它位置为0的向量），则EI值就会更大。也就是说，如果在状态转移的过程中，从某一时刻到下一时刻的跳转越确定，则EI值就会倾向于越高。


根据公式2，我们发现，EI在马尔科夫链的情景下可以被分解为两项，即：

[math]\displaystyle{ \begin{aligned} EI &= \frac{1}{N}\cdot \sum_{i=1}^N\left(P_i\cdot \log P_i - P_i\cdot\log \bar{P}\right)\\ &=\underbrace{-\langle H(P_i)\rangle}_{确定性项}+\underbrace{H(\bar{P})}_{非简并性项} \end{aligned} }[/math]

 

 

 

 

(2)

其中，第一项：[math]-\langle H(P_i)\rangle\equiv \frac{1}{N}\sum_{i=1}^N H(P_i)[/math]为每个行向量[math]P_i[/math]的负熵的平均值，它刻画了整个马尔科夫转移矩阵的确定性（determinism）；第二项：[math]H(\bar{P})[/math]为平均行向量的熵，其中[math]\bar{P}\equiv \frac{1}{N}\sum_{i=1}^N P_i [/math]为所有N个行向量的平均行向量，它刻画了整个马尔科夫转移矩阵的非简并性或非退化性(non-degeneracy)。

我们定义一个马尔科夫链转移矩阵P的简并性为：[math]\displaystyle{ Degeneracy \equiv \log N - H(\bar{P})=\log N + \sum_{j=1}^N \bar{P}_{\cdot j}\log \bar{P}_{\cdot j}=\sum_{j=1}^N \frac{\sum_{i=1}^Np_{ij}}{N}\log \left(\sum_{i=1}^Np_{ij}\right) }[/math]

这一项为简并性或叫退化性，为了防止其为负数，所以加上了[math]\log N[/math][17]。这里的“简并性”的含义是：如果知道了系统的当前状态，能不能反推系统在上一时刻的状态的能力，如果可以推断，则这个马尔科夫矩阵的简并性就会比较低，也就是非简并的；而如果很难推断，则马尔科夫矩阵就是简并的，也即退化的。


最大化有效信息的划分

开头提到，对马尔科夫链做粗粒化的其中一个动机就是，想抹除一些微观动力学上的噪音，而得到在宏观上更具规律和因果性的动力学。所以，我们可以利用EI这一马尔科夫链的因果效应度量，来衡量聚类后的宏观动力学的因果强度。

但是，EI的大小和状态空间大小有关，这一性质在我们比较不同尺度的马尔科夫链的时候非常不方便，我们需要一个尽可能不受尺度效应影响的因果效应度量。因此，我们需要对有效信息EI做一个归一化处理，得到和系统尺寸无关的一个量化指标。根据Erik Hoel和Tononi等人的工作，要用均匀分布即最大熵分布下的熵值，即[math]\displaystyle{ \log N }[/math]来做分母对EI进行归一化，这里的[math]N[/math]为状态空间[math]\mathcal{X}[/math]中的状态的数量[17]。那么归一化后的EI便等于：

[math]\displaystyle{ Eff=\frac{EI}{\log N} }[/math]

进一步定义归一化指标也称为有效性（effectiveness）。

有了有效信息这一度量指标后，我们便可以讨论马尔科夫链上的因果涌现了。对于一个马尔科夫链，观察者可以建立多尺度视角去观测，区分出微观和宏观。首先，原始的马尔科夫概率转移矩阵[math]\displaystyle{ P_N }[/math]即定义了微观动力学；其次，经过对微观态的粗粒化映射（coarse-graining）后（通常反映为对微观状态的分组），观察者可以得到对应的宏观状态（即分组），这些宏观状态彼此之间的概率转移可以用宏观的马尔科夫概率转移矩阵[math]\displaystyle{ P_M }[/math]来刻画。对两个动力学分别可以计算EI，如果宏观的EI大于微观EI，我们说该系统发生了因果涌现。

因果涌现示意图

可以定义一个新的指标直接度量因果涌现的程度：

[math]\displaystyle{ ce = Eff(P_M) - Eff(P_N) = \frac{EI(P_M)}{\log M} - \frac{EI(P_N)}{\log N} }[/math]

这里[math]P[/math]为微观状态的马尔科夫概率转移矩阵，维度为：[math]N\times N[/math]，这里N为微观的状态数；而[math]P'[/math]为对[math]P[/math]做粗粒化操作之后得到的宏观态的马尔科夫概率转移矩阵，维度为[math]M\times M[/math]，其中[math]M<N[/math]为宏观状态数。如果计算得出的[math]\displaystyle{ ce \gt 0 }[/math]，则称该系统发生了因果涌现，否则没有发生。

下面，我们展示一个具体的因果涌现的例子：

马尔科夫链示例

[math]\displaystyle{ P_m=\begin{pmatrix} &1/3 &1/3 &1/3 &0& \\ &1/3 &1/3 &1/3 &0& \\ &1/3 &1/3 &1/3 &0& \\ &0 &0 &0 &1& \\ \end{pmatrix} }[/math],

[math]\displaystyle{ P_M=\begin{pmatrix} &1 &0 & \\ &0 &1 & \\ \end{pmatrix} }[/math].

[math]\begin{aligned}&Det(P_m)=0.81\ bits,\\&Deg(P_m)=0\ bits,\\&EI(P_m)=0.81\ bits\end{aligned}[/math] [math]\begin{aligned}&Det(P_M)=1\ bits,\\&Deg(P_M)=0\ bits,\\&EI(P_M)=1\ bits\end{aligned}[/math]

在这个例子中，微观态的转移矩阵是一个4*4的矩阵，其中前三个状态彼此以1/3的概率相互转移，这导致该转移矩阵具有较小的确定性，因此EI也不是很大为0.81。然而，当我们对该矩阵进行粗粒化，也就是把前三个状态合并为一个状态a，而最后一个状态转变为一个宏观态b。这样所有的原本三个微观态彼此之间的转移就变成了宏观态a到a内部的转移了。因此，转移概率矩阵也就变成了[math]P_M[/math]，它的EI为1。在这个例子中，可以计算它的因果涌现度量为：

[math]\displaystyle{ ce=Eff(P_M)-Eff(P_m)=1-0.405=0.595 }[/math]


当[math]\displaystyle{ ce }[/math]越大，说明宏观动力学的因果效应更强，也说明宏观动力学能够更好地表示系统的因果涌现。


因果涌现：动力学可逆性 Dynamical Reversibility

基于可逆性的因果涌现理论是一种量化因果涌现的新框架，提出了近似动力学可逆性（[math]\displaystyle{ \Gamma_{\alpha} }[/math]）的概念，用于量化马尔科夫动力学接近可逆动力学的程度[18]。 该理论的核心思想是指出所谓的因果涌现其实等价于动力学可逆性的涌现。

直观地说，如果将一个马尔科夫动力学看作是一个通信信道，那么，每一状态下的概率转移就看作是一个信息通路（information pathway），于是该系统发生因果涌现就等价于系统中存在着冗余的信息通路——存在着一些状态转移与其它状态转移相似或线性相关，而消除这些冗余通路的最好办法就是利用SVD找到等于零或接近于零的奇异值以及对应的奇异向量，这实际上正是最大化有效信息的精髓所在。事实上，我们可以证明，最大化宏观动力学有效信息的粗粒化方案近似等价于将系统的动力学向较大奇异值对应的奇异向量方向做投影的粗粒化方案。而且，系统所能达到的最大信息传输效率的提升（归一化的近似动力学可逆性），就可以定义为因果涌现的强度。

所以，给定一个系统的马尔科夫转移矩阵，通过对它进行奇异值分解，将奇异值[math]\sigma_i[/math]的[math]\displaystyle{ \alpha }[/math]次方之和定义为马尔科夫动力学的近似可逆性度量（[math]\displaystyle{ \Gamma_{\alpha}\equiv \sum_{i=1}^N\sigma_i^{\alpha} }[/math]），用以刻画最终指标偏向于衡量确定性或简并性的程度。通常我们取[math]\alpha=1[/math]，表示平权考虑确定性或简并性对最终可逆性度量的影响。

该指标与有效信息具有高度的相关性，因而也可以用于刻画动力学的因果效应强度。基于可逆性的因果涌现理论与已有的其他因果涌现理论最大的区别就是不需要指定粗粒化策略，仅从马尔科夫转移矩阵的奇异值谱就能判断因果涌现的发生，直接定义所谓清晰因果涌现（clear causal emergence）和模糊因果涌现（vague causal emergence）的概念，在一定程度上规避了因果涌现结果依赖于粗粒化策略的问题和最大化有效信息的难以优化求解的问题。其次，这套理论指出一个动力系统的因果涌现特性仅仅是系统线性化算符的奇异值谱决定的，因而是一套更加反映系统动力学特性本身的因果涌现理论。相关详情请参照基于可逆性的因果涌现理论词条。


定义

对于给定的马尔可夫链[math]\displaystyle{ \chi }[/math]和对应的转移概率矩阵(Transitional Probability Matrix，简称TPM) [math]\displaystyle{ P }[/math] ，如果 [math]\displaystyle{ P }[/math] 同时满足：

  1. P是可逆矩阵，即存在矩阵 [math]\displaystyle{ P^{-1} }[/math]，使得 [math]\displaystyle{ PP^{-1}=I }[/math]，[math]\displaystyle{ I }[/math] 为单位矩阵；
  2. [math]\displaystyle{ P^{-1} }[/math] 也是另一个马尔可夫链 [math]\displaystyle{ \chi^{-1} }[/math] 的有效TPM，

则 [math]\displaystyle{ \chi }[/math] 和 [math]\displaystyle{ P }[/math] 可以称为严格动力学可逆的。


关于动力学可逆性，我们有如下定理：对于一个给定的马尔科夫链[math]\displaystyle{ \chi }[/math]和对应的TPM [math]\displaystyle{ P }[/math]，当且仅当 [math]\displaystyle{ P }[/math] 是置换矩阵(Permutation Matrix)的时候，[math]\displaystyle{ P }[/math] 是严格动力学可逆的。所谓的置换矩阵是指每一个行向量都是独热向量（one-hot vector，即只有一个元素是1，其余元素均为零的向量），每两个行向量都不相同。

纯粹的置换矩阵在所有可能的TPM中非常稀少，所以大多数的TPM并不是严格动力学可逆的。因此，需要一个指标来刻画任意一个TPM接近动力学可逆的程度。

考虑 [math]\displaystyle{ P }[/math] 的秩为 [math]\displaystyle{ r }[/math] ，当且仅当 [math]\displaystyle{ r\lt N }[/math] （ [math]\displaystyle{ N }[/math] 为矩阵的维数）的时候，[math]\displaystyle{ P }[/math] 是不可逆的；且 [math]\displaystyle{ P }[/math] 越退化对应着越小的 [math]\displaystyle{ r }[/math]。然而，非退化（满秩）的矩阵 [math]\displaystyle{ P }[/math] 并不总是动力学可逆的，因为：

  1. 尽管[math]\displaystyle{ P^{-1} }[/math]存在，[math]\displaystyle{ P^{-1} }[/math]并不一定是满足归一化条件的合法TPM（即矩阵中每个元素都大于等于0，小于等于1，且每一行满足加和为1的条件）；
  2. 如前所述，若 [math]\displaystyle{ P }[/math] 满足动力学可逆性，则 [math]\displaystyle{ P }[/math] 必为置换矩阵。


所有置换矩阵的行向量都是独热向量。这一特性可以被矩阵 [math]\displaystyle{ P }[/math] 的弗罗贝尼乌斯范数（Frobenius norm，F-norm）刻画。事实上，当且仅当P的行向量是独热向量的时候，矩阵P的弗罗贝尼乌斯范数取最大值。因此，我们可以由矩阵P的秩r和矩阵的弗罗贝尼乌斯范数找到P的近似动力学可逆性与矩阵奇异值之间的联系。

首先，矩阵的秩可以被写作：[math]\displaystyle{ r=\sum_{i=1}^{N}\sigma_{i}^{0} }[/math]，其中[math]\displaystyle{ \sigma_{i} }[/math]是矩阵 [math]\displaystyle{ P }[/math] 的第 [math]\displaystyle{ i }[/math] 个奇异值，[math]\displaystyle{ N }[/math] 是马尔科夫链的状态数。

紧接着，矩阵的F-norm可以被写作：[math]\displaystyle{ {||P||}_{F}^{2}=\sum_{i=1}^{N}\sigma_{i}^{2} }[/math]。这也是所有奇异值的平方和。可以看出矩阵的秩和F-norm都与奇异值相联系。

由于绝大部分马尔科夫链的状态转移矩阵都不是置换矩阵，所以我们需要一个指标来刻画其靠近（可逆的）置换矩阵的程度。这就是矩阵 [math]\displaystyle{ P }[/math] 的 [math]\displaystyle{ \alpha }[/math] 阶近似动力学可逆性。定义为：

[math]\displaystyle{ \begin{aligned} \Gamma_{\alpha}=\sum_{i=1}^{N}\sigma_{i}^{\alpha}\end{aligned} }[/math]

其中[math]\displaystyle{ \alpha\in(0,2) }[/math]是参数。实际上，当[math]\displaystyle{ \alpha\ge1 }[/math]时，[math]\displaystyle{ \Gamma_{\alpha} }[/math]是 [math]\displaystyle{ P }[/math] 的沙滕范数（Schatten norm）；当[math]\displaystyle{ 0\lt \alpha\lt 1 }[/math]时，[math]\displaystyle{ \Gamma_{\alpha} }[/math]是P的准范数（quasinorm）[19][20][21][22]。使用这个定义来刻画近似动力学可逆性是合理的，因为完全动力学可逆性可以通过最大化[math]\displaystyle{ \Gamma_{\alpha} }[/math]来得到。我们可以证明如下结论：对于任意[math]\displaystyle{ \alpha\in(0,2) }[/math]，[math]\displaystyle{ \Gamma_{\alpha} }[/math]的最大值是N，当且仅当P是置换矩阵的时候能取到该最大值。


最大化动力学可逆性的划分

文章[18]提供了一种基于奇异值分解的粗粒化方法，以获得宏观层面的简化TPM。其基本思想是将 [math]\displaystyle{ P }[/math] 中的行向量 [math]\displaystyle{ P_{i},\forall i \in [1,N] }[/math]投影到 [math]\displaystyle{ P }[/math] 的几个较大奇异值对应的奇异向量所张成的子空间上，从而保留P的主要信息，并保持 [math]\displaystyle{ \Gamma }[/math] 不变。

最大化动力学可逆性的合理性来源于这套SVD方案定义的因果涌现与最大化有效信息的因果涌现定义的密切关系。在Erik Hoel因果涌现理论的框架内，给定马尔可夫动力学系统的因果涌现的大小依赖于粗粒化策略的选择，可以通过宏观动力学的EI最大化寻找最优的粗粒化策略。

我们可以证明，要想获得最大化有效信息的粗粒化策略也近似等价于是向最大奇异值对应的奇异向量方向投影的策略。严格推导可见基于有效信息的因果涌现理论。由于对于离散的马尔可夫链来说，粗粒化策略往往都是比较严格的0或1构成的向量，因此并不能简单地根据奇异向量来直接投影。但是，如果粗粒化策略尽可能地与最大的几个奇异向量平行，那么就可以让粗粒化后动力学的EI尽可能更大，因此最大化EI的必要条件近似可以认为就是往最大奇异值对应的奇异向量空间做投影。在实际应用中，粗粒化策略 [math]\displaystyle{ \phi }[/math] 的约束（例如，宏观动力学TPM的分组和归一化）会阻止 [math]\displaystyle{ \phi_j }[/math] 与奇异向量完美对齐。因此，基于SVD的粗粒化策略只是最大化EI的近似必要条件。


划分方法

当我们想把 [math]\displaystyle{ N }[/math] 个状态的马尔科夫链划分成 [math]\displaystyle{ M }[/math] 组时，待选的划分数量可以由 第二类斯特林数（Stirling numbers of the second kind）表示，[math]\displaystyle{ S(N, M) = \frac{1}{M!} \sum_{k=0}^{M} (-1)^{M-k} \binom{M}{k} k^N }[/math]

待选的划分数量极大，比如[math]\displaystyle{ S(20,10)=1.107e32 }[/math]，几乎不可能通过遍历所有待选划分来寻找最优划分，所以，我们需要使用一些方法来有效的寻找到近似最优解。


逐层贪婪算法

在这个问题中，我们可以采取 自底向上合并（Bottom-Up）的层次聚类（Hierarchical Clustering）方法，即从最细粒度（每个状态单独成组，一共 [math]\displaystyle{ N }[/math] 组）开始，每一步合并两个最相似的组，直到达到目标分组数 [math]\displaystyle{ M }[/math] 。

具体步骤如下：

  1. 初始化：每个状态单独成组，[math]\displaystyle{ A_0={{1},{2},…,{N}} }[/math]
  2. 在每一层 [math]\displaystyle{ k }[/math] 内：
    1. 遍历所有组对 [math]\displaystyle{ G_p }[/math] 和 [math]\displaystyle{ G_q }[/math]：
      1. 合并为新组 [math]\displaystyle{ G_{new} = G_p \cup G_q }[/math]
      2. 获得新划分 [math]\displaystyle{ =A_k \setminus \{ G_p , G_q \} \cup G_{new} }[/math]
      3. 根据此划分获得新动力学 [math]\displaystyle{ P' }[/math]
      4. 计算 [math]\displaystyle{ P' }[/math] 对应的划分指标
    2. 更新划分 [math]\displaystyle{ A_{k} }[/math] 为该层最优的划分
  3. 当剩余分组数 [math]\displaystyle{ |A_k| = M }[/math] 时停止

该算法的复杂度为 [math]\displaystyle{ O \left (\sum_{k=M}^N k^2 \right ) }[/math]，介乎[math]\displaystyle{ O \left ( N^2 \right ) }[/math]和[math]\displaystyle{ O \left ( N^3 \right ) }[/math]之间。

该算法的缺点在于其容易陷入局部最优。一旦进入了某一分支，就无法跳出。为了找到更好的解，我们可以在其基础上使用遗传算法等方法帮助跳出局部最优。


谱分解

谱分解（Spectral Clustering）的思想是利用矩阵的谱，即特征向量，所展示的物理含义进行聚类和划分。


最大化有效信息的谱分解方法

该方法是将原始网络对应的转移矩阵做特征值分解，并使用这些特征向量来对节点进行聚类，之后再归并网络。

输入：原始包含 [math]\displaystyle{ N }[/math] 个状态的马尔科夫链 [math]\displaystyle{ P }[/math]和距离超参[math]\displaystyle{ \epsilon }[/math]；

  1. 对 [math]\displaystyle{ P }[/math] 进行矩阵的特征值分解，得到特征值集合 [math]\displaystyle{ \Lambda=\{\lambda_i\}^N_{i=1} }[/math] 与特征向量集合 [math]\displaystyle{ E=\{e_i\}^N_{i=1} }[/math]；
  2. 构建新的有偏的特征向量集合 [math]\displaystyle{ E’=\{\lambda_ie_i|\lambda_i≠0\}^N_{i=1} }[/math]，该集合中的每个向量被其对应的特征值所加权，同时去除了零特征值和对应的特征向量（新的状态数量为[math]\displaystyle{ N'=rank(A) }[/math]）。直观地说，忽略特征值为0的特征向量是有意义的，因为它对应了动力学中的简并性，并且非零特征值和相应的特征向量包含了丰富的网络拓扑结构信息；
  3. 依据[math]\displaystyle{ E' }[/math]计算状态间的相似矩阵[math]\displaystyle{ D_{N'×N'} }[/math]：
    1. 如果状态 [math]\displaystyle{ v_i }[/math] 和 [math]\displaystyle{ v_j }[/math] 分别在对方的马尔可夫毯中，则通过新的特征向量计算两个节点的cosine相似性作为两者间的相似度 [math]\displaystyle{ d_{ij} }[/math] 和 [math]\displaystyle{ d_{ji} }[/math]；
    2. 否则将两个节点间的相似度设为无穷大 [math]\displaystyle{ \infty }[/math]（可以设个比较大的值，如10000）；
  4. 基于相似度矩阵[math]\displaystyle{ D_{N'×N'} }[/math]和一个超参[math]\displaystyle{ \epsilon }[/math]（需要线性搜索，可以选择EI最大的参数），使用OPTICS算法（是一种基于密度的划分算法，旨在识别数据集中不同密度的划分结构）进行划分，输出对应超参[math]\displaystyle{ \epsilon }[/math]下的划分方式。

该方案的时间复杂度为[math]\displaystyle{ O(N^3) }[/math]。

缺点：

  1. 对大规模动力学进行粗粒化时，对转移矩阵进行特征值分解的计算复杂度也较高
  2. 聚类算法所需的距离超参[math]\displaystyle{ \epsilon }[/math]需要线性搜索，需要加入人为的先验知识


寻找成块划分的谱分解方法

当马尔科夫矩阵存在块状结构，或者状态明显可以被划分为几类时，根据这样的划分，该矩阵就是成块的。对于简单的[math]\displaystyle{ P }[/math]，如图6中的[math]\displaystyle{ \bar{P} }[/math]所示，我们能够很快的找出它的成块划分[math]\displaystyle{ \Phi }[/math]。然而，在实际应用中，我们经常遇到的矩阵可能并不完美。例如，成块的矩阵的状态排序可能被打乱（如图（6）中的[math]\displaystyle{ P_1 }[/math]），或者矩阵包含了如[math]\displaystyle{ P_2 }[/math]的噪声（如图6中的[math]\displaystyle{ P }[/math]，[math]\displaystyle{ P = P_1 + P_2 }[/math]，[math]\displaystyle{ P_1^TP_2 = 0 }[/math]）。在这种情况下，我们通常面临的是像[math]\displaystyle{ P }[/math]这样的矩阵。我们既无法确定它是否成块，也无法决定它的划分，甚至不知道它的马尔科夫秩。

图6：Zhang[4] 文章中的示意图。图中左面四个矩阵都是成块的马尔科夫矩阵，而右面的P_2是一个噪声矩阵，(P_1)^T P_2 = 0


针对这一问题，Anru Zhang[4][4]提出了一种寻找最优划分的方法。通过这种方法，我们可以找到最优的划分，并基于这个划分判断一个马尔科夫矩阵是否成块。具体步骤如下：

  1. 获取一个n*n维的马尔科夫矩阵P，或者是马尔科夫链的频率采样；
  2. 获取[math]\displaystyle{ r\lt n }[/math] 的马尔科夫秩，或指定一个[math]\displaystyle{ r }[/math]；
  3. 对P进行SVD分解，[math]\displaystyle{ P = U \Sigma V^T }[/math]，其中 [math]\displaystyle{ U }[/math] 为左奇异向量，[math]\displaystyle{ V }[/math] 为右奇异向量，[math]\displaystyle{ \Sigma = diag(\sigma_{1},\sigma_{2},...,\sigma_{N}) }[/math] 是一个对角矩阵，包含所有有序奇异值。
  4. 通过下列公式得到最优划分 [math]\displaystyle{ \Omega_1, \Omega_2, ... , \Omega_r }[/math]

[math]\displaystyle{ \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]\displaystyle{ v_s }[/math]是第[math]\displaystyle{ s }[/math]类[math]\displaystyle{ \Omega_s }[/math]的特征向量。这个算法的核心思想是：在最优的划分中，[math]\displaystyle{ \Omega_s }[/math]中的状态[math]\displaystyle{ i }[/math]的右奇异向量会和[math]\displaystyle{ v_s }[/math]距离最小，也就是说，让每个微观态的右特征向量和它对应的宏观态的右特征向量尽可能接近。

这个算法看起来非常熟悉，实际上它与 k-Means 聚类算法非常相似。让[math]\displaystyle{ v_s }[/math]和组里的[math]\displaystyle{ V[i:r] }[/math]距离最小的方法，就是让[math]\displaystyle{ v_s }[/math]成为[math]\displaystyle{ V[i:r] }[/math]这若干个点的中心。回顾一下k-Means算法把n个点聚成r类的目标函数，其中第[math]\displaystyle{ i }[/math]个点被归到第[math]\displaystyle{ c^i }[/math]类，而[math]\displaystyle{ \mu_j }[/math]为第[math]\displaystyle{ j }[/math]类点的中心点：

[math]\displaystyle{ J(c^1, ... , c^n, \mu_1, ... , \mu_r) = \sum_{i=1}^n || x_i - \mu_{c^i} ||^2 }[/math]

我们看到这两个目标函数非常相似。因此，我们可以通过上述算法或 k-Means 聚类来获取最优的状态划分，进而判断马尔科夫矩阵的成块性。

论文[4]中并未提供使用右奇异向量向量进行聚类的物理原因，但提供了该算法的解与最优解的差距的上界。


最大化可逆性的谱分解方法

动力学可逆性的文章[18]提供了另一种基于奇异值分解的粗粒化方法，以获得宏观层面的简化TPM。其基本思想是将 [math]\displaystyle{ P }[/math] 中的行向量 [math]\displaystyle{ P_{i},\forall i \in [1,N] }[/math] 投影到 [math]\displaystyle{ P }[/math] 的几个较大奇异值对应的奇异向量所张成的子空间上，从而保留 [math]\displaystyle{ P }[/math] 的主要信息，并保持 [math]\displaystyle{ \Gamma }[/math] 不变。

跟上面的方法类似，该方法的基本思路是将P中的所有行向量 [math]\displaystyle{ P_{i} }[/math] 视为维数为 [math]\displaystyle{ N }[/math] 的数据向量。首先对这些行向量进行PCA降维，其次将其聚类为 [math]\displaystyle{ r }[/math] 个簇，其中 [math]\displaystyle{ r }[/math] 是根据奇异值谱的阈值[math]\displaystyle{ \epsilon }[/math]截取的：

  1. 对P进行SVD分解：[math]\displaystyle{ P=U\cdot \Sigma \cdot V^{T}, }[/math]
  2. 选择一个[math]\displaystyle{ \epsilon }[/math]作为阈值来截断奇异值谱，并得到[math]\displaystyle{ r_{\epsilon} }[/math]作为保留状态的个数；
  3. 通过计算[math]\displaystyle{ \tilde{P}\equiv P\cdot V_{N\times r_{\epsilon}} }[/math]对P中的所有[math]\displaystyle{ P_{i} }[/math]进行降维，其中[math]\displaystyle{ V_{N\times r_{\epsilon}} }[/math]由[math]\displaystyle{ P\cdot P^{T} }[/math]的前[math]\displaystyle{ r_{\epsilon} }[/math]个特征向量构成；
  4. 通过 K-means 算法将 [math]\displaystyle{ \tilde{P} }[/math] 中的所有行向量聚类为 [math]\displaystyle{ r }[/math] 组，得到投影矩阵[math]\displaystyle{ \Phi }[/math]，其定义为：

[math]\displaystyle{ \Phi_{ij} =\begin{cases} 1, & \text{如果}\tilde{P_{i}}\text{属于第r组}\\ 0, & \text{其他情况} \end{cases} }[/math]


其他方法

寻找最优划分本质上是个组合优化问题，所以我们也能使用其他合适的方法来解决该问题。

首先是遗传算法。它是一种受生物进化启发的全局优化算法，通过模拟自然选择（选择适应度高的个体作为父代）、交叉（通过父代基因重组生成子代）和变异（以一定概率随机改变子代基因）等机制，在解空间中搜索最优解。其核心流程包括：

  1. 种群初始化：随机生成候选解的群体。
  2. 在每一步迭代过程中：
    1. 选择：选择最大化目标（划分标准）的个体作为父代。
    2. 交叉：通过父代基因重组生成子代。
    3. 变异：以一定概率随机改变子代基因。

所以，我们可以把上述的逐层贪婪算法和谱分解方法的解作为初始解，通过遗传算法，帮助跳出局部最优解，并寻找到新的全局更优解。


另外，我们也能使用梯度下降方法，搭配深度学习方法，比如图神经网络，来寻找最优划分。比如，在Klein等人的论文[23]，以及Griebenow 等人[24]的文献中，把粗粒化方案写成一个分组矩阵，同时将这个分组矩阵参数化，这样把最终的粗粒化网络表达为分组矩阵与原始网络邻接矩阵的乘积，并由此计算EI。之后，使用自动微分方法执行梯度下降算法从而优化分组矩阵，以使得EI最大化。

输入：包含[math]N[/math]个节点的原始网络[math]G[/math]，其对应邻接矩阵为：[math]\displaystyle{ A }[/math]，粗粒化后的网络所包含的节点数：[math]\displaystyle{ K }[/math]；输出：宏观网络[math]G'[/math]，对应的邻接矩阵为：[math]\displaystyle{ B }[/math]，以及对应的从[math]A[/math]到[math]B[/math]的粗粒化矩阵：[math]\displaystyle{ M }[/math]

  1. 针对一个含有[math]\displaystyle{ N }[/math]个节点的网络[math]\displaystyle{ A }[/math]，随机初始化一个分组矩阵[math]\displaystyle{ M\in \mathbb{R}^{N×K} }[/math]，[math]\displaystyle{ K }[/math]表示分组的大小，其中矩阵里面的每个元素[math]\displaystyle{ m_{iμ}=Pr⁡(v_i\in v_{\mu}) }[/math]，表示微观节点[math]\displaystyle{ v_i }[/math]属于宏观节点[math]\displaystyle{ v_{\mu} }[/math]的概率
  2. 根据微观网络和分组矩阵构建宏观网络 [math]\displaystyle{ B }[/math] ，具体计算方法分为两个步骤：1）[math]\displaystyle{ M^T A M }[/math]；2）归一化前一步骤得到的矩阵。
  3. 使用带动量的梯度下降方法优化 [math]\displaystyle{ M }[/math] ，优化目标是最大化宏观网络的有效信息EI。

时间复杂度：[math]\displaystyle{ O(N^3) }[/math]

缺点：

  1. 初始化分组矩阵的选择，多次重复实验会得到不一样的结果。
  2. 依赖神经网络的超参，如学习率、迭代次数等。


对划分的限制

在上述的算法介绍中，并未强调各方法对划分节点的限制。但在实际问题场景里，我们可能需要在寻找划分时加入限制来满足问题的特性。

比如，在社区检测问题中，具有只允许相邻的节点聚合在一起的约束。这要求状态划分必须遵循状态转移图的拓扑结构，同一聚合状态内的节点必须在状态转移图中直接相连或通过短路径可达（如K-邻域）。这种限制确保聚合状态具有空间连贯性，同时大幅减少可行划分数量。

在因果涌现分组问题[23]中，要求同属一组的状态需同属一个马尔科夫毯（Markov blanket）。在一个概率图中，一个节点的马尔科夫毯包括其父节点、子节点以及子节点的其他父节点（配偶节点）。这些节点共同作用，使得该节点与图中的其他节点条件独立，起到一个“信息防火墙”、确保组内状态动态仅通过固定接口与外部交互的作用。该约束确保：

  1. 因果独立性：当我们把整个组看作一个宏观状态时，由于组内所有状态的马尔科夫毯相同，意味着这个宏观状态与外部状态的交互（输入和输出）在统计上是一致的；
  2. 信息压缩与模式识别：具有相同马尔科夫毯的状态在系统中扮演相似的角色（即它们与其他状态之间的依赖关系相同）。将这些状态聚合成一个宏观状态，相当于识别出了系统中的一个功能模块或模式；
  3. 因果涌现：在复杂系统中，这种聚合可能导致宏观层面出现新的因果特性（即涌现）。如果聚合后的宏观状态能产生比单个状态更明确、更强大的因果作用，则称为因果涌现。此时，宏观状态的马尔可夫转移模型可以更有效地预测系统行为。

但也有其他的划分问题不需遵从这类限制。比如在学术网络、网页流网络的场景上的模式识别时，我们可能会判断某些从未合作过的学者具有相似的行为模式，而把他们划分在一组内。


聚合方法

当我们决定了最优的划分之后，就能根据该划分得到聚合后的矩阵了。

假设划分是把若干个状态合成一组，而其他状态各为一组，我们要做的是对原本的矩阵内，待聚合状态对应的若干行和若干列聚合成一行和一列，而其他元素则保持不动。

矩阵聚合.png


把一个 [math]\displaystyle{ N \times N }[/math]的矩阵聚合成 [math]\displaystyle{ N' \times N' }[/math]的矩阵 可以分成两步：

  1. 把待聚合的列直接进行相加，得到 [math]\displaystyle{ N \times N' }[/math] 的矩阵；
  2. 把待聚合的行进行加权平均，得到 [math]\displaystyle{ N \times N }[/math] 的矩阵


入流的聚合方法：直接加和

我们可以把马尔科夫链看做一个状态为节点，转移概率为边。

我们把待聚合节点都放进一个大节点的框 [math]\displaystyle{ M }[/math] 里，会发现每个外面的节点转移进大节点的概率就等于转移进 [math]\displaystyle{ M }[/math] 中所有节点的概率的和，即[math]\displaystyle{ p_{iM} = \sum_{j \in M} p_{ij} }[/math]。如此就能满足每行加和为一的条件。

入流聚合.png

所以，对入流（若干列）的聚类方法就是把这若干列直接加和，得到的列等于各节点转移到大节点的概率。


出流的聚合方法：加权平均

计算大节点的入流比较简单，而计算大节点的出流会相对较复杂，需要考虑待聚合节点之间的随机游走和状态分布，来计算流出到每个节点的概率。我们可以把出流的聚合看作是对若干行的线性加权平均。


权重一致的平均

在现有的方法中，大多采用最简单的思路：对状态进行权重一致的平均。当我们只看矩阵本身，认为所有的状态都是平等的，应该分配同等的权重，即[math]\displaystyle{ p_{Mj} = \frac{1}{n} \sum_{j \in M} p_{ji} }[/math]。

但是，这种情况只适用于成块的状态划分：由于成块性的同属于一组的微观状态转移到另一组的概率相同的性质，这样的平权平均是没有问题的。但是，当我们对不成块的状态划分进行平权平均时，就会导致宏观和微观的马尔科夫动力学不一致的情况。


权重为平稳分布的平均

当一个随机游走粒子进入到大节点时，它下一步游走到大节点里或大节点外的其他节点的概率时，应考虑其当前处于大节点内的哪个小节点的概率。

设游走粒子在 [math]\displaystyle{ t }[/math] 时刻处于小节点 [math]\displaystyle{ j }[/math] 的概率为[math]\displaystyle{ s_j(t) }[/math]，其下一步游走到其他粒子的概率应为[math]\displaystyle{ p_{Mj}(t) = \frac{1}{\sum_{j \in M} s_j(t)} \sum_{j \in M} s_j(t)\ p_{ji} }[/math] 。

假设游走已达到了平稳分布[math]\displaystyle{ \pi }[/math]，游走粒子处于不同小节点 [math]\displaystyle{ j }[/math] 的概率为[math]\displaystyle{ \pi_j }[/math]，所以其下一步游走到其他粒子的概率应为[math]\displaystyle{ p_{Mj} = \frac{1}{\sum_{j \in M} \pi_j} \sum_{j \in M} \pi_j\ p_{ji} }[/math]，而留在大节点的概率为[math]\displaystyle{ p_{MM} = \frac{1}{\sum_{j \in M} \pi_j} \sum_{j \in M} \pi_j \ p_{jj} }[/math] 。按照权重为平稳分布对若干行进行加权平均，得到的行仍会满足加和为一的条件。


为了说明如何获得简化的TPM，首先定义一个矩阵，称为稳态流矩阵：[math]\displaystyle{ F_{ij} \equiv \pi_i \cdot P_{ij}, \, \forall i,j \in [1, N], }[/math]。其中，[math]\displaystyle{ \pi }[/math] 是 [math]\displaystyle{ P }[/math] 的稳态分布，即满足[math]\displaystyle{ \pi P = \pi }[/math]。

其次，我们将根据 [math]\displaystyle{ \Phi }[/math] 和 [math]\displaystyle{ F }[/math] 得到粗粒化后的稳态流矩阵：

[math]\displaystyle{ F' = \Phi^T \cdot F \cdot \Phi, }[/math]

其中，[math]\displaystyle{ F' }[/math] 是粗粒化后的平稳流矩阵。在归并的过程中，保持所有矩阵上的流量总量是不变的。

最后，粗粒化后的 [math]\displaystyle{ P' }[/math] 可直接通过归一化每一行得到：

[math]\displaystyle{ P'_i = \frac{F'_i}{ \sum_{j=1}^{N} (F'_i)_j }, \, \forall i \in [1, N]. }[/math]

将 [math]\displaystyle{ N }[/math] 个 [math]\displaystyle{ P'_i }[/math] 拼在一起就能得到最终的 [math]\displaystyle{ P' }[/math] 。可以证明，这样的粗粒化方案可以保证让粗粒化操作和时间演化操作具有可交换的性质（证明可以参考[18]）。

按照该方法进行加和平均，得到的宏观和微观的马尔科夫动力学就会满足一致的条件。虽然该方法直观明了，该方法的问题在于其需要在系统演变到平稳分布前，计算每个 [math]\displaystyle{ t }[/math] 时刻的状态分布概率 [math]\displaystyle{ s_j(t) }[/math]。在达到平稳分布前直接使用[math]\displaystyle{ \pi }[/math]来聚合还是不会满足一致性，所以比较麻烦，而HOM的计算方式则可以解决这个问题。


高阶网络聚类 Higher-Order Macronodes (HOM)

在Klein等人的论文[23]，以及Griebenow 等人[24]的文献中，作者们采用了高阶依赖项建模（HOMs）[25]来对网络节点进行聚合。其目的是为了保证聚合后的宏观网络与原始的微观网络具有相同的随机游走动力学特性。该方法也可以用于对马尔科夫链状态进行合并。

Higher-order Macronodes，是由Jian Xu等人[25]提出，根据给定的微观马尔科夫链和节点划分方案，得出一个和微观转移矩阵保持一致的宏观转移矩阵的方法。下面将简单介绍构建宏观 HOM 转移概率矩阵的计算方法。具体方法请参考网络归并。

假设我们把若干个节点视为待合并的节点[math]\displaystyle{ A }[/math]。

  1. 计算合并节点的入流：对于其他节点到待合并节点[math]\displaystyle{ A }[/math]的连边，直接加和即可。
  2. 计算合并节点的出流：对于待合并节点[math]\displaystyle{ A }[/math]到其他节点的连边，需要考虑三种情况：
    1. 当待合并节点[math]\displaystyle{ A }[/math]间存在连边时；
    2. 当待合并节点[math]\displaystyle{ A }[/math]指向同一个输出节点时；
    3. 当待合并节点[math]\displaystyle{ A }[/math]指向不同输出节点时；
  3. 待合并节点[math]\displaystyle{ A }[/math]的出流的计算方式都不一样，如图（7）所示：
图7：HOM节点边权合并示意图

论文[23]中提到，原始的微观动力学和使用该方法聚合后的宏观动力学和能满足其定义的一致性。


不同粗粒化方法之间的异同

高阶网络聚类 Higher-order Macronodes v.s. 成块性

我们通过一个简单的成块划分例子来展示，成块性的粗粒化与HOM是如何对应的。该例子的状态转移图如下图所示：

图6：例子示意图

与之相应的状态转移概率矩阵为： [math]\displaystyle{ \begin{aligned} P = P_{s_i \rightarrow s_j} = \left [ \begin{array}{c:c:cc:c:c} 0 & 0 & 0.3 & 0.1 & 0.6 & 0\\ \hdashline 0 & 0 & 0.2 & 0.3 & 0 & 0.5\\ \hdashline 0 & 0 & 0 & 0.7 & 0.2 & 0.1\\ 0 & 0 & 0.7 & 0 & 0.2 & 0.1\\ \hdashline 1.0 & 0 & 0 & 0 & 0 & 0\\ \hdashline 0 & 1.0 & 0 & 0 & 0 & 0\\ \end{array} \right ] \end{aligned} }[/math] 。对于这个转移概率矩阵，[math]\displaystyle{ \{1, 2, \{34\}, 5, 6\} }[/math]是一个成块的划分。其对应的宏观转移概率： [math]\displaystyle{ \begin{aligned} \left [ \begin{array}{c:c:c:c:c} 0 & 0 & 0.4 & 0.6 & 0\\ \hdashline 0 & 0 & 0.5 & 0 & 0.5\\ \hdashline 0 & 0 & 0.7 & 0.2 & 0.1\\ \hdashline 1.0 & 0 & 0 & 0 & 0\\ \hdashline 0 & 1.0 & 0 & 0 & 0\\ \end{array} \right ] \end{aligned} }[/math]


根据上一部分对HOM方法介绍，我们注意到，HOM中各情况的计算方式本质上都是对待合并节点[math]\displaystyle{ A }[/math]中各节点的出流[math]\displaystyle{ W_i^{out} }[/math]（即图中的红色连边）的加权求和。而成块性的特性决定了，群组里的各节点到其他群组的出流[math]\displaystyle{ W_i^{out} }[/math]都是一样的。所以，对相同的[math]\displaystyle{ W_i^{out} }[/math]进行任意加权平均的结果仍然是[math]\displaystyle{ W_i^{out} }[/math]。这意味着，无论根据哪种情况采用哪种加权方法，对于成块划分来说，最终的结果都是一致的，且与成块性的粗粒化结果相同。

由此，我们总结，成块性的聚合方式和HOM的计算方式，无论入流还是出流都完全一致。


成块性 v.s. 有效信息/动力学可逆性

当一个马尔科夫矩阵满足成块性时，同一组的行会受到约束而相似，所以这种成块的矩阵的EI或动力学可逆性也会很高。这就是成块性和因果涌现之间的物理联系。

但是，这两个概念并不完全等价。成块性的目标是想在保证宏观和微观动力学的一致性的情况下去对状态做分组，而因果涌现的主要目标是想最大化粗粒化之后宏观动力学的因果强度。

在普通的成块性例子里，它们粗粒化之后能得到一个EI或者可逆性很高的宏观动力学。但是，当我们考虑一些极端的例子，如下图所示，我们强行让每一组里的状态转移概率都尽可能的不一样，但是加起来又是一样的，满足了成块性的划分。按照这种成块的划分来去做粗粒化，它并不一定得到高的EI或可逆性。相反，我们用另一个不成块的分组来做粗粒化，得到的宏观动力学的EI和可逆性反而会更高。这现象说明：

  1. 基于成块性的最优划分不一定等于基于可逆性的最优划分。
  1. 成块性的分组数量是基于马尔科夫秩的，而可逆性的分组数量是对奇异值谱的截断位置。





  1. ↑ Hafner, Danijar, et al. "Mastering atari with discrete world models." arXiv preprint arXiv:2010.02193 (2020).
  2. ↑ 2.0 2.1 2.2 2.3 Coarse graining. Encyclopedia of Mathematics. URL: http://encyclopediaofmath.org/index.php?title=Coarse_graining&oldid=16170
  3. ↑ 3.0 3.1 Buchholz, Peter. "Exact and ordinary lumpability in finite Markov chains." Journal of applied probability 31.1 (1994): 59-75.
  4. ↑ 4.0 4.1 4.2 4.3 4.4 Zhang, Anru, and Mengdi Wang. "Spectral state compression of markov processes." IEEE transactions on information theory 66.5 (2019): 3202-3231.
  5. ↑ Singh, Satinder, Tommi Jaakkola, and Michael Jordan. "Reinforcement learning with soft state aggregation." Advances in neural information processing systems 7 (1994).
  6. ↑ Duan, Yaqi, Tracy Ke, and Mengdi Wang. "State aggregation learning from markov transition data." Advances in Neural Information Processing Systems 32 (2019).
  7. ↑ Coarse graining. Encyclopedia of Mathematics. URL: http://encyclopediaofmath.org/index.php?title=Coarse_graining&oldid=16170
  8. ↑ 8.0 8.1 概率论与统计学5——马尔科夫链（Markov Chain） - 做大饼馅儿的韭菜的文章 - 知乎 https://zhuanlan.zhihu.com/p/274775796
  9. ↑ 9.0 9.1 9.2 Kemeny, John G., and J. Laurie Snell. Finite markov chains. Vol. 26. Princeton, NJ: van Nostrand, 1969. https://www.math.pku.edu.cn/teachers/yaoy/Fall2011/Kemeny-Snell_Chapter6.3-4.pdf
  10. ↑ Crutchfield, James P. (1994). "The Calculi of Emergence: Computation, Dynamics, and Induction". Physica D: Nonlinear Phenomena. 75: 11–54.
  11. ↑ Shalizi, C. R.; Crutchfield, J. P. (2001). "Computational Mechanics: Pattern and Prediction, Structure and Simplicity". Journal of Statistical Physics. 104 (3/4): 817–879.
  12. ↑ Barnett, L., & Seth, A. K. (2023). Dynamical independence: discovering emergent macroscopic processes in complex dynamical systems. Physical Review E, 108(1), 014304
  13. ↑ 13.0 13.1 13.2 13.3 13.4 Rosas, F. E., Geiger, B. C., Luppi, A. I., Seth, A. K., Polani, D., Gastpar, M., & Mediano, P. A. (2024). Software in the natural world: A computational approach to hierarchical emergence. arXiv preprint arXiv:2402.09090.
  14. ↑ B. C. Geiger and C. Temmel, “Lumpings of Markov chains, entropy rate preservation, and higher-order lumpability,” J. Appl. Probab., vol. 51, no. 4, pp. 1114–1132, Dec. 2014, extended version: arXiv:1212.4375.
  15. ↑ Hoel E P, Albantakis L, Tononi G. Quantifying causal emergence shows that macro can beat micro[J]. Proceedings of the National Academy of Sciences, 2013, 110(49): 19790-19795.
  16. ↑ Hoel E P. When the map is better than the territory[J]. Entropy, 2017, 19(5): 188.
  17. ↑ 17.0 17.1 Hoel, Erik P.; Albantakis, L.; Tononi, G. (2013). "Quantifying causal emergence shows that macro can beat micro". Proceedings of the National Academy of Sciences. 110 (49): 19790–19795.
  18. ↑ 18.0 18.1 18.2 18.3 Zhang, J., Tao, R., Leong, K.H. et al. Dynamical reversibility and a new theory of causal emergence based on SVD. npj Complex 2, 3 (2025).
  19. ↑ Schatten norm from Wikipedia. https://en.wikipedia.org/wiki/Schatten norm
  20. ↑ Recht, B., Fazel, M., Parrilo, P.A.: Guaranteed minimum-rank solutions of linear matrix equations via nuclear norm minimization. SIAM review 52(3), 471–501 (2010)
  21. ↑ Chi, Y., Lu, Y.M., Chen, Y.: Nonconvex optimization meets low-rank matrix factorization: An overview. IEEE Transactions on Signal Processing 67(20), 52395269 (2019)
  22. ↑ Cui, S., Wang, S., Zhuo, J., Li, L., Huang, Q., Tian, Q.: Towards discriminability and diversity: Batch nuclear-norm maximization under label insufficient situations. In: Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pp. 3941–3950 (2020)
  23. ↑ 23.0 23.1 23.2 23.3 Klein B, Hoel E. The emergence of informative higher scales in complex networks[J]. Complexity, 2020, 20201-12.
  24. ↑ 24.0 24.1 Griebenow R, Klein B, Hoel E. Finding the right scale of a network: efficient identification of causal emergence through spectral clustering[J]. arXiv preprint arXiv:190807565, 2019.
  25. ↑ 25.0 25.1 Xu, J., Wickramarathne, T. L., & Chawla, N. V. Representing higher-order dependencies in networks[J]. Science advances, 2016, 2(5), e1600028.
取自“https://wiki.swarma.org/index.php?title=马尔科夫链的粗粒化&oldid=45247”