图数据挖掘：基于概率的流行病模型

1 导引

在上一篇博客《图数据挖掘：网络中的级联行为》中介绍了用基于决策的模型来对级联行为进行建模，该模型是基于效用(Utility)的且是是确定性的，主要关注于单个节点如何根据其邻居的情况来做决策，需要大量和数据相关的先验信息。这篇博客就让我们来介绍基于概率的传播模型，这种模型基于对数据的观测来构建，不过不能对因果性进行建模。

2 基于随机树的流行病模型

接下来我们介绍一种基于随机树的传染病模型，它是分支过程（branching processes）的一种变种。在这种模型中，一个病人可能接触d" role="presentation">d$d$个其他人，对他们中的每一个都有概率q>0" role="presentation">q>0$q>0$将其传染，如下图所示：

接下来我们来看当d" role="presentation">d$d$和q" role="presentation">q$q$取何值时，流行病最终会消失（die out），也即满足

limh→∞ph=0" role="presentation">limh→∞ph=0$$\underset{h\to \mathrm{\infty }}{lim}{p}_h=0$$

这里ph" role="presentation">ph$p_h$为在深度h" role="presentation">h$h$处存在感染节点的概率（是关于q" role="presentation">q$q$和d" role="presentation">d$d$的函数）。如果流行病会永远流行下去，则上述极限应该>0" role="presentation">>0$>0$。

ph" role="presentation">ph$p_h$满足递归式：

ph=1−(1−q⋅ph−1)d" role="presentation">ph=1−(1−q⋅ph−1)d$$p_h=1-{(1-q\cdot p_{h-1})}^d$$

这里(1−q⋅ph−1)d" role="presentation">(1−q⋅ph−1)d${(1-q\cdot p_{h-1})}^d$表示在距离根节点h" role="presentation">h$h$深度处没有感染节点的概率。

接下来我们通过对函数

f(x)=1−(1−q⋅x)d" role="presentation">f(x)=1−(1−q⋅x)d$$f(x)=1-(1-q\cdot x{)}^d$$

进行迭代来得到limh→∞ph" role="presentation">limh→∞ph$\underset{h\to \mathrm{\infty }}{lim}{p}_h$。我们从根节点x=1" role="presentation">x=1$x=1$(因为p1=1" role="presentation">p1=1$p_1=1$)开始，依次迭代得到x1=f(1),x2=f(x1),x3=f(x2)" role="presentation">x1=f(1),x2=f(x1),x3=f(x2)$x_1=f(1),x_2=f(x_1),x_3=f(x_2)$。事实上，该迭代最终会收敛到不动点f(x)=x" role="presentation">f(x)=x$f(x)=x$，如下图所示：

这里x" role="presentation">x$x$是在深度h−1" role="presentation">h−1$h-1$处存在感染节点的概率，f(x)" role="presentation">f(x)$f(x)$是在深度为h" role="presentation">h$h$处存在感染节点的概率，q" role="presentation">q$q$为感染概率，d" role="presentation">d$d$为节点的度。
如果我们想要传染病最终消失，那么迭代f(x)" role="presentation">f(x)$f(x)$的结果必须要趋向于0" role="presentation">0$0$，也即不动点需要为0。而这也就意味着f(x)" role="presentation">f(x)$f(x)$必须要在y=x" role="presentation">y=x$y=x$下方，如下所示：

如何控制f(x)" role="presentation">f(x)$f(x)$必须要在y=x" role="presentation">y=x$y=x$下方呢？我们先来分析下f(x)" role="presentation">f(x)$f(x)$的图像形状，我们有以下结论：

f(x)" role="presentation">f(x)$f(x)$是单调的：对0≤x,q≤1,d>1" role="presentation">0≤x,q≤1,d>1$0\le x,q\le 1,d>1$，f′(x)=q⋅d(1−qx)d−1>0" role="presentation">f′(x)=q⋅d(1−qx)d−1>0$f^{\prime }(x)=q\cdot d(1-qx{)}^{d-1}>0$，故f(x)" role="presentation">f(x)$f(x)$是单调的。

f′(x)" role="presentation">f′(x)$f^{\prime }(x)$是非增的：f′(x)=q⋅d(1−qx)d−1" role="presentation">f′(x)=q⋅d(1−qx)d−1$f^{\prime }(x)=q\cdot d(1-qx{)}^{d-1}$会着x" role="presentation">x$x$减小而减小。

而f(x)" role="presentation">f(x)$f(x)$低于y=x" role="presentation">y=x$y=x$，则需要满足

f′(0)=q⋅d<1" role="presentation">f′(0)=q⋅d<1$$f^{\prime }(0)=q\cdot d<1$$

综上所述，我们有结论：

limh→∞ph=0 when q⋅d<1" role="presentation">limh→∞ph=0 when q⋅d<1$$\underset{h\to \mathrm{\infty }}{lim}{p}_h=0\text{ when }q\cdot d<\mathbf{1}$$

这里R0=q⋅d" role="presentation">R0=q⋅d$R_0=q\cdot d$表示每个被感染的个体在期望意义上所产生的新的病体数，我们将其称为基本再生数（reproductive number），它决定了传染病病是否会流行：

• 若R0≥1" role="presentation">R0≥1$R_0\ge 1$: 流行病永远不会消失且感染人数会以指数速度上升。
• 若R0≤1" role="presentation">R0≤1$R_0\le 1$: 流行病会以指数速度快速消失。

3 SIR与SIS流行病模型

3.1 模型范式

在病毒的传播中，有两个最基本的参数：

• 出生率β" role="presentation">β$\beta$ 被已感染邻居攻击的概率
• 死亡率δ" role="presentation">δ$\delta$ 已感染节点治愈的概率

网络中的节点可以在以下四个状态（S+E+I+R）之间做转移：

• 易感期(susceptible)： 节点患病之前，处于容易被邻居传染的时期，也称敏感期。
• 潜伏期(exposed)：节点已被感染，但是还没具备能力去传染别人。
• 传染期(infectious)：节点已被感染，且能够以一定的概率把疾病传染给那些处于易感期的邻居，也称感染期。
• 移除期(removed)：当一个节点经历了完整的传染期，就不再被考虑了，因为它不会再受感染，也不会对其它节点构成威胁，也称隔离期。

状态转移图如下图所示（图中的Z" role="presentation">Z$Z$表示人工免疫）：

其中状态转移的概率由我们前面提到的模型参数β" role="presentation">β$\beta$和δ" role="presentation">δ$\delta$控制。

3.2 SIR模型

在SIR模型中，节点经历S-I-R三个阶段：

事实上，该模型可用于对水痘和鼠疫的建模，也即一旦我治愈了，那我就永远不会再被感染了。

假设模型满足完美混合（即网络是完全图），则模型的动力方程为：

dSdt=−βSIdRdt=δIdIdt=βSI−δI" role="presentation">dSdt=−βSIdRdt=δIdIdt=βSI−δI$$\begin{array}{rl}& \frac{dS}{dt}=-\beta SI\\ & \frac{dR}{dt}=\delta I\\ & \frac{dI}{dt}=\beta SI-\delta I\end{array}$$

处于S" role="presentation">S$S$、I" role="presentation">I$I$、R" role="presentation">R$R$状态的节点数量随着时间变化曲线如下图所示：

3.3 SIS模型

SIS模型中节点只有S-I两个阶段，它假设已经治愈的节点会立即变为易感节点。节点的状态转移图如下：

这里我们把s=βδ" role="presentation">s=βδ$s=\frac{\beta }{\delta }$定义为病毒的“力量”（strength）。

该模型可用于对流感的建模，也即已被感染的节点经过治愈后会重新回到易感状态。

同样我们假设模型满足完美混合（即网络是完全图），则模型的动力方程为：

dSdt=−βSI+δIdIdt=βSI−δI" role="presentation">dSdt=−βSI+δIdIdt=βSI−δI$$\begin{array}{rl}& \frac{dS}{dt}=-\beta SI+\delta I\\ & \frac{dI}{dt}=\beta SI-\delta I\end{array}$$

处于S" role="presentation">S$S$、I" role="presentation">I$I$状态的节点数量随着时间变化曲线如下图所示：

3.4 传染阈值

接下来我们考虑SIS模型中的传染阈值(epidemic threshold)τ" role="presentation">τ$\tau$。对于SIS" role="presentation">SIS$SIS$模型而言，流行阈值可以是任意的。

我们设图G" role="presentation">G$G$的传染阈值为τ" role="presentation">τ$\tau$。如果病毒的“力量”s=βδ<τ" role="presentation">s=βδ<τ$s=\frac{\beta }{\delta }<\tau$（这里β" role="presentation">β$\beta$指病毒的死亡率，δ" role="presentation">δ$\delta$指病毒的出生率），则疾病的流行就不会发生（它最终会消失）。事实上，图G" role="presentation">G$G$的传染阈值τ" role="presentation">τ$\tau$可以表示为

τ=1λ1,A" role="presentation">τ=1λ1,A$$\tau =\frac{1}{{\lambda }_{1,A}}$$

这里λ1,A" role="presentation">λ1,A${\lambda }_{1,A}$为图G" role="presentation">G$G$的邻接矩阵最大的特征值。这个定理看起来非常神奇，因为我们只用λ1,A" role="presentation">λ1,A${\lambda }_{1,A}$就捕捉到了整个图的属性！

以下是在AS图上，当s" role="presentation">s$s$大于、小于或等于传染阈值τ" role="presentation">τ$\tau$时的感染节点数量随时间变化图：

如果我们再考虑不同的初始感染人数，则会得到以下的感染人数变化图像：

3.5 一个埃博拉的例子

在一个埃博拉的例子^[1]中，设置如下的转换状态：

当设置R0=1.5-2.0" role="presentation">R0=1.5-2.0$R_0=1.5\text{-}2.0$时，总死亡人数随时间变化如下：

参考

[1] Gomes M F C, y Piontti A P, Rossi L, et al. Assessing the international spreading risk associated with the 2014 West African Ebola outbreak[J]. PLoS currents, 2014, 6.
[2] http://web.stanford.edu/class/cs224w/
[3] Easley D, Kleinberg J. Networks, crowds, and markets: Reasoning about a highly connected world[M]. Cambridge university press, 2010.
[4] Barabási A L. Network science[J]. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2013, 371(1987): 20120375.

• 1 导引
• 2 基于随机树的流行病模型
• 3 SIR与SIS流行病模型
•     3.1 模型范式
•     3.2 SIR模型
•     3.3 SIS模型
•     3.4 传染阈值
•     3.5 一个埃博拉的例子
• 参考

__EOF__