具有Logistic增长和年龄结构的SEIQR传染病模型

SEIQR infectious disease models with Logistic growth and age structure

摘要

本文建立了具有\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}增长和年龄结构的\mathrm{S}\mathrm{E}\mathrm{I}\mathrm{Q}\mathrm{R}传染病模型.首先，定义了基本再生数R_{\mathrm{0}}，并证明R_{\mathrm{0}}是决定疾病灭绝或存在的阈值.当时，E_{\mathrm{0}}是全局渐近稳定的；当R_{\mathrm{0}}>\mathrm{1}时，E_{\mathrm{0}}是不稳定的，疾病将持续存在.其次，通过考察相应特征方程的根分布，研究了地方病平衡点E_{\mathrm{*}}的稳定性，证明了当R_{\mathrm{0}}>\mathrm{1},\tau =\mathrm{0}时，系统的地方病平衡点E_{\mathrm{*}}是局部渐近稳定的.接着，当地方病平衡点E_{\mathrm{*}}不稳定时，证明了在其周围存在\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支.最后，利用\mathrm{M}\mathrm{a}\mathrm{t}\mathrm{l}\mathrm{a}\mathrm{b}对传染病模型进行数值模拟，验证所获结论的有效性.

Abstract

In this paper，an infectious disease model with \mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c} growth and age structure was established. First， the basic regeneration number R_{\mathrm{0}} was defined and proved to be the threshold that determines the extinction or survival of the disease. When ，E_{\mathrm{0}} was globally asymptotically stable. When R_{\mathrm{0}}>\mathrm{1}，E_{\mathrm{0}} was unstable， the disease will persist. Second， the stability of the endemic equilibrium point E_{\mathrm{*}} was studied by examining the root distribution of the corresponding characteristic equation， it was proved that if R_{\mathrm{0}}>\mathrm{1},\tau =\mathrm{0}， the endemic equilibrium point of the system was locally asymptotically stable. Then， when the endemic equilibrium point was unstable， it was proved that there were \mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f} branches around it. Finally， the \mathrm{M}\mathrm{a}\mathrm{t}\mathrm{l}\mathrm{a}\mathrm{b} software was used to simulate the infectious disease model and the validity of the conclusions was verified.

Graphical abstract

关键词

\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}

增长 / 年龄结构 / 基本再生数 /

\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}

分支

Key words

\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}

growth / age structure / basic regeneration number /

\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}

branches

数学模型是了解传染病传播和控制的基本工具.自\mathrm{K}\mathrm{e}\mathrm{r}\mathrm{m}\mathrm{a}\mathrm{c}\mathrm{k}等^［1］提出仓室模型以来，易感染-感染-恢复（susceptible-infectious-recovered，SIR）模型已成为传染病在人群中传播的基础数学理论，并应用于特定疾病，如麻疹、疟疾、霍乱、季节性流感、\mathrm{C}\mathrm{O}\mathrm{V}\mathrm{I}\mathrm{D}-\mathrm{19}等^［2］.然而，大多数易感者的输入率通常被假定为常数，对于死亡率高或持续时间相对较长的疾病，在某些实际情况下，假设易感者的输入为\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}增长输入可能更合理^［3］.\mathrm{L}\mathrm{i}等^［4］考虑了具有\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}增长率和饱和治疗以及双线性发病率的\mathrm{S}\mathrm{I}\mathrm{R}传染病模型，并研究了模型的稳定性以及分支.\mathrm{W}\mathrm{a}\mathrm{n}\mathrm{g}等^［5］建立了一种具有病毒突变和易感人群\mathrm{l}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}\mathrm{s}增长的随机\mathrm{S}\mathrm{E}\mathrm{I}\mathrm{R}（susceptible-exposed-infectious-recovered）传染病模型，研究了模型的稳定性并获得了传染病生存和灭绝的充分条件.

在研究传染病模型的动力学行为时，时滞是一个不可忽视的因素^［6-8］.例如，感染手足口病的个体在接触病毒后的一段时间 （一般为\mathrm{3}~\mathrm{7} d） 后才会具有传染性.另外，年龄也是传染病传播、控制、预防和建模的重要因素之一，因为不同年龄的人经历相同疾病可能有不同死亡率和感染率，即使同一个个体，在病程的不同时期传染能力也不尽相同，所以许多学者在模型中加入了年龄结构来研究某些传染病模型传播的动力学.早期建立和研究年龄结构流行病模型的是\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{p}\mathrm{e}\mathrm{n}\mathrm{s}\mathrm{t}\mathrm{e}\mathrm{n}\mathrm{d}^［9］，此后，研究年龄结构流行病模型的工作者不断出现，大多数作者讨论的是SIS（susceptible-infectious-susceptible）、SIR、SEIR模型^［10-12］.此外，隔离被认为是治疗患者的有效方法，也是清除过程之一.通过自我检查、自我隔离等手段可以有效地提高恢复率以及降低死亡率.许多学者在模型中考虑了隔离，从而让模型更加符合实际情况^［13-15］.

基于以上分析，本文建立了一个具有易感群体\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}增长和年龄结构的\mathrm{S}\mathrm{E}\mathrm{I}\mathrm{Q}\mathrm{R}（susceptible-exposed-infectious-quarantined-recovered）传染病模型，将人群分为易感者类、潜伏者类、感染者类、隔离者类和康复者类，并且个体从潜伏者到感染者的状态转移率以及潜伏者的因病死亡率都与年龄有关.首先，本研究定义了模型的基本再生数；其次，计算出无病平衡点和地方病平衡点，并且研究了无病平衡点和地方病平衡点的稳定性；最后，通过数值模拟来验证理论的有效性.

1 建立模型

总人口在流行病学上被划分为易感者、潜伏者、感染者、隔离者和康复者.S\left(t\right)表示时间t时的易感者的个体数量，E\left(t,a\right)表示时间为t、感染年龄为a的潜伏者的个体数量，I\left(t\right)表示时间t时的感染者的个体数量，Q\left(t\right)表示时间t时的隔离者的个体数量，R\left(t\right)表示时间t时的康复者的个体数量.因此，建立具有\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}\mathrm{s}增长和年龄结构的\mathrm{S}\mathrm{E}\mathrm{I}\mathrm{Q}\mathrm{R}传染病模型，模型的流程图如图\mathrm{1}所示.

建立模型如下所示：

\left\{\begin{array}{l}\frac{\mathrm{d}S\left(t\right)}{\mathrm{d}t}=bS\left(t\right)\left(\mathrm{1}-\frac{S\left(t\right)}{K}\right)-\beta S\left(t\right)I\left(t\right),\\ \frac{\partial E\left(t,a\right)}{\partial t}+\frac{\partial E\left(t,a\right)}{\partial a}=-\left(d_{\mathrm{1}}+\sigma \left(a\right)+\alpha \left(a\right)\right)E\left(t,a\right),\\ \frac{\mathrm{d}I\left(t\right)}{\mathrm{d}t}={\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)E\left(t,a\right)\mathrm{d}a-\left(d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)I\left(t\right),\\ \frac{\mathrm{d}Q\left(t\right)}{\mathrm{d}t}=mI\left(t\right)-\left(d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}\right)Q\left(t\right),\\ \frac{\mathrm{d}R\left(t\right)}{\mathrm{d}t}={\gamma }_{\mathrm{1}}I\left(t\right)+{\gamma }_{\mathrm{2}}Q\left(t\right)-dR\left(t\right),\\ E\left(t,\mathrm{0}\right)=\beta S\left(t\right)I\left(t\right),\end{array}\right.（1）

初始条件：

\begin{array}{c}S\left(\mathrm{0}\right)=S_{\mathrm{0}}>\mathrm{0},E\left(\mathrm{0},a\right)=E_{\mathrm{0}}\left(a\right)>\mathrm{0},\\ I\left(\mathrm{0}\right)=I_{\mathrm{0}}>\mathrm{0},Q\left(\mathrm{0}\right)=Q_{\mathrm{0}}\ge \mathrm{0},R\left(\mathrm{0}\right)=R_{\mathrm{0}}\ge \mathrm{0}\end{array}.

2 基本再生数

首先讨论系统\left(\mathrm{1}\right)的基本再生数R_{\mathrm{0}}.显然，系统\left(\mathrm{1}\right)存在无病平衡点E_{\mathrm{0}}=\left(S^{\mathrm{0}},\mathrm{0,0},\mathrm{0,0}\right)，其中，S^{\mathrm{0}}=K.R相对于系统\left(\mathrm{1}\right)是独立的，因此只研究系统\left(\mathrm{1}\right)的前四个方程，令

\left\{\begin{array}{l}b{S}_{\mathrm{*}}\left(\mathrm{1}-\frac{S_{\mathrm{*}}}{K}\right)-\beta S_{\mathrm{*}}{I}_{\mathrm{*}}=\mathrm{0},\\ \frac{\mathrm{d}{E}_{\mathrm{*}}\left(a\right)}{\mathrm{d}a}=-\chi \left(a\right)E_{\mathrm{*}}\left(a\right),\\ {\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)E_{\mathrm{*}}\left(a\right)\mathrm{d}a-\left(d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)I_{\mathrm{*}}=\mathrm{0},\\ m{I}_{\mathrm{*}}-\left(d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}\right)Q_{\mathrm{*}}=\mathrm{0},\end{array}\right.（2）

其中，\chi \left(a\right)=d_{\mathrm{1}}+\sigma \left(a\right)+\alpha \left(a\right),E_{\mathrm{*}}\left(\mathrm{0}\right)=\beta S_{\mathrm{*}}{I}_{\mathrm{*}}.

求解上述方程组可得

\begin{array}{l}{E}_{\mathrm{*}}\left(a\right)=E_{\mathrm{*}}\left(\mathrm{0}\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\chi \left(\theta \right)\mathrm{d}\theta },\\ S_{\mathrm{*}}=\frac{d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}{\beta {\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\chi \left(\theta \right)\mathrm{d}\theta }\mathrm{d}a},\\ I_{\mathrm{*}}=\frac{E_{\mathrm{*}}\left(\mathrm{0}\right){\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\chi \left(\theta \right)\mathrm{d}\theta }\mathrm{d}a}{d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m},\\ Q_{\mathrm{*}}=\frac{m}{d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}}{I}_{\mathrm{*}}\mathrm{。}\end{array}

由系统\left(\mathrm{2}\right)的第一个方程可得E_{\mathrm{*}}\left(\mathrm{0}\right)=b{S}_{\mathrm{*}}\left(\mathrm{1}-\frac{S_{\mathrm{*}}}{K}\right),b\left(\mathrm{1}-\frac{S_{\mathrm{*}}}{K}\right)=\beta I_{\mathrm{*}}，则有I_{\mathrm{*}}=\frac{b}{\beta }\left(\mathrm{1}-\frac{S_{\mathrm{*}}}{K}\right)，即

I_{\mathrm{*}}=\frac{b}{\beta }\left(\mathrm{1}-\frac{d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}{K\beta {\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\chi \left(\theta \right)\mathrm{d}\theta }\mathrm{d}a}\right)，

因此可以得到地方病平衡点E_{\mathrm{*}}=\left(S_{\mathrm{*}},E_{\mathrm{*}}\left(a\right),I_{\mathrm{*}},Q_{\mathrm{*}}\right).

定义基本再生数

R_{\mathrm{0}}=\frac{K}{S_{\mathrm{*}}}=\frac{K\beta {\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\chi \left(\theta \right)\mathrm{d}\theta }\mathrm{d}a}{d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}.（3）

3 无病平衡点的稳定性

在上述讨论的基础上，研究无病平衡点E_{\mathrm{0}}的稳定性.

定理1 当时，E_{\mathrm{0}}是局部渐近稳定的；当R_{\mathrm{0}}>\mathrm{1}时，E_{\mathrm{0}}是不稳定的.

证明 由于无病平衡点E_{\mathrm{0}}=\left(S^{\mathrm{0}},\mathrm{0,0},\mathrm{0}\right)，令

\begin{array}{c}s\left(t\right)=S\left(t\right)-S^{\mathrm{0}},e\left(t,a\right)=E\left(t,a\right),\\ i\left(t\right)=I\left(t\right),q\left(t\right)=Q\left(t\right)\end{array}.

则系统\left(\mathrm{1}\right)在E_{\mathrm{0}}处的线性化系统为

\left\{\begin{array}{l}\frac{\mathrm{d}s\left(t\right)}{\mathrm{d}t}=-bs\left(t\right)-\beta S^{\mathrm{0}}i\left(t\right),\\ \frac{\partial e\left(t,a\right)}{\partial t}+\frac{\partial e\left(t,a\right)}{\partial a}=-\chi \left(a\right)e\left(t,a\right),\\ \frac{\mathrm{d}i\left(t\right)}{\mathrm{d}t}={\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)e\left(t,a\right)\mathrm{d}a-\left(d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)i\left(t\right),\\ \frac{\mathrm{d}q\left(t\right)}{\mathrm{d}t}=mi\left(t\right)-\left(d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}\right)q\left(t\right),\end{array}\right.（4）

其中，e\left(t,\mathrm{0}\right)=\beta S^{\mathrm{0}}i\left(t\right).

将s\left(t\right)=s_{\mathrm{0}}{\mathrm{e}}^{\lambda t},e\left(t,a\right)=e_{\mathrm{0}}\left(a\right){\mathrm{e}}^{\lambda t},i\left(t\right)=i_{\mathrm{0}}{\mathrm{e}}^{\lambda t},q\left(t\right)=q_{\mathrm{0}}{\mathrm{e}}^{\lambda t}代入上述方程组可得

\left\{\begin{array}{l}\left(\lambda +b\right)s_{\mathrm{0}}=-\beta S^{\mathrm{0}}{i}_{\mathrm{0}},\\ \frac{\mathrm{d}{e}_{\mathrm{0}}\left(a\right)}{\mathrm{d}a}=-\left(\lambda +\chi \left(a\right)\right)e_{\mathrm{0}}\left(a\right),\\ \left(\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)i_{\mathrm{0}}={\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)e_{\mathrm{0}}\left(a\right)\mathrm{d}a,\\ \left(\lambda +d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}+m\right)q_{\mathrm{0}}=m{i}_{\mathrm{0}},\end{array}\right.（5）

其中，e_{\mathrm{0}}\left(\mathrm{0}\right)=\beta S^{\mathrm{0}}{i}_{\mathrm{0}}.

求解上述方程组可得

\begin{array}{c}{e}_{\mathrm{0}}\left(a\right)=e_{\mathrm{0}}\left(\mathrm{0}\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta },\\ i_{\mathrm{0}}=\frac{e_{\mathrm{0}}\left(\mathrm{0}\right){\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta }\mathrm{d}a}{\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m},\\ q_{\mathrm{0}}=\frac{m}{\lambda +d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}}{i}_{\mathrm{0}}\end{array}.

将i_{\mathrm{0}}代入e_{\mathrm{0}}\left(\mathrm{0}\right)=\beta S^{\mathrm{0}}{i}_{\mathrm{0}}，得

e_{\mathrm{0}}\left(\mathrm{0}\right)=\beta S^{\mathrm{0}}\frac{e_{\mathrm{0}}\left(\mathrm{0}\right)}{\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}{\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta }\mathrm{d}a，

则有

F\left(\lambda \right)=\beta S^{\mathrm{0}}\frac{\mathrm{1}}{\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}{\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta }\mathrm{d}a-\mathrm{1}=\mathrm{0}，（6）

即\left(\mathrm{6}\right)是系统\left(\mathrm{4}\right)所对应的特征方程.

若\lambda \in \mathbb{R}，则F\left(\lambda \right)是单调递减的实值连续函数，所以F\left(\mathrm{0}\right)=R_{\mathrm{0}}-\mathrm{1}.因此，当R_{\mathrm{0}}>\mathrm{1}时，特征方程至少有一个正实根，即无病平衡点E_{\mathrm{0}}是不稳定的.否则，当时，F\left(\lambda \right)=\mathrm{0}没有非负实根.此时，讨论\lambda是复根的情形.不妨设\lambda =k+\omega \mathrm{i},k\ge \mathrm{0}是F\left(\lambda \right)=\mathrm{0}的任意复根，则

这里可以看出矛盾，特征方程F\left(\lambda \right)=\mathrm{0}的根只有负实部，即，E_{\mathrm{0}}是局部渐近稳定的.

定理2 若，则无病平衡点E_{\mathrm{0}}是全局渐近稳定的.

证明 根据定理1，只需证明当时，无病平衡点E_{\mathrm{0}}是全局吸引的，即证明

\underset{t\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\left(S\left(t\right),E\left(t,a\right),I\left(t\right),Q\left(t\right)\right)=\left(S^{\mathrm{0}},\mathrm{0,0},\mathrm{0}\right).

根据系统\left(\mathrm{4}\right)的第一个方程

\underset{t\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\mathrm{s}\mathrm{u}\mathrm{p}\frac{\mathrm{d}\stackrel{⌣}{S}\left(t\right)}{\mathrm{d}t}\le b\stackrel{⌣}{S}\left(t\right)\left(\mathrm{1}-\frac{\stackrel{⌣}{S}\left(t\right)}{K}\right)，

求解\frac{\mathrm{d}\stackrel{⌣}{S}\left(t\right)}{\mathrm{d}t}=b\stackrel{⌣}{S}\left(t\right)\left(\mathrm{1}-\frac{\stackrel{⌣}{S}\left(t\right)}{K}\right)，可以获得\stackrel{⌣}{S}\left(t\right)=\left(\mathrm{1}+\frac{\mathrm{1}}{{\mathrm{e}}^{bt+Kc}}\right)K，其中，c是一个常数，即\underset{t\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\stackrel{⌣}{S}\left(t\right)=K，所以\underset{t\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}S\left(t\right)=K.因此，对于\forall \epsilon >\mathrm{0},\exists t_{\mathrm{0}}>\mathrm{0},使得S\left(t\right)\le S^{\mathrm{0}}+\epsilon对于所有t>t_{\mathrm{0}}成立，则有E\left(t,\mathrm{0}\right)\le \beta \left(S^{\mathrm{0}}+\epsilon \right)I\left(t\right).

接下来考虑如下线性系统，

\left\{\begin{array}{l}\frac{\partial \stackrel{⌣}{E}\left(t,a\right)}{\partial t}+\frac{\partial \stackrel{⌣}{E}\left(t,a\right)}{\partial a}=-\chi \left(a\right)\stackrel{⌣}{E}\left(t,a\right),\\ \frac{\mathrm{d}\stackrel{⌣}{I}\left(t\right)}{\mathrm{d}t}={\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)\stackrel{⌣}{E}\left(t,a\right)\mathrm{d}a-\left(d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)\stackrel{⌣}{I}\left(t\right),\\ \frac{\mathrm{d}\stackrel{⌣}{Q}\left(t\right)}{\mathrm{d}t}=m\stackrel{⌣}{I}\left(t\right)-\left(d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}\right)\stackrel{⌣}{Q}\left(t\right),\end{array}\right.（7）

其中，\stackrel{⌣}{E}\left(t,\mathrm{0}\right)=\beta \left(S^{\mathrm{0}}+\epsilon \right)\stackrel{⌣}{I}\left(t\right).

与定理\mathrm{1}的证明类似，可以获得\left(\mathrm{7}\right)解的形式为

\begin{array}{c}\left(\stackrel{⌣}{E}\left(t,a\right),\stackrel{⌣}{I}\left(t\right),\stackrel{⌣}{Q}\left(t\right)\right)=\\ \left({\stackrel{⌣}{E}}_{\mathrm{0}}\left(a\right){\mathrm{e}}^{{\lambda }_{\mathrm{0}}t},{\stackrel{⌣}{I}}_{\mathrm{0}}{\mathrm{e}}^{{\lambda }_{\mathrm{0}}t},{\stackrel{⌣}{Q}}_{\mathrm{0}}{\mathrm{e}}^{{\lambda }_{\mathrm{0}}t}\right),\end{array}

其中，{\stackrel{⌣}{E}}_{\mathrm{0}}\left(a\right)\mathrm{、}{\stackrel{⌣}{I}}_{\mathrm{0}}\mathrm{、}{\stackrel{⌣}{Q}}_{\mathrm{0}}都是正的，且{\lambda }_{\mathrm{0}}为所对应特征方程的特征根.

利用特征线法来求解\stackrel{⌣}{E}\left(t,a\right)，

\stackrel{⌣}{E}\left(t,a\right)=\left\{\begin{array}{l}\beta S^{\mathrm{0}}\left(t-a\right)\stackrel{⌣}{I}\left(t-a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\chi \left(\theta \right)\mathrm{d}\theta },a\le t,\\ {\stackrel{⌣}{E}}_{\mathrm{0}}\left(a-t\right){\mathrm{e}}^{-{\int }_{a-t}^t\chi \left(\theta \right)\mathrm{d}\theta },a>t,\end{array}\right.

可以看出，当t\ge t_{\mathrm{0}}时，E\left(t,a\right)\le \stackrel{⌣}{E}\left(t,a\right),I\left(t\right)\le \stackrel{⌣}{I}\left(t\right),Q\left(t\right)\le \stackrel{⌣}{Q}\left(t\right).

因为时，系统\left(\mathrm{4}\right)的所有特征值都具有负实部，所以存在\epsilon >\mathrm{0}，使得，故

\underset{t\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\left(E\left(t,a\right),I\left(t\right),Q\left(t\right)\right)=\left(\mathrm{0,0},\mathrm{0}\right)，

因此，若，则\underset{t\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\left(S\left(t\right),E\left(t,a\right)\right.,I\left(t\right),\left.Q\left(t\right)\right)=\left(S^{\mathrm{0}},\mathrm{0,0},\mathrm{0}\right)，即E_{\mathrm{0}}是全局渐近稳定的.

4 地方病平衡点的稳定性及\mathbf{H}\mathbf{o}\mathbf{p}\mathbf{f}分支

本节研究地方病平衡点E_{\mathrm{*}}=\left(S_{\mathrm{*}},E_{\mathrm{*}}\left(a\right),I_{\mathrm{*}},Q_{\mathrm{*}}\right)的稳定性和它附近\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支的存在性.

为了证明E_{\mathrm{*}}的局部渐近稳定性，将系统\left(\mathrm{1}\right)在E_{\mathrm{*}}附近进行线性化.

令x\left(t\right)=S\left(t\right)-S_{\mathrm{*}},y\left(t,a\right)=E\left(t,a\right)-E_{\mathrm{*}}\left(a\right),g\left(t\right)=I\left(t\right)-I_{\mathrm{*}},h\left(t\right)=Q\left(t\right)-Q_{\mathrm{*}}，可以得到

\left\{\begin{array}{l}\frac{\mathrm{d}x\left(t\right)}{\mathrm{d}t}=bx\left(t\right)\left(\mathrm{1}-\frac{\mathrm{2}{S}_{\mathrm{*}}}{K}\right)-\beta \left(I_{\mathrm{*}}x\left(t\right)+S_{\mathrm{*}}g\left(t\right)\right),\\ \frac{\partial y\left(t,a\right)}{\partial t}+\frac{\partial y\left(t,a\right)}{\partial a}=-\chi \left(a\right)E\left(t,a\right),\\ \frac{\mathrm{d}g\left(t\right)}{\mathrm{d}t}={\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)y\left(t,a\right)\mathrm{d}a-\left(d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)g\left(t\right),\\ \frac{\mathrm{d}h\left(t\right)}{\mathrm{d}t}=mg\left(t\right)-\left(d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}\right)h\left(t\right),\end{array}\right.（8）

其中，y\left(t,\mathrm{0}\right)=\beta \left(I_{\mathrm{*}}x\left(t\right)+S_{\mathrm{*}}g\left(t\right)\right).

为了分析E{}_{\mathrm{*}}{}^{}的渐近行为，令x\left(t\right)=x_{\mathrm{0}}{\mathrm{e}}^{\lambda t},y\left(t,a\right)=y_{\mathrm{0}}\left(a\right){\mathrm{e}}^{\lambda t},g\left(t\right)=g_{\mathrm{0}}{\mathrm{e}}^{\lambda t},h\left(t\right)=h_{\mathrm{0}}{\mathrm{e}}^{\lambda t}，则有如下特征问题：

\left\{\begin{array}{l}\left[\lambda -b\left(\mathrm{1}-\frac{S_{\mathrm{*}}}{K}\right)\right]x_{\mathrm{0}}=-\beta \left(I_{\mathrm{*}}{x}_{\mathrm{0}}+S_{\mathrm{*}}{g}_{\mathrm{0}}\right),\\ \frac{\mathrm{d}{y}_{\mathrm{0}}\left(a\right)}{\mathrm{d}a}=-\left(\lambda +\chi \left(a\right)\right)y_{\mathrm{0}}\left(a\right),\\ \left(\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)g_{\mathrm{0}}={\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right)y_{\mathrm{0}}\left(a\right)\mathrm{d}a,\\ \left(\lambda +d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}\right)h_{\mathrm{0}}=m{g}_{\mathrm{0}},\end{array}\right.（9）

其中，y_{\mathrm{0}}\left(\mathrm{0}\right)=\beta \left(I_{\mathrm{*}}{x}_{\mathrm{0}}+S_{\mathrm{*}}{g}_{\mathrm{0}}\right).

求解系统\left(\mathrm{9}\right)可以获得

\begin{array}{c}{x}_{\mathrm{0}}=-\frac{y_{\mathrm{0}}\left(\mathrm{0}\right)}{\lambda -b\left(\mathrm{1}-\frac{\mathrm{2}{S}_{\mathrm{*}}}{K}\right)},\\ y_{\mathrm{0}}\left(a\right)=y_{\mathrm{0}}\left(\mathrm{0}\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta },\\ g_{\mathrm{0}}=\frac{y_{\mathrm{0}}\left(\mathrm{0}\right)}{\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}{\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta }\mathrm{d}a,\\ h_{\mathrm{0}}=\frac{m}{\lambda +d_{\mathrm{3}}+{\mu }_{\mathrm{2}}+{\gamma }_{\mathrm{2}}}{g}_{\mathrm{0}}\mathrm{。}\end{array}

结合上面四个式子，可以获得地方病平衡点E_{\mathrm{*}}所对应的特征方程

\begin{array}{l}\lambda -b\left(\mathrm{1}-\frac{\mathrm{2}{S}_{\mathrm{*}}}{K}\right)=-\beta I_{\mathrm{*}}+\\ \beta S_{\mathrm{*}}\frac{\lambda -b\left(\mathrm{1}-\frac{\mathrm{2}{S}_{\mathrm{*}}}{K}\right)}{\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}{\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta }\mathrm{d}a.\end{array}

考虑到特征方程的复杂性和\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支，假设潜伏者的疾病死亡率\sigma \left(a\right)={\sigma }_{\mathrm{*}}，同时假设当被感染者的感染年龄a小于感染时滞\tau时，被感染者没有感染性，因此\alpha \left(a\right)=\mathrm{0}，而当被感染者的感染年龄a大于或等于感染时滞\tau时，就有可能被感染，则有\alpha \left(a\right)=c_{\mathrm{*}}\left(>\mathrm{0}\right)，因此，\alpha \left(a\right)的表达式为

其中，c_{\mathrm{*}}为一个正常数，可视为潜伏个体到感染个体的转化率.因此，

{\int }_{\mathrm{0}}^{+\mathrm{\infty }}\alpha \left(a\right){\mathrm{e}}^{-{\int }_{\mathrm{0}}^a\left(\lambda +\chi \left(\theta \right)\right)\mathrm{d}\theta }\mathrm{d}a=\frac{c_{\mathrm{*}}{\mathrm{e}}^{-\left(\lambda +d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}\right)\tau }}{\lambda +d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}+c_{\mathrm{*}}}.

通过直接计算，可以得到对应的特征方程

\begin{array}{l}\lambda -\phi =-\beta I_{\mathrm{*}}+\\ \beta S_{\mathrm{*}}\frac{\lambda -\phi }{\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m}\frac{c_{\mathrm{*}}{\mathrm{e}}^{-\left(\lambda +d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}\right)\tau }}{\lambda +d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}+c_{\mathrm{*}}},\end{array}

其中，\phi =b\left(\mathrm{1}-\frac{\mathrm{2}{S}_{\mathrm{*}}}{K}\right).

合并同类项可得

\mathrm{\Pi }\left(\lambda ,\tau \right)=\frac{f\left(\lambda ,\tau \right)}{g\left(\lambda \right)}=\frac{{\lambda }^{\mathrm{3}}+a_{\mathrm{2}}{\lambda }^{\mathrm{2}}+a_{\mathrm{1}}\lambda +a_{\mathrm{0}}+\left(b_{\mathrm{1}}\lambda +b_{\mathrm{0}}\right){\mathrm{e}}^{-\lambda \tau }}{\left(\lambda +d_{\mathrm{2}}+{\mu }_{\mathrm{1}}+{\gamma }_{\mathrm{1}}+m\right)\left(\lambda +d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}+c_{\mathrm{*}}\right)}=\mathrm{0},

其中，

\begin{array}{l}{a}_{\mathrm{2}}=c_{\mathrm{*}}+d_{\mathrm{1}}+d_{\mathrm{2}}+m+\beta I_{\mathrm{*}}+{\gamma }_{\mathrm{1}}+{\mu }_{\mathrm{1}}+{\sigma }_{\mathrm{*}}-\phi ,\\ a_{\mathrm{1}}=\left(d_{\mathrm{2}}+m+{\gamma }_{\mathrm{1}}+{\mu }_{\mathrm{1}}\right)\left(\beta I_{\mathrm{*}}+{\sigma }_{\mathrm{*}}+d_{\mathrm{1}}+c_{\mathrm{*}}-\phi \right)+\beta I_{\mathrm{*}}\left({\sigma }_{\mathrm{*}}+c_{\mathrm{*}}\right)-\phi \left({\sigma }_{\mathrm{*}}+c_{\mathrm{*}}+d_{\mathrm{1}}\right),\\ a_{\mathrm{0}}=\left({\sigma }_{\mathrm{*}}+\mathrm{2}{c}_{\mathrm{*}}\right)\left(d_{\mathrm{2}}+m+{\gamma }_{\mathrm{1}}+{\mu }_{\mathrm{1}}\right)\left(\beta I_{\mathrm{*}}-\phi \right),\\ b_{\mathrm{1}}=\beta c_{\mathrm{*}}{S}_{\mathrm{*}}{\mathrm{e}}^{-\left(d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}\right)\tau },\\ b_{\mathrm{0}}=\beta c_{\mathrm{*}}{S}_{\mathrm{*}}\phi {\mathrm{e}}^{-\left(d_{\mathrm{1}}+{\sigma }_{\mathrm{*}}\right)\tau }.\end{array}

进一步探究不同时滞情况下，E_{\mathrm{*}}的稳定性以及\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支的存在性.

定理3 若R_{\mathrm{0}}>\mathrm{1},\tau =\mathrm{0}，并且a_{\mathrm{2}}\left(a_{\mathrm{1}}+b_{\mathrm{1}}\right)-\left(a_{\mathrm{0}}+b_{\mathrm{0}}\right)>\mathrm{0}，则系统的地方病平衡点E_{\mathrm{*}}是局部渐近稳定的.

证明 若\tau =\mathrm{0}，则f\left(\lambda ,\tau \right)={\lambda }^{\mathrm{3}}+a_{\mathrm{2}}{\lambda }^{\mathrm{2}}+\left(a_{\mathrm{1}}+b_{\mathrm{1}}\right)\lambda +\left(a_{\mathrm{0}}+b_{\mathrm{0}}\right).

经过简单的计算，

a_{\mathrm{2}}=c_{\mathrm{*}}+d_{\mathrm{1}}+d_{\mathrm{2}}+m+\beta I_{\mathrm{*}}+{\gamma }_{\mathrm{1}}+{\mu }_{\mathrm{1}}+{\sigma }_{\mathrm{*}}-\phi =\frac{b}{R_{\mathrm{0}}}>\mathrm{0}.

若a_{\mathrm{2}}\left(a_{\mathrm{1}}+b_{\mathrm{1}}\right)-\left(a_{\mathrm{0}}+b_{\mathrm{0}}\right)>\mathrm{0}，由Routh-Hurwitz判据可知，特征方程的所有根都有负实部，即系统的地方病平衡点E_{\mathrm{*}}是局部渐近稳定的.

设\lambda =\omega \mathrm{i}\left(\omega >\mathrm{0}\right)是f\left(\lambda ,\tau \right)的纯虚根，则有

-{\omega }^{\mathrm{3}}\mathrm{i}-a_{\mathrm{2}}{\omega }^{\mathrm{2}}+a_{\mathrm{1}}\omega \mathrm{i}+a_{\mathrm{0}}+\left(b_{\mathrm{1}}\omega \mathrm{i}+b_{\mathrm{0}}\right)\left(\mathrm{c}\mathrm{o}\mathrm{s}\omega \tau -\mathrm{i}\mathrm{s}\mathrm{i}\mathrm{n}\omega \tau \right)=\mathrm{0}.

将实部和虚部分离，

\left\{\begin{array}{l}-{\omega }^{\mathrm{3}}+a_{\mathrm{1}}\omega =b_{\mathrm{1}}\omega \mathrm{c}\mathrm{o}\mathrm{s}\omega \tau -b_{\mathrm{0}}\mathrm{s}\mathrm{i}\mathrm{n}\omega \tau ,\\ -a_{\mathrm{2}}{\omega }^{\mathrm{2}}+a_{\mathrm{0}}=-b_{\mathrm{0}}\mathrm{c}\mathrm{o}\mathrm{s}\omega \tau -b_{\mathrm{1}}\omega \mathrm{s}\mathrm{i}\mathrm{n}\omega \tau ,\end{array}\right.

则有

\begin{array}{l}{\omega }^{\mathrm{6}}+\left(a_{\mathrm{2}}^{\mathrm{2}}-\mathrm{2}{a}_{\mathrm{1}}\right){\omega }^{\mathrm{4}}+\left(a_{\mathrm{1}}^{\mathrm{2}}-b_{\mathrm{1}}^{\mathrm{2}}-\mathrm{2}{a}_{\mathrm{0}}{a}_{\mathrm{2}}\right){\omega }^{\mathrm{2}}+\\ \left(a_{\mathrm{0}}^{\mathrm{2}}-b_{\mathrm{0}}^{\mathrm{2}}\right)=\mathrm{0}.\end{array}

令

\begin{array}{l}z={\omega }^{\mathrm{2}},p=a_{\mathrm{2}}^{\mathrm{2}}-\mathrm{2}{a}_{\mathrm{1}},\\ q=a_{\mathrm{1}}^{\mathrm{2}}-b_{\mathrm{1}}^{\mathrm{2}}-\mathrm{2}{a}_{\mathrm{0}}{a}_{\mathrm{2}},r=a_{\mathrm{0}}^{\mathrm{2}}-b_{\mathrm{0}}^{\mathrm{2}},\end{array}

则上述方程可以转化为

H\left(z\right)=z^{\mathrm{3}}+p{z}^{\mathrm{2}}+qz+r.

对H\left(z\right)关于z求导\frac{\mathrm{d}H\left(z\right)}{\mathrm{d}z}=\mathrm{3}{z}^{\mathrm{2}}+\mathrm{2}pz+q，记\mathrm{\Delta }=\mathrm{4}{p}^{\mathrm{2}}-\mathrm{12}q.

引理1^［16］ 对于方程H\left(z\right)=\mathrm{0}，

1） 若，则方程至少有一个正根；

2） 若r\ge \mathrm{0},\mathrm{\Delta }\le \mathrm{0}，则方程没有正根；

3） 若r\ge \mathrm{0},\mathrm{\Delta }>\mathrm{0}，则方程有正根当且仅当z_{\mathrm{*}}>\mathrm{0}，且H\left(z_{\mathrm{*}}\right)\le \mathrm{0}.

引理2 假设H^{\text{'}}\left(z_{\mathrm{*}}\right)\ne \mathrm{0}，若\tau ={\tau }_k,k=\mathrm{0,1},\mathrm{2},\dots，则方程f\left(\lambda ,\tau \right)=\mathrm{0}有一对纯虚根\pm {\omega }_{\mathrm{*}}\mathrm{i}，并且

\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left\{\frac{\mathrm{d}\mathrm{R}\mathrm{e}\left(\lambda \right)}{\mathrm{d}\tau }{|}_{\lambda ={\omega }_{\mathrm{*}}\mathrm{i}}\right\}\ne \mathrm{0}.

证明 设z_{\mathrm{*}}是H\left(z\right)=\mathrm{0}的正实数解，则{\omega }_{\mathrm{*}}=\sqrt[]{z_{\mathrm{*}}}是方程的唯一正实数解，所以当\tau ={\tau }_k,k=\mathrm{0,1},\mathrm{2}\dots时，f\left(\lambda ,\tau \right)=\mathrm{0}有两个纯虚根\pm {\omega }_{\mathrm{*}}\mathrm{i}，其中，

这里，k=\mathrm{0,1},\mathrm{2},\dots，且

\rho =\frac{-b_{\mathrm{1}}{\omega }_{\mathrm{*}}\left(-{{\omega }_{\mathrm{*}}}^{\mathrm{3}}+a_{\mathrm{1}}{\omega }_{\mathrm{*}}\right)+a_{\mathrm{2}}{b}_{\mathrm{0}}{{\omega }_{\mathrm{*}}}^{\mathrm{2}}-a_{\mathrm{0}}{b}_{\mathrm{0}}}{b_{\mathrm{0}}^{\mathrm{2}}-b_{\mathrm{1}}{\omega }_{\mathrm{*}}}.

将f\left(\lambda ,\tau \right)=\mathrm{0}两边关于\tau进行微分，可得

{\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau }\right)}^{-\mathrm{1}}=\frac{\left(\mathrm{3}{\lambda }^{\mathrm{2}}+\mathrm{2}{a}_{\mathrm{2}}\lambda +a_{\mathrm{1}}\right){\mathrm{e}}^{\lambda \tau }}{\lambda \left(b_{\mathrm{1}}\lambda +b_{\mathrm{0}}\right)}+\frac{b_{\mathrm{1}}}{\lambda \left(b_{\mathrm{1}}\lambda +b_{\mathrm{0}}\right)}-\frac{\tau }{\lambda }，

又因为

\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left\{\mathrm{R}\mathrm{e}\left[{\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau }\right)}^{-\mathrm{1}}{|}_{\lambda ={\omega }_{\mathrm{*}}\mathrm{i}}\right]\right\}=\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left\{\mathrm{R}\mathrm{e}\left[\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau }\right){|}_{\lambda ={\omega }_{\mathrm{*}}\mathrm{i}}\right]\right\}，

所以

\begin{array}{l}\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left\{\mathrm{R}\mathrm{e}\left[\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau }\right){|}_{\lambda ={\omega }_{\mathrm{*}}\mathrm{i}}\right]\right\}=\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left\{\mathrm{R}\mathrm{e}\left[\frac{\left(\mathrm{3}{\lambda }^{\mathrm{2}}+\mathrm{2}{a}_{\mathrm{2}}\lambda +a_{\mathrm{1}}\right){\mathrm{e}}^{\lambda \tau }}{\lambda \left(b_{\mathrm{1}}\lambda +b_{\mathrm{0}}\right)}+\frac{b_{\mathrm{1}}}{\lambda \left(b_{\mathrm{1}}\lambda +b_{\mathrm{0}}\right)}-\frac{\tau }{\lambda }{|}_{\lambda ={\omega }_{\mathrm{*}}\mathrm{i}}\right]\right\}=\\ \mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\left\{\frac{\mathrm{3}{\omega }_{\mathrm{*}}^{\mathrm{6}}+\mathrm{2}p{\omega }_{\mathrm{*}}^{\mathrm{4}}+q{\omega }_{\mathrm{*}}^{\mathrm{2}}}{b_{\mathrm{1}}^{\mathrm{2}}{\omega }_{\mathrm{*}}^{\mathrm{4}}+b_{\mathrm{0}}{\omega }_{\mathrm{*}}^{\mathrm{2}}}\right\}=\frac{z_{\mathrm{*}}}{b_{\mathrm{1}}^{\mathrm{2}}{\omega }_{\mathrm{*}}^{\mathrm{4}}+b_{\mathrm{0}}{\omega }_{\mathrm{*}}^{\mathrm{2}}}{H}^{\text{'}}\left(z_{\mathrm{*}}\right)\ne \mathrm{0}.\end{array}

基于以上分析，可总结出如下定理.

定理4 1） 若r\ge \mathrm{0},\mathrm{\Delta }\le \mathrm{0}，则对任意的\tau \ge \mathrm{0}时，地方病平衡点E_{\mathrm{*}}是局部渐近稳定的；2） 若或r\ge \mathrm{0}，且\mathrm{\Delta }>\mathrm{0},z_{\mathrm{*}}>\mathrm{0}，H\left(z_{\mathrm{*}}\right)\le \mathrm{0}，则对于任意的\tau \in \left.\left[\mathrm{0},{\tau }_k\right.\right)，地方病平衡点E_{\mathrm{*}}是局部渐近稳定的；3） 若条件2）满足，且H^{\text{'}}\left(z_{\mathrm{*}}\right)\ne \mathrm{0}，则\tau ={\tau }_k时，系统\left(\mathrm{1}\right)在地方病平衡点E_{\mathrm{*}}处出现\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支.

5 数值模拟

利用数值模拟，以说明系统\left(\mathrm{1}\right)无病平衡点和地方病平衡点的稳定性，及其周期性行为和非周期性行为的存在性.

1） 为了说明系统\left(\mathrm{1}\right)无病平衡点的全局稳定性，假设系统\left(\mathrm{1}\right)的参数如下：

\begin{array}{c}b=\mathrm{0.6},K=\mathrm{100},\beta =\mathrm{0.0164},c_{\mathrm{*}}=\mathrm{0.22},\\ d_{\mathrm{1}}=\mathrm{0.05},m=\mathrm{0.43},\sigma \left(a\right)=\mathrm{0.4},\end{array}

\begin{array}{c}{d}_{\mathrm{2}}=\mathrm{0.10},{\mu }_{\mathrm{1}}=\mathrm{0.05},{\gamma }_{\mathrm{1}}=\mathrm{0.25},{\gamma }_{\mathrm{2}}=\mathrm{0.714},\\ d_{\mathrm{3}}=\mathrm{0.02},{\mu }_{\mathrm{2}}=\mathrm{0.03}.\end{array}

在这种情况下，.如图\mathrm{2}所示，无病平衡点E_{\mathrm{0}}=\left(\mathrm{100,0},\mathrm{0,0}\right)是全局渐近稳定的.

进一步研究在不同情况下，地方病平衡点E_{\mathrm{*}}的渐近行为.

2） 为了说明系统\left(\mathrm{1}\right)地方病平衡点的局部稳定性，假设系统\left(\mathrm{1}\right)的参数如下：

\begin{array}{c}b=\mathrm{0.6},K=\mathrm{1000},\beta =\mathrm{0.0164},\alpha \left(a\right)=\mathrm{0.8},\\ d_{\mathrm{1}}=\mathrm{0.05},m=\mathrm{0.15},\sigma \left(a\right)=\mathrm{0.4},\end{array}

\begin{array}{c}{d}_{\mathrm{2}}=\mathrm{0.05},{\mu }_{\mathrm{1}}=\mathrm{0.05},{\gamma }_{\mathrm{1}}=\mathrm{0.12},{\gamma }_{\mathrm{2}}=\mathrm{0.55},\\ d_{\mathrm{3}}=\mathrm{0.02},{\mu }_{\mathrm{2}}=\mathrm{0.03}.\end{array}

在这种情况下，R_{\mathrm{0}}=\mathrm{2.8367}>\mathrm{1}.如图\mathrm{3}所示，地方病平衡点E_{\mathrm{*}}是局部渐近稳定的.

3） 为了说明不同时滞下E_{\mathrm{*}}的稳定性，假设系统\left(\mathrm{1}\right)的参数如下：

\begin{array}{l}b=\mathrm{0.5},K=\mathrm{1000},\beta =\mathrm{0.00164},\\ d_{\mathrm{1}}=\mathrm{0.05},m=\mathrm{0.06},\sigma \left(a\right)=\mathrm{0.2},\end{array}

\begin{array}{c}{d}_{\mathrm{2}}=\mathrm{0.05},{\mu }_{\mathrm{1}}=\mathrm{0.05},{\gamma }_{\mathrm{1}}=\mathrm{0.12},\\ {\gamma }_{\mathrm{2}}=\mathrm{0.30},d_{\mathrm{3}}=\mathrm{0.10},{\mu }_{\mathrm{2}}=\mathrm{0.05}.\end{array}

分别取\tau =\mathrm{0.8},\tau =\mathrm{0.5},\tau =\mathrm{0.1}，图\mathrm{4}显示了在不同时滞下E_{\mathrm{*}}的稳定性，比较图3中不同\tau的曲线，可以发现随着时滞\tau的减小，易感者和潜伏者的数量会减少，感染者和隔离者的数量会增加.然而，如果\tau >{\tau }_k，E_{\mathrm{*}}将变得不稳定，并将导致在E_{\mathrm{*}}附近出现\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支.

4） 为了说明不同时滞下E_{\mathrm{*}}的稳定性，假设系统的参数如下：

\begin{array}{c}b=\mathrm{0.5},K=\mathrm{1000},\beta =\mathrm{0.00164},\\ d_{\mathrm{1}}=\mathrm{0.05},m=\mathrm{0.35},\sigma \left(a\right)=\mathrm{0.2},\end{array}

\begin{array}{c}{d}_{\mathrm{2}}=\mathrm{0.05},{\mu }_{\mathrm{1}}=\mathrm{0.05},{\gamma }_{\mathrm{1}}=\mathrm{0.12},\\ {\gamma }_{\mathrm{2}}=\mathrm{0.30},d_{\mathrm{3}}=\mathrm{0.10},{\mu }_{\mathrm{2}}=\mathrm{0.05}.\end{array}

如图\mathrm{5}所示，分别取\tau =\mathrm{3},\tau =\mathrm{2.5},\tau =\mathrm{1.5}，地方病平衡点E_{\mathrm{*}}失去稳定性，并且出现了周期解，可以看出感染者的人数随着时滞的增大而增大，潜伏者随着时滞的增大而减少.

6 结论

本文建立了具有\mathrm{L}\mathrm{o}\mathrm{g}\mathrm{i}\mathrm{s}\mathrm{t}\mathrm{i}\mathrm{c}\mathrm{s}增长和年龄结构的\mathrm{S}\mathrm{E}\mathrm{I}\mathrm{Q}\mathrm{R}传染病模型.研究得到，当时，E_{\mathrm{0}}是全局渐近稳定的；当R_{\mathrm{0}}>\mathrm{1}时，E_{\mathrm{0}}是不稳定的，疾病将持续存在.证明了当R_{\mathrm{0}}>\mathrm{1},\tau =\mathrm{0}时，系统的地方病平衡点E_{\mathrm{*}}是局部渐近稳定的.当E_{\mathrm{*}}不稳定时，证明了在其周围存在\mathrm{H}\mathrm{o}\mathrm{p}\mathrm{f}分支.最后利用\mathrm{M}\mathrm{a}\mathrm{t}\mathrm{l}\mathrm{a}\mathrm{b}软件对传染病模型进行数值模拟，验证了所获结论的有效性.

参考文献

[1] KERMACK W O， MCKENDRICK A G. A contribution to the mathematical theory of epidemics［J］. Proceedings of the Royal Society of London. Series A， Containing Papers of a Mathematical and Physical Character， 1927， 115（772）： 700-721.

[2] LI W， JI J， HUANG L， et al. Global dynamics of a controlled discontinuous diffusive SIR epidemic system［J/OL］. Applied Mathematics Letters， 2021， 121［2024-12-25］.

[3] TIAN C， LIU Z， RUAN S G. Asymptotic and transient dynamics of SEIR epidemic models on weighted networks［J］. European Journal of Applied Mathematics， 2023， 34（2）： 238-261.

[4] LI J， TENG Z， WANG G， et al. Stability and bifurcation analysis of an SIR epidemic model with logistic growth and saturated treatment［J］. Chaos， Solitons and Fractals， 2017， 99： 63-71.

[5] WANG Q， XIANG K N， ZHU C H， et al. Stochastic SEIR epidemic models with virus mutation and logistic growth of susceptible populations［J］. Mathematics and Computers in Simulation， 2023， 212： 289-309.

[6] LYU Y， CHEN L， CHEN F， et al. Stability and bifurcation in an SI epidemic model with additive Allee effect and time delay［J/OL］. International Journal of Bifurcation and Chaos， 2021， 31（4）［2024-12-25］.

[7] MA Y， CUI Y， WANG M. Global stability and control strategies of a SIQRS epidemic model with time delay［J］. Mathematical Methods in the Applied Sciences， 2022， 45（13）： 8269-8293.

[8] ZHANG Z， KUNDU S， TRIPATHI J P， et al. Stability and Hopf bifurcation analysis of an SVEIR epidemic model with vaccination and multiple time delays［J/OL］. Chaos， Solitons and Fractals， 2020， 131［2024-12-25］.

[9] HOPPENSTEND F. An age-structured epidemic model［J］. Journal of the Franklin Institute， 1994， 197： 325-333.

[10] LIU P Y， LI H X. Global behavior of a multi-group SEIR epidemic model with age structure and spatial diffusion［J］. Mathematical Biosciences and Engineering， 2020， 17（6）： 7248-7273.

[11] HUANG J C， KANG H， LU M， et al. Stability analysis of an age structured epidemic model with vaccination and standard incidence rate［J/OL］. Nonlinear Analysis： Real World Applications， 2022， 66［2024-12-25］.

[12] YU Y， TAN Y， TANG S Y. Stability analysis of the COVID-19 model with age structure under media effect［J/OL］. Computational and Applied Mathematics， 2023， 42（5）［2024-12-25］.

[13] LIU Q， JIANG D， SHI N. Threshold behavior in a stochastic SIQR epidemic model with standard incidence and regime switching［J］. Applied Mathematics and Computation， 2018， 316： 310-325.

[14] WANG S， DING Y， LU H， et al. Stability and bifurcation analysis of SIQR for the COVID-19 epidemic model with time delay［J］. Mathematical Biosciences and Engineering， 2021， 18（5）： 5505-5524.

[15] YANG J， SHI X， SONG X， et al. Threshold dynamics of a stochastic SIQR epidemic model with imperfect quarantine［J/OL］. Applied Mathematics Letters， 2023， 136［2024-12-25］.

[16] ZHAO H， LIN Y， DAI Y. An SIRS epidemic model incorporating media coverage with time delay［J/OL］. Computational and Mathematical Methods in Medicine， 2014， 2014（1）［2024-12-25］.

基金资助

国家自然科学基金项目(11801243)

甘肃省高等学校青年博士项目(2022QB-056)