一类具有隔离项的传染病模型稳定性分析

Stability Analysis of a Kind of Epidemic Disease Model with Isolation Term

摘要: 研究了一类具有带常规预防和隔离防护措施的传染病模型，运用下一代生成矩阵法求得系统的基本再生数并证明了平衡点的存在性，通过Hurwitz判别法得到无病平衡点与地方病平衡点的局部渐进稳定性，构造Lyapunov函数并利用Lasalle不变集原理得到平衡点全局渐近稳定的充分条件。

Abstract: A kind of infectious disease model with routine prevention and isolation protection measures is studied. The basic reproductive number of the system is obtained by the next generation matrix method and the existence of the equilibrium point is proved. The local asymptotic stability of the disease-free equilibrium point and the endemic equilibrium point was obtained by the Hurwitz dis-criminant method, the Lyapunov function is constructed and the sufficient condition for the global asymptotic stability of the equilibrium point is obtained by the Lasalle invariant set principle.

1. 引言

自古以来，传染病都威胁着人类的生存与生活，研究传染病模型是国内外研究者的重要课题。经过国内外学者多年对传染病模型的研究，已取得了不少成果 [1] - [9] 。在早期的数学模型中，多只考虑影响传染病传播的部分直接因素，如感染期，传播方式，易感性等。随着传染病模型的研究的更加深入，人们在建立模型时，考虑了更多能影响传染病传播的因素，如人口条件，地理环境，经济等因素。因此，有必要加入更多的参数和变量，使建立的数学模型更加真实。2002年底到2003年初，SARS迅速在国际上蔓延开来，超过28个国家遭受SARS疫情，在世卫组织的努力下，确定了SARS-CoV的全基因组序列，许多专家学者也提出了预测SARS病毒的数学模型。2019年底爆发的COVID-19，传染性极强且传播速度快，为成功解决突如其来的疫情，必须实施预防措施，在这些预防措施中，建立数学模型并分析对怎样实施预防措施有着重要的作用 [10] [11] [12] 。在实际生活中，许多传染病有潜伏期，控制传染病传播最有效的措施就是隔离。Anwarud Din等人 [13] 为控制COVID-19这类疫情的传播，使用两种控制措施建立模型并进行分析，得到COVID-19在种群中灭绝或持续存在的条件。胡瑞等人 [14] 结合COVID-19的传播性质建立了数学模型，并得出模型的平衡点的稳定性和不同参数对疾病传播的敏感性分析。梁桂珍等人 [15] [16] 考虑了具有连续接种措施的传染病模型，得到疾病流行与否的阈值并证明平衡点的全局稳定性。Swati Tyagi等人 [17] 考虑了COVID-19的两种感染形式，结合隔离预防措施建立数学模型，得出系统的基本再生数并对系统平衡点做稳定性分析。在以上文献的研究下，本文考虑加入常规预防措施以及隔离措施对传染病传播进行控制，在第2节提出具有两种隔离措施并考虑种群有两种感染形式的数学模型。在第3节，运用下一代生成矩阵法确定疾病是否流行的阈值以及平衡点的存在性。第4节通过Hurwitz判别法以及Lasalle不变集原理得到系统的局部稳定和全局稳定性。第5节通过数值模拟验证结果的正确性。

2. 模型建立

假设 $N\left(t\right)$ 为某时刻的人口总数，将 $N\left(t\right)$ 分为六部分为 $S\left(t\right)$ ， $E\left(t\right)$ ， $I\left(t\right)$ ， $A\left(t\right)$ ， $Q\left(t\right)$ ， $R\left(t\right)$ ，分别表示t时刻易感染者，潜伏者，早期无症状感染者，后期有症状感染者，隔离者和恢复者的数量。t时刻人口总数为

$N\left(t\right)=S\left(t\right)+E\left(t\right)+I\left(t\right)+A\left(t\right)+Q\left(t\right)+R\left(t\right).$

基于假设可建立如下模型：

$\{\begin{array}{l}\frac{\text{d}S\left(t\right)}{\text{d}t}=\Lambda -\left(1-\mu \right){\beta }_{1}SI-\left(1-\mu \right){\beta }_{2}SA-dS\\ \frac{\text{d}E\left(t\right)}{\text{d}t}=\left(1-\mu \right){\beta }_{1}SI+\left(1-\mu \right){\beta }_{2}SA-{\delta }_{1}E-{\delta }_{2}E-dE\\ \frac{\text{d}I\left(t\right)}{\text{d}t}={\delta }_{1}E+\upsilon A-{\alpha }_{1}I-{\gamma }_{1}I-{\epsilon }_{1}I-dI\\ \frac{\text{d}A\left(t\right)}{\text{d}t}={\delta }_{2}E-\upsilon A-{\gamma }_{2}A-{\epsilon }_{2}A-dA\\ \frac{\text{d}Q\left(t\right)}{\text{d}t}={\epsilon }_{1}I+{\epsilon }_{2}A-{\alpha }_{2}Q-{\gamma }_{3}Q-dQ\\ \frac{\text{d}R\left(t\right)}{\text{d}t}={\gamma }_{1}I+{\gamma }_{2}A+{\gamma }_{3}Q-dR\end{array}$ (1)

在系统(1)中， $\Lambda$ 表示人群的输入率， ${\beta }_1$ 表示有症状感染者与易感者的传染率， ${\beta }_2$ 表示无症状感染者与易感者的传染率， ${\delta }_1$ 表示潜伏者转为有症状感染者的比率， ${\delta }_2$ 表示潜伏者转为无症状感染者的比率， $\upsilon$ 表示无症状感染者进入有症状阶段的比率， ${\alpha }_1$ ， ${\alpha }_2$ 分别表示有症状感染者与隔离者的因病死亡率， ${\epsilon }_1$ ， ${\epsilon }_2$ 分别表示有症状感染者与无症状感染者的隔离率， ${\gamma }_1$ ， ${\gamma }_2$ ， ${\gamma }_3$ 分别表示后期有症状感染，无症状感染和隔离者的恢复率， $\mu$ 表示日常常规预防传染病的有效率，d表示人群的自然死亡率，考虑模型的生物意义，所有参数均为非负数。

由系统(1)可知

$\frac{\text{d}N}{\text{d}t}=\Lambda -{\alpha }_{1}I-{\alpha }_{2}Q-dN.$

利用比较原理，得

$\frac{\text{d}N}{\text{d}t}\le \Lambda -dN.$

从而，有

$0<N\left(t\right)\le \frac{\Lambda }{d}-\left[\frac{\Lambda }{d}-N\left(0\right)\right]{\text{e}}^{-dt}.$

因此

$\underset{t\to +\infty }{\lim }\sup N\left(t\right)\le \frac{\Lambda }{d}.$

显然可行域为 $\Omega =\left\{\left(S,E,I,A,Q,R\right)\in R_{+}^6:S+E+I+A+Q+R\le \frac{\Lambda }{d}\right\}$ ，集合 $\Omega$ 是系统(1)的正向不变集。

3. 基本再生数与平衡点的存在性

通过下一代生成矩阵法得到系统(1)的基本再生数 $R_0$ ，构造矩阵 $\mathcal{F}$ 和 $\mathcal{V}$ ，

$\mathcal{F}=\left(\begin{array}{c}\left(1-\mu \right){\beta }_{1}SI+\left(1-\mu \right){\beta }_{2}SA\\ 0\\ 0\end{array}\right),\mathcal{V}=\left(\begin{array}{c}\left({\delta }_1+{\delta }_2+d\right)E\\ -{\delta }_{1}E-vA+\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)I\\ -{\delta }_{2}E+\left(v+{\epsilon }_2+{\gamma }_2+d\right)A\end{array}\right).$

$\mathcal{F}$ 和 $\mathcal{V}$ 在无病平衡点 $P_0$ 处的Jacobi矩阵为

$V=\left[\begin{array}{ccc}{\delta }_1+{\delta }_2+d& 0& 0\\ -{\delta }_1& {\alpha }_1+{\epsilon }_1+{\gamma }_1+d& -v\\ -{\delta }_2& 0& v+{\epsilon }_2+{\gamma }_2+d\end{array}\right],$

$V^{-1}=\left[\begin{array}{ccc}\frac{1}{{\delta }_1+{\delta }_2+d}& 0& 0\\ \frac{{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v}{\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)}& \frac{1}{{\alpha }_1+{\epsilon }_1+{\gamma }_1+d}& \frac{v}{\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)}\\ \frac{{\delta }_2}{\left({\delta }_1+{\delta }_2+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)}& 0& \frac{1}{v+{\epsilon }_2+{\gamma }_2+d}\end{array}\right].$

因此，

$F{V}^{-1}=\left[\begin{array}{ccc}{m}_1& m_2& m_3\\ 0& 0& 0\\ 0& 0& 0\end{array}\right].$

其中

$m_1=\frac{\left[{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v\right]\left(1-\mu \right){\beta }_1{S}^0+{\delta }_2\left(1-\mu \right){\beta }_2{S}^0\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)}{\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)},$

$m_2=\frac{\left(1-\mu \right){\beta }_1{S}^0}{{\alpha }_1+{\epsilon }_1+{\gamma }_1+d},$

$m_3=\frac{v\left(1-\mu \right){\beta }_1{S}^0+\left(1-\mu \right){\beta }_2{S}^0\left(v+{\epsilon }_2+{\gamma }_2+d\right)}{\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)}.$

因此系统(1)的基本再生数为

$R_0=\rho \left(F{V}^{-1}\right)=\frac{\left[{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v\right]\left(1-\mu \right){\beta }_1{S}^0+{\delta }_2\left(1-\mu \right){\beta }_2{S}^0\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)}{\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)}$ .

定理1当 $R_0<1$ 时，系统(1)总存在一个无病平衡点 $P_0\left(S_0,0,0,0,0,0\right)$ ，当 $R_0>1$ 时，系统(1)存在一个地方病平衡点 $P^{*}\left(S^{*},E^{*},I^{*},A^{*},Q^{*},R^{*}\right)$ 。

证明 令系统(1)中 $E\left(t\right)=0,I\left(t\right)=0,A\left(t\right)=0,Q\left(t\right)=0,R\left(t\right)=0$ ，则系统(1)有无病平衡点为 $P_0\left(S_0,0,0,0,0,0\right)$ ，其中 $S_0=\frac{\Lambda }{d}$ 。

$\{\begin{array}{l}\Lambda -\left(1-\mu \right){\beta }_1{S}^{*}{I}^{*}-\left(1-\mu \right){\beta }_2{S}^{*}{A}^{*}-d{S}^{*}=0\\ \left(1-\mu \right){\beta }_1{S}^{*}{I}^{*}+\left(1-\mu \right){\beta }_2{S}^{*}{A}^{*}-\left({\delta }_1+{\delta }_2+d\right)E^{*}=0\\ {\delta }_1{E}^{*}+\upsilon A^{*}-\left({\alpha }_1+{\gamma }_1+{\epsilon }_1+d\right)I^{*}=0\\ {\delta }_2{E}^{*}-\left(\upsilon +{\gamma }_2+{\epsilon }_2+d\right)A^{*}=0\\ {\epsilon }_1{I}^{*}+{\epsilon }_2{A}^{*}-\left({\alpha }_2+{\gamma }_3+d\right)Q^{*}=0\\ {\gamma }_1{I}^{*}+{\gamma }_2{A}^{*}+{\gamma }_3{Q}^{*}-d{R}^{*}=0\end{array}$ (2)

$S^{*}=\frac{1}{d}\left[\Lambda -\frac{\left({\delta }_1+{\delta }_2+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)A^{*}}{{\delta }_2}\right],$ (3)

$E^{*}=\frac{1}{{\delta }_2}\left(v+{\epsilon }_2+{\gamma }_2+d\right)A^{*},$ (4)

$I^{*}=\frac{{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v}{{\delta }_2\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)}{A}^{*},$ (5)

$Q^{*}=\left[\frac{{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v}{{\delta }_2\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left({\alpha }_2+{\gamma }_3+d\right)}+\frac{{\epsilon }_2}{{\alpha }_2+{\gamma }_3+d}\right]A^{*},$ (6)

$R^{*}=\frac{1}{d}\left[\frac{{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v}{{\delta }_2\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)}\left({\gamma }_1+\frac{{\gamma }_3}{{\alpha }_2+{\gamma }_3+d}\right)+{\gamma }_2+{\epsilon }_2\right]A^{*},$ (7)

$A^{*}=\frac{\Lambda {\delta }_2}{\left({\delta }_1+{\delta }_2+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)}-\frac{d{\delta }_2\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)}{\left(1-\mu \right)\left[{\beta }_1{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\beta }_1{\delta }_{2}v-{\beta }_2{\delta }_2\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\right]}.$ (8)

当 $R_0>1$ 时，有 $A^{*}>0$ ，且 $S^{*},E^{*},I^{*},A^{*},Q^{*},R^{*}$ 均大于零，所以系统(1)存在地方病平衡点 $P^{*}\left(S^{*},E^{*},I^{*},A^{*},Q^{*},R^{*}\right)$ 。

4. 平衡点的稳定性分析

定理2 若 $R_0<1$ ，则系统(1)的无病平衡点 $P_0$ 局部渐近稳定；若 $R_0>1$ ，则系统(1)的无病平衡点 $P_0$ 不稳定。

证明 考虑子系统

$\{\begin{array}{l}{S}^{\prime }\left(t\right)=\Lambda -\left(1-\mu \right){\beta }_{1}SI-\left(1-\mu \right){\beta }_{2}SA-dS\hfill \\ E^{\prime }\left(t\right)=\left(1-\mu \right){\beta }_{1}SI+\left(1-\mu \right){\beta }_{2}SA-{\delta }_{1}E-{\delta }_{2}E-dE\hfill \\ I^{\prime }\left(t\right)={\delta }_{1}E+\upsilon A-{\alpha }_{1}I-{\gamma }_{1}I-{\epsilon }_{1}I-dI\hfill \\ A^{\prime }\left(t\right)={\delta }_{2}E-\upsilon A-{\gamma }_{2}A-{\epsilon }_{2}A-dA\hfill \end{array}$ (9)

由系统(9)第二个方程可得

$\left(1-\mu \right){\beta }_{1}SI+\left(1-\mu \right){\beta }_{2}SA-{\delta }_{1}E-{\delta }_{2}E-dE\le \left(1-\mu \right){\beta }_1{S}^{0}I+\left(1-\mu \right){\beta }_2{S}^{0}A-{\delta }_{1}E-{\delta }_{2}E-dE.$

有线性系统

$\{\begin{array}{l}{X}^{\prime }\left(t\right)=\left(1-\mu \right){\beta }_1{S}^{0}Y+\left(1-\mu \right){\beta }_2{S}^{0}Z-{\delta }_{1}X-{\delta }_{2}X-dX\hfill \\ Y^{\prime }\left(t\right)={\delta }_{1}X+\upsilon Z-{\alpha }_{1}Y-{\gamma }_{1}Y-{\epsilon }_{1}Y-dY\hfill \\ Z^{\prime }\left(t\right)={\delta }_{2}X-\upsilon Z-{\gamma }_{2}Z-{\epsilon }_{2}Z-dZ\hfill \end{array}$ (10)

其系数矩阵为

$J=\left(\begin{array}{ccc}-\left({\delta }_1+{\delta }_2+d\right)& \left(1-\mu \right){\beta }_1{S}^0& \left(1-\mu \right){\beta }_2{S}^0\\ {\delta }_1& -\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)& v\\ {\delta }_2& 0& -\left(v+{\epsilon }_2+{\gamma }_2+d\right)\end{array}\right),$

特征方程为

$x^3+a_1{x}^2+a_{2}x+a_3=0,$ (11)

其中

$a_1={\delta }_1+{\delta }_2+d+{\alpha }_1+{\epsilon }_1+{\gamma }_1+d+v+{\epsilon }_2+{\gamma }_2+d,$ (12)

$\begin{array}{l}{a}_2=\left({\delta }_1+{\delta }_2+d+{\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)\\ +\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left({\delta }_1+{\delta }_2+d\right)-{\delta }_1\left(1-\mu \right){\beta }_1{S}^0-{\delta }_2\left(1-\mu \right){\beta }_2{S}^0,\end{array}$ (13)

$\begin{array}{c}{a}_3=\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)\\ -\left[{\delta }_1\left(v+{\epsilon }_2+{\gamma }_2+d\right)+{\delta }_{2}v\right]\left(1-\mu \right){\beta }_1{S}^0-{\delta }_2\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(1-\mu \right){\beta }_2{S}^0\\ =\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)\left(1-R_0\right).\end{array}$ (14)

显然，当 $R_0<1$ 时， $a_1,a_2,a_3>0$ ，此外

${\Delta }_1=\left|a_1\right|>0,{\Delta }_2=\left|\begin{array}{cc}{a}_1& 1\\ a_3& a_2\end{array}\right|=a_1{a}_2-a_3>0.$

其中

$\begin{array}{c}{a}_1{a}_2-a_3=a_1\left({\delta }_1+{\delta }_2+d+{\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)\\ +\left({\delta }_1+{\delta }_2+d+{\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\\ -{\delta }_1\left(1-\mu \right){\beta }_1{S}^0\left({\delta }_1+{\delta }_2+d+{\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\\ +{\delta }_{2}v\left(1-\mu \right){\beta }_1{S}^0-{\delta }_2\left(1-\mu \right){\beta }_2{S}^0\left({\delta }_1+{\delta }_2+d+v+{\epsilon }_2+{\gamma }_2+d\right)\end{array}$

$\begin{array}{c}>\left({\delta }_1+{\delta }_2+d+{\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)\left[\left({\delta }_1+{\delta }_2+d\right)\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)-{\delta }_1\left(1-\mu \right){\beta }_1{S}^0\right]\\ +\left({\delta }_1+{\delta }_2+d+v+{\epsilon }_2+{\gamma }_2+d\right)\left[\left({\delta }_1+{\delta }_2+d\right)\left(v+{\epsilon }_2+{\gamma }_2+d\right)-{\delta }_2\left(1-\mu \right){\beta }_2{S}^0\right]\\ +{\delta }_{2}v\left(1-\mu \right){\beta }_1{S}^0.\end{array}$

利用Hurwitz判别法可知，特征方程 $x^3+a_1{x}^2+a_{2}x+a_3=0$ 的根均有负实部，则当 $R_0<1$ 时，无病平衡点 $P_0$ 局部渐进稳定。

定理3 [18] 若给定的模型(1)写成：

$\frac{\text{d}X}{\text{d}t}=F\left(X,Y\right),\frac{\text{d}Y}{\text{d}t}=G\left(X,Y\right),G\left(X,0\right)=0,$

其中 $X=S,Y={\left(E,I,A\right)}^{\text{T}}$ 表示未感染者和感染者。无病平衡点表示为

$P_0=\left(X_0,0\right)=\left(\frac{\Lambda }{d},0\right),$

要使 $P_0$ 全局渐近稳定，则满足以下条件：

(I) $\frac{\text{d}X}{\text{d}t}=F\left(X_0,0\right)$ ， $X_0$ 全局渐进稳定。

(II) $G\left(X,Y\right)=MY-\hat{G}\left(X,Y\right),\hat{G}\left(X,Y\right)\ge 0$ 。这里 $M=D_{Y}G\left(X_0,0\right)$ 是M矩阵(在M矩阵中，所有非对

角线元素都是非负的)。

如果模型(1)中给定微分方程满足定理3的条件，则在 $R_0<1$ 的条件下，无病平衡点 $P_0=\left(X_0,0\right)$ 在系统(1)中全局渐近稳定。

定理4 若 $R_0<1$ ，则系统(1)的无病平衡点 $P_0$ 全局渐近稳定。

证明 将定理3运用到系统(1)，得到

$F\left(X_0,0\right)=\Lambda -\mu S^0,G\left(X,Y\right)=AY-\hat{G}\left(X,Y\right),$

$A=\left(\begin{array}{ccc}-\left({\delta }_1+{\delta }_2+d\right)& \left(1-\mu \right){\beta }_1{S}^0& \left(1-\mu \right){\beta }_2{S}^0\\ {\delta }_1& -\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)& v\\ {\delta }_2& 0& -\left(v+{\epsilon }_2+{\gamma }_2+d\right)\end{array}\right),$

$\hat{G}\left(X,Y\right)=AY-G\left(X,Y\right)=\left[\begin{array}{c}\left(1-\mu \right)\left(S^0-S\right)\left({\beta }_{1}I+{\beta }_{2}A\right)\\ 0\\ 0\end{array}\right].$

由于 $S^0-S\ge 0$ ，故 $\hat{G}\left(X,Y\right)\ge 0$ ，且A矩阵是M矩阵，根据文献 [18] ，故 $P_0$ 全局渐近稳定。

定理5 若 $R_0>1$ ，地方病平衡点 $P^{*}\left(S^{*},E^{*},I^{*},A^{*},Q^{*},R^{*}\right)$ 在 $\Omega \backslash \left\{P_0\right\}$ 上全局渐进稳定。

证明 构造Lyapunov函数

$V=K_1\left(S-S^{*}-S^{*}\ln \frac{S}{S^{*}}\right)+K_2\left(E-E^{*}-E^{*}\ln \frac{E}{E^{*}}\right)+K_3\left(I-I^{*}-I^{*}\ln \frac{I}{I^{*}}\right)+K_4\left(A-A^{*}-A^{*}\ln \frac{A}{A^{*}}\right).$

沿着系统(1)的轨迹对函数V求导

$\begin{array}{c}\frac{\text{d}V}{\text{d}t}=K_1\left(1-\frac{S^{*}}{S}\right)\left[\Lambda -\left(1-\mu \right){\beta }_{1}SI-\left(1-\mu \right){\beta }_{2}SA-dS\right]\\ +K_2\left(1-\frac{E^{*}}{E}\right)\left[\left(1-\mu \right){\beta }_{1}SI+\left(1-\mu \right){\beta }_{2}SA-\left({\delta }_1+{\delta }_2+d\right)E\right]\\ +K_3\left(1-\frac{I^{*}}{I}\right)\left[{\delta }_{1}E+vA-\left({\alpha }_1+{\epsilon }_1+{\gamma }_1+d\right)I\right]\\ +K_4\left(1-\frac{A^{*}}{A}\right)\left[{\delta }_{2}E-\left(v+{\epsilon }_2+{\gamma }_2+d\right)A\right].\end{array}$

在平衡点 $P^{*}\left(S^{*},E^{*},I^{*},A^{*},Q^{*},R^{*}\right)$ 处，系统(1)满足以下条件

$\Lambda =\left(1-\mu \right){\beta }_1{S}^{*}{I}^{*}-\left(1-\mu \right){\beta }_2{S}^{*}{A}^{*}-d{S}^{*},$

${\delta }_1+{\delta }_2+d=\frac{\left(1-\mu \right){\beta }_1{S}^{*}{I}^{*}+\left(1-\mu \right){\beta }_2{S}^{*}{A}^{*}}{E^{*}},$

${\alpha }_1+{\epsilon }_1+{\gamma }_1+d=\frac{{\delta }_1{E}^{*}+v{A}^{*}}{I^{*}},$

$v+{\epsilon }_2+{\gamma }_2+d=\frac{{\delta }_2{E}^{*}}{A^{*}}.$

代入 $\frac{\text{d}V}{\text{d}t}$ ，原方程可写为

$\frac{\text{d}V}{\text{d}t}=-\frac{K_1\mu {\left(S-S^{*}\right)}^2}{S}+g\left(y_1,y_2,y_3,y_4\right).$

令

$\frac{S}{S^{*}}=y_1,\frac{E}{E^{*}}=y_2,\frac{I}{I^{*}}=y_3,\frac{A}{A^{*}}=y_4,S^{*}{I}^{*}=a,S^{*}{A}^{*}=b,{\delta }_1{E}^{*}=c,v{A}^{*}=d,{\delta }_2{E}^{*}=e.$

$\begin{array}{c}g\left(y_1,y_2,y_3,y_4\right)=y_1{y}_3\left[-K_1\left(1-\mu \right){\beta }_{1}a+K_2\left(1-\mu \right){\beta }_{1}a\right]\\ +y_1{y}_4\left[-K_1\left(1-\mu \right){\beta }_{2}b+K_2\left(1-\mu \right){\beta }_{2}b\right]\\ +y_3\left[K_1\left(1-\mu \right){\beta }_{1}a-K_{3}c-K_{3}d\right]\\ +y_4\left[K_1\left(1-\mu \right){\beta }_{2}b+K_{3}d-K_{4}e\right]\\ +y_2\left[-K_2\left(1-\mu \right){\beta }_{1}a-K_2\left(1-\mu \right){\beta }_{2}b+K_{3}c+K_{4}e\right]\end{array}$

$\begin{array}{c}+K_3\left(c+d\right)+K_{4}e+\frac{1}{y_1}\left[-K_1\left(1-\mu \right){\beta }_{1}a-K_1\left(1-\mu \right){\beta }_{2}b\right]\\ +\frac{1}{y_2}\left[-K_2\left(1-\mu \right){\beta }_{1}a{y}_1{y}_3-K_2\left(1-\mu \right){\beta }_{2}b{y}_1{y}_4\right]\\ +\frac{1}{y_3}\left[-K_{3}c{y}_2-K_{3}d{y}_4\right]-\frac{1}{y_4}{K}_{4}e{y}_2\\ +K_1\left(1-\mu \right)\left({\beta }_{1}a+{\beta }_{2}b\right)+K_2\left(1-\mu \right)\left({\beta }_{1}a+{\beta }_{2}b\right).\end{array}$

令 $y_1{y}_3$ ， $y_1{y}_4$ ， $y_4$ ， $y_3$ ， $y_2$ 为零，可以得到

$\begin{array}{l}{K}_1=K_2,\\ K_3=\frac{\left(1-\mu \right){\beta }_{1}a}{c+d}{K}_1,\\ K_4=\left[\frac{\left(1-\mu \right){\beta }_{2}b}{e}+\frac{\left(1-\mu \right){\beta }_{1}ad}{\left(c+d\right)e}\right]K_1.\end{array}$

令 $K_1=K_2=c+d$ ，得到

$\begin{array}{l}{K}_3=\left(1-\mu \right){\beta }_{1}a{K}_1,\\ K_4=\frac{\left(1-\mu \right){\beta }_{2}b\left(c+d\right)+\left(1-\mu \right){\beta }_{1}ad}{e}.\end{array}$

则

$\begin{array}{c}g\left(y_1,y_2,y_3,y_4\right)=\left(1-\mu \right){\beta }_{1}ac\left[3-\frac{1}{y_1}-\frac{y_1{y}_3}{y_2}-\frac{y_2}{y_3}\right]\\ +\left(1-\mu \right){\beta }_{1}ad\left[4-\frac{1}{y_1}-\frac{y_1{y}_3}{y_2}-\frac{y_4}{y_3}-\frac{y_2}{y_4}\right]\\ +\left(1-\mu \right){\beta }_{2}b\left(c+d\right)\left[3-\frac{1}{y_1}-\frac{y_1{y}_4}{y_2}-\frac{y_2}{y_4}\right].\end{array}$

由算术平均值与几何平均值可知

$\begin{array}{l}\frac{1}{y_1}+\frac{y_1{y}_3}{y_2}+\frac{y_2}{y_3}\ge 3,\\ \frac{1}{y_1}+\frac{y_1{y}_3}{y_2}+\frac{y_4}{y_3}+\frac{y_2}{y_4}\ge 4,\\ \frac{1}{y_1}+\frac{y_1{y}_4}{y_2}+\frac{y_2}{y_4}\ge 3.\end{array}$

经计算可知

$\begin{array}{c}\frac{\text{d}V}{\text{d}t}=-\frac{K_1\mu {\left(S-S^{*}\right)}^2}{S}+\left(1-\mu \right){\beta }_{1}ac\left[3-\frac{1}{y_1}-\frac{y_1{y}_3}{y_2}-\frac{y_2}{y_3}\right]\\ +\left(1-\mu \right){\beta }_{1}ad\left[4-\frac{1}{y_1}-\frac{y_1{y}_3}{y_2}-\frac{y_4}{y_3}-\frac{y_2}{y_4}\right]\\ +\left(1-\mu \right){\beta }_{2}b\left(c+d\right)\left[3-\frac{1}{y_1}-\frac{y_1{y}_4}{y_2}-\frac{y_2}{y_4}\right].\end{array}$

故 $\frac{\text{d}V}{\text{d}t}\le 0$ ，当且仅当 $y_1=y_2=y_3=y_4=1$ 时， $S=S^{*},E=E^{*},I=I^{*},A=A^{*}$ 等号成立，且 $\left\{\left(S,E,I,A\right)\in R_{+}^4|S+E+I+A\le \frac{\Lambda }{\mu },\stackrel{\dot{}}{V}=0\right\}$ 的最大单点集为 $\left\{\left(S^{*},E^{*},I^{*},A^{*}\right)\right\}$ 。由LaSalle不变集原理 [19] 可得，故当 $R_0>1$ 时，系统(1)的平衡点 $P^{*}\left(S^{*},E^{*},I^{*},A^{*},Q^{*},R^{*}\right)$ 全局渐近稳定。

5. 数值模拟

取模型中各参数 $\Lambda =1$ ， $\mu =0.8$ ， $d=0.05$ ， ${\beta }_1=0.05$ ， ${\beta }_2=0.05$ ， ${\delta }_1=0.3$ ， ${\delta }_2=0.2$ ， $\upsilon =0.1$ ， ${\alpha }_1=0.07$ ， ${\alpha }_2=0.05$ ， ${\gamma }_1=0.05$ ， ${\gamma }_2=0.1$ ， ${\gamma }_3=0.05$ ， ${\epsilon }_1=0.1$ ， ${\epsilon }_2=0.2$ 。其中S，E，I，A，Q，R的初值为 $S=20$ ， $E=0$ ， $I=10$ ， $A=5$ ， $Q=5$ ， $R=0$ 。得到 $R_0<1$ ，将初值中 $\Lambda =1$ 改为 $\Lambda =8$ ，得到 $R_0>1$ 。如图1和图2：

Figure 1. Model 1

图1. 模型1

Figure 2. Model 2

图2. 模型2

从图1看出，各仓室的人数很快趋于稳定，在传染病模型中加入常规预防和隔离措施后，传染病传播能在短时间内得到有效控制，所以当 $R_0<1$ 时，无病平衡点全局渐近稳定。从图2看出，当 $R_0>1$ 时，各仓室人数在短暂波动后，最终趋于稳定，地方病平衡点的稳定性得到验证。