数学物理学报  2017, Vol. 37 Issue (3): 577-592   PDF    
扩展功能
加入收藏夹
复制引文信息
加入引用管理器
Email Alert
RSS
本文作者相关文章
韦爱举
张新建
王俊义
李科赞
一类埃博拉传染病模型的动力学分析
韦爱举1, 张新建1, 王俊义2, 李科赞1,3     
1. 桂林电子科技大学数学与计算科学学院 广西桂林 541004;
2. 广西无线宽带通信与信号处理重点实验室 广西桂林 541004;
3. 广西密码学与信息安全重点实验室 广西桂林 541004
摘要:该文考虑了人类和动物耦合传播的情况,提出了一类具有隔离仓室和动物仓室的埃博拉传染病模型.然后运用第二代生成矩阵的方法得到了基本再生数的表达式,利用稳定性和动力系统的相关理论证明了无病平衡点、边界平衡点、共存平衡点的存在性及全局稳定性,并对模型中的参数进行敏感性指数分析.最后通过数值模拟验证了这些结果的正确性.该文工作对于如何有效预防和控制埃博拉病毒的传播具有重要的意义.
关键词埃博拉病毒    基本再生数    全局稳定性    敏感性分析    
Dynamics Analysis of an Ebola Epidemic Model
Wei Aiju1, Zhang Xinjian1, Wang Junyi2, Li Kezan1,3     
1. School of Mathematics and Computing Science, Guilin University of Electronic Technology, Guangxi Guilin 541004;
2. Guangxi Key Laboratory of Wireless Wideband Communication and Signal Processing, Guangxi Guilin 541004;
3. Guangxi Key Laboratory of Cryptography and Information Security, Guangxi Guilin 541004
Abstract: In this paper, we consider the human-animal coupling epidemic where an Ebola epidemic model with isolation compartment and animal compartment is proposed. Then, the basic reproductive number is obtained by using the method of the second generation matrix. Applying the stability theory and dynamical system theory, we prove both the existence and global stability of disease-free equilibrium, boundary equilibrium and coexistence equilibrium. Moreover, sensitivity analysis of model's parameters is also performed. Finally, the theoretical results are all verified by numerical simulations. This work has important significance for guiding us to prevent and control the Ebola spread.
Key words: Ebola virus     The basic reproduction number     Global stability     Sensitivity analysis    
1 引言

纵观人类历史, 传染病一直是影响人类健康的宿敌.长久以来, 人们都想尽一切办法预防和控制传染病带来的影响, 对传染病动力学的分析一直是人们研究的热点.在文献[1]中, Kermack和McKendrick构造了著名的SIR仓室模型, 随后又建立了SIS仓室模型, 提出了研究疾病流行的“阈值理论”, 为今后人们研究传染病动力学奠定了基础.

随着人们对传染病动力学研究的深入, 仓室模型得到了许多的应用, 例如用于媒介传染、垂直传播、接触传播等不同传播方式的研究以及用于一些具体疾病的研究[2-9].

1976年, 在非洲同时爆发的两起疫情中(其中一起疫情出现在埃博拉河附近的一处村庄, 因此该病称为埃博拉病毒病), 埃博拉病毒首次出现在了人们的视野中.据世界卫生组织报道, 最近一次疫情发生于2014年3月, 西非出现首批病例, 而此次疫情是最大最复杂的埃博拉疫情, 死亡人数和病例均超过了以往疫情总和.西非这次出现的埃博拉病毒疫情的持续攀升现象强调我们需要进一步了解埃博拉病毒的流行病学特征, 以及控制干预埃博拉病毒传播的影响.在1996年, 为了模拟扎伊尔地区的埃博拉爆发情况, Fauci等建立了SIR和SEIR模型[10].在文献[11]中, 作者主要通过收集相关数据进行实验, 对早期发生的埃博拉传播动力学进行评估比较, 评估了基本公共卫生措施对疾病传播的影响以及分析了隔离控制策略对基本再生数的影响.之后也有一些学者对埃博拉病毒进行了相关研究[12-13].

埃博拉出血热是由埃博拉病毒引起的一种急性、出血性人兽共患传染病, 病死率可高达90%, 是人类目前已知最为烈性的传染病之一[14-15].埃博拉病毒具有人兽共患的特点, 该病毒既可以在动物之间传播, 又可以在人群中传播, 并且它也可以从动物传播到人.埃博拉是通过密切接触到感染动物的血液、分泌物、器官或其它体液而传到人, 随后通过人际间传播加以蔓延, 人群易感者既可以通过破损皮肤或粘膜直接接触感染者的血液、分泌物、器官或其它体液的方式感染该疾病, 也可以通过接触受到这些液体污染的表面和材料(如床上用品、衣物)等方式感染[16].不过在以往的研究中, 人们只分析该疾病在人群中流行的情况, 并未考虑动物群体中疾病的流行情况以及人和动物共存状态下疾病传播的情况.另外, 由于目前针对埃博拉病毒还没有研制出完全有效的疫苗, 因此预防和控制的一项重要措施就是控制传染源[17-18], 其中最有效可行的办法是对患者进行隔离, 不过这方面的理论研究还相当少见.

为此, 本文根据埃博拉病毒的传播特点, 对埃博拉病毒的传播特征进行了综合考虑, 建立带有隔离仓室和动物仓室的埃博拉传染病模型, 证明了解的相关性质, 利用稳定性理论研究了多种平衡点的稳定性条件, 并对模型参数作了敏感性分析, 从而把握埃博拉流行的稳态情况.此工作对于了解和研究埃博拉病毒的发生和发展情况原理, 以及有效地预防和控制此类疾病有着重要的意义.

文章结构安排如下:在第二节中, 给出了关于埃博拉传染病模型的一些定义及引理.在第三节中, 建立了带有隔离仓室和动物仓室的${S_h}{I_h}{Q_h}{R_h}{S_{\rm{b}}}{I_b}$模型, 并证明了解的正性及有界性.第四节分别证明了无病平衡点、共存平衡点以及边界平衡点的存在性和稳定性.第五节我们对模型中的参数进行了敏感性分析.在第六节中, 通过数值模拟验证了理论的正确性.第七节是本文的结论.

2 预备知识

为了后续理论分析的需要, 首先介绍本文将要用到的一些预备知识.

引理2.1  (Lasalle不变原理[19])考虑自治系统

$ \begin{equation} {x'}(t) = f(x), \end{equation} $ (2.1)

其中$f\in C(D, {{\boldsymbol{\text{R}}}^{n}}), D\subset {{\boldsymbol{\text{R}}}^{n}}$.

$V$是系统(2.1) 的定义在开子集$\Omega \subset D$内的一个Lyapunov函数, $V$$\Omega$上连续, 令

$ E = \{ x \in \Omega |V'(x) = 0\}, $

$M$是系统(2.1) 在$E$中的最大不变子集, 从$\Omega $内出发的任一正半轨${\Gamma _ + }({x_0})({x_0} \in \Omega )$, 恒在$\Omega $中且有界, 则轨线${\Gamma _ + }({x_0})$$\omega $极限集$\omega ({\Gamma _ + }({x_0})) \subset M$, 且有

$ \mathop {\lim }\limits_{t \to + \infty } {\rm dist}(x(t, {x_0}), M) = 0, $

其中dist$( \cdot, \cdot )$表示距离.

推论2.1[19] 在引理2.1的条件下, 若$M = \{ {x^*}\}$, 这里$f({x^*}) = 0$, 则系统(2.1) 的平衡点${x^*}$$\Omega $内是全局吸引的.

推论2.2[19] 设$V(x)$是系统(2.1) 的Lyapunov函数且有下界, 当$\left\| x \right\| \to + \infty$时, $V(x) \to + \infty$, 则系统(2.1) 的任一正半轨${\Gamma _ + }({x_0})$是有界的, 且当$t \to + \infty $时该轨线趋向于$E=\{x\in {{\boldsymbol{\text{R}}}^{n}}|{V}'(x)=0\}$中的最大不变集$M$.特别地, 若$M = \{ {x^*}\}, f({x^*}) = 0$, 则平衡点${x^*}$全局渐近稳定.

引理2.2 (Bendixson-Dulac判别法[20])考虑平面自治系统

$ {{x'}(t) = P(x, y)}, {{y'}(t) = Q(x, y)}. $

若在单连通域$G$内存在函数$B(x, y)\in {{C}^{1}}(G, \boldsymbol{\text{R}})$, 对任意$(x, y) \in G$, 都有

$ \begin{eqnarray*} \frac{{\partial (BP)}}{{\partial x}} + \frac{{\partial (BQ)}}{{\partial y}} \ge 0~~( \le 0), \end{eqnarray*} $

且不在$G$的任一子区域内恒为零, 则该系统不存在全部位于$G$内的闭轨线和具有有限个奇点的奇异闭轨线.

定义2.1[21]  变量$h$关于参数$l$的敏感性指数依赖于$h$$l$的微分, 定义为

$ \begin{eqnarray*} \Upsilon _l^h: = \frac{{\partial h}}{{\partial l}} \times \frac{l}{h}. \end{eqnarray*} $
3 模型的建立与解的相关性质

根据埃博拉病毒的传播特点, 我们把人群分为四类:易感者, 感染者, 隔离者以及恢复者, 并分别用${S_h}(t), {I_h}(t), {Q_h}(t)$${R_h}(t)$表示$t$时刻相应人群的数量; 把动物分为两类:易感者和染病者, 分别用$ {S_b}(t), {I_b}(t)$表示$t$时刻相应动物的数量, $t$时刻总人口数和动物总数分别为${N_h}(t) = {S_h}(t) + {I_h}(t) + {Q_h}(t) + {R_h}(t), {N_b}(t) = {S_b}(t) + {I_b}(t)$.

现作如下假设:潜伏期很短, 忽略不计; 患者治愈后具有终身免疫能力; 被隔离的人群完全断绝与外界的接触, 不具有传染性; 人口总数是变动的, 考虑出生人口和外来人口, 并且一定时期内, 一个国家或地区的人口出生率和迁移率大体是稳定的, 所以假定输入率为常数.考虑人口的自然死亡和因病死亡, 不考虑人口的迁出和输出.有关参数定义见表 1.

表 1 模型中各参数的意义及范围

基于埃博拉病毒的实际传播特点, 具有隔离仓室和动物仓室的该传染病模型可表示如下

$ \begin{equation} \left\{ {\begin{array}{*{20}{l}} {S_h'(t) = {\Lambda _h} - {\beta _h}{S_h}{I_h} - {\beta _{bh}}{S_h}{I_b} - {\mu _h}{S_h}}, \\ {I_h'(t) = {\beta _h}{S_h}{I_h} + {\beta _{bh}}{S_h}{I_b} - {\mu _h}{I_h} - {\alpha _h}{I_h} - {\delta _h}{I_h} - {\gamma _h}{I_h}}, \\ {Q_h'(t) = {\delta _h}{I_h} - {\mu _h}{Q_h} - {d_h}{Q_h} - {\varepsilon _h}{Q_h}}, \\ {R_h'(t) = {\gamma _h}{I_h} + {\varepsilon _h}{Q_h} - {\mu _h}{R_h}}, \\ {S_b'(t) = {\Lambda _b} - {\beta _b}{S_b}{I_b} - {\mu _b}{S_b}}, \\ {I_b'(t) = {\beta _b}{S_b}{I_b} - {\mu _b}{I_b} - {\alpha _b}{I_b}}. \end{array}} \right. \end{equation} $ (3.1)

在以往对埃博拉传染病的研究中, 人们只分析了人群中该疾病的流行情况, 并未考虑其在人群和动物群体中共同存在时的流行稳态问题, 并且没有考虑隔离人群的情况.因此本文根据埃博拉病毒具有人兽共患的特点以及目前最有效的方式是将已感染该疾病的人群进行隔离的实际情况, 建立了如(3.1) 式所示的模型.此外, 根据各个变量间的转移关系, 本文所建模型(3.1) 的仓室图见图 1.

图 1 埃博拉传染病模型的传播机制图

由系统(3.1) 的前面四个方程可知

$ \begin{eqnarray*} {({S_h} + {I_h} + {Q_h} + {R_h})'} \le {\Lambda _h} - {\mu _h}({S_h} + {I_h} + {Q_h} + {R_h}), \end{eqnarray*} $

解得

$ \begin{eqnarray*} {S_h}(t) + {I_h}(t) + {Q_h}(t) + {R_h}(t) \le \frac{{{\Lambda _h}}}{{{\mu _h}}} + ({S_h}(0) + {I_h}(0) + {Q_h}(0) + {R_h}(0) - \frac{{{\Lambda _h}}}{{{\mu _h}}}){{\rm e}^{ - {\mu _h}t}}. \end{eqnarray*} $

由系统(3.1) 的后面两个方程可知

$ \begin{eqnarray*} {({S_b} + {I_b})'} \le {\Lambda _b} - {\mu _b}({S_b} + {I_b}), \end{eqnarray*} $

解得

$ \begin{eqnarray*} {S_b}(t) + {I_b}(t) \le \frac{{{\Lambda _b}}}{{{\mu _b}}} + ({S_b}(0) + {I_b}(0) - \frac{{{\Lambda _b}}}{{{\mu _b}}}){{\rm e}^{ - {\mu _b}t}}. \end{eqnarray*} $

$ \begin{eqnarray*} \Omega = \{ ({S_h}, {I_h}, {Q_h}, {R_h}, {S_b}, {I_b}) \in R_ + ^6|{S_h} + {I_h} + {Q_h} + {R_h} \le \frac{{{\Lambda _h}}}{{{\mu _h}}}, {S_b} + {I_b} \le \frac{{{\Lambda _b}}}{{{\mu _b}}}\}, \end{eqnarray*} $

故系统(3.1) 的所有正半解只要从$\Omega$中出发, 必定含在$\Omega$中, 因此可得到如下引理.

引理3.1 集合$\Omega$是系统(3.1) 的正向不变集.

4 平衡点的存在性及稳定性

为了弄清楚埃博拉病毒在两个种群(即人群和动物群体)中流行的稳态情况, 下面我们将运用相关的理论知识分别研究无病平衡点, 共存平衡点以及边界平衡点的存在性和稳定性问题.

4.1 无病平衡点的存在性及稳定性

我们利用第二代生成矩阵的方法计算基本再生数.由文献[22]的相关定义可知

$ \begin{eqnarray*} F = \left[{\begin{array}{*{20}{c}} {{\beta _h}{S_h}}~~&{{\beta _{bh}}{S_h}}\\ 0~~&{{\beta _b}{S_b}} \end{array}} \right], \end{eqnarray*} $
$ \begin{eqnarray*} V = \left[{\begin{array}{*{20}{c}} {{\mu _h} + {\alpha _h} + {\delta _h} + {r_h}}~~&0\\ 0~~&{{\mu _b} + {\alpha _b}} \end{array}} \right], \end{eqnarray*} $

因此

$ \begin{eqnarray*} F{V^{ - 1}} = \left[{\begin{array}{*{20}{c}} {\displaystyle\frac{{{\beta _h}{S_h}}}{{{\mu _h} + {\alpha _h} + {\delta _h} + {r_h}}}}~~&{\displaystyle\frac{{{\beta _{bh}}{S_h}}}{{{\mu _b} + {\alpha _b}}}}\\[3mm] 0~~&{\displaystyle\frac{{{\beta _b}{S_b}}}{{{\mu _b} + {\alpha _b}}}} \end{array}} \right]. \end{eqnarray*} $

从而基本再生数为

$ \begin{eqnarray*} {R_0}= \rho (F{V^{ - 1}}) = \max \{ {R_{01}}, {R_{02}}\}, \end{eqnarray*} $

其中$\rho $表示矩阵$F{V^{ - 1}}$的谱半径, ${R_{01}}, {R_{02}}$分别表示人群患者和动物患者的基本再生数, 且

$ {R_{01}} =\frac{{{\beta _h}{\Lambda _h}}}{{{\mu _h}{{({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})}}}}, {R_{02}}= \frac{{{\beta _b}{\Lambda _b}}}{{{\mu _b}({\mu _b} + {\alpha _b})}}. $

定理4.1.1 当${R_0} <1$时, 系统(3.1) 存在无病平衡点${E_0}(\frac{{{\Lambda _h}}}{{{\mu _h}}}, 0, 0, 0, \frac{{{\Lambda _b}}}{{{\mu _b}}}, 0)$, 且无病平衡点在$\Omega$内是局部渐近稳定的; 当${R_0}>1$时, 无病平衡点${E_0}$是不稳定的.

 令(3.1) 式右端全部为零且${I_h} = 0, {I_b} = 0$

$ {S_h} = \frac{{{\Lambda _h}}}{{{\mu _h}}}, ~{Q_h} = 0, ~{R_h} = 0, ~{S_b} = \frac{{{\Lambda _b}}}{{{\mu _b}}}, $

因此存在无病平衡点${E_0}(\frac{{{\Lambda _h}}}{{{\mu _h}}}, 0, 0, 0, \frac{{{\Lambda _b}}}{{{\mu _b}}}, 0)$.根据第二代生成矩阵分析方法, 由文献[22]中定理2可知, 当${R_0} <1$时, 系统(3.1) 的无病平衡点${E_0}$是局部渐近稳定的.当${R_0}>1$时, ${E_0}$是不稳定的.定理得证.

定理4.1.2 当${R_0} <1$时, 无病平衡点${E_0}(\frac{{{\Lambda _h}}}{{{\mu _h}}}, 0, 0, 0, \frac{{{\Lambda _b}}}{{{\mu _b}}}, 0)$$\Omega$内是全局渐近稳定的.

 由于系统(3.1) 的其它方程均不含$Q$$R$, 所以只需要考虑如下子系统

$ \begin{equation} \left\{ {\begin{array}{*{20}{l}} {S_h'(t) = {\Lambda _h} - {\beta _h}{S_h}{I_h} - {\beta _{bh}}{S_h}{I_b} - {\mu _h}{S_h}}, \\ {I_h'(t) = {\beta _h}{S_h}{I_h} + {\beta _{bh}}{S_h}{I_b} - {\mu _h}{I_h} - {\alpha _h}{I_h} - {\delta _h}{I_h} - {\gamma _h}{I_h}}, \\ {S_b'(t) = {\Lambda _b} - {\beta _b}{S_b}{I_b} - {\mu _b}{S_b}}, \\ {I_b'(t) = {\beta _b}{S_b}{I_b} - {\mu _b}{I_b} - {\alpha _b}{I_b}}. \end{array}} \right. \end{equation} $ (4.1)

由定理4.1.1可知, 当${R_0} <1$时, 无病平衡点$ {E_0}$是局部渐近稳定的, 接下来我们只需要证明$(S_h^0, 0, S_b^0, 0)$是系统(4.1) 的全局吸引子.

我们考虑系统(4.1) 的一个解$({S_h}(t), {I_h}(t), {S_b}(t), {I_b}(t))$, 根据系统(4.1) 的第二和第四个方程可得

$ \begin{equation} \left\{ {\begin{array}{*{20}{l}} {I_h'(t) \le {\beta _h}S_h^0{I_h} + {\beta _{bh}}S_h^0{I_b} - {\mu _h}{I_h} - {\alpha _h}{I_h} - {\delta _h}{I_h} - {\gamma _h}{I_h}}, \\ {I_b'(t) \le {\beta _b}S_b^0{I_b} - {\mu _b}{I_b} - {\alpha _b}{I_b}}. \end{array}} \right. \end{equation} $ (4.2)

我们先考虑以下微分方程

$ \begin{equation} \left\{ {\begin{array}{*{20}{l}} { {X_h'} = {\beta _h}S_h^0{X_h} + {\beta _{bh}}S_h^0{Y_b} - {\mu _h}{X_h} - {\alpha _h}{X_h} - {\delta _h}{X_h} - {\gamma _h}{X_h}}, \\ { {Y_b'}= {\beta _b}S_b^0{Y_b} - {\mu _b}{Y_b} - {\alpha _b}{Y_b}}, \end{array}} \right. \end{equation} $ (4.3)

因为系统(4.3) 是一个线性系统, 所以(4.3) 的全局稳定性可由以下雅可比矩阵$J$决定

$ \begin{eqnarray*} J = \left[{\begin{array}{*{20}{c}} {{\beta _h}S_h^0-({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})}&{{\beta _{bh}}S_h^0}\\ 0&{{\beta _b}S_b^0-({\mu _b} + {\alpha _b})} \end{array}} \right]. \end{eqnarray*} $

又因为

$ \begin{eqnarray*} \begin{array}{c} {\beta _b}S_b^0 - ({\mu _b} + {\alpha _b}) = ({\mu _b} + {\alpha _b})({R_{02}} - 1), \\ {\beta _h}S_h^0 - ({\mu _h} + {\alpha _h} + {\delta _h} + {r_h}) = ({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})({R_{01}} - 1), \end{array} \end{eqnarray*} $

所以当${R_{02}} <1, {R_{01}} < 1$${R_0} <1$时, $J$的所有特征值均有负实部, 所以系统(4.3) 满足

$ \begin{eqnarray*} \mathop {\lim }\limits_{t \to +\infty }{{{X_h}} }(t) = 0, \mathop {\lim }\limits_{t \to +\infty } { {{Y_b}}}(t) = 0. \end{eqnarray*} $

由比较原理可知, 当$t \to + \infty$时, 系统(4.1) 满足$I_h(t) \to 0, I_b(t) \to 0$.因此(4.1) 的极限系统为

$ \begin{equation} \left\{ {\begin{array}{*{20}{c}} {S_h'(t) = {\Lambda _h} - {\mu _h}{S_h}, }\\ \begin{array}{l} S_b'(t) = {\Lambda _b} - {\mu _b}{S_b}.\\ \end{array} \end{array}} \right. \end{equation} $ (4.4)

显然系统(4.4) 的平衡点$(S_h^0, S_b^0)$是全局渐近稳定的.由极限系统理论[19]可知, $(S_h^0, 0, S_b^0, 0)$是系统(4.1) 的全局吸引子.因此无病平衡点${E_0}(\frac{{{\Lambda _h}}}{{{\mu _h}}}, 0, 0, 0, \frac{{{\Lambda _b}}}{{{\mu _b}}}, 0)$$\Omega$内是全局渐近稳定的.定理得证.

4.2 共存平衡点的存在性及稳定性

定理4.2 当${R_{02}}> 1$时, 系统(3.1) 存在共存平衡点$\mathop E\limits^ - (\mathop {{S_h}}\limits^ -, \mathop {{I_h}}\limits^ -, \mathop {{Q_h}}\limits^ -, \mathop {{R_h}}\limits^ -, \mathop {{S_b}}\limits^ -, \mathop {{I_b}}\limits^ - )$, 且共存平衡点在$\Omega $内是全局渐近稳定的.

 令(3.1) 式右端全部为零得

$ \begin{eqnarray*} % to remove numbering (before each equation) \mathop {{S_h}}\limits^ - &=& \frac{{{\Lambda _h} - ({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})\mathop {{I_h}}\limits^ - }}{{{\mu _h}}}, \\ \mathop {{R_h}}\limits^ - &=& \frac{{{r_h}({\mu _h} + {d_h} + {\varepsilon _h}) + {\varepsilon _h}{\delta _h}}}{{{\mu _h}({\mu _h} + {d_h} + {\varepsilon _h})}}\mathop {{I_h}}\limits^ -, \\ \mathop {{Q_h}}\limits^ - &=& \frac{{{\delta _h}}}{{{\mu _h} + {d_h} + {\varepsilon _h}}}\mathop {{I_h}}\limits^ -, \\ \mathop {{S_b}}\limits^ - &=& \frac{{{\mu _b} + {\alpha _b}}}{{{\beta _b}}}, \\ \mathop {{I_b}}\limits^ - &=&\frac{{{\mu _b}}}{{{\beta _b}}}({R_{02}} - 1), \end{eqnarray*} $

其中$\mathop {{I_h}}\limits^ - $满足

$ \begin{equation} a_1{{\overline I _h}^2} + a_2{{\overline I _h}} +a_3 = 0. \end{equation} $ (4.5)

其中$a_1={\beta _h}({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})$, $a_2=({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})({\beta _{bh}}\mathop {{I_b}}\limits^ - + {\mu _h}) - {\Lambda _h}{\beta _h}$, $a_3=- {\Lambda _h}{\beta _{bh}}\mathop {{I_b}}\limits^ -$.由于$a_1>0, a_3 <0$, 故方程(4.5) 存在唯一的正解$\mathop {{I_h}}\limits^ - $.

为了研究共存平衡点的稳定性问题, 我们先研究以下动物系统

$ \begin{equation} \begin{array}{l} S_b'(t) = {\Lambda _b} - {\beta _b}{S_b}{I_b} - {\mu _b}{S_b} \buildrel \Delta \over = H({S_b}, {I_b}), \\ I_b'(t) = {\beta _b}{S_b}{I_b} - {\mu _b}{I_b} - {\alpha _b}{I_b} \buildrel \Delta \over = G({S_b}, {I_b}). \end{array} \end{equation} $ (4.6)

系统(4.6) 在平衡点$(\mathop {{S_b}}\limits^ -, \mathop {{I_b}}\limits^ - )$的雅可比矩阵为

$ \begin{eqnarray*} {J_1} = \left[{\begin{array}{*{20}{c}} {-\displaystyle\frac{{{\beta _b}{\Lambda _b}}}{{{\mu _b} + {\alpha _b}}}}&{-({\mu _b} + {\alpha _b})}\\ {{\mu _b}({R_{02}}-1)}&0 \end{array}} \right]. \end{eqnarray*} $

雅可比矩阵${J_1}$对应的特征方程为

$ \begin{eqnarray*} {\lambda ^2} + \frac{{{\beta _b}{\Lambda _b}}}{{{\mu _b} + {\alpha _b}}}\lambda + ({\mu _b} + {\alpha _b}){\mu _b}({R_{02}} - 1) = 0. \end{eqnarray*} $

因为当${R_{02}} > 1$时, 有

$ \begin{eqnarray*} \begin{array}{l} {\lambda _{\rm{1}}}{\rm{ + }}{\lambda _{\rm{2}}}{\rm{ = - }}\displaystyle\frac{{{\beta _b}{\Lambda _b}}}{{{\mu _b} + {\alpha _b}}} < 0, \\ {\lambda _{\rm{1}}}{\lambda _{\rm{2}}} = ({\mu _b} + {\alpha _b}){\mu _b}({R_{02}} - 1) > 0. \end{array} \end{eqnarray*} $

所以${\lambda _{\rm{1}}} < 0, {\lambda _2} < 0$, 系统(4.6) 的平衡点$(\mathop {{S_b}}\limits^ -, \mathop {{I_b}}\limits^ - )$是局部渐近稳定的.

取Dulac函数为$B = \frac{1}{{{I_b}}}$, 则

$ \begin{eqnarray*} \frac{{\partial (BH)}}{{\partial {S_b}}} + \frac{{\partial (BG)}}{{\partial {I_b}}} = - {\beta _b} - \frac{{{\mu _b}}}{{{I_b}}} < 0. \end{eqnarray*} $

故系统(4.6) 在$\{({{S}_{b}}, {{I}_{b}})\in \boldsymbol{\text{R}}_{+}^{2}|{{S}_{b}}+{{I}_{b}}\le \frac{{{\Lambda }_{b}}}{{{\mu }_{b}}}\}$内无闭轨, 结合平衡点$(\mathop {{S_b}}\limits^ -, \mathop {{I_b}}\limits^ - )$的局部稳定性可知, 当${R_{02}} > 1$时, 动物系统(4.6) 的平衡点$(\mathop {{S_b}}\limits^ -, \mathop {{I_b}}\limits^ - )$是全局渐近稳定的.

接下来我们证明系统(4.1) 的共存平衡点的全局稳定性.因为当${R_{02}} > 1$时, 动物系统(4.6) 的平衡点是全局渐近稳定的, 即$\mathop {\lim }\limits_{t \to +\infty } {S_b} = \mathop {{S_b}}\limits^ -, \mathop {\lim }\limits_{t \to +\infty } {I_b} = \mathop {{I_b}}\limits^ - $, 所以只需考虑如下极限子系统

$ \begin{eqnarray*} \begin{array}{*{20}{l}} {S_h'(t) = {\Lambda _h} - {\beta _h}{S_h}{I_h} - {\beta _{bh}}{S_h}\mathop {{I_b}}\limits^ - - {\mu _h}{S_h}, }\\ {I_h'(t) = {\beta _h}{S_h}{I_h} + {\beta _{bh}}{S_h}\mathop {{I_b}}\limits^ - - {\mu _h}{I_h} - {\alpha _h}{I_h} - {\delta _h}{I_h} - {\gamma _h}{I_h}.} \end{array} \end{eqnarray*} $

构造如下形式的Lyapunov函数

$ \begin{eqnarray*} U({S_h}, {I_h}) =\bigg( {{S_h} - \mathop {{S_h}}\limits^ - - \mathop {{S_h}}\limits^ - \ln \frac{{{S_h}}} {{\mathop {{S_h}}\limits^ - }}} \bigg) + \bigg( {{I_h} - \mathop {{I_h}}\limits^ - - \mathop {{I_h}}\limits^ - \ln \frac{{{I_h}}}{{\mathop {{I_h}} \limits^ - }}} \bigg). \end{eqnarray*} $

显然$U({S_h}, {I_h})\geq0$, 并且当$\left\| {({S_h}, {I_h})} \right\| \to + \infty $时, $U({S_h}, {I_h}) \to + \infty$.因为系统在共存平衡点处满足以下关系

$ \begin{eqnarray*} \begin{array}{l} {\Lambda _h} = {\beta _h}\mathop {{S_h}}\limits^ - \mathop {{I_h}}\limits^ - + {\beta _{bh}}\mathop {{S_h}}\limits^ - \mathop {{I_b}}\limits^ - + {\mu _h}\mathop {{S_h}}\limits^ -, \\ {\mu _h} + {\alpha _h} + {\delta _h} + {\gamma _h} = \displaystyle\frac{{{\beta _h}\mathop {{S_h}}\limits^ - \mathop {{I_h}}\limits^ - + {\beta _{bh}}\mathop {{S_h}}\limits^ - \mathop {{I_b}}\limits^ - }}{{\mathop {{I_h}}\limits^ - }}, \end{array} \end{eqnarray*} $

因此可得

$ \begin{eqnarray*} % to remove numbering (before each equation) {U'} &=& \bigg( {1 - \frac{{\mathop {{S_h}}\limits^ - }}{{{S_h}}}}\bigg)S_h' + \bigg( {1 - \frac{{\mathop {{I_h}}\limits^ - }}{{{I_h}}}} \bigg)I_h' \\ &=&\bigg( {1 - \frac{{\mathop {{S_h}}\limits^ - }}{{{S_h}}}} \bigg) \bigg[{{\beta _h}\mathop {{S_h}}\limits^-\mathop {{I_h}}\limits^-+ {\beta _{bh}}\mathop {{S_h}}\limits^-\mathop {{I_b}}\limits^ - + {\mu _h}\mathop {{S_h}}\limits^ - - {\beta _h}{S_h}{I_h} - {\beta _{bh}}{S_h}{I_b} - {\mu _h}{S_h}}\bigg] \\ & & +\bigg( {1 - \frac{{\mathop {{I_h}}\limits^ - }}{{{I_h}}}} \bigg) \bigg[{{\beta _h}{S_h}{I_h} + {\beta _{bh}}{S_h}{I_b}-\frac{{{\beta _h}\mathop {{S_h}}\limits^-\mathop {{I_h}}\limits^-+ {\beta _{bh}}\mathop {{S_h}}\limits^ - \mathop {{I_b}}\limits^ - }}{{\mathop {{I_h}}\limits^ - }}{I_h}} \bigg]\\ &=&\Big( {{\mu _h} + {\beta _h}\mathop {{I_h}}\limits^ - }\Big) \mathop {{S_h}}\limits^ - \bigg( {2 - \frac{{{S_h}}}{{\mathop {{S_h}}\limits^ - }} - \frac{{\mathop {{S_h}}\limits^ - }}{{\mathop {{S_h}}\limits^{} }}} \bigg) + {\beta _{bh}}\mathop {{S_h}}\limits^ - \mathop {{I_b}}\limits^ - \bigg( {3 - \frac{{\mathop {{S_h}}\limits^ - }}{{\mathop {{S_h}}\limits^{} }} - \frac{{{I_h}}}{{\mathop {{I_h}}\limits^ - }} - \frac{{{S_h}}}{{\mathop {{S_h}}\limits^ - }}\frac{{\mathop {{I_h}} \limits^ - }}{{\mathop {{I_h}}\limits^{} }}} \bigg) . \end{eqnarray*} $

根据均值不等式可得

$ {2 - \frac{{{S_h}}}{{\mathop {{S_h}}\limits^ - }} - \frac{{\mathop {{S_h}}\limits^ - }}{{\mathop {{S_h}}\limits^{} }}} \le 0, {3 - \frac{{\mathop {{S_h}}\limits^ - }}{{\mathop {{S_h}}\limits^{} }} - \frac{{{I_h}}}{{\mathop {{I_h}}\limits^ - }} - \frac{{{S_h}}}{{\mathop {{S_h}}\limits^ - }}\frac{{\mathop {{I_h}}\limits^ - }}{{\mathop {{I_h}}\limits^{} }}}\le 0, $

$U' \le 0$.当且仅当${S_h} = \mathop {{S_h}}\limits^ -, {I_h} = \mathop {{I_h}}\limits^ - $时, $U' = 0$.所以由Lasalle不变原理(见引理2.1) 及推论2.2可知, 当${R_{02}}> 1$时, 共存平衡点是全局渐近稳定的.定理得证.

4.3 边界平衡点的存在性及稳定性

由(4.5) 式可知, 当${\mathop {{I_h}}\limits^ - } = 0$时, 有${\mathop {{I_b}}\limits^ - } = 0$.所以当人群中感染者灭绝时, 除了无病平衡点以外, 不存在如$(S_h^*, 0, Q_h^*, R_h^*, S_b^*, I_b^*)$这样的边界平衡点.下面我们讨论另一类边界平衡点(即动物感染者灭绝时)的存在性和稳定性.

当动物感染者灭绝时, 即${I_b} = 0$, 设其对应的边界平衡点为${E_h}(S_h^*, I_h^*, Q_h^*, R_h^*, S_b^*, 0)$, 其中

$ \begin{eqnarray*} % to remove numbering (before each equation) S_b^*&=&\frac{{{\Lambda _b}}}{{{\mu _b}}}, \\ S_h^*&=&\frac{{{\Lambda _h} - ({\mu _h} + {\alpha _h} + {\delta _h} + {r_h}){I_h}^ * }}{{{\mu _h}}}, \\ Q_h^* &=& \frac{{{\delta _h}}}{{{\mu _h} + {d_h} + {\varepsilon _h}}}{I_h}^ *, \\ R_h^* &=& \frac{{{r_h}({\mu _h} + {d_h} + {\varepsilon _h}) + {\varepsilon _h}{\delta _h}}}{{{\mu _h}({\mu _h} + {d_h} + {\varepsilon _h})}}{I_h}^ *, \end{eqnarray*} $

$I_h^*$满足

$ \begin{eqnarray*} {\beta _h}({\mu _h} + {\alpha _h} + {\delta _h} + {r_h}){I{_h^*}^2} + [{\mu _h}({\mu _h} + {\alpha _h} + {\delta _h} + {r_h})-{\Lambda _h}{\beta _h}]{I_h^*} = 0. \end{eqnarray*} $

解得

$ \begin{eqnarray*} \displaystyle I_h^* = \frac{{{\Lambda _h}}}{{{\mu _h} + {\alpha _h} + {\delta _h} + {r_h}}} - \frac{{{\mu _h}}}{{{\beta _h}}} = \frac{{{\mu _h}}}{{{\beta _h}}}({R_{01}} - 1). \end{eqnarray*} $

因此当${R_{01}}> 1$时, 有$I_h^* > 0$, 即模型存在边界平衡点${E_h}$.令${\Omega _b} = \{ ({S_h}, {I_h}, {Q_h}, {R_h}, {S_b}, {I_b}) \in \Omega |{I_b} = 0\}$, 我们可得到以下结论.

定理4.3 当${R_{02}} < 1 < {R_{01}}$时, ${E_h}$$\Omega /{\Omega _b}$内是全局渐近稳定的.

 系统(3.1) 在边界平衡点${E_h}(S_h^*, I_h^*, Q_h^*, R_h^*, S_b^*, 0)$处的雅可比矩阵为

$ \begin{eqnarray*} {J_2} = \left[{\begin{array}{*{20}{c}} {-{\beta _h}I_h^*-{\mu _h}}& ~~{-{\beta _h}S_h^*}~~&0&0&0&{ - {\beta _{bh}}S_h^*}\\ {{\beta _h}I_h^*}&0&0&0&0&{{\beta _{bh}}S_h^*}\\ 0&{{\delta _h}}&{ - {\mu _h} - {d_h} - {\varepsilon _h}}&0&0&0\\ 0&{{r_h}}&{{\varepsilon _h}}&~~{ - {\mu _h}}~~&0&0\\ 0&0&0&0&{ - {\mu _b}}&{ - {\beta _b}S_b^*}\\ 0&0&0&0&0&{{\beta _b}S_b^* - {\mu _b} - {\alpha _b}} \end{array}} \right]. \end{eqnarray*} $

雅可比矩阵${J_2}$对应的特征方程为

$ \begin{eqnarray*} (\lambda + {\mu _h} + {d_h} + {\varepsilon _h})(\lambda + {\mu _h})(\lambda + {\mu _b})[\lambda-({\beta _b}\frac{{{\Lambda _b}}}{{{\mu _b}}}-{\mu _b}-{\alpha _b})][{\lambda ^2} + ({\beta _h}I_h^* + {\mu _h})\lambda + \beta _h^2S_h^*I_h^*] = 0. \end{eqnarray*} $

显然${\lambda _3} = - {\mu _h} - {d_h} - {\varepsilon _h}, {\lambda _4} = - {\mu _h}, {\lambda _5} = - {\mu _b}, {\lambda _6} = {\beta _b}\frac{{{\Lambda _b}}}{{{\mu _b}}} - {\mu _b} - {\alpha _b} < 0$, 且其它两个特征值${\lambda _1}$${\lambda _2}$满足

$ \begin{eqnarray*} P(\lambda ) = {\lambda ^2} + ({\beta _h}I_h^* + {\mu _h})\lambda + \beta _h^2S_h^*I_h^*{\rm{ = 0}}. \end{eqnarray*} $

因为

$ {\lambda _{\rm{1}}}{\rm{ + }}{\lambda _{\rm{2}}}{\rm{ = - (}}{\beta _h}I_h^* + {\mu _h}) < 0, {\lambda _{\rm{1}}}{\lambda _{\rm{2}}} = \beta _h^2S_h^*I_h^* > 0, $

所以${\lambda _{\rm{1}}} < 0, {\lambda _2} < 0, $故系统(3.1) 的边界平衡点${E_h}(S_h^*, I_h^*, Q_h^*, R_h^*, S_b^*, 0)$是局部渐近稳定的.

接下来证明边界平衡点的全局吸引性, 以此证明它的全局稳定性.构造Lyapunov函数$V = {I_b}$, 则

$ \frac{{dV}}{{dt}}{|_{({\rm{3.1}})}} = {\beta _b}{S_b}{I_b} - {\mu _b}{I_b} - {\alpha _b}{I_b} \le ({\mu _b} + {\alpha _b})({R_{02}} - 1){I_b}. $

因此当${R_{02}} < 1$时, 就有$V' \le 0$.当且仅当${I_b}(t) = 0$$V' = 0$.由Lasalle不变原理知$\mathop {\lim }\limits_{t \to + \infty } {I_b}(t) = 0$.考虑极限系统

$ \begin{eqnarray*} {S_b'(t) = {\Lambda _b} - {\mu _b}{S_b}.} \end{eqnarray*} $

易知当${R_{02}}< 1$时, 该子系统的平衡点是全局渐近稳定的, 从而

$ \begin{eqnarray*} \left\{ {\begin{array}{*{20}{l}} {\mathop {\lim }\limits_{t \to +\infty } {S_b}(t) = S_b^*, }\\ {\mathop {\lim }\limits_{t \to +\infty } {I_b}(t) = 0.} \end{array}} \right. \end{eqnarray*} $

考虑系统(4.1) 的极限子系统为

$ \begin{equation} \begin{array}{*{20}{l}} {S_h'(t) = {\Lambda _h} - {\beta _h}{S_h}{I_h} - {\mu _h}{S_h} \buildrel \Delta \over = H({S_h}, {I_h}), }\\ {I_h'(t) = {\beta _h}{S_h}{I_h} - {\mu _h}{I_h} - {\alpha _h}{I_h} - {\delta _h}{I_h} - {\gamma _h}{I_h} \buildrel \Delta \over = G({S_h}, {I_h}).} \end{array} \end{equation} $

取Dulac函数为$B = \frac{1}{{{I_h}}}$, 则

$ \begin{eqnarray*} \frac{{\partial (BH)}}{{\partial {S_h}}} + \frac{{\partial (BG)}}{{\partial {I_h}}} = - {\beta _h} - \frac{{{\mu _h}}}{{{I_h}}} < 0. \end{eqnarray*} $

由引理2.2可知, 系统(4.7) 在$\{({{S}_{h}},{{I}_{h}})\in \boldsymbol{\text{R}}_{+}^{2}|{{S}_{h}}+{{I}_{h}}\le \frac{{{\Lambda }_{h}}}{{{\mu }_{h}}}\}$内无闭轨, 平衡点$(S_h^*, I_h^*)$是全局渐近稳定的.故当$t \to + \infty$时, $({S_h}, {I_h}) \to (S_h^*, I_h^*)$.

又因为

$ \begin{eqnarray*} \begin{array}{l} Q_h'(t) = {\delta _h}{I_h} - {\mu _h}{Q_h} - {d_h}{Q_h} - {\varepsilon _h}{Q_h}, \\ R_h'(t) = {\gamma _h}{I_h} + {\varepsilon _h}{Q_h} - {\mu _h}{R_h}, \end{array} \end{eqnarray*} $

所以当$t \to + \infty$时, 有

$ \begin{eqnarray*} && \mathop {R_h} \to \displaystyle\frac{{{r_h}({\mu _h} + {d_h} + {\varepsilon _h}) + {\varepsilon _h}{\delta _h}}}{{{\mu _h}({\mu _h} + {d_h} + {\varepsilon _h})}}\mathop {{I_h}^ *}, \\ &&\mathop {Q_h} \to \displaystyle\frac{{{\delta _h}}}{{{\mu _h} + {d_h} + {\varepsilon _h}}}{I_h}^ * . \end{eqnarray*} $

${I_b}=0$时, 虽然同时存在无病平衡点和边界平衡点, 但当${R_{02}} < 1 < {R_{01}}$时, 即${R_0} = \max \{ {R_{01}}, {R_{02}}\} > 1$时, 无病平衡点是不稳定的, 再结合边界平衡点的局部稳定性和全局吸引性, 系统(3.1) 的边界平衡点是全局渐近稳定的.定理得证.

5 参数敏感性分析

基本再生数${R_0}$的敏感性影响着系统参数对于模型(3.1) 的鲁棒性, 通过敏感性分析可以得到对于${R_0}$影响较大的参数, 从而提出更好的防控措施.下面对基本再生数${R_0}$中的参数进行分析.

由定义2.1, 可计算得到基本再生数${R_0}$关于模型中各个参数的敏感性指数.例如${\beta _h}$, ${\mu _h}$, ${\alpha _h}$这三个参数的敏感性指数可推导如下

$ \begin{eqnarray*} && \Upsilon _{{\beta _h}}^{R_0}: = \frac{{\partial {R_0}}}{{\partial {\beta _h}}} \times \frac{{{\beta _h}}}{{R_0}}{\rm{ = }}\left\{ {\begin{array}{*{20}{c}} {1, ~{R_{01}} > {R_{02}}}, \\ {0, ~{R_{01}} < {R_{02}}}, \end{array}} \right.\\ && \Upsilon _{{\mu _h}}^{R_0}: = \frac{{\partial {R_0}}}{{\partial {\mu _h}}} \times \frac{{{\mu _h}}}{{R_0}}{\rm{ = }}\left\{ {\begin{array}{{ll}} {\rm{ - }}\frac{{{\mu _h}}}{{{\mu _h}{\rm{ + }}{\alpha _h}{\rm{ + }}{\delta _h}{\rm{ + }}{r_h}}} - 1, ~&{R_{01}} > {R_{02}}, \\ 0, ~&{R_{01}} < {R_{02}}, \end{array}} \right.\\ && \Upsilon _{{\alpha _h}}^{R_0}: = \frac{{\partial {R_0}}}{{\partial {\alpha _h}}} \times \frac{{{\alpha _h}}}{{R_0}}{\rm{ = }}\left\{ {\begin{array}{{ll}} {\rm{ - }}\frac{{{\alpha _h}}}{{{\mu _h}{\rm{ + }}{\alpha _h}{\rm{ + }} {\delta _h}{\rm{ + }}{r_h}}}, ~&{R_{01}} > {R_{02}}, \\ 0, ~&{R_{01}} < {R_{02}}, \end{array}} \right. \end{eqnarray*} $

其它参数的敏感性指数可类推得到.

${R_{01}} > {R_{02}}$时, ${R_0} = \max \{ {R_{01}}, {R_{02}}\} = {R_{01}}$, 所得结果如表 2所示.

表 2 基本再生数${R_0}$关于模型中各参数的敏感性分析

同理, 当${R_{01}} < {R_{02}}$时, ${R_0} = \max \{ {R_{01}}, {R_{02}}\} = {R_{02}}$, 则各参数的敏感性指数如表 3所示.

表 3 基本再生数${R_0}$关于模型中各参数的敏感性分析

通过参数的敏感性分析可知, 当${R_{01}}>{R_{02}}$时, ${R_0}=\max \{{R_{01}}, {R_{02}}\}={R_{01}}$, 最敏感的参数之一是${\beta _h}$, 即人与人之间的传染率.由于$\Upsilon _{{\beta _h}}^{R_0}: = \frac{{\partial {R_0}}}{{\partial {\beta _h}}} \times \frac{{{\beta _h}}} {{R_0}}{\rm{ = 1, }}$所以人与人之间的传染率${\beta _h}$增大(或减小)10%时, 基本再生数${R_0}$也随着增大(或减小)10%.类似的, 当${R_{02}} > {R_{01}}$时, ${R_0}= \max \{ {R_{01}}, {R_{02}}\} = {R_{02}}$, 动物与动物之间的传染率${\beta _b}$增大(或减小)10%时, 基本再生数${R_0}$也随着增大(或减小)10%.因此, 为了更为有效地预防和控制埃博拉病毒, 我们应该尽量不要去埃博拉爆发的区域, 避免直接接触埃博拉人群患者和动物患者.医护人员或从事屠宰动物工作的相关人员应该做好充分的保护措施, 避免直接接触感染者的血液、分泌物、器官或其它体液等.对于受到这些液体污染的表面和材料(如床上用品、衣物等), 应该及时进行有效的处理, 以免再次受到感染.

6 数值模拟

下面通过数值模拟验证系统无病平衡点、共存平衡点和边界平衡点的稳定性.不失一般性, 我们考虑如下3种情形.

情形1  取参数${\Lambda _h} = 10, {\Lambda _b} = 50, {\mu _h} = 0.1, {\mu _b} = 0.8, {r_h} = 0.32, {\alpha _h} = 0.35, $${\delta _h} = 0.4, {d_h} = 0.25, {\varepsilon _h} = 0.25, {\alpha _b} = 0.4, $${\beta _h} = 0.01, {\beta _{bh}} = 0.02, {\beta _b} = 0.012, {S_h}(0) = 30, $${I_h}(0) = 78, {Q_h}(0) = 70, {R_h}(0) = 30, {S_b}(0) = 60, {I_b}(0) = 10$, 此时${R_{01}} = 0.85, {R_{02}} = 0.625, ~$满足${R_0} = \max \{ {R_{01}}, {R_{02}}\} < 1$, 则无病平衡点是全局渐近稳定的(见图 2).

图 2 系统无病平衡点的稳定性

情形2  取参数${\Lambda _h} = 10, ~{\Lambda _b} = 100, ~{\mu _h} = 0.1, ~{\mu _b} = 0.8, ~{r_h} = 0.3, ~{\alpha _h} = 0.4, ~{\delta _h} = 0.45, ~{d_h} = 0.3, ~{\varepsilon _h} = 0.3, ~{\alpha _b} = 0.7, $$~{\beta _h} = 0.03, ~{\beta _{bh}} = 0.03, ~{\beta _b} = 0.002, ~{S_h}(0) = 30, ~{I_h}(0) = 78, ~{Q_h}(0) = 70, ~{R_h}(0) = 30, ~{S_b}(0) = 90, ~{I_b}(0) = 20$, 此时${R_{01}} = 2.4, {R_{02}} = 1.67$, 满足$R_{02}>1, $则共存平衡点是全局渐近稳定的(见图 3).

图 3 系统共存平衡点的稳定性

情形3  取参数${\Lambda _h} = 10, {\Lambda _b} = 50, {\mu _h} = 0.1, {\mu _b} = 0.8, {r_h} = 0.3, {\alpha _h} = 0.4, $${\delta _h} = 0.45, {d_h} = 0.3, {\varepsilon _h} = 0.3, {\alpha _b} = 0.4, {\beta _h} = 0.03, {\beta _{bh}} = 0.03, $${\beta _b} = 0.012, {S_h}(0) = 30, {I_h}(0) = 78, {Q_h}(0) = 70, {R_h}(0) = 30, {S_b}(0) = 60, {I_b}(0) = 10, $此时${R_{01}} = 2.4 > 1, {R_{02}} = 0.625 <1, $则边界平衡点是全局渐近稳定的(见图 4).

图 4 系统边界平衡点的稳定性

在第5节中, 我们利用定义2.1分析了各个参数的敏感性指数, 下面我们通过Matlab软件作图模拟基本再生数${R_0}$对于模型中部分参数的变化情况(如图 5, 6所示).

图 5 基本再生数${R_0}$与模型中参数${\beta _h}$, ${\mu _b}$的关系. (a) ${R_0}$与模型中参数${\beta _h}$的关系; (b) ${R_0}$与模型中参数${\mu _b}$的关系

图 6 参数${\beta _h}$, ${\mu _b}$${R_0}$的综合影响

图 5可以看出, 当${R_{01}} > {R_{02}}$时, 基本再生数${R_0}$${\beta _h}$的增加而增大, 当${\beta _h} <0.0075$时, ${R_0}$几乎趋于一个常数, 此时${R_0} = \max \{{R_{01}}, {R_{02}}\} = {R_{02}}$.当${R_{01}} < {R_{02}}$时, ${R_0}$${\mu _b}$的增加而减小, 当$0<{\mu _b}<0.2$时, ${R_0}$的速率减小得很快.当${\mu _b}>0.2$时, ${R_0}$几乎趋于一个常数, 此时${R_0} = \max \{ {R_{01}}, {R_{02}}\}={R_{01}}$.而从图 6可知, 当${R_{01}} > {R_{02}}$时, ${\beta _h}$${R_0}$的影响较大.当${R_{01}} < {R_{02}}$时, ${\mu _b}$${R_0}$的影响较大.

对于埃博拉病毒目前尚无获批的疫苗, 隔离是控制疾病蔓延的有效方式, 接下来我们模拟分析隔离措施对基本再生数${R_0}$和人群中染病个体${I_h}$的影响效果(如图 7, 8所示).

图 7 ${R_0}$与模型中参数${\delta _h}$的关系

图 8 人群染病个体${I_h}$与模型中参数${\delta _h}$的关系

图 7图 8可以看出, 基本再生数${R_0}$随着隔离参数${\delta _h}$的增大而减小, 并且随着${\delta _h}$的增大, 人群中染病个体达到平衡状态所用的时间越短, 这说明采取隔离措施能有效的控制疾病的爆发.因此在埃博拉病毒的传播过程中, 我们应该及时地对染病者或疑似病例进行隔离观察治疗, 以便有效地控制疾病的蔓延.

7 结论

本文根据埃博拉病毒的传播特点, 考虑了人类和动物耦合传播的情况, 研究了一类带有隔离仓室和动物仓室的埃博拉传染病模型, 通过第二代生成矩阵的方法得到了决定疾病是否爆发的阈值条件, 也就是基本再生数${R_{01}}$${R_{02}}$.当人群患者和动物患者的基本再生数都小于1时, 仅存在无病平衡点, 此时疾病消亡.当动物患者的基本再生数大于1时, 存在共存平衡点, 此时疾病同时存在于两个群体中.当只有人群患者的基本再生数大于1时, 存在边界平衡点, 疾病只在人类群体中蔓延.并利用敏感性指数分析研究了模型中参数对基本再生数的影响, 最敏感的参数是人与人之间的传染率${\beta _h}$, 人群的输入率${\Lambda _h}$, 动物之间的传染率${\beta _b}$以及动物的输入率${\Lambda _b}$, 由此可见, 人们在平时要注意个人卫生, 避免直接接触埃博拉患者.医护人员或从事屠宰动物工作的相关人员应该做好充分的保护措施, 避免直接接触感染者的血液、分泌物、器官或其它体液等.对于受到这些液体污染的表面和材料, 应该及时进行有效的处理, 以免再次受到感染.在人流变动比较频繁或动物输入量比较大的地区, 更应该严格做好相关疾病的宣传和检查工作, 避免隐性患者的输入.对于感染者或疑似病例应该及时进行隔离治疗, 从而及时有效地控制疾病的传播.

参考文献
[1] Kermack W O, Mckendrick A G. Contributions to the mathematical theory of epidemics. B Math Biol, 1991, 53: 33–118.
[2] Busenberg S, Cooke K L. The effect of integral conditions in certain equations modeling epidemics and population growth. J Math Biol, 1980, 10: 13–32. DOI:10.1007/BF00276393
[3] 娄洁, 马知恩. 一类有被动免疫的流行病模型研究. 数学物理学报, 2003, 23: 357–368.
Lou J, Ma Z E. A class of epidemic model with passive immunity research. Acta Math Sci, 2003, 23: 357–368. DOI:10.3321/j.issn:1003-3998.2003.03.015
[4] 凌琳, 刘苏雨, 蒋贵荣. 具有饱和接触率和垂直传染的SIRS传染病模型分岔分析. 数学物理学报, 2014, 34: 1415–1425.
Ling L, Liu S Y, Jiang G R. Bifurcation analysis of SIRS epidemic model with saturated contact rate and vertical infection. Acta Math Sci, 2014, 34: 1415–1425.
[5] Anderson R M, Grenfell B T. Quantitative investigations of differentions of different vaccination polies for the control of congenital rubella syndrome(CRS) in the United Kingdom. J Hyg, 1986, 96: 305–333. DOI:10.1017/S0022172400066079
[6] Esteva L, Vargas C. Analysis of a dengue disease transmission mode. Math Biosci, 1998, 150: 131–151. DOI:10.1016/S0025-5564(98)10003-2
[7] Arino O, Parra RBDL. Alcala 1st International Conference on Mathematical Ecology. Madrid Spain:Univ Alcala, 1998
[8] Cooke K L, Yorke J A. Some equations modelling growth processes and gonorrhea epidemics. Math Biosci, 1973, 16: 75–101. DOI:10.1016/0025-5564(73)90046-1
[9] Feng Z L, Huang W Z. On the role of variable latent periods in mathematical model for tuberculosis. J Dyn Differ Equ, 2001, 13: 425–452. DOI:10.1023/A:1016688209771
[10] Fauci A S, Vargas C. Ebola-underscoring the global disparities in health care resources. New Engl J Med, 2014, 371: 1084–1086. DOI:10.1056/NEJMp1409494
[11] Gerardo C, Hiroshi N. Transmission dynamics and control of Ebola virus disease (EVD):a review. Chowell and Nishiura BMC Medicine, 2014, 12: 1–16. DOI:10.1186/1741-7015-12-1
[12] Joseph A L, Martial L N M, Jorge A A. Dynamics and control of Ebola virus transmission in Montserrado, Liberia:a mathematical modeling analysis. Lancet Infect Dis, 2014, 14: 1189–1195. DOI:10.1016/S1473-3099(14)70995-8
[13] A mathematical look at the Ebola Virus. http://www.docin.com/p-947410939.html
[14] Dixon M G, Schafer I J. Ebola viral disease outbreak-West Africa, 2014. Morb Mortal Wkly Rep, 2014, 63: 548–551.
[15] Ebola virus disease outbreak[2016-11-29]. http://www.who.int/csr/disease/ebola/en/
[16] 埃博拉病毒病实况报道[2016-11-29]. http://www.who.int/mediacentre/factsheets/fs103/zh/
The Ebola virus disease commentary[2016-11-29]. http://www.who.int/mediacentre/factsheets/fs103/zh/
[17] 张雅为. 埃博拉出血热-人兽共患传染病. 现代畜牧兽医, 2014, 10: 37–40.
Zhang Y W. Ebola haemorrhagic fever-zoonotic infectious disease. Modern Journal of Animal Husbandry and Veterinary Medicine, 2014, 10: 37–40. DOI:10.3969/j.issn.1672-9692.2014.05.012
[18] 郑学星, 杨松涛, 王化磊, 夏咸柱. 埃博拉疫苗研究新进展. 中国病原生物学杂志, 2013, 7: 652–655.
Zheng X X, Yang S T, Wang H L, Xia X Z. Ebola vaccine research progress. Journal of Pathogen Biology, 2013, 7: 652–655.
[19] 马知恩, 周义仓, 王稳地, 靳祯. 传染病动力学的数学建模与研究(第1版). 北京:科学出版社, 2004: 58–59.
Ma Z E, Zhou Y C, Wang W D, Jin Z. The Mathematical Modeling and Research of Infectious Diseases (Version 1). Beijing:Science Press, 2004: 58–59.
[20] 马知恩, 周义仓. 常微分方程定性与稳定性方法(第1版). 北京: 科学出版社, 2001: 158
Ma Z E, Zhou Y C. Qualitative and Stability of Ordinary Differential Qquation (Version 1). Beijing:Science Press, 2001:158
[21] Nakul C, James M H, Jim M C. Determining important parameters in the spread of malaria through the sensitivity analysis of a mathematical model. Bulletin Math Biol, 2008, 70: 1272–1296. DOI:10.1007/s11538-008-9299-0
[22] Van den Driessche P, Watmough J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Math Biosci, 2002, 180: 29–48. DOI:10.1016/S0025-5564(02)00108-6