猴痘病毒传播的传染病动力学建模与定性分析

Infectious Disease Dynamics Modeling and Qualitative Analysis of Monkeypox Virus Transmission

摘要: 本文通过对猴痘病毒的传播机制的分析，并考虑了天花疫苗对猴痘病毒的交叉防护，构建了包含隔离措施和疫苗接种的SEIQVR猴痘病毒模型。通过应用下一代矩阵法，计算了模型的基本再生数$R_0=max\left\{R_{01},R_{02}\right\}$ ，其中人类种群的基本再生数为$R_{01}$ ，动物种群的基本再生数为$R_{02}$ ，并验证了相关的阈值定理。根据$R_{01}$ 和$R_{02}$ 的值，我们确定了三类平衡点的存在性，包括无病平衡点、无猴痘地方病平衡点及共存的地方病平衡点。利用Hurwitz判据，我们证明了这些平衡点的局部稳定性，并通过构造Lyapunov函数结合LaSalle不变性原理评估了其全局稳定性。数值模拟结果进一步验证了所得到结论的合理性。

Abstract: This study analyzes the transmission mechanism of the monkeypox virus and considers the cross-protection of the smallpox vaccine against monkeypox. We construct an SEIQVR monkeypox model that incorporates isolation measures and vaccination. By applying the next generation matrix method, we calculate the basic reproduction number $R_0=max\left\{R_{01},R_{02}\right\}$ , where $R_{01}$ is the basic reproduction number for the human population and $R_{02}$ is the basic reproduction number for the animal population, and we verify the related threshold theorem. Based on the values of $R_{01}$ and $R_{02}$ , we determine the existence of three types of equilibrium points: The disease-free equilibrium, the monkeypox endemic equilibrium, and the coexistence endemic equilibrium. Using the Hurwitz criterion, we prove the local stability of these equilibrium points and assess their global stability by constructing a Lyapunov function in conjunction with the LaSalle invariance principle. Numerical simulation results further validate the reasonableness of the conclusions obtained.

1. 引言

猴痘是一种病毒性人畜共患病，其临床症状范围从中度到重度不等。该病主要通过与非洲中部和南部热带雨林地区的农村和城市动物直接接触传播给人类。作为一种典型的动物源性病原体，猴痘病毒的自然宿主范围广泛，包括多种啮齿目、土拨鼠科、松鼠科以及类人猿科的黑猩猩等。鉴于该病毒在上述宿主动物中的持续循环，其在自然界中非常难以消除。该病毒具备跨物种传播的能力，存在从动物宿主向人类传播的潜在风险。传播途径主要包括直接接触感染者的体液、吸入含病毒的飞沫核，以及接触被病毒污染的环境或物品表面。这些途径均有可能引发病毒在人群中的感染和传播[1] [2]。

自2021年起至2022年9月2日，全球多地报告了2017例猴痘新发病例，其中导致500人死亡。此次疫情并非该病毒首次肆虐。猴痘病毒最初于2000年在刚果被记录，该地区曾是天花疫情的重灾区。猴痘病例主要位于非洲中西部农村地区。刚果疫情之后，贝宁、喀麦隆、加蓬和尼日利亚等国也报告了病例。尼日利亚在过去三年中病例最多，报告了2007例疑似病例，70例确诊，死亡率达3%。2007年，美国也爆发了一起与宠物感染相关的疫情，确诊70例。随着病例增加，深入研究猴痘的流行特性、传播方式和源头，以及开发治疗和疫苗变得迫切。科学防控策略的制定对阻止病毒传播极为关键。

近期研究致力于揭示猴痘病毒的传播机制。Bhunu和Mushayabasa [3]首次运用SIR模型来探讨猴痘在跨物种交互群体中的传播规律，并通过Lyapunov稳定性理论证实了模型中地方性平衡状态下的局部稳定性。为持续保护易感人群免受猴痘病毒感染，Usman [4]提出了一个猴痘病毒动力学模型，研究了平衡点的稳定性。Somma等人[5]建立了一个包含人类和啮齿动物两个相互作用宿主群体的猴痘病毒传播数学模型，考虑了隔离措施和公众教育活动对疾病传播的影响。随后，Peter等人[6]提出了一个确定性传染病动力学模型，以探究猴痘病毒的传播机制，该研究不仅确定了无病平衡点和地方平衡点的局部和全局稳定性条件，还揭示了模型中存在一种后向分岔现象。通过引入针对猴痘疫情的疫苗接种策略，Lasisi等人[7]通过建立猴痘病毒传播的数学模型，深入剖析了其传播机制，特别强调了在疫情暴发时为高风险群体接种疫苗的关键作用。此外，Khan等人[8]利用随机微分方程模型，分析了一个广义的随机猴痘流行模型，基于疾病的复杂性，将猴痘病毒作为一个分区系统，采取随机微分方程组进一步描述了猴痘病毒的传播特性。Madubueze等[9]考虑了环境传播途径、检疫措施以及疫苗接种等多种因素，提出了一个猴痘病毒的动力学模型，并确定了三个重要的平衡点，通过基本再生数的分析探讨了这些平衡状态的局部与全局稳定性。Majee和Suvankar [10]构建了一个包含七个隔室的分数阶数学模型，并详细研究了该模型解的存在性和唯一性。Yapkan等人[11]提出了一类具有治疗仓室的SITR新型猴痘传播模型，确定减少传染人群的优化控制策略，进行稳定性和灵敏度分析，并考虑了最优控制策略。Ahmad [12]建立了包括暴露后接种、暴露前接种和隔离的确定性SVIR传染病模型，来描述猴痘病毒的人际传播动态。

本文在文献[7]和[9]-[12]的基础上，对猴痘病毒传播进行了深入研究，提出了包含接种疫苗人群和隔离人群的SEIQVR模型，以更精确地描述病毒的传播动态。本文的结构安排如下：第二节中，我们详细构建了SEIQVR仓室模型，并验证了模型的正不变性及解的有界性。第三节探讨了无病平衡点、无猴痘地方病平衡点及共存的地方病平衡点的存在性，并分析了它们的局部和全局稳定性。第四节通过数值模拟验证了结论的合理性。最后一部分对本文的工作进行了总结。

2. SEIQVR传染病动力学模型

2.1. 模型建立

考虑到猴痘病毒在人群与动物群体间的传播特性，将人群划分为如下六个仓室：易感者(S)、潜伏者类(E)、猴痘感染者类(I)、隔离者类(Q)、治愈类者(R)、接种疫苗者类(V)，并用$S\left(t\right),E\left(t\right),I\left(t\right),Q\left(t\right),R\left(t\right),V\left(t\right)$ 分别表示t时刻相应仓室人群的数量。将动物种群划分为四个仓室，即易感动物者类($S_A$ )，潜伏动物者类($E_A$ )，患病动物者类($I_A$ )，恢复动物者类($R_A$ )，分别用$S_A\left(t\right),E_A\left(t\right),I_A\left(t\right),R_A\left(t\right)$ 分别表示t时刻相应动物种群的数量。此外，设t时刻总人口数$N_H\left(t\right)$ 和动物种群总数$N_A\left(t\right)$ 分别定义如下：

$N_H\left(t\right)=S\left(t\right)+E\left(t\right)+I\left(t\right)+Q\left(t\right)+R\left(t\right)+V\left(t\right)$ ，$N_A\left(t\right)=S_A\left(t\right)+E_A\left(t\right)+I_A\left(t\right)+R_A\left(t\right)$

为考虑猴痘病毒在人类与动物之间的相互作用，我们做出如下假设：

1) 考虑到没有针对猴痘病毒的专门疫苗，而天花疫苗对猴痘的有效性仅为85%左右。接种假设免疫人群以$\epsilon$ 的速率向易感人群转移。

2) 被隔离患者与外界没有接触无法传播疾病。

3) 仅考虑人群的出生率，自然死亡率与因病致死率，不考虑人口迁入与迁出率。

4) 对于啮齿类动物，假设它们能够依靠自身的恢复力来战胜疾病，在进入恢复者类后，不考虑二次感染的可能性。

考虑到生物学意义，所有模型参数为正数，具体含义见表1。本研究围绕猴痘病毒引发的⼈畜共患病特点，结合现行防控策略，深入分析了关键因素间的传递规律。通过对确诊感染者的隔离措施以及对易感人群的疫苗接种方案，构建了一个详尽的仓室流程图。该流程图清晰展示了猴痘病毒在人类与动物宿主之间的传播路径，如图1所示。

Table 1. The meaning parameters of SEIQVR model

表1. 参数意义

参数	意义
$\Lambda ,{\Lambda }_1$	人类动物自然出生率
$\beta ,{\beta }_1$	人与人之间，动物之间的有效接率
$\xi$	表示S与$I_A$ 之间的交叉感染率
$\mu$	人类自然死亡率
${\mu }_1$	动物自然死亡率
$\phi$	感染者的隔离率
$\omega$	动物自我治愈率
d	人类因感染猴痘死亡率
q	潜伏者的隔离率
$\epsilon$	接种者免疫丧失率
$\delta$	疫苗有效保护率
$r_1,r_2$	治愈率，个体自我修护率
$\theta$	动物潜伏向感染的转化率
$d_1$	动物因感染猴痘死亡率

Figure 1. Flow chart of the transformation between different populations and rodents of the SEIQVR model

图1. SEIQVR模型的各类人群及啮齿动物之间转化流程图

鉴于猴痘病传播特点，根据仓室流程图1，我们可以得到包含隔离仓室、疫苗接种仓室和动物室的SEAIQR传染病模型。其微分方程形式如下：

$\{\begin{array}{l}{S}^{\prime }\left(t\right)=\Lambda -\left(\beta SI+\xi S{I}_A\right)-\mu S+\epsilon V-\delta S\hfill \\ E^{\prime }\left(t\right)=\beta SI+\xi S{I}_A-M_{1}E\hfill \\ I^{\prime }\left(t\right)=\sigma E-M_{2}I\hfill \\ Q^{\prime }\left(t\right)=\phi I-M_{3}Q+qE\hfill \\ R^{\prime }\left(t\right)=\tau Q+M_{4}I-\mu R\hfill \\ V^{\prime }\left(t\right)=\delta S-\epsilon V-\mu V\hfill \\ {S^{\prime }}_A\left(t\right)={\Lambda }_1-{\beta }_1{S}_A{I}_A-{\mu }_1{S}_A\hfill \\ {E^{\prime }}_A\left(t\right)={\beta }_1{S}_A{I}_A-M_5{E}_A\hfill \\ {I^{\prime }}_A\left(t\right)=\theta E_A-M_6{I}_A\hfill \\ {R^{\prime }}_A\left(t\right)=\omega I_A-{\mu }_1{R}_A\hfill \end{array}$ (2.1)

其中$M_1=\mu +\sigma +q$ ，$M_2=\mu +d+\phi +r_1+r_2$ ，$M_3=\mu +d+\tau$ ，$M_4=r_1+r_2$ ，$M_5=\theta +{\mu }_1$ ，$M_6=\omega +{\mu }_1+d_1$

2.2. 正不变性与解的有界性

令$x=\left(S,E,I,Q,R,V,S_A,E_A,I_A,R_A\right)$ 及$f\left(x\right)=\left(f_1\left(x\right),f_2\left(x\right),\cdots ,f_{10}\left(x\right)\right)$ ，其中

$\{\begin{array}{l}{f}_1\left(x\right)=\Lambda -\left(\beta SI+\xi S{I}_A\right)-\mu S+\epsilon V-\delta S\hfill \\ f_2\left(x\right)=\beta SI+\xi S{I}_A-M_{1}E\hfill \\ f_3\left(x\right)=\sigma E-M_{2}I\hfill \\ f_4\left(x\right)=\phi I-M_{3}Q+qE\hfill \\ f_5\left(x\right)=\tau Q+M_{4}I-\mu R\hfill \\ f_6\left(x\right)=\delta S-\epsilon V-\mu V\hfill \\ f_7\left(x\right)={\Lambda }_1-{\beta }_1{S}_A{I}_A-{\mu }_1{S}_A\hfill \\ f_8\left(x\right)={\beta }_1{S}_A{I}_A-M_5{E}_A\hfill \\ f_9\left(x\right)=\theta E_A-M_6{I}_A\hfill \\ f_{10}\left(x\right)=\omega I_A-{\mu }_1{R}_A\hfill \end{array}$ (2.2)

则系统(2.1)可写成如下向量模式：

$x^{\prime }=f\left(x\right)$ (2.3)

定理1 1) 集合$R_{+}^{10}$ 关于系统(2.1)是正不变的。

2) 考虑两个区域${\Omega }_1$ 和${\Omega }_2$ ，其中$\Omega ={\Omega }_1\times {\Omega }_2$ ，

${\Omega }_1=\left\{\left(S\left(t\right),E\left(t\right),I\left(t\right),Q\left(t\right),R\left(t\right),V\left(t\right)\right)\in R_6^{+}:0\le S\left(t\right)+E\left(t\right)+I\left(t\right)+Q\left(t\right)+R\left(t\right)+V\left(t\right)\le \frac{\Lambda }{\mu }\right\}$ ，

${\Omega }_2=\left\{\left(S_A\left(t\right),E_A\left(t\right),I_A\left(t\right),R_A\left(t\right)\right)\in R_4^{+}:0\le S_A\left(t\right)+E_A\left(t\right)+I_A\left(t\right)+R_A\left(t\right)\le \frac{{\Lambda }_1}{{\mu }_1}\right\}$

则区域Ω关于系统(2.1)是正不变的。

证明：1) 对任意的$x\left(t\right)\in R_{+}^{10}$ ，

$\{\begin{array}{l}{f_1\left(x\right)|}_{S=0}=\Lambda +\epsilon V>0\\ {f_2\left(x\right)|}_{E=0}=\beta SI+\xi S{I}_A\ge 0\\ {f_3\left(x\right)|}_{I=0}=\sigma E\ge 0\\ {f_4\left(x\right)|}_{Q=0}=\phi I\ge 0\\ {f_5\left(x\right)|}_{R=0}=\tau Q+\left(r_1+r_2\right)I\ge 0\\ {f_6\left(x\right)|}_{V=0}=\delta S\ge 0\\ {f_7\left(x\right)|}_{S_A=0}={\Lambda }_1>0\\ {f_8\left(x\right)|}_{E_A=0}={\beta }_1{S}_A{I}_A\ge 0\\ {f_9\left(x\right)|}_{I_A=0}=\theta E_A\ge 0\\ {f_{10}\left(x\right)|}_{R_A=0}=\omega I_A\ge 0\end{array}$ (2.4)

故对任意的$x\left(0\right)\ge 0$ ，都有$x\left(t\right)\ge 0$ ，根据不变性原理[12] below知集合$R_{+}^{10}$ 关于系统(2.1)是正不变的。

2) 根据系统(2.1)，总人口数为$N_H\left(t\right)$ 和啮齿类动物种群总数$N_A\left(t\right)$ 满足：

$\frac{\text{d}{N}_H\left(t\right)}{\text{d}t}=\Lambda -\mu N_H\left(t\right)-d\left(Q+I\right)\le \Lambda -\mu N_H\left(t\right)$ ，$\frac{\text{d}{N}_A\left(t\right)}{\text{d}t}={\Lambda }_1-{\mu }_1{N}_A\left(t\right)-d_1{I}_A\le {\Lambda }_1-{\mu }_1{N}_A\left(t\right)$ ，

由比较定理[12]可知，对任意的$t\ge 0$ ，

$N_H\left(t\right)\le \frac{\Lambda }{\mu }-{\text{e}}^{-\mu t}\left(\frac{\Lambda }{\mu }-N_H\left(0\right)\right)\le \frac{\Lambda }{\mu }$ ；

$N_A\left(t\right)\le \frac{{\Lambda }_1}{{\mu }_1}-{\text{e}}^{-{\mu }_{1}t}\left(\frac{{\Lambda }_1}{{\mu }_1}-N_A\left(0\right)\right)\le \frac{{\Lambda }_1}{{\mu }_1}$

所以区域Ω关于系统(2.1)是正不变的。

3. 动力学分析

3.1. 基本再生数

当疾病消亡时，即$E=I=Q=R=E_A=I_A=R_A=0$ 易知模型(2.1)始终存在无病平衡点：

$E_0\left(S^0,0,0,0,0,V^0,\frac{{\Lambda }_1}{{\mu }_1},0,0,0\right)$

其中$S^0=\frac{\Lambda \left(\epsilon +\mu \right)}{\mu \left(\epsilon +\mu +\delta \right)}$ ，$V^0=\frac{\Lambda \delta }{\mu \left(\epsilon +\mu +\delta \right)}$ 。令$u=\left(E,I,Q,E_A,I_A\right)$ 表示染病仓室，且

$ℱ=\left(\begin{array}{c}\beta SI+\xi S{I}_A\\ 0\\ 0\\ {\beta }_1{S}_A{I}_A\\ 0\end{array}\right)$ ，$\mathbb{V}=\left(\begin{array}{c}{M}_{1}E\\ M_{2}I-\sigma E\\ M_{3}Q-qE-\varphi I\\ M_5{E}_A\\ M_6{I}_A-\theta E_A\end{array}\right)$

其中$M_1=\mu +\sigma +q$ ，$M_4=r_1+r_2$ ，$M_2=\mu +d+\varphi +r_1+r_2$ ，$M_5=\theta +{\mu }_1$ ，$M_3=\mu +d+\tau$ ，$M_6=\omega +{\mu }_1+d_1$ 。$ℱ$ 和$\mathbb{V}$ 在$E_0$ 处雅可比矩阵分别为：$F=Dℱ\left(E_0\right)$ ，$V=D\mathbb{V}\left(E_0\right)$ ，计算得：

$F=\left(\begin{array}{ccccc}0& \frac{\Lambda \left(\epsilon +\mu \right)}{\mu \left(\epsilon +\mu +\delta \right)}\beta & 0& 0& \frac{\Lambda \left(\epsilon +\mu \right)}{\mu \left(\epsilon +\mu +\delta \right)}\xi \\ 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0\\ 0& 0& 0& 0& {\beta }_1\frac{{\Lambda }_1}{{\mu }_1}\\ 0& 0& 0& 0& 0\end{array}\right)$ ，$V=\left(\begin{array}{ccccc}{M}_1& 0& 0& 0& 0\\ -\sigma & M_2& 0& 0& 0\\ -q& -\phi & M_3& 0& 0\\ 0& 0& 0& M_5& 0\\ 0& 0& 0& -\theta & M_6\end{array}\right)$

通过再生矩阵方法[13]，知模型(2.1)的基本再生数为：

$R_0=\rho \left(F{V}^{-1}\right)=\max \left\{R_{01},R_{02}\right\}$ (3.1)

其中$R_{01}=\frac{\beta \Lambda \left(\epsilon +\mu \right)\sigma }{\mu \left(\delta +\epsilon +\mu \right)M_1{M}_2}$ ，$R_{02}=\frac{{\Lambda }_1{\beta }_1\theta }{{\mu }_1{M}_5{M}_6}$ ($R_{01,02}$ 分别表示人类的基本再生数与啮齿类动物的基本再生数)。

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

定理3.2.1 当$R_0<1$ 时，模型(2.1)的无病平衡点$E_0$ 是局部渐近稳定的；当$R_0>1$ 时，$E_0$ 是不稳定的。

证明：模型(2.1)在$E_0$ 处的雅可比矩阵为

$J_0={\left[\begin{array}{cc}{A}_{6\times 6}& B_{4\times 6}\\ 0_{6\times 4}& D_{4\times 4}\end{array}\right]}_{10\times 10}$ (3.2)

$A=\left[\begin{array}{cccccc}-\mu -\delta & 0& -\beta x_1& 0& 0& \delta \\ 0& -M_1& \beta x_1& 0& 0& 0\\ 0& \sigma & -M_2& 0& 0& 0\\ 0& q& \phi & -M_3& 0& 0\\ 0& 0& M_4& \tau & -\mu & 0\\ \delta & 0& 0& 0& 0& -\epsilon -\mu \end{array}\right]$ ，$B=\left[\begin{array}{cccc}0& 0& x_1\xi & 0\\ 0& 0& -x_1\xi & 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right]$ ，$D=\left[\begin{array}{cccc}-{\mu }_1& 0& -x_7{\beta }_1& 0\\ 0& -M_5& x_7{\beta }_1& 0\\ 0& \theta & -M_6& 0\\ 0& 0& \omega & -{\mu }_1\end{array}\right]$

简单计算，可得矩阵$J_0$ 的特征值：

${\lambda }_{1,2}=-{\mu }_1$ ，${\lambda }_{3,4}=-\mu \omega$ ，${\lambda }_5=-\left(\mu +\delta \right)$ ，${\lambda }_6=-M_3$ ，${\lambda }_7=-M_5$ ，${\lambda }_8=-M_6\left(1-R_{02}\right)$ ，${\lambda }_9,{\lambda }_{10}$ 由${\lambda }^2+\left(M_1+M_2\right)\lambda +M_1{M}_2\left(1-R_{01}\right)=0$ 给出。

由韦达定理可知只有当$R_{01}<1$ 时，${\lambda }_9,{\lambda }_{10}<0$ 。又因为当$R_{02}<1$ 时，${\lambda }_8=-M_6\left(1-R_{02}\right)<0$ 。因此当$R_{01},R_{02}<1$ 时，即$R_0<1$ 时，矩阵$J_0$ 所有特征根具有负实部，此时模型(2.1)的无病平衡点$E_0$ 是局部渐近稳定的。而当$R_{01},R_{02}$ 有一个大于1时，矩阵的特征值至少有一个正实部，此时无病平衡点$E_0$ 是不稳定的。

定理3.2.2 当$R_0<1$ 时，模型(2.1)的无病平衡点$E_0$ 是全局渐近稳定的；当$R_0>1$ 时，$E_0$ 是不稳定的。

证明：由于系统(2.1)最后四个方程中不含人类种群，考虑下述子系统：

$\{\begin{array}{l}{S^{\prime }}_A\left(t\right)={\Lambda }_1-{\beta }_1{S}_A{I}_A-{\mu }_1{S}_A\\ {E^{\prime }}_A\left(t\right)={\beta }_1{S}_A{I}_A-M_5{E}_A\\ {I^{\prime }}_A\left(t\right)=\theta E_A-M_6{I}_A\\ {R^{\prime }}_A\left(t\right)=\omega I_A-{\mu }_1{R}_A\end{array}$ (3.3)

令$x_A=\left(E_A,I_A\right)$ 表示染病仓室，$y_A=\left(S_A,R_A\right)$ 表示未染病仓室，则模型(3.3)可写

$\{\begin{array}{l}{x^{\prime }}_A=\left(F_A-V_A\right)x_A-f_A\left(x_A,y_A\right)\\ {y^{\prime }}_A=g_A\left(x_A,y_A\right)\end{array}$

其中

$F_A=\left[\begin{array}{cc}0& \frac{{\beta }_1{\Lambda }_1}{{\mu }_1}\\ 0& 0\end{array}\right]$ ，$V_A=\left[\begin{array}{cc}{M}_5& 0\\ -\theta & M_6\end{array}\right]$ ，$g_A\left(x_A,y_A\right)=\left(\begin{array}{c}{\Lambda }_1-{\beta }_1{S}_A{I}_A-{\mu }_1{S}_A\\ \omega I_A-{\mu }_1{R}_A\end{array}\right)$ ，$f_A\left(x_A,y_A\right)=\left(\begin{array}{c}{\beta }_1\left(\frac{{\Lambda }_1}{{\mu }_1}-S_A\right)I_A\\ 0\end{array}\right)$

根据定理1，可得$f_A\left(x_A,y_A\right)\ge 0$ 。构造Lyapunov函数

${\Phi }_A\left(x_A\right)={\omega }_A^{\text{T}}{V}_A^{-1}{x}_A$

其中${\omega }_A^{\text{T}}$ 为矩阵$V_A^{-1}{F}_A$ 关于特征值$R_{02}$ 的左特征向量。取${\omega }_A^{\text{T}}=\left(0,1\right)$ ，计算得到：

${\Phi }_A\left(x_A\right)={\omega }^{\text{T}}{V}_A^{-1}{x}_A=\frac{{\Lambda }_1{\beta }_1\theta }{{\mu }_1{M}_5{M}_6}{I}_A$

则当$R_{02}<1$ 时，${\Phi }_A\left(x_A\right)$ 沿着模型(3.3)的全导数为：

$\frac{\text{d}}{\text{d}t}{\Phi }_A\left(x_A\right)={\omega }^{\text{T}}{V}_A^{-1}\left(F-V\right)x_A-{\omega }^{\text{T}}{V}_A^{-1}{f}_A\left(x_A,y_A\right)=\left(R_{02}-1\right){\omega }^{\text{T}}{x}_A-{\omega }^{\text{T}}{V}_A^{-1}{f}_A\left(x_A,y_A\right)\le 0$ 。

进而得到

$\frac{\text{d}}{\text{d}t}{\Phi }_A\left(x_A\right)=\left(R_{02}-1\right)I_A-\frac{\theta }{M_5{M}_6}\left(\frac{{\Lambda }_1}{{\mu }_1}-S_A\right)I_A\le 0$ (3.4)

易知${{\Phi }^{\prime }}_A\left(x_A\right)=0$ 当且仅当$S_A=\frac{{\Lambda }_1}{{\mu }_1}$ ，$I_A=0$ ，$E_A=0$ ，$R_A=0$ 。由LaSalle不变原理[14] [15]知，当$R_{02}<1$ 时，模型(3.3)的无病平衡点是全局渐近稳定的，即

$\underset{t\to +\infty }{\lim }{S}_A\left(t\right)=\frac{{\Lambda }_1}{{\mu }_1}$ ，$\underset{t\to +\infty }{\lim }{E}_A\left(t\right)=0$ ，$\underset{t\to +\infty }{\lim }{I}_A\left(t\right)=0$ ，$\underset{t\to +\infty }{\lim }{R}_A\left(t\right)=0$ ，

由于当$R_{02}<1$ 时，$\underset{t\to +\infty }{\lim }{I}_A\left(t\right)=0$ ，此时由前六个方程控制的人类种群系统对应的极限系统为：

$\{\begin{array}{l}{S}^{\prime }\left(t\right)=\Lambda -\beta SI-\mu S+cV-\delta S\\ E^{\prime }\left(t\right)=\beta SI-M_{1}E\\ I^{\prime }\left(t\right)=\sigma E-M_{2}I\\ Q^{\prime }\left(t\right)=\phi I-M_{3}Q+qE\\ R^{\prime }\left(t\right)=\tau Q+M_{4}I-\mu R\\ V^{\prime }\left(t\right)=\delta S-cV-\mu V\end{array}$ (3.5)

令$x_H=\left(E,I,Q\right)$ ，$y_H=\left(S,R,V\right)$ 分别表示人类染病与未染病仓室，则模型(3.5)可写成

$\{\begin{array}{l}{x^{\prime }}_H=\left(F_H-V_H\right)x_H-f_H\left(x_H,y_H\right)\\ {y^{\prime }}_H=g_H\left(x_H,y_H\right)\end{array}$

其中

$F_H=\left[\begin{array}{ccc}\beta S^0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\right]$ ，$V_H=\left[\begin{array}{ccc}{M}_1& 0& 0\\ -\theta & M_2& 0\\ -q& -\phi & M_3\end{array}\right]$

$g_H\left(x_H,y_H\right)=\left(\begin{array}{c}\Lambda -\beta SI-\mu S+cV-\delta S\\ \tau Q+M_{4}I-\mu R\\ \delta S-\epsilon V-\mu V\end{array}\right)$ ，$f_H\left(x_H,y_H\right)=\left(\begin{array}{c}\beta \left(S^0-S\right)I\\ 0\\ \delta \left(S^0-S\right)-\left(\epsilon +\mu \right)\left(V^0-V\right)\end{array}\right)$

根据定理1，可得$f_H\left(x_H,y_H\right)\ge 0$ ，构造Lyapunov函数

$V_H\left(x_H\right)={\omega }^{\text{T}}{V}_H^{-1}{x}_H$

其中${\omega }^{\text{T}}$ 为矩阵$V_H^{-1}{F}_H$ 关于特征值$R_0$ 的左特征向量。经计算得：

$V_H^{-1}F=\left[\begin{array}{ccc}0& \frac{\beta S^0}{M_1}& 0\\ 0& \frac{\beta S^0\sigma }{M_1{M}_2}& 0\\ 0& \frac{\beta S^0\left(M_{2}q+\phi \sigma \right)}{M_1{M}_2{M}_3}& 0\end{array}\right]$

可取${\omega }^{\text{T}}=\left(0,1,0\right)$ ，计算得到：

${\Phi }_H\left(x_H\right)={\omega }^{\text{T}}{V}_H^{-1}{x}_H=\frac{\beta S^0\sigma }{M_1{M}_2}I$

则当$R_{01}<1$ 时，${\Phi }_H\left(x_H\right)$ 关于模型(3.5)的全导数为：

$\begin{array}{c}{{\Phi }^{\prime }}_H\left(x_H\right)={\omega }^{\text{T}}{V}_H^{-1}\left(F-V\right)x_H-{\omega }^{\text{T}}{V}_H^{-1}{f}_H\left(x_H,y_H\right)\\ =\left(R_{01}-1\right){\omega }^{\text{T}}{x}_H-{\omega }^{\text{T}}{V}_H^{-1}{f}_H\left(x_H,y_H\right)\\ =\left(R_{01}-1\right)I-\frac{\sigma }{M_1{M}_2}\beta \left(S^0-S\right)I\le 0\end{array}$

易知${{\Phi }^{\prime }}_H\left(x_H\right)=0$ 当且仅当$S=S^0$ ，$V=V^0$ ，$E=I=Q=0$ ，由LaSalle不变原理知，当$R_{01}<1$ 时，模型(3.5)的无病平衡点是全局渐近稳定的。

综上可知，根据极限系统理论和定理3.1，对于系统(2.1)当$R_{02}<1,R_{01}<1$ 时，即$R_0<1$ 时模型(2.1)的无病平衡点$E_0$ 是全局渐近稳定的。当$R_0>1$ 时，模型(3.3)的无病平衡点不稳定。

3.3. 无动物地方病平衡点的存在性及稳定性

若当子系统(3.3)存在无病平衡点，而子系统(3.5)存在一个地方病平衡点，该平衡点称为模型(2.1)无动物平衡点。

定理3.3.1 当$R_{01}>1$ 时，系统(2.1)存在无动物地方病平衡点$E_1$ 。

证明：通过系统(2.1)，令$E_A=I_A=R_A=0$ 可得到一个平衡点$E_1\left(S^{\text{*}},E^{\text{*}},I^{\text{*}},Q^{\text{*}},R^{\text{*}},V^{\text{*}},S_A^{*},0,0,0\right)$ 。

$S_A^{*}=\frac{{\Lambda }_1}{{\mu }_1}$ ，$S^{*}=\frac{M_1{M}_2}{\beta \sigma }$ ，$E^{*}=\frac{\Lambda \beta \sigma \left(\epsilon +\mu \right)-M_1{M}_2\mu \left(\delta +\epsilon +\mu \right)}{M_1\beta \left(\epsilon +\mu \right)}=\frac{\Lambda }{M_1}\left(1-\frac{1}{R_{01}}\right)$ ，

$I^{*}=\frac{\Lambda \beta \sigma \left(\epsilon +\mu \right)-M_1{M}_2\mu \left(\delta +\epsilon +\mu \right)}{M_1{M}_2\beta \left(\epsilon +\mu \right)}=\frac{\Lambda \sigma }{M_1{M}_2}\left(1-\frac{1}{R_{01}}\right)$ ，

$Q^{*}=\frac{AB}{M_1{M}_2{M}_3\beta \left(\epsilon +\mu \right)\sigma } =\frac{\Lambda \left(M_{2}q+\phi \sigma \right)}{M_1{M}_2{M}_3}\left(1-\frac{1}{R_{01}}\right)$ ，

$R^{*} =\frac{AC}{M_1{M}_2{M}_3\beta \left(\epsilon +\mu \right)\sigma \mu }=\frac{\Lambda \mu \left(M_{2}q+\phi \sigma +M_3{M}_4\sigma \right)}{{\mu }_1{M}_1{M}_2{M}_3}\left(1-\frac{1}{R_{01}}\right)$ ，$V^{*}=\frac{M_1{M}_2\delta }{\beta \sigma \left(\epsilon +\mu \right)}$ ，

其中$A=\left[\Lambda \beta \sigma \left(\epsilon +\mu \right)-M_1{M}_2\mu \left(\delta +\epsilon +\mu \right)\right]$ ，$B=\left(M_{2}q+\phi \sigma \right)$ ，$C=\left(M_{2}q\tau +\phi \tau \sigma +M_3{M}_4\sigma \right)\mu$ 。

综上可知，当$R_{01}>1$ 时，$E_1$ 的每个元素都是非负的，这与种群动态一致，故定理3.3.1得证。

定理3.3.2 若$R_{01}>1$ ，则模型(2.1)的无动物地方病平衡点当$R_{02}<1$ 时局部渐近稳定，当$R_{02}>1$ 时不稳定。

证明：系统(2.1)在$E_1$ 处的雅可比矩阵为：

$J_0={\left[\begin{array}{cc}{A}_{6\times 6}& B_{4\times 6}\\ 0_{6\times 4}& D_{4\times 4}\end{array}\right]}_{10\times 10}$

$A=\left[\begin{array}{cccccc}-\beta I^{*}-\mu -\delta & 0& -\beta S^{*}& 0& 0& c\\ \beta I^{**}& -M_1& \beta S^{*}& 0& 0& 0\\ 0& \sigma & -M_2& 0& 0& 0\\ 0& q& \phi & -M_3& 0& 0\\ 0& 0& M_4& \tau & -\mu & 0\\ \delta & 0& 0& 0& 0& -\epsilon -\mu \end{array}\right]$ ，$B={\left[\begin{array}{cccc}0& 0& -\xi S^{*}& 0\\ 0& 0& \xi S^{*}& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right]}_{6\times 4}$ ，

$D={\left[\begin{array}{cccc}-{\mu }_1& 0& -S_A^{*}{\beta }_1& 0\\ 0& -M_5& S_A^{*}{\beta }_1& 0\\ 0& \theta & -M_6& 0\\ 0& 0& \omega & -{\mu }_1\end{array}\right]}_{4\times 4}$ ，$M=\left(\begin{array}{ccc}-\mu & -M_1& 0\\ \beta I^{*}& -M_1& -\beta S^{*}\\ 0& -\sigma & -M_2\end{array}\right)$

矩阵$J_0$ 的特征根为：${\lambda }_{1,2}=-{\mu }_1$ ，${\lambda }_3=-M_3$ ，${\lambda }_4=-\mu$ ，${\lambda }_6=-M_5$ ，${\lambda }_5=M_6\left(R_{02}-1\right)$ ，${\lambda }_7=-\left(\epsilon +\delta +\mu \right)$ ，${\lambda }_{8,9,10}$ 是上述矩阵M的特征根。它们满足多项式方程：${\lambda }^3+a_1{\lambda }^2+a_2{\lambda }^1+a_3=0$ ，其中$a_1=M_1+M_2$ ，

$a_2=\mu \left(M_1+M_2\right)+\frac{\Lambda \sigma \beta }{M_2}\left(1-\frac{1}{R_{01}}\right)$ ，$a_3=\Lambda \sigma \beta \left(1-\frac{1}{R_{01}}\right)$ ，根据$M_i\left(i=1,2\right)$ 的定义，知当$R_{01}>1$ ，$a_1>0$ ，$a_3>0$ ，且

$a_1{a}_2-a_3=\mu {\left(M_1+M_2\right)}^2+\frac{M_1\beta \sigma \Lambda }{M_2}\left(1-\frac{1}{R_{01}}\right)>0$

由Routh-Hurwitz准则可知当$R_{01}>1$ 时，${\lambda }_{8,9,10}$ 具有有负实部，又当$R_{02}<1$ 时，故而${\lambda }_5=M_6\left(R_{02}-1\right)<0$ ，则矩阵$J_0$ 的特征根均具有负实部，所以系统(2.1)的无动物地方病平衡点$E_1$ 在$R_{01}>1,R_{02}<1$ 时局部渐近稳定。而当$R_{01}>1,R_{02}>1$ 时，$J_0$ 特征根至少有一个具有正实部，此时$E_1$ 不稳定。

定理3.3.3 当$R_{01}>1,R_{02}<1$ 时，系统(2.1)的无动物地方病平衡点$E_1$ 是全局渐近稳定。

证明：在无动物地方病平衡点$E_1$ 处，$R_{02}<1$ 对于动物处于无病状态，由定理3.2.2可知，此时动物子

系统(3.3)是全局稳定的，即$\underset{t\to +\infty }{\lim }{I}_A\left(t\right)=0$ ，此时人类子系统对应的极限系统化为(3.5).对系统(3.5)，依据文献[16]，构造Lyapunov函数：

$L\left(x\right)=\left(S-S^{*}-S^{*}\ln \frac{S}{S^{*}}\right)+a_2\left(E-E^{*}-E^{*}\ln \frac{E}{E^{*}}\right) +a_3\left(I-I^{*}-I^{*}\ln \frac{I}{I^{*}}\right)+\left(V-V^{*}-V^{*}\ln \frac{V}{V^{*}}\right)$

其中$a_2=1$ ，$a_3=\frac{S^{*}\beta }{M_2}$ 。则函数$L\left(x\right)$ 关于模型(2.1)的全导数为：

$\begin{array}{c}{L^{\prime }\left(x\right)|}_{\left(2.1\right)}=C-\left(\beta I+\mu +\delta \right)S+\epsilon V-\frac{S^{*}\Lambda }{S}+a_2\beta S^{*}I-\frac{S^{*}\epsilon \Lambda }{S}+a_2\beta IS-a_2{M}_{1}E\\ -\frac{E^{*}\beta SI}{E}+\frac{S^{*}\beta }{M_2}\sigma E-a_3{M}_2-a_3\frac{I^{*}\sigma E}{I}+\delta S-\epsilon V-\mu V-\frac{V^{*}\delta S}{V}\\ \triangleq G\left(x\right)\end{array}$

其中，$C=\text{Λ}+\mu S^{*}+\delta S^{*}+E^{*}{M}_1+a_3{I}^{*}{M}_2+V^{*}\left(\epsilon +\mu \right)$ 。

取$\frac{S}{S^{*}}=x$ ，$\frac{E}{E^{*}}=y$ ，$\frac{I}{I^{*}}=z$ ，$\frac{Q}{Q^{*}}=u$ ，$\frac{V}{V^{*}}=n$ ，则有：

$\begin{array}{c}G\left(x\right)=C+\left(a_2\beta -\beta \right)S^{*}{I}^{*}zx-\mu S^{*}-\frac{\delta S^{*}x}{n}-a_2\frac{\beta S^{*}{I}^{*}zx}{y}+\left(-a_2{M}_1+a_3\sigma \right)E^{*}y\\ -a_3\frac{I^{*}\sigma E^{*}y}{z}+\left(-a_3{M}_2+\beta S^{*}\right)I^{*}z-\mu V^{*}n-\frac{V^{*}\Lambda \epsilon }{x}-\frac{\Lambda }{x}\end{array}$

化简得：$G\left(x\right)=\mu S^{*}\left(2-x-\frac{1}{x}\right)+c{V}^{*}\left(2-\frac{n}{x}-\frac{x}{n}\right)+\beta I^{*}{S}^{*}\left(3-\frac{1}{x}-\frac{y}{z}-\frac{xz}{y}\right)+V^{*}\mu \left(3-n-\frac{x}{n}-\frac{1}{x}\right)$

易知$G\left(x\right)=0$ 当且仅当$x=y=z=n=1$ ，此时$S=S^{*}$ ，$E=E^{*}$ ，$I=I^{*}$ ，$V=V^{*}$ ，记

$D=\left\{x|{L^{\prime }\left(x\right)|}_{\left(2.1\right)}\equiv 0,x\in {\Omega }_1\right\}$ ，

则系统(3.5)的最大不变集K满足$K=\left\{\left(S^{*},E^{*},I^{*},Q^{*},R^{*},V^{*}\right)\right\}\subset D$ ，故由LaSalle不变性原理知：

$\underset{t\to +\infty }{\lim }S\left(t\right)=S^{*}$ ，$\underset{t\to +\infty }{\lim }E\left(t\right)=E^{*}$ ，$\underset{t\to +\infty }{\lim }I\left(t\right)=I^{*}$ ，$\underset{t\to +\infty }{\lim }V\left(t\right)=V^{*}$ ，$\underset{t\to +\infty }{\lim }R\left(t\right)=R^{*}$ ，$\underset{t\to +\infty }{\text{lim}}Q=Q^{*}$

由极限系统理论知，当$R_{02}<1,R_{01}>1$ 时，原系统(2.1)无动物地方病平衡点$E_1$ 是全局渐近稳定的。

3.4. 共存状态下地方病平衡点的存在性及稳定性

除了无病平衡点和无动物地方病平衡点外，系统(2.1)还可能存在一个共存地方病平衡点$E_2$ ，此时疾病在人类和动物中持续存在。

定理3.4.1 当$R_{02}>1$ 时，系统(2.1)具有唯一的地方平衡点$E_2$ 。

证明：由系统可解的系统的另一个地方病$E_2\left(S^{**},E^{**},I^{**},Q^{**},R^{**},V^{**},S_A^{**},E_A^{**},I_A^{**},R_A^{**}\right)$ ，其中

$S^{**}=\frac{\Lambda }{\beta I^{**}+\xi I_A^{**}+\mu +\delta -\frac{\epsilon \delta }{\epsilon +\mu }}$ ，$E^{\text{**}}=\frac{M_2}{\sigma }{I}^{**}$ ，$Q^{**}=\frac{\phi \sigma +q{M}_2}{\sigma M_3}{I}^{**}$ ，$R^{**}=\frac{\phi \sigma \tau +q{M}_2\tau +\sigma M_3{M}_4}{\sigma M_3\mu }{I}^{**}$

$V^{**}=\frac{\Lambda \delta }{\left(\epsilon +\mu \right)\left(\beta I^{**}+\xi I_A^{**}+\mu \right)}$ ，$S_A^{**}=\frac{M_5{M}_6}{\theta {\beta }_1}$ ，$E_A^{**}=\frac{{\Lambda }_1}{M_5}\left(1-\frac{1}{R_{02}}\right)$ ，$I_A^{**}=\frac{{\Lambda }_1\theta }{M_5{M}_6}\left(1-\frac{1}{R_{02}}\right)$ ，

$R_A^{**}=\frac{{\Lambda }_1\theta }{M_5{M}_6\omega }\left(1-\frac{1}{R_{02}}\right)$ ，

其中，$I^{\text{**}}$ 满足

$F\left(I^{**}\right)=I^{**}{}^2{M}_1{M}_2\beta +I^{**}\left[\left(\xi I_A^{**}+\mu +\delta -\frac{\epsilon \delta }{\epsilon +\mu }\right)M_1{M}_2-\sigma \beta \Lambda \right]-\sigma \xi \Lambda I_A^{**}=0$ ，

考虑到$F\left(0\right)=-\sigma \beta \Lambda T{M}_1{M}_2<0$ 且$\Delta ={\left[\left(\xi I_A^{**}+\mu +\delta -\frac{\epsilon \delta }{\epsilon +\mu }\right)M_1{M}_2-\sigma \beta \Lambda \right]}^2+4M_1{M}_2\beta \sigma \xi \Lambda I_A^{**}>0$ 则由

韦达定理知，$F\left(I^{\text{**}}\right)=0$ 存在唯一正解。故当$R_{02}>1$ 时，$E_2$ 存在且唯一。

定理3.4.2 当$R_{02}>1$ 时，系统(2.1)的地方病平衡点$E_2$ 是局部渐近稳定的。

证明：由定理3.4.1可知，对任意的$R_{01}$ ，当且仅当$R_{02}>1$ 时系统(2.1)具有地方平衡点。系统(2.1)在$E_2$ 处的雅可比矩阵为：

$J_2={\left[\begin{array}{cc}{A}_{6\times 6}& B_{4\times 6}\\ 0_{6\times 4}& D_{4\times 4}\end{array}\right]}_{10\times 10}$ ，$A={\left[\begin{array}{cccccc}{A}_{11}& 0& -\beta S^{*}& 0& 0& c\\ \beta I^{**}+\xi I_A^{**}& -M_1& \beta S^{*}& 0& 0& 0\\ 0& \sigma & -M_2& 0& 0& 0\\ 0& q& \phi & -M_3& 0& 0\\ 0& 0& M_4& \tau & -\mu & 0\\ \delta & 0& 0& 0& 0& -\epsilon -\mu \end{array}\right]}_{6\times 6}$

$B={\left[\begin{array}{cccc}0& 0& -\xi S^{**}& 0\\ 0& 0& \xi S^{**}& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right]}_{6\times 4}$ ，$D={\left[\begin{array}{cccc}-{\beta }_1{I}_A^{**}-{\mu }_1& 0& -S_A^{**}{\beta }_1& 0\\ {\beta }_1{I}_A^{**}& -M_5& S_A^{**}{\beta }_1& 0\\ 0& \theta & -M_6& 0\\ 0& 0& \omega & -{\mu }_1\end{array}\right]}_{4\times 4}$

其中$A_{11}=-\beta I^{**}-\xi I_A^{**}-\mu -\delta$ 。矩阵$J_2$ 的特征值分别为：${\lambda }_1=-{\mu }_1$ ，${\lambda }_2=-\mu$ ，${\lambda }_3=-M_3$ ，${\lambda }_4=-M_6$ ，

${\lambda }_5=-M_5$ ，${\lambda }_6=-\beta I_A^{**}$ ，${\lambda }_7=-M_2$ ，${\lambda }_8=-\mu$ ，${\lambda }_9=-\epsilon -\mu$ ，${\lambda }_{10} =-\beta I^{**}-\xi I_A^{**}-\mu -\frac{\beta S^{**}\sigma \left(\mu +\delta \right)}{M_1{M}_2}-\frac{\mu \delta }{\epsilon +\mu }$ 。

对任意的$R_{01}$ ，当且仅当$R_{02}>1$ 时，即$E_2$ 存在。对任意的$R_{01}$ ，$E_2$ 对应的雅可比行列式的特征值均具有负实部，故在$R_{02}>1$ 时存在，对任意$R_{01}$ ，$E_2$ 均局部稳定。

定理3.4.3 当且仅当$R_{02}>1$ 时，系统(2.1)共存状态下的地方病平衡点$E_2$ 当$R_{01}>1$ 时全局渐近稳定。

证明：注意到动物子系统(3.3)是独立的，且当$R_{02}>1$ 时，系统(3.3)存在正平衡点$E_2^1=\left(S_A^{**},E_A^{**},I_A^{**},R_A^{**}\right)$

构造Lyapunov函数：

$L_A\left(x\right)=\left(S_A-S_A^{**}-S_A^{**}\ln \frac{S_A}{S_A^{**}}\right)+\left(E_A-E_A^{**}-E_A^{**}\ln \frac{E_A}{E_A^{**}}\right)+\frac{M_5}{\theta }\left(I_A-I_A^{**}-I_A^{**}\ln \frac{I_A}{I_A^{**}}\right)$

沿着系统关于$L_A\left(x\right)$ 求全导数为：

${{L^{\prime }}_A\left(x\right)|}_{\left(3.3\right)}=C_A-{\mu }_1{S}_A+\left(-{\beta }_1+\left(1-\frac{E_A^{**}}{E_A}\right){\beta }_1\right)S_A{I}_A+\left(S_A^{**}{\beta }_1-\frac{M_5{M}_6}{\theta }\right)I_A-\frac{I_A^{**}{M}_5{E}_A^{**}}{I_A}-\frac{{\Lambda }_1{S}_A^{**}}{S_A}=A\left(x\right)\le 0$

其中$C_A=S_A^{**}{\mu }_1+{\Lambda }_1+\frac{M_5{M}_6{I}_A^{**}}{\theta }+M_5{E}_A^{**}$

取$\frac{S_A}{S_A^{**}}=x_A$ ，$\frac{E_A}{E_A^{**}}=y_A$ ，$\frac{I_A}{I_A^{**}}=z_A$ ，则有：

$\begin{array}{c}A\left(x\right)=C_A-{\mu }_1{S}_A^{**}{x}_A+\left(-{\beta }_1+\left(1-\frac{1}{y_A}\right){\beta }_1\right)S_A^{**}{I}_A^{**}{x}_A{z}_A+\left(S_A^{**}{\beta }_1-\frac{M_5{M}_6}{\theta }\right)I_A^{**}{z}_A-\frac{y_A{M}_5{E}_A^{**}}{z_A}-\frac{{\Lambda }_1}{x_A}\\ ={\mu }_1{S}_A^{**}\left(2-\frac{S_A^{**}}{S_A}-\frac{S_A}{S_A^{**}}\right)+S_A^{**}{\beta }_1{I}_A^{**}\left(3-\frac{S_A^{**}}{S_A}-\frac{S_A{I}_A{E}_A^{**}}{S_A^{**}{I}_A^{**}{E}_A}-\frac{I_A^{**}{E}_A}{I_A{E}_A^{**}}\right)\end{array}$

易知$A\left(x\right)=0$ 当且仅当$x_A=y_A=z_A=1$ ，此时，$S_A=S_A^{**}$ ，$E_A=E_A^{**}$ ，$I_A=I_A^{**}$ 。记

$D_A=\left\{x|{L^{\prime }\left(x\right)|}_{\left(3.3\right)}\equiv 0,x\in {\Omega }_2\right\}$ ，

则最大不变集K为$K=\left\{E_2^1\right\}\subset D_A$ 。故由LaSalle不变性原理知：

$\underset{t\to +\infty }{\lim }{S}_A\left(t\right)=S_A^{**}$ ，$\underset{t\to +\infty }{\text{lim}}{E}_A\left(t\right)=E_A^{**}$ ，$\underset{t\to +\infty }{\text{lim}}{I}_A\left(t\right)=I_A^{**}$ ，$\underset{t\to +\infty }{\text{lim}}{R}_A\left(t\right)=R_A^{**}$

即当$R_{02}>1$ 时，系统(3.3)的平衡点$E_2^1$ 是全局吸引的。

取$\xi I_A^{**}=\pi$ ，由于$\underset{t\to +\infty }{\text{lim}}{I}_A\left(t\right)=I_A^{**}$ ，则人类子系统对应的极限系统

$\{\begin{array}{l}{S}^{\prime }\left(t\right)=\Lambda -\beta SI-\mu S+\epsilon V-\delta S-\pi S\\ E^{\prime }\left(t\right)=\beta SI+\pi S-M_{1}E\\ I^{\prime }\left(t\right)=\sigma E-M_{2}I\\ Q^{\prime }\left(t\right)=\phi I-M_{3}Q+qE\\ V^{\prime }\left(t\right)=\delta S-\epsilon V-\mu V\end{array}$ (3.6)

易知，系统(3.6)在$R_{01}>1$ 时存在正平衡点$E_2^2=\left(S^{**},E^{**},I^{**},Q^{**},R^{**},V^{**}\right)$ 。

构造Lyapunov函数：

$L_H\left(x\right)=\left(S-S^{**}-S^{**}\ln \frac{S}{S^{**}}\right)+\left(E-E^{**}-E^{**}\ln \frac{E}{E^{**}}\right)+\frac{S^{**}\beta }{M_2}\left(I-I^{**}-I^{**}\ln \frac{I}{I^{**}}\right)+\left(V-V^{**}-V^{**}\ln \frac{V}{V^{**}}\right)$

则$L_H\left(x\right)$ 关于系统(3.6)的全导数为：

$\begin{array}{c}{{L^{\prime }}_H\left(x\right)|}_{\left(3.6\right)}=C_H-\mu S-\frac{V^{**}\delta S}{V}-\frac{E^{**}\beta SI}{E}-\frac{E^{**}\pi S}{E}+\left(-M_1+\frac{S^{**}\beta \sigma }{M_2}\right)E\\ -\frac{S^{**}{I}^{**}\beta \sigma E}{M_{2}I}-\mu V-\frac{S^{**}\epsilon V}{S}-\frac{S^{**}\Lambda }{S}\\ =H\left(x\right)\end{array}$

其中$C_H=\Lambda +M_1{E}^{**}+S^{**}{I}^{**}\beta +\left(\epsilon +\mu \right)V^{**}+\left(\delta +\mu +\pi \right)S^{**}$ 。

令$\frac{S}{S^{**}}=x_H$ ，$\frac{E}{E^{**}}=y_H$ ，$\frac{I}{I^{**}}=z_H$ ，$\frac{V}{V^{**}}=u_H$ ，则

$\begin{array}{c}H\left(x\right)=C_H-\mu S^{**}{x}_H-\frac{\delta S^{**}{x}_H}{u_H}-\frac{\beta S^{**}{I}^{**}{x}_H{z}_H}{y_H}-\frac{\pi S^{**}{x}_H}{y_H}\\ +\left(-M_1+\frac{S^{**}\beta \sigma }{M_2}\right)E^{**}{y}_H-\frac{S^{**}\beta \sigma E^{**}{y}_H}{M_2{z}_H}-\mu V^{**}{u}_H-\frac{\epsilon V^{**}{u}_H}{x_H}-\frac{\Lambda }{x_H}\\ = \mu S^{**}\left(2-x_H-\frac{1}{x_H}\right)+c{V}^{**}\left(2-\frac{x_H}{u_H}-\frac{u_H}{x_H}\right)+\pi S^{**}\left(3-y_H-\frac{x_H}{y_H}-\frac{1}{x_H}\right)\\ +\mu V^{**}\left(3-\frac{u_H}{x_H}-u_H-\frac{1}{x_H}\right)+\beta S^{**}{I}^{**}\left(3-\frac{z_H{x}_H}{y_H}-\frac{y_H}{z_H}-\frac{1}{x_H}\right)\end{array}$

易知$H\left(x\right)=0$ 当且仅当$x_H=y_H=z_H=u_H=1$ ，此时$S=S^{**}$ ，$E=E^{**}$ ，$I=I^{**}$ ，$V=V^{**}$ ，进而系统(3.6)的最大不变集合$K_H$ 为此时有：

$K_H=\left\{E_2^2\right\}\subset D_H=\left\{x|{{L^{\prime }}_H\left(x\right)|}_{\left(3.6\right)}\equiv 0,x\in {\Omega }_1\right\}$ ，

根据LaSalle不变集原理得到

$\underset{t\to +\infty }{\text{lim}}S=S^{**}$ ，$\underset{t\to +\infty }{\text{lim}}E=E^{**}$ ，$\underset{t\to +\infty }{\text{lim}}I=I^{**}$ ，$\underset{t\to +\infty }{\text{lim}}V=V^{**}$ ，$\underset{t\to +\infty }{\text{lim}}Q=Q^{**}$ ，$\underset{t\to +\infty }{\text{lim}}R=R^{**}$

当$R_{02}>1$ 时，系统(3.6)的平衡点$E_2^2$ 是全局吸引的。

根据定理3.3.2和极限系统理论[17]知，当$R_{02}>1$ 时，系统(2.1)的地方病平衡点$E_2$ 当$R_{01}>1$ 时是全局渐近稳定的。

4. 数值模拟

为验证理论分析结果的正确性，本节将运用数值模拟方法进行实证验证。表2中对应的参数满足$R_{01}\approx 0.3323<1$ ，$R_{02}\approx 0.2176<1$ ，对应的解曲线为图2。若将表2中参数$\beta$ 取0.022，${\beta }_1$ 取0.012。其他参数值不变，则此时$R_{01}\approx 3.3226>1$ ，$R_{02}\approx 2.176>1$ ，对应的解曲线如图2所示。

Table 2. The parameter values of SEIQVR model

表2. $R_{01}<1,R_{02}<1$ 参数取值

参数	$\Lambda$	$\beta$	$\mu$	$\delta$	$\epsilon$	$\sigma$	$\xi$	d	$r_1+r_2$	$\theta$
取值	0.5	0.0022	0.003	0.12	0.1	0.2	0.5	0.2	0.83	0.64
参数	${\beta }_1$	${\mu }_1$	$\phi$	${\text{Λ}}_1$	$d_1$	q	$\omega$
取值	0.0012	0.002	0.85	0.31	0.2	0.75	0.7

在图2~图4中，图中的纵坐标表示模型中各仓室量、动物数量，横坐标表示时间。若取模型的参数如表2所示，对应的基本再生数$R_{01}\approx 0.3323<1$ ，$R_{02}\approx 0.2176<1$ ，我们可以得到图2，可以发现随着时间的演化，染病动物数量及类数量终近于零，疾病将不会发展称为地方病而最终消亡。这验证了定理3.2.1的结论。

Figure 2. Solution curves of the model at $R_{01}\approx 0.3323<1,R_{02}\approx 0.2176<1$

图2. $R_{01}\approx 0.3323<1,R_{02}\approx 0.2176<1$ 时模型的解曲线

若替换表2中的参数$\beta$ = 022，${\beta }_1$ = 0012，其余参数不变，此时$R_{01}\approx 3.3226>1$ ，$R_{02}\approx 0.2176<1$ ，系统存在无动物平衡点，且随着时间的演化，动物种群疾病最终消亡，而人类种群疾病持续存在，患病人类数量最终趋于无动物地方病平衡点，见图3这与定理3.3.1的结论一致。

Figure 3. Solution curves of the model at $R_{01}\approx 0.3323>1,R_{02}\approx 0.2176<1$

图3. $R_{01}\approx 0.3323>1,R_{02}\approx 0.2176<1$ 时模型的解曲线

若将表2中的参数$\beta$ 替换为0.022，${\beta }_1$ 替换为0.012，其余参数不变，计算可得$R_{01}\approx 3.3226>1$ ，$R_{02}\approx 2.176>1$ ，模拟得到图4此时疾病在动物及人类种群持续存在而形成了地方病，染病动物和人类数量最终趋于共存的地方病平衡点，这与定理3.4.1的结论一致。

Figure 4. Solution curves of the model at $R_{01}\approx 0.3323>1,R_{02}\approx 0.2176>1$

图4. $R_{01}\approx 0.3323>1,R_{02}\approx 0.2176>1$ 时模型的解曲线

5. 结论

为了分析猴痘病毒的传播特征，我们构建了一个包含六个人类仓室物仓室的耦合传染病动力学模型。通过计算二代再生矩阵，我们获得了模型的人类种群基本再生数($R_{01}$ )和动物种群基本再生数($R_{02}$ )，进而确定了耦合系统的基本再生数$R_0=\max \left\{R_{01},R_{02}\right\}$ ，验证了阈值定理和地方病平衡点的存在性。在此基础上，利用李雅普诺夫稳定性定理，我们分析了在特定条件下无动物地方病平衡点和共存地方病平衡点的稳定性，具体结果如表3所示。此外，通过数值模拟验证了我们所获得结论的合理性与有效性。这些结果为理解猴痘病毒在不同宿主之间的传播机制提供了重要的理论依据。然而，本模型尚未充分考虑到疫情后期可能出现的疾病复发、季节性感染以及合并感染等因素的影响，这在一定程度上限制了其应用的广泛性。因此，未来的研究应综合考虑这些关键因素，以对猴痘模型进行更为深入的探讨与分析。

Table 3. Balance stenced stability conditions

表3. 平衡点存在性及稳定性条件

平衡点	存在条件	局部稳定	全局稳定
无病平衡点	一直存在	$R_0<1$	$R_0<1$
无动物地方病平衡点	$R_{01}>1$	$R_{01}>1,R_{02}<1$	$R_{01}>1,R_{02}<1$
共存地方病平衡点	$R_{02}>1$	$R_{02}>1,\forall R_{01}$	$R_{02}>1,R_{01}>1$

基金项目

重庆市教委科学技术研究项目(KJQN201801136)；重庆理工大学研究生教育高质量发展行动计划资助成果(gzlcx20243274)。

NOTES

^*通讯作者。