下载掌阅APP，畅读海量书库
立即打开
本节分 3 个小节来详细地介绍经典的传播动力学模型、计算机模拟方法，以及常用理论分析方法。
在真实社会中，计算机病毒、信息、谣言、行为传播和金融风险扩散等现象都可以描述为复杂网络上的传播动力学 [1] 。根据描述对象的不同，网络传播动力学可分为生物传播、社会传播和社会—生物传播。生物传播主要刻画的是一类简单传播，即节点两次连续接触导致感染的概率是相同的，包括计算机病毒传播和传染性疾病传播等。相比之下，社会传播主要描述的是一类复杂传播，即再次接触时的感染概率依赖于先前的接触，包括革新采纳、健康活动和金融行为等。社会—生物传播则主要关注信息和疾病共同演化传播过程。生物传播和社会传播已经有一些经典的模型。例如，描述生物传播的易感态—感染态（ SI）模型、易感态—感染态—易感态（SIS）模型和易感态—感染态—恢复态（ SIR）模型，刻画社会传播的线性阈值模型（简称为阈值模型）。图 2-3 给出了 4 个经典传播模型的示意图。
（1）易感态—感染态（SI）模型[图 2-3（a）]。 SI模型主要用于描述一些致命的疾病传播，例如艾滋病。在任意时刻，节点只能处于易感态或感染态。在 t 时刻，每个感染态节点 i 以概率 β 将疾病传递给每个易感态邻居节点。节点一旦被感染，将无法康复。
图2-3 经典传播模型示意图
（a）SI模型，（b）SIS模型，（c）SIR模型，（d）阈值模型
（2）易感态—感染态—易感态（SIS）模型[图 2-3（b）]。 SIS模型是一个经典的可逆传播过程，可用于描述感冒这类容易康复，但可能被再次感染的疾病。与SI模型类似，在每个时刻，每个节点只能处于易感态或感染态。在 t 时刻，每个感染态节点 i 以概率 β 将疾病传递给每个易感态邻居节点，并且以概率 γ 恢复成易感态。疾病有效传播概率为 λ ＝ β/γ 。
（3）易感态—感染态—恢复态（ SIR）模型[图 2-3（ c）]。 SIR模型是经典的非可逆传播过程，它可用于刻画麻疹、水痘等这类康复后不再被感染的疾病。在任意时刻，每个节点只能处于易感态、感染态或恢复态。在 t 时刻，每个感染态节点 i 以概率 β 将疾病传递给每个易感态邻居节点，并且以概率 γ 变为恢复态。疾病有效传播概率为 λ ＝ β/γ 。节点一旦进入恢复态，将不再参与后续传播过程。
（4）阈值模型[图 2-3（ d）]。社会传播往往存在加强效应，即节点采纳某行为时，需同时考虑它的所有邻居节点状态。每个时刻，节点只能处于活跃态或非活跃态。每个非活跃态节点变为活跃态，当且仅当它的活跃邻居比例超过某一指定阈值。图 2-3（ d）中的采纳阈值 T ＝ 0.2。节点一旦变为活跃态，它的状态不再发生改变。
最终的传播范围和爆发阈值是研究网络传播动力学稳态时的两个重点关注问题。传播范围是指在稳态或终态时，网络中最终被感染的节点比例，它是统计物理所观察的序参量。爆发阈值是一个临界传播率，当传播率低于它时，最终传播范围为零，即系统处于非活跃态（吸收态）；否则，一定比例的节点被感染，即系统处于活跃态。从统计物理角度来看，爆发阈值为系统的相变点。图2-4 展示了传播范围随传播率 λ 的变化。系统存在一个临界点 λ c ，当 λ < λ c 时，传播范围为零（吸收态）；否则，最终传播范围非零（活跃态）。
图2-4 传播范围随传播率的变化
计算机模拟是分析网络传播动力学的一个重要研究手段。一方面，在计算模拟过程中不需要过多的假设；另一方面，计算机模拟能描述复杂的接触形式，并且能准确刻画各种因素对传播动力学的影响。由于真实的计算机病毒传播、流行病传播和行为传播往往需要付出惨痛的代价，因此只能借助计算机模拟，来探究它们的传播过程，最终控制传播。在第 2.2.2 节中，首先介绍同步更新和异步更新两个常用更新节点状态的方法，从而得到各类节点比例随时间的变化及传播范围。然后，介绍易感性、可变性和相对波动方差三种判断疾病模拟爆发阈值的方法。
传播动力学的两种常见模拟方法：同步更新和异步更新 [144，145] 。从不同的角度看待真实现象的演化过程，可用两种不同的方法更新个体状态。从长期来看，所有个体状态都在同时更新，即同步更新；从连续时间来看，系统在每个时刻仅发生一个事件，即个体状态是异步更新。下面以SIS模型为例来阐述两种不同的更新方法。
（1）同步更新。根据上一步状态，每个节点同时更新它的当前状态。运用同步更新方法模拟SIS模型，其过程如下：
（i）初始时刻，随机选择 ρ 0 比例的种子节点（感染态），剩余节点都处于易感态。用队列 Q 1 和 Q 2 分别存放当前感染节点和新感染的节点；
（ii）遍历队列 Q 1 中的每个感染态节点，并以概率 β 尝试感染它的所有易感态邻居。若成功感染易感态邻居，则将其存放到队列 Q 2 中；
（iii）遍历队列 Q 1 中的每个感染态节点，以概率 γ 恢复成易感态。若成功恢复，则从队列 Q 1 中移除；
（iv）将队列 Q 2 中的所有元素移至队列 Q 1 中，并清空队列 Q 2 中的所有元素；
（v）更新时间 t → t ＋1，重复第（ii）步到第（iv）步，直到 t ＝ t max 或者没有感染态节点为止。
同步更新的空间复杂度为 O （ N ＋ E ）。每个时间步的时间复杂度为 O （‹ k › N ）。
（2）异步更新。异步更新最典型的方法为Gillespie算法 [146] ，广泛地运用于博弈、传播和阈值模型等。运用异步更新来模拟SIS传播步骤如下：
（i）初始时刻随机选择
ρ
0
比例节点处于感染态，其余节点都处于易感态。用队列
和
分别存放感染态节点和活跃边。活跃边是指易感态节点和感染态节点之间的连边；
（ii）从队列
和
中选择一个事件的发生。对于事件
i
，它的发生概率为∏
i
＝ o
i
/ （
γN
A
（
t
）
＋ βN
E
（
t
））。其中，
o
i
表示事件
i
发生的概率，
N
A
（
t
）和
N
E
（
t
）分别表示在
t
时刻，系统中活跃节点的数量和活跃边的数量；
（iii）若发生事件为恢复节点
i
，则移除队列
中的节点
i
；移除队列
中节点
i
所连到的活跃边；添加节点
i
现在连接的活跃边到队列
中。若发生事件为感染节点
i
，则添加节点
i
到队列
中；移除队列
中的节点
i
所连到的活跃边；添加节点
i
现在的活跃边到队列
中；
（iv）更新当前系统的时间为 t → t ＋ τ 。其中， τ ＝ 1 /[ γN A （ t ） ＋ βN E （ t ）] [27] ；
（v）重复第（ ii）—（ iv）步，直到 t ＝ t max 或者没有感染态节点为止。
异步更新的空间复杂度为 O （ E ），每个时间步的时间复杂度为 O [（ N ＋ E ） 2 ]。
运用同步更新或异步更新方法，能得到稳态时最终的感染比例。从上述描述不难发现，两种更新方法存在一定的差异性，但它们之间也存在一定的联系。在每个时间步，更新 f 比例的节点。若 f ＝ 1 / N ，则为异步更新；若 f ＝ 1.0，则为同步更新 [21] 。对于SIR模型或其他传播动力学过程，可根据上述方法进行模拟。值得注意的是，运用同步更新和异步更新所得的传播范围有一定的差异（图2-5），并且相对应的理论方法也有所不同 [144，147] 。
图2-5 同步更新和异步更新模拟SIS模型
爆发阈值是传播动力学所关注的另一个重要问题。通常，模拟阈值也视为疾病爆发的真实阈值，可用于检验理论阈值的准确性 [28，147] 。因此，准确地判断模拟阈值至关重要。下面将介绍三种常用判断模拟阈值的判断方法。
（1）易感性。易感性广泛地运用于确定SIS传播和渗流等临界点，其表达式 [27] 为
其中， ρ 表示传播范围，即序参量。在爆发阈值处， x 呈现出峰值。
（2）可变性。可变性最早运用于确定磁化系统中的临界点 [148] 。 Shu等利用可变性方法来确定SIS模型和SIR模型的模拟阈值，其表达式 [28] 为
在临界点处， Δ 呈现出峰值。
（3）相对波动方差。相对波动方差可运用于确定SIR模型的模拟阈值和渗流临界点，其表达式 [149] 为
在临界点处， v R 呈现出峰值。
理论分析网络传播动力学过程受到了学者们的青睐。在运用理论方法描述复杂网络上的传播动力学时，我们会面临两个挑战。由于真实网络度分布异质性强、权重分布异质强、高集群系数、模块和社区等复杂的结构 [48，131，134，150] ，第一个挑战就是如何描述错综复杂的接触网络，第二个挑战是如何刻画相邻节点间的动力学关联性。动力学关联性是指一个节点被两个邻居感染存在关联性 [151] 。现已有众多理论方法描述网络上的传播动力学，常见的有异质平均场、淬火平均场、动态信息传递方法和渗流理论。
（1）异质平均场方法。 Pastor-Satorras和Vespignani最早运用异质平均场方法求解复杂网络上的传播动力学
[26]
。该理论方法假设度相同的节点不存在任何差异性。下面以SIS模型为例，介绍如何运用异质平均场求解传播范围和爆发阈值。在
t
时刻，度为
k
的节点处于感染态的概率记为
ρ
k
（
t
）。系统中处于感染态的节点比例为
ρ
（
t
）
＝
。当
t
→∞ ，
ρ
（
t
）即为最终的感染节点比例。在无关联网络中，易感态节点连接到一个感染态邻居的概率
Θ
（
t
），其表达式为
因此， ρ k （ t ）的演化方程为
在等式（2-9）的右边，第一项表示在 t 时刻度为 k 的感染态节点的恢复概率，第二项表示易感态节点被邻居感染的概率。
在 ρ k （0）→0 附近线性展开等式（2-9），可得
其中，
。 Jacobian矩阵
C
＝
，其元素为
其中，
为Dirac-Delta函数。在初始时刻，若
ρ
（
t
）指数增长，则全局疾病爆发，这意味着矩阵
C
的最大特征值‹
k
2
›/‹
k
› －1 大于零。因此，疾病爆发阈值为
其中，‹
k
›和‹
k
2
›分别表示度分布的一阶矩和二阶矩。在均匀网络（如ER网络）上，SIS疾病的爆发阈值为
＝ 1 /（‹
k
› ＋1）。对于度分布为幂率分布的无标度网络[即
，在热力学极限下且
γ
D
≤3 时，不存在疾病爆发阈值。因为，当
γ
D
>3 时，‹
k
2
›发散。异质平均场方法的广泛运用源于两个原因：①它仅需要度分布信息来预测传播范围和爆发阈值，②它能揭示网络结构对爆发阈值的影响。学者们已将异质平均场用于描述权重
[152]
、度关联
[153]
、结构多重性对传播
[154
-
156]
的影响。
（2）淬火平均场方法。由于异质平均场不能描述网络结构的所有信息，学者们利用邻接矩阵
A
来刻画接触网络所有特性。其他利用邻接矩阵的方法可归类为淬火平均场方法，例如离散马尔科夫链方法
[157]
和
N
缠绕方法
[158，159]
。在
t
时刻，易感态节点
i
被邻居感染的概率为
。其中，
ρ
j
（
t
） ＝ 1－
s
j
（
t
）表示在
t
时刻节点
j
处于感染态的概率。不难得出
ρ
i
（
t
）的演化方程为
等式（2-13）右边的第一项表示在
t
时刻节点
i
恢复的概率，第二项表示节点
i
在
t
时刻被感染的概率。在
t
时刻，处于感染态节点的比例为
ρ
（
t
）
＝
1
/ N
。初始时刻，只有少许节点处于感染态（即
ρ
i
（0）→0），对等式（2-13）在
ρ
i
（
t
）→0 处线性展开，写成矩阵形式为
其中，
ρ
i
（
t
）是向量
的第
i
个元素。利用与得到等式（2-12）类似的方法，疾病爆发阈值为
其中，∧
A
是矩阵
A
的最大特征值。利用淬火平均场得到的疾病爆发阈值仅依赖网络结构。在无关联的无标度网络中，当
γ
D
时，
∝‹
k
› /‹
k
2
›，与等式（2-12）相同；当
γ
D
>2.5 时，
，意味着在热力学极限下不存在爆发阈值
[160]
，与异质平均场所预测的值恰恰相反。有学者在文献
[33，37，161]
中讨论了异质平均场和淬火平均场所预测的理论阈值的差异性。
Goltsev等
[162]
定义了参与反比率（ IPR），用于量化特征值∧
A
所对应的特征向量的局域性。特征值∧
A
的参与反比率定义为
。其中，
f
i
（∧
A
）为∧
A
对应的特征向量
的第
i
个元素。若
为非局域特征向量，
v
（∧
A
）∝
O
（0）；若
v
（∧
A
）∝
O
（1），
为局域特征向量。 Goltsev等
[162]
指出：当
γ
D
<2.5 时，
为非局域特征向量，意味着当
时就有一定比例的节点被疾病感染；当
γ
D
>2.5 时，
为局域特征向量，意味着
时仅有中心节点和它的邻居被感染，由于波动性疾病最终缓慢消亡。因此，这种局域态并不能构成真正的活跃态，真实阈值更接近等式（2-12）。更多相关讨论见第3.1节。
Pastor-Satorras和Castellano
[161]
指出∧
A
的特征向量
具有不同的局域特性：当
γ
D
>2.5 时，
局域于中心节点；当
γ
D
<2.5 时，局域于
K
核值
[163]
高的节点。为了进一步理解淬火平均场失效的原因，需要理解
的物理意义，即特征向量中心性
[20]
。每个节点的特征向量中心性正比于它邻居特征向量中心性之和。中心节点具有高特征向量中心性，从而导致它们的邻居也具有高特征向量中心性。然而，这些中心节点的特征向量中心性又取决于它的邻居，导致中心节点的特征向量中心性被高估。类似地，可知易感态节点被感染的概率被高估了。根据等式（2-13），
ρ
i
（
t
）随
ρ
j
（
t
）增长，
ρ
j
（
t
）随
ρ
i
（
t
）增长。因此，同一条边来回地传递疾病导致“回音室（ echo chamber）”现象，最终使得感染态节点比例被高估
[164]
。
（3）动态信息传递方法。为了弥补淬火平均场方法的不足，同时又保留它的优点，即考虑网络的淬火结构，学者们发展了动态信息传递方法。 Karrer和Newman首次利用动态信息传递方法来求解复杂网络上的SIR模型 [165] 。最近，Shrestha等拓展后用于求解SIS模型 [166] 。在动态信息传递方法中，假设处于“空穴”态的节点不能传递信息给邻居，但允许邻居将信息传递给它。这消除了等式（2-13）所导致的“回音室”现象，并且考虑了一些动力学关联性。Karrer和Newman发现，动态信息传递方法能准确描述树形网络上的传播动力学 [165] 。基于动态信息传递方法， ρ i （ t ）的演化方程为
其中，
θ
j
→
i
（
t
）表示在
t
时刻，节点
j
被除
i
之外的其他邻居感染的概率。等式（2-16）右边的第一项表示节点
i
恢复的概率，第二项表示被邻居感染的概率。若节点
j
从感染态恢复，
θ
j
→
i
（
t
）将减小。若节点
j
被除
i
之外的其他邻居感染，
θ
j
→
i
（
t
）将增加，其概率为
λ
[1 －
ρ
j
（
t
）]
。其中，
N
j
表示节点
j
的邻居集合。考虑上述两个因素，
θ
j
→
i
（
t
）可写为
根据等式（2-16）和等式（2-17），可得每个节点的状态演化方程。系统需要 2 E ＋ N 个微分方程才能来描述。初始时刻， θ j → i （0）→0。在 θ j → i （0）＝ 0 处线性化等式（2-17），可得
其中，
B
为非回溯矩阵
[167]
。
θ
j
→
i
（
t
）为向量
的元素。矩阵
B
的元素为
其中，
δ
il
为Dirac-Delta函数。
的物理意义为：当
i
≠
l
，连边
l
→
h
可以影响连边
j
→
i
。与得到等式（2-12）的方法类似，可得疾病爆发阈值
，其表达式为
其中，∧ B 为矩阵 B 的最大特征值。
动态信息传递方法广泛地运用于传播动力学 [165，167] 、渗流 [168 - 170] 和级联失效 [171 - 173] 。动态信息传递方法广泛运用源于：①它用非回溯矩阵来描述网络结构，②它假设“空穴”态节点不能将信息传递给邻居，从而能描述部分动力学关联性。大量的实验模拟表示，动态信息传递方法能比较准确地描述SIS和SIR模型 [30] 。
然而，动态信息传递方法有两个缺陷：①方程复杂度高，②只有在树形网络上是精确的。为了弥补第一个缺陷，可以假设每条边没有传递信息给邻居的概率相同。因此，易感态节点连接到一个感染态邻居的概率为
若节点按照度大小分类，等式（2-21）可写为
将等式（2-22）代入等式（2-9），Barthélemy等发现，在无标度网络上，疾病传播呈现出层次性。此时，疾病爆发阈值
，其表达式为
为了弥补动态信息传递方法的第二个缺陷，需要减少有限环导致的“回音室”效应。 Radicchi和Castellano拓展了一套更为准确的动态信息传递方法，他们消除了三角形所导致的冗余感染，使得理论方法能更准确地描述真实网络 [164] 。对于真实网络中的度关联、模体和社区等复杂的网络拓扑结构，需要进一步发展准确的理论方法。
（4）边渗流理论。由于SIS是可逆传播过程，而SIR是非可逆传播过程，这导致理论求解SIR模型有一些特殊的理论方法。其中，边渗流理论是运用最为广泛的方法 [31，174，175] 。 SIR模型的传播范围可映射为求解边渗流理论中极大连通子图大小。假设感染态节点以概率 λ 将疾病传递给邻居，并在 τ 步后变为恢复态。那么，一条边的有效传递率为 ϕ ，即感染态节点在恢复前将疾病传递给易感态邻居的概率。若疾病传播随时间连续演化，可得
有效传递率 ϕ 可映射成边渗流中的边占有概率。对于树形网络结，图大小 g ，其值为
其中，
G
0
（
u
）
＝
为度分布
p
（
k
）的生成函数，
u
表示随机选择一条边没有连向极大连通子图的概率。为了求解
u
，需考虑两种不同情况：①边没有被占有，其概率为 1－
ϕ
，②被占有的边没有连到极大连通子图，其概率为
ϕ G
1
（
u
）。其中，
G
1
（
u
）
＝
为剩余度分布
Q
（
k
）＝
P
（
k
＋1）（
k
＋1）/‹
k
›的生成函数。联立上述两个方面，可得
数值求解等式（2-25）和等式（2-26），可得极大连通子图大小，即SIR疾病传播范围。若流行病爆发，等式（2-26）存在一个非平凡根 u <1。若没有流行病爆发，等式（2-26）仅有平凡根 u ＝ 1。当等式（2-26）的左右两端在 u ＝ 1 处相切时，可得临界传播率 ϕ c ，即
联立等式（2-24），疾病爆发阈值
，其值为
基于边渗流理论，很容易得到疾病传播范围和爆发阈值。经典的边渗流理论仅适合于无限大小的树形网络 [176，177] 。为了克服经典边渗流理论的不足，Noёl等发展了针对有限无关联树形网络的边渗流理论 [178] 。 Marder发展了边渗流理论来得到SIR传播模型的演化过程 [179] 。最近，学者们拓展了渗流理论来研究集群系数 [180] 、度关联 [181] 、社区结构 [182] 和关系多重性 [183] 对SIR疾病传播的影响。最近，Newman将边渗流理论拓展后用于分析共演化动力学 [106，110] ，发现两种疾病共存的临界条件。 Parshani等将边渗流理论用于求解SIS传播模型 [36] 。