具有饱和发生率和饱和治疗函数的离散SIS传染病模型的分岔分析

Bifurcation Analysis of Discrete SIS Model with Saturated Incidence and Saturated Treatment Function

摘要: 本文探讨了一类具有饱和发生率和饱和治疗函数的离散SIS模型，得到模型无病平衡点的类型及双曲性利用中心流形定理与分岔理论得到无病平衡点产生跨临界分岔和flip分岔的条件。最后给出分岔的生物学解释。

Abstract: In this paper, we discuss the discrete SIS models with saturation incidence and saturation treat-ment function. We obtain the type and hyperbolicity of the disease-free equilibrium point. We also use the central manifold theorem and the bifurcation theory to obtain the conditions of Tran-scritical bifurcation and flip bifurcation for disease-free equilibrium point. Finally, the biological explanation of bifurcation is given.

1. 引言

近几十年来，为了预防和控制传染病的传播，许多传染病模型被数学家们提出 [1] [2] [3] 。传染病模型主要分为连续模型和离散模型两大类，由于公共卫生部门收集的传染病数据都是离散数据，相比于连续模型，离散模型更容易理解和应用，且离散模型具有更加丰富的动力学行为，离散模型也可以通过连续模型离散化得到，因此，近年来离散传染病模型的研究越来越多 [4] [5] 。关于离散传染病模型的研究主要集中于模型平衡点的分岔，稳定性和持久性。研究传染病模型的分岔性质可以得到传染病模型的解受外界干扰时的性态变化和其对应模型拓扑结构的变化的影响，从而得到传染病模型相关的动力学性质，这对于传染病的防控措施与流行规律的理论研究是非常重要的。

在传染病模型中，许多因素影响着模型的动力学行为，如治疗函数、发生率等。Zhang和Liu于文献 [6] 研究了具有饱和治疗函数与饱和发生率的SIR传染病模型的倒向分岔，Wang和Ruan于文献 [7] 中研究了具有常数治疗函数与双线性发生率的SIR传染病模型的鞍结点分岔、Hopf分岔和同宿分岔。在已有研究的基础上，综合考虑当医疗条件有限且感染者数量较大时延迟治疗所产生的影响，以及在特定时间里，每个染病者接触易感者的数量有限，Zhou和Zhang在文献 [8] 中研究了如下具有饱和发生率和饱和治疗函数的连续SIS模型：

$\{\begin{array}{l}\frac{\text{d}S}{\text{d}t}=A-dS-\frac{\lambda SI}{1+kI}+vI+\frac{rI}{1+\alpha I},\\ \frac{\text{d}I}{\text{d}t}=\frac{\lambda SI}{1+kI}-\left(d+v+u\right)I-\frac{rI}{1+\alpha I},\end{array}$ (1)

这里的S和I分别代表了t时刻的易感者和染病者数量；A代表人口的补充率；d是自然死亡率；u是疾

病致死率；v指的是感染性个体的自然恢复率；饱和发生率为 $\frac{\lambda SI}{1+kI}$ ，表示单位时间内易感者接触感染者的数量；饱和治疗函数被定义为 $\frac{rI}{1+\alpha I}$ ，r指的是治愈率； $\alpha$ 描述的是延迟治疗所产生的影响。在文献 [8] 中，作者得到模型(1)的基本再生数 $R_0=\frac{\lambda A}{d\left(d+u+v+r\right)}$ ，如 $R_0>1$ ，无病平衡点 $$ 是不稳定的，

如果 $$ ，无病平衡点 $E_0$ 是稳定的。将模型(1)离散化并进行泰勒展开，我们得到如下的离散SIS传染病模型

$\{\begin{array}{l}{S}_{n+1}=\left(1-d\right)S_n-\lambda I_n\left(S_n+\frac{A}{d}\right)+k\lambda I_n^2\left(S_n+\frac{A}{d}\right)+v{I}_n+r{I}_n-\alpha r{I}_n^2,\\ I_{n+1}=\lambda I_n\left(S_n+\frac{A}{d}\right)-k\lambda I_n^2\left(S_n+\frac{A}{d}\right)-\left(d+u+v-1\right)I_n-r{I}_n+\alpha r{I}_n^2,\end{array}$ (2)

$S_n$ 和 $I_n$ 指的是n时刻易感者和染病者的数量，假定 $S\left(0\right)>0,I\left(0\right)>0$ ，系数参数 $A,d,u,v,\lambda ,r$ 都是正的， $k,\alpha$ 非负。

本文的安排如下：在第二节讨论无病平衡点 $E_0$ 的双曲性，在第三节，应用中心流形定理研究并得到非双曲无病平衡点产生跨临界分岔和flip分岔的条件，在第四节进行了总结。

2. 平衡点分析

模型(2)的无病平衡点为 $$ ，其雅克比矩阵为

$DF\left(E_0\right)=\left[\begin{array}{cc}1-d& -\frac{\lambda A}{d}+v+r\\ 0& \text{}\frac{\lambda A}{d}-\left(d+u+v+r-1\right)\end{array}\right],$

特征值为

${\lambda }_1=\left(1-d\right),{\lambda }_2=\frac{\lambda A}{d}-\left(d+u+v+r-1\right).$

定义如下符号

$B_1=\frac{d\left(d+u+v+r\right)}{\lambda },B_2=\frac{d\left(d+u+v+r-2\right)}{\lambda }.$

定理2.1 无病平衡点 $E_0$ 有以下性质：

i) 当或 $A=B_2$ 时， $E_0$ 非双曲；

ii) 当 $B_1<A<B_2$ 时， $E_0$ 全局渐近稳定；当 $A>B_2$ 或 $A<B_1$ 时， $E_0$ 是鞍点。

3. 跨临界分岔和Flip分岔

本节考虑无病平衡点 $E_0$ 的跨临界分岔和flip分岔，选取A为分岔参数。

从第2节可知，当 $A=B_1$ 时，有，当 $A=B_2$ 有 ${\lambda }_2=-1,{\lambda }_1=1-d$ ，定义如下符号

$M_1=d^2,M_2=-d\left(d+u\right),M_3=d\left(d-2\right),M_4=d\left(2-d-u\right),$

$\omega =\frac{d\left(v+r\right)-\lambda A}{\lambda A-d\left(u+v+r\right)},{\omega }_1=\frac{M_2-\lambda A_1}{\lambda A_1-M_1},{\omega }_2=\frac{M_4-\lambda A_2}{\lambda A_2-M_3}.$

定理3.1 当 $A=B_1$ 时，系统(2)在无病平衡点 $E_0$ 产生跨临界分岔。具体地讲，对A的微小扰动使得 $A<B_1$ 时，无病平衡点 $E_0$ 是双曲且稳定的；当 $A>B_1$ 时无病平衡点 $E_0$ 是双曲且不稳定的。

证明：把A看作分岔参数且把F写成 $F_A$ 来体现对A的依赖性，易得

$D{F}_A\left(E_0\right)=\left[\begin{array}{cc}1-d& -\frac{\lambda A}{d}+v+r\\ 0& \text{}\frac{\lambda A}{d}-\left(d+u+v+r-1\right)\end{array}\right].$

特征值 ${\lambda }_1,{\lambda }_2$ 对应的特征向量为

${\left(1,0\right)}^{\tau },{\left(\omega ,1\right)}^{\tau },$ (3)

利用特征向量组(3)得到如下变换

$\left[\begin{array}{c}S\\ I\end{array}\right]=\left[\begin{array}{cc}1& \omega \\ 0& 1\end{array}\right]\left[\begin{array}{c}X\\ Y\end{array}\right],$ (4)

由变换(4)可以把系统(1)变为

$\left[\begin{array}{c}X\\ Y\\ A\end{array}\right]\mapsto \left[\begin{array}{ccc}1-d& 0& 0\\ 0& 1& 0\\ 0& 0& 1\end{array}\right]\left[\begin{array}{c}X\\ Y\\ A\end{array}\right]+\left[\begin{array}{c}-\left(1+\omega \right)\left[\lambda Y\left(Xu+\omega Y\right)-k\lambda Y^2\left(X+\omega Y+\frac{A}{d}\right)+\alpha r{Y}^2\right]\\ \lambda Y\left(X+\omega Y\right)-k\lambda Y^2\left(X+\omega Y+\frac{A}{d}\right)+\alpha r{Y}^2\\ 0\end{array}\right].$ (5)

当时，有 ${\lambda }_2=1,{\lambda }_1=1-d$ ，由文献 [9] 的定理2.1.4可知，无病平衡点 $\left(X,Y\right)=\left(0,0\right)$ 的稳定性与分岔可通过研究在中心流形上一参数族映射来确定，令 $\xi =B_1,A_1=A-\xi$ ，在 $A_1$ 附近的中心流形有如下形式：

$W_{loc}\left(0,0\right)=\left\{\left(X,Y,A_1\right)\in R^3|X=h\left(Y,A_1\right),h\left(0,0\right)=0,Dh\left(0,0\right)=0,\left|Y\right|<\delta ,\left|A_1\right|<\delta \right\}.$

假设中心流形如下

其中， $\mathcal{O}\left(3\right)$ 表示次数大于或等于3的项，由中心流形的不变性可得

$\begin{array}{l}\mathcal{N}\left(h\left(Y,A_1\right)\right)=h\left(Y,A_1\right)-\left(1-d\right)h\left(Y,A_1\right)\\ +\left\{\left(1+\omega \right)\lambda Y\left[h\left(Y,A_1\right)+\omega Y\right]-k\lambda Y^2\left[h\left(X,A_1\right)+\omega Y+\frac{A}{d}\right]+\alpha r{Y}^2\right\}\\ =0.\end{array}$

对比 $Y^2,Y{A}_1,A_1^2$ 的系数解得

$a=\frac{\left(M_1+M_2\right)\left[-\lambda d{M}_2+k\lambda \xi M_1-\alpha rd{M}_1\right]}{{\left(d{M}_1\right)}^2},b=c=0,$

因此系统(1)在 $\left(X,Y\right)=\left(0,0\right)$ 附近的中心流形为

$X=h\left(Y,A_1\right)=a{Y}^2+\mathcal{O}\left(3\right).$

将其代入(5)中的第二个等式可得一维映射

$Y\mapsto {\Phi }_{A_1}\left(Y\right)=Y+\left[\lambda Y\left(a{Y}^2+{\omega }_{1}Y\right)-k\lambda Y^2\left(a{Y}^2+{\omega }_{1}Y+\frac{A_1+\xi }{d}\right)+\alpha r{Y}^2\right]+\mathcal{O}\left(5\right).$

当 $A=B_1$ 时，验算如下式子成立：

$\frac{\partial {\Phi }_{A_1}}{\partial A_1}\left(0,0\right)=0,\frac{\partial {\Phi }_{A_1}}{\partial Y\partial A_1}\left(0,0\right)\ne 0,\frac{{\partial }^2\Phi }{\partial Y^2}\left(0,0\right)\ne 0.$ (6)

根据文献 [9] 可知系统(2)在 $E_0$ 处产生跨临界分岔，在 $\left(X,Y\right)=\left(0,0\right)$ 附近系统(2)的拓扑结构与式子(7)的拓扑结构是局部拓扑等价的。

$Y\mapsto {\Phi }_{A_1}\left(Y\right)=Y+\left[\lambda Y\left(a{Y}^2+{\omega }_{1}Y\right)-k\lambda Y^2\left(a{Y}^2+{\omega }_{1}Y+\frac{A_1+\xi }{d}\right)+\alpha r{Y}^2\right]$ (7)

定理3.2 当 $A=B_2$ 时，系统(2)在无病平衡点 $E_0$ 处经历flip分岔。具体地讲，当 $$ 时，系统在无病平衡点 $E_0$ 处分岔出稳定的2周期轨；当 ${\alpha }_2<0$ 时，系统(2)在无病平衡点 $E_0$ 处分岔出的2周期轨是不稳定的。

证明：令 $\zeta =B_2,A_2=A-\zeta$ ，在 $A_2$ 附近的中心流形为

$X=h\left(Y,A_2\right)=a^{*}{Y}^2+\mathcal{O}\left(3\right).$

这里

$a^{*}=\frac{\left(M_3+M_4\right)\left[\lambda d{M}_4+k\lambda \zeta M_3-\alpha rd{M}_3\right]}{{\left(d{M}_3\right)}^2}$

类似定理3.1的计算可得到如下一维映射：

$Y\mapsto {\Phi }_{A_2}\left(Y\right)=-Y+\left[\lambda Y\left(a^{*}{Y}^2+{\omega }_{2}Y\right)-k\lambda Y^2\left(a^{*}{Y}^2+{\omega }_{2}Y+\frac{A_2+\zeta }{d}\right)+\alpha r{Y}^2\right]+\mathcal{O}\left(5\right).$

当 $A=\frac{d\left(d+u+v+r-2\right)}{\lambda }$ 时，有 ${\lambda }_2=-1,{\lambda }_1=1-d$ ，系统(2)在无病平衡点有

这里

${\alpha }_1=-\frac{2{\lambda }^2\left(M_3+M_4\right)}{M_3^2}-\frac{2k\lambda }{d},{\alpha }_2=\frac{\lambda M_3}{M_4}\left(1-2k\right)-\frac{\lambda k\zeta }{d}+2\lambda a^{*}+\alpha r.$

由文献 [10] 中的定理3.5.1可知，当 $A=B_2,{\alpha }_1\ne 0,{\alpha }_2\ne 0$ 时，系统(2)在无病平衡点 $E_0$ 处经历flip分岔。

4. 分岔生物学解释

由定理3.1和定理3.2可知，在无病平衡点 $E_0$ 处，系统(2)的局部拓扑结构对参数A非常敏感，A穿过 $B_1$ 时，系统(2)在的双曲性和稳定性都发生了改变；当A穿过 $B_2$ 时，系统(2)在 $E_0$ 经历倍周期分岔。这表明传染病的流形状态也发生相应改变。通过研究传染病模型的分岔性质，人们可以掌握传染病传播的规律，从而采取有效的措施防止传染病的传播。比如控制外来人口的数量、降低接触率、研发高效疫苗等，调整模型的系数参数，最终使模型的无病平衡点趋于稳定的状态。

基金项目

国家自然科学基金(No: 11561019)；广东省攀登计划专项资金项目(No: pdjha0304/No: pdjh2019b0297)；全国大学生创新创业训练计划项目(201910579725)。