一种基于SEIS的封闭系统传染病模型

A Closed System Epidemic Model Based on SEIS

摘要: 本文针对一个封闭系统中的某种传染病进行研究，建立基于SEIS的感染病模型，讨论了系统中染病人数的变化情况。分析了潜伏者对传染病传播的影响，对于初始潜伏者是工作人员的情况和初始潜伏者是其他人员的情况分别讨论，调整模型参数，得到系统中人数变化图。接着分析了模型平衡点的稳定性，结果表明当R₀<1时，λ<0，平衡点A⁰是局部渐近稳定的；当R₀＞1时，λ＞0 故A⁰是不稳定的，存在地方性平衡点。我们还讨论了主要参数对感染人数的变化情况，结果表明传染率为影响染病人数的主要影响因素，其次是治愈率，因此发现传染病采取隔离，戴口罩等措施是十分必要的。

Abstract: This paper studies an infectious disease in a closed system, establishes an infectious disease model based on SEIS, and discusses the change of the number of infected people in the system. This paper analyzes the influence of the lurk on the spread of infectious diseases, discusses the situation that the initial lurk is a worker and the initial lurk is other people respectively, adjusts the model parameters, and obtains the change chart of the number of people in the system. Then the stability of the equilibrium point is analyzed. The results show that when R₀<1, λ<0, equilibrium point A⁰ is locally asymptotically stable; when  R₀＞1,  λ＞0, so A⁰ is unstable, and there is a local equilibrium point. We also discussed the changes of the main parameters on the number of infected people. The results show that the infection rate is the main factor affecting the number of infected people, followed by the cure rate. Therefore, it is necessary to isolate infectious diseases and wear masks.

1. 引言

传染病(Infectious Diseases)是由各种病原体引起的能在人与人、动物与动物或人与动物之间相互传播的一类疾病。传染病是一种可以从一个人或其他物种，经过各种途径传染给另一个人或物种的感染病。在现实生活中，传染病广泛流传。传染病的特点是有病原体，有传染性和流行性，感染后常有免疫性。

本文分别讨论初始潜伏者是公共场所的工作人员和其它人员时，系统中得病人数的变化情况。本文以沈阳市某一商城为例，假设该城目前处于封城状态，即为一封闭系统，商城为其中的一个公共场所进行研究。当初始潜伏者是工作人员时，初始潜伏者只在固定的地方接触随机的人员；当初始潜伏者是其他人员时，初始潜伏者在随机的地方接触随机的人员。我们运用再生矩阵的方法计算了模型的控制再生数，讨论了平衡点的稳定性，当 $R_{\text{0}}<1$ 时，可得 $\lambda <0$，故平衡点 $A^{\text{0}}$ 是局部渐近稳定的；当 $R_{\text{0}}>1$ 时 $\lambda >0$， $A^{\text{0}}$ 是不稳定的，存在地方病平衡点。本文还重点分析了传染率，治愈率对传染病传播的影响，利用数值模拟了参数取值不同时染病人数的变化曲线，结果表明传染率对传染病传播的影响最大，因此隔离是传染病发生时的首要措施。接着我们讨论了疫苗对传染病传播的影响，数值模拟结果表明疫苗的出现大大降低了感染人数的峰值，延长了峰值出现的时间，因此研制疫苗等手段也必不可少。

2. 模型建立

2.1. 传染病的微分动力模型

利用动力学方法建立传染病模型是应用数学的一个重要分支，为流行病传播数学模型的研究做出了很大贡献。常见的传染病模型按照传染病类型分为SI、SIR、SIRS、SEIR、SEIS模型等，按照传播机理又分为基于常微分方程、偏微分方程、网络动力学的不同类型。对于带有潜伏期的传染病，假设潜伏者按指数分布发展成染病者。

通常感染人数为 $\frac{\beta SI}{S+I}$ [1]，考虑到工作人员和其他人员在特定场所的接触人数和活动范围不同，本文引入 $\frac{\beta S{I}_1}{1+S+I_1}+\frac{\beta S{I}_2}{1+S+I_2}$ 作为感染人数，得到其动力学性态的完整分析，如图1所示：

Figure 1. Schematic diagram of SEIS infectious disease model

图1. SEIS传染病模型示意图

2.2. 基于SEIS传染病的模型建立

初始潜伏者分别为工作人员和其他人员时，其接触人员数量和活动范围不同。当初始潜伏者是工作人员时，初始潜伏者只在固定的地方接触随机的人员；当初始潜伏者是其他人员时，初始潜伏者在随机的地方接触随机的人员。综合考虑封闭系统中的情况，我们将总人口分为四个不同仓室：易感者(S)，潜伏期患者(E)，工作人员感染者(I₁)，其他人员感染者(I₂)，做出如下假设：

1) 假设工作人员和其他人员被感染的概率分别与工作人员和其他人员在人群中占的比例相等。

2) 康复者会二次感染

根据上述假设和图1的传播示意图，我们建立如下模型：

$\frac{\text{d}S}{\text{d}t}=\gamma \left(I_1+I_2\right)-\frac{\beta S{I}_1}{1+S+I_1}-\frac{\beta S{I}_2}{1+S+I_2}-\alpha S$ (1)

$\frac{\text{d}E}{\text{d}t}=\frac{\beta S{I}_1}{1+S+I_1}+\frac{\beta S{I}_2}{1+S+I_2}-\epsilon E-\alpha E$ (2)

$\frac{\text{d}{I}_1}{\text{d}t}=\epsilon E_1-\left(\alpha +\gamma \right)I_1$ (3)

$\frac{\text{d}{I}_{\text{2}}}{\text{d}t}=\epsilon E_{\text{2}}-\left(\alpha +\gamma \right)I_{\text{2}}$ (4)

其中，S表示封闭系统健康总人数，即易感染人数；S是健康总人数中工作人员数量A和健康总人数中其他人员数量B，其表达式为：

$S=A+B$ (5)

E表示潜伏者，E₁表示工作人员潜伏者，E₂表示其他人员潜伏者，其表达式为：

$E_1=\frac{A}{S}E$ (6)

$E_2=\frac{B}{S}E$ (7)

$I_1$ 表示工作人员染病者数量； $I_2$ 表示其他人员染病者数量；N染病者接触的平均人数； $\gamma$ 表示染病者到健康人的恢复系数； $\beta$ 表示传染率， ${\beta }_{\text{1}}$ 表示工作人员传染率， ${\beta }_{\text{2}}$ 表示其他人员传染率； $\alpha$ 表示死亡率； $\epsilon$ 表示潜伏者转为染病者的转化率。

3. 实验结果

3.1. 初始潜伏者是公共场所的工作人员

当初始潜伏者是工作人员时，初始潜伏者只在固定的地方接触随机的人员，但其在该公共场所的时间较其他人员的时间长，所以接触的人更多，即 ${\beta }_1>{\beta }_2$。

设置 $E_1=1$，通过Python模拟上述模型，可视化所得数据，如图2所示：

Figure 2. The change chart of latent and infected persons with time when the initiators were people in public places

图2. 初始者为公共场所人员时潜伏者和染病者随时间变化图

其中纵坐标分别表示总人数、潜伏者、其他人员染病和工作人员染病的人数，横坐标表示时间。通过上图可以看出该病在25~50天左右达到感染高峰，在200天左右，疾病基本灭绝。

3.2. 初始潜伏者是公共场所的其他人员

当初始潜伏者是其他人员时，初始潜伏者在随机的地方接触随机的人员，但在该公共场所的时间较其他人员的时间短。

设置 $E_2=1$，通过Python模拟上述模型，可视化所得数据，如图3所示：

Figure 3. When the initial lurker is other people, the change chart of lurker and infected person with time

图3. 初始潜伏者为其他人员时潜伏者和染病者随时间变化图

其中纵坐标分别表示总人数、潜伏者、其他人员染病和工作人员染病的人数，横坐标表示时间。通过上图可以看出该病在50~75天左右达到感染高峰，在175~200天作用，疾病基本灭绝。

通过比较图2和图3，可以看出当初始潜伏者是工作人员时较初始潜伏者是其他人员时的图像斜率较陡。可以反映出当初始潜伏者是工作人员接触的人更多，疾病传播的人数更多。

4. 模型分析与检验

4.1. 平衡点的存在性

我们把上述模型中的工作人员感染人数和其他人员感染人数变为一个变量，如下所示：

$\{\begin{array}{l}\frac{\text{d}S}{\text{d}t}=\gamma I-\frac{\beta SI}{N}-\alpha S\\ \frac{\text{d}E}{\text{d}t}=\frac{\beta SI}{N}-\epsilon E-\alpha E\\ \frac{\text{d}I}{\text{d}t}=\epsilon E-\left(\alpha +\gamma \right)I\end{array}$ (8)

其中， $D=\left\{\left(S,E,I\right)\in R_{+}^{\text{3}}:N=S+E+I\le \frac{1}{\alpha }\right\}$ 是模型(8)的不变集。

设无病平衡点为 $A^{\text{0}}=\left(S^0,E^0,I^0\right)$，由于 $I=\text{0}$，则 $E=\text{0}$，若 $\gamma I^{\text{0}}-\frac{\beta S^0{I}^0}{N}-\alpha S^0=0$ 成立，解得 $S^{\text{0}}=\frac{\text{1}}{\alpha }$，因此模型(8)恒存在无病平衡点为 $\left(\frac{1}{\alpha },0,0\right)$。

设 $A^{\text{*}}=\left(S^{\text{*}},E^{\text{*}},I^{\text{*}}\right)$ 为模型(8)的任意正平衡点，当 $I\ne \text{0}$ 时，由 $\frac{\text{d}I}{\text{d}t}=0$ 得 $E^{\text{*}}=\frac{\left(\alpha +\gamma \right)I^{\text{*}}}{\epsilon }$

代入 $\frac{\text{d}E}{\text{d}t}=0$，得 $S^{\text{*}}=\frac{\left(\alpha +\epsilon \right)\left(\alpha +\gamma \right)}{\epsilon \beta }$，代入 $\frac{\text{d}S}{\text{d}t}=0$，得 $I^{\text{*}}=\frac{\alpha \left(\epsilon +\alpha \right)\left(\alpha +\gamma \right)}{\epsilon \beta \gamma -\beta \left(\epsilon +\alpha \right)\left(\alpha +\gamma \right)}$，因此正平衡点为

$A^{*}=\left(\frac{\left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)}{\epsilon \beta },\frac{\left(\alpha +\gamma \right)I^{*}}{\epsilon },\frac{\alpha \left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)}{\epsilon \beta \gamma -\beta \left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)}\right)$。

4.2. 平衡点的稳定性分析

模型(8)总存在一个无病平衡点 $A^0=\left(\frac{1}{\alpha },0,0\right)$，我们应用再生矩阵的方法 [2] 来计算模型(8)的控制再生数并验证 $A^{\text{0}}$ 的局部稳定性 [3]。

此模型有两个仓室分别为： $I\left(t\right),E\left(t\right)$，划分X,Y向量为 $X=\left[\begin{array}{c}E\left(t\right)\\ I\left(t\right)\end{array}\right],Y=\left[S\left(t\right)\right]$

$X=\left[\begin{array}{c}E\left(t\right)\\ I\left(t\right)\end{array}\right]=\left[\begin{array}{c}\frac{\beta SI}{N}\\ 0\end{array}\right]+\left[\begin{array}{c}\epsilon E+\alpha E\\ -\epsilon E+\alpha I+\gamma I\end{array}\right]=F_{1,2}\left(E,I\right)+V_{1,2}\left(E,I\right)$

因为 $S+I+E=N$，因此我们取 $N=\text{1}$

$F,V$ 的雅可比矩阵如下：

$F=Jacobian\left(F_{1,2}\left(E,I\right)\right)=\left[\begin{array}{cc}0& \beta S\\ 0& 0\end{array}\right]$

$V=Jacobian\left(V_{1,2}\left(E,I\right)\right)=\left[\begin{array}{cc}\epsilon +\alpha & \text{0}\\ -\epsilon & \alpha +\gamma \end{array}\right]$

则可得模型(8)的控制再生数：

$R_{\text{0}}={\rho \left(F{V}^{-1}\right)|}_{A_0}=\frac{\epsilon \beta }{\left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)\alpha }$

下面利用特征方程理论来证明平衡点 $A_{\text{0}}$ 的稳定性 [4]。

系统在平衡点处的雅可比矩阵如下：

$G_0=\left[\begin{array}{ccc}-\alpha & 0& \gamma -\frac{\beta }{\alpha }\\ 0& -\epsilon -\alpha & \frac{\beta }{\alpha }\\ 0& \epsilon & -\alpha -\gamma \end{array}\right]$ (9)

得到 ${\lambda }^{\text{2}}+\lambda \left(\text{2}\alpha +\gamma +\epsilon \right)+\left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)-\frac{\beta \epsilon }{\alpha }=\text{0}$

因此 $\frac{{\lambda }^{\text{2}}}{\left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)}+\lambda \frac{\left(\text{2}\alpha +\gamma +\epsilon \right)}{\left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)}=R_{\text{0}}-1$，其中 ${\lambda }^{\text{2}}>0,\alpha >0,\beta >0,\gamma \ge 0$。

因此当 $R_{\text{0}}<1$ 时，可得 $\lambda <0$，故平衡点 $A^{\text{0}}$ 是局部渐近稳定的；当 $R_{\text{0}}>1$ 时， $\lambda >0$ 故 $A^{\text{0}}$ 是不稳定的。

当 $R_{\text{0}}>1$ 时， $\frac{\epsilon \beta }{\left(\alpha +\gamma \right)\left(\alpha +\epsilon \right)\alpha }>1$，则 $S^{\text{*}},E^{*},I^{\text{*}}$ 均为正数，故模型(8)存在地方病平衡点 $A^{\text{*}}$ [5]。

5. 数值模拟各参数对传染病模型的影响

5.1. N对传染病模型的影响

我们分别取N = 100, 1000, 10000, 100000来模拟，如图4~7所示：

从下图可以看出N值基本只影响疫情规模的时间推移情况，对其他数据不存在影响。

Figure 4. When n = 100, the change chart of three kinds of people with time

图4. n = 100时三种人数随时间变化图

Figure 5. When n = 1000, the change chart of three kinds of people with time

图5. n = 1000时三种人数随时间变化图

Figure 6. When n = 10,000, the change chart of three kinds of people with time

图6. n = 10,000时三种人数随时间变化图

Figure 7. When n = 100,000, the change chart of three kinds of people with time

图7. n = 100,000时三种人数随时间变化图

5.2. 传染率对传染病模型的影响

我们分别取 $\beta =0.4,0.8,1.2,1.6$ 进行模拟，如图8~11所示。

Figure 8. When $\beta =0.4$, the change chart of three kinds of people with time

图8. $\beta =0.4$ 时三种人数随时间变化图

Figure 9. $\beta =0.8$, the change chart of three kinds of people with time

图9. $\beta =0.8$ 时三种人数随时间变化图

Figure 10. When $\beta =1.2$, the change chart of three kinds of people with time

图10. $\beta =1.2$ 时三种人数随时间变化图

Figure 11. When $\beta =1.6$, the change chart of three kinds of people with time

图11. $\beta =1.6$ 时三种人数随时间变化图

由图可知， $\beta$ 的大小直接影响到的是E的数值大小及变化趋势，进而影响到了S与I值的大小及变化趋势。且从图中可以直观了解到 $\beta$ 的值对整个模型的趋势是有很大影响的，因此戴口罩，减少外出时间居家隔离等措施都可以行之有效的减少 $\beta$ 的大小，对阻碍病毒传播有一定效果。

5.3. 治愈率对传染病传播的影响

当存在有效药物时，治愈率 $\gamma$ 将提高，我们分别设置 $\gamma =0.05,0.2,0.4,0.6$ 对比如图12~15所示。

Figure 12. When $\gamma =\text{0}\text{.05}$, the change chart of three kinds of people with time

图12. $\gamma =\text{0}\text{.05}$ 时三种人数随时间变化图

Figure 13. When $\gamma =\text{0}\text{.2}$, the change chart of three kinds of people with time

图13. $\gamma =\text{0}\text{.2}$ 时三种人数随时间变化图

Figure 14. When $\gamma =\text{0}\text{.4}$, the change chart of three kinds of people with time

图14. $\gamma =\text{0}\text{.4}$ 时三种人数随时间变化图

Figure 15. When $\gamma =\text{0}\text{.6}$, the change chart of three kinds of people with time

图15. $\gamma =\text{0}\text{.6}$ 时三种人数随时间变化图

由上图可以清晰的看出， $\gamma$ 值的大小影响了I和E值的极值大小。因此提高医疗救治水平可以有效减少同时间累积的感染人数，提高效率。

5.4. 疫苗对传染病传播的影响

当存在疫苗时，部分易感人群S会直接进入免疫人群，因此我们增加免疫人群这一参数R重新建立模型，假设每天新增的免疫者数量与易感人群数量成正比，比例系数为k，建立模型为：

$\{\begin{array}{l}\frac{\text{d}S}{\text{d}t}=-\frac{\beta SI}{N}-k\frac{S}{N}-\alpha S\\ \frac{\text{d}E}{\text{d}t}=\frac{\beta SI}{N}-\epsilon E-\alpha E\\ \frac{\text{d}I}{\text{d}t}=\epsilon E-\left(\alpha +\gamma \right)I\\ \frac{\text{d}R}{\text{d}t}=\gamma I+k\frac{S}{N}-\alpha R\end{array}$ (10)

打疫苗前后对比如图16、图17所示：

Figure 16. Time chart of the number of the first three vaccinations

图16. 打疫苗前三种人数随时间变化图

Figure 17. Time chart of the three kinds of people after vaccination

图17. 打疫苗后三种人数随时间变化图

可以看到，疫苗使开始感染的人数有所下降，感染到达峰值也有所延长。

6. 结论

我们构建了一类在封闭系统中员工和其他人员分开的SEIS模型，探究了潜伏者对传染病传播的影响。结果表明潜伏者为工作人员时染病人数变化的曲线斜率更大，因此在生活中像商场等封闭环境中工作人员实施戴口罩等防护措施是十分必要的。对模型的定性分析表明控制再生数 $R_{\text{0}}$ 是一个临界阈值，决定了传染病传播规模和发展趋势。接着对平衡点进行稳定性分析，结果表明当 $R_{\text{0}}<1$ 时，可得 $\lambda <0$，故平衡点 $A^{\text{0}}$ 是局部渐近稳定的；当 $R_{\text{0}}>1$ 时， $\lambda >0$ 故 $A^{\text{0}}$ 是不稳定的，存在地方性平衡点。接着我们通过对调整主要参数来数值模拟传染病人数变化情况，可以直观的看出传染病的传染率和治愈率为影响传染病传播的主要因素，因此在传染病传播期间隔离应该为首要措施，民众应该加强防护意思，采取勤洗手，出门佩戴口罩等来阻断传播，除此之外还应该提高治愈率，例如加强疫苗研制等。

NOTES

^*通讯作者。