一类斑块环境下口蹄疫传染病模型的动力学分析
Dynamical Analysis of a Foot-and-Mouth Disease Epidemic Model in Patchy Environment
DOI: 10.12677/aam.2025.148366, PDF, HTML, XML,    科研立项经费支持
作者: 盛俊杰:天津职业技术师范大学理学院,天津
关键词: 口蹄疫病毒斑块环境基本再生数李雅普诺夫函数 Foot-and-Mouth Disease Virus Patchy Environment Basic Reproduction Number Lyapunov Function
摘要: 针对口蹄疫病毒在斑块环境中的传播现象,建立了一个斑块口蹄疫传染病模型。分别就单斑块和两个斑块中病毒动力学行为展开研究。通过下一代矩阵理论,定义了模型的基本再生数,分析了单斑块模型无病平衡点及地方病平衡点的存在性和局部稳定性;通过构造合适的李雅普诺夫函数证明了单斑块模型平衡点的全局渐进稳定性。在两个对称的斑块环境中,通过李亚普诺夫函数法进一步确定了多斑块模型无病平衡点和地方病平衡点的全局渐进稳定性。
Abstract: A patchy foot-and-mouth disease epidemic model was established to study the transmission dynamics of the virus in patch environment. The viral transmission behaviors were investigated in both single-patch and two-patch systems. Using the next-generation matrix theory, the basic reproduction number of the model was defined, and the existence as well as local stability of the disease-free equilibrium and endemic equilibrium in the single-patch model were analyzed. By constructing an appropriate Lyapunov function, the global asymptotic stability of the equilibria in the single-patch model was demonstrated. In a symmetric two-patch environment, the global asymptotic stability of the disease-free equilibrium and endemic equilibrium in the multi-patch model was further confirmed through the Lyapunov function method.
文章引用:盛俊杰. 一类斑块环境下口蹄疫传染病模型的动力学分析[J]. 应用数学进展, 2025, 14(8): 1-16. https://doi.org/10.12677/aam.2025.148366

1. 引言

口蹄疫(Foot-and-Mouth Disease, FMD),是由一种名为口蹄疫病毒的小核糖核酸科病毒引起的,主要包括了SAT (1-3),O,A,C和Asial七种血清型[1]。对于口蹄疫病毒在连续空间上的传播,国内外学者已做了一定的研究。Wang等[2] [3]研究了延迟和非局部的口蹄疫模型的空间传播。Yang等[4]对空间异质环境下带年龄结构的非局部扩散口蹄疫模型及其动力学做了详尽的研究。Shuai等[5]主要研究的是霍乱病毒在异质宿主种群中的传播现象,其中异质宿主种群由共享共同水源的几个斑块的同质宿主种群组成。本文就口蹄疫病毒在离散斑块环境下的传播机理,建立了一个多斑块口蹄疫传染病模型并分析了模型的长时间动力学行为。利用Van den Driessche P等[6]和Diekmann O等[7]中研究传染病模型的下一代矩阵理论,定义并计算了疾病传播的基本再生数,分析了单斑块模型无病平衡点和地方病平衡点的存在和局部稳定性,并通过构造合适的李雅普诺夫函数,证明了模型平衡点的全局渐进稳定性。针对两个斑块模型,在对称的假设下进一步利用李雅普诺夫法证明了模型无病平衡点和地方病平衡点的全局吸引性。

2. 模型建立

对于两个斑块内口蹄疫病毒的传播方式,我们建立了如下的模型进行研究。其中 S i , I i 分别表示在第 i 个斑块内易感动物和患病动物的种群密度,而 V i 则表示的是在第 i 个斑块内口蹄疫病毒的密度。

{ d S 1 dt = Λ 1 μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1 d I 1 dt = μ 1 I 1 + β 11 S 1 I 1 + β 12 S 1 V 1 d V 1 dt = r 1 I 1 k 1 V 1 d 21 V 1 + d 12 V 2 d S 2 dt = Λ 2 μ 2 S 2 β 21 S 2 I 2 β 22 S 2 V 2 d I 2 dt = μ 2 I 2 + β 21 S 2 I 2 + β 22 S 2 V 2 d V 2 dt = r 2 I 2 k 2 V 2 + d 21 V 1 d 12 V 2 (1)

其他字母表示如下表:

字母

生物学意义

Λ i

i个斑块内动物的输入率

μ i

i个斑块内动物的自然死亡率

r i

i个斑块由患病动物释放病毒率

k i

i个斑块内病毒的自然降解率

β i1

i个斑块内由易感动物转变为患病动物的感染率

β i2

i个斑块内通过病毒转变为患病动物的感染率

d 21

病毒从第一个斑块到第二个斑块的转移率

d 12

病毒从第二个斑块到第一个斑块的转移率

3. 一个斑块中口蹄疫病毒的动力学行为

首先,我们考虑如下单斑块模型的动力学行为

{ d S 1 dt = Λ 1 μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1   d I 1 dt = μ 1 I 1 + β 11 S 1 I 1 + β 12 S 1 V 1 d V 1 dt = r 1 I 1  k 1 V 1 (2)

3.1. 基本再生数的定义

在模型(2)中,令 I 1 = V 1 =0 得系统的无病平衡点 E 0 =( Λ 1 μ 1 ,0,0 ) ,将模型在 E 0 处线性化得下面线性化方程组

{ d S 1 dt = Λ 1 μ 1 Λ 1 μ 1 + β 11 Λ 1 μ 1 I 1 + β 12 Λ 1 μ 1 V 1 d I 1 dt = μ 1 I 1 + β 11 Λ 1 μ 1 I 1 + β 12 Λ 1 μ 1 V 1 d V 1 dt = r 1 I 1  k 1 V 1 (3)

注意系统(3)中变量 I 1 V 1 的方程不依赖于变量 S 1 ,故可利用 I 1 V 1 的方程来定义模型的基本再生数。运用文献中的下一代矩阵理论[7],得其相关仓室为

F=( β 11 Λ 1 μ 1 β 12 Λ 1 μ 1 0 0 ), V=( μ 1 0 r 1 k 1 )

因此,模型(2)的基本再生数 R 10 可由 F V 1 的最大特征值 ρ( F V 1 ) 来表示,经计算得

R 10 = Λ 1 ( β 11 k 1 + β 12 r 1 ) u 1 2 k 1

类似第二个斑块上的基本再生数为

R 20 = Λ 2 ( β 21 k 2 + β 22 r 2 ) u 2 2 k 2

3.2. 地方性平衡点的计算

为了计算得到地方性平衡点 E 1 * ,令模型(2)的右边都等于0,即

{ Λ 1 μ 1 S 1 * β 11 S 1 * I 1 * β 12 S 1 * V 1 * =0 μ 1 I 1 * + β 11 S 1 * I 1 * + β 12 S 1 * V 1 * =0 r 1 I 1 * k 1 V 1 * =0 (4)

由此计算得出 E 1 * =( S 1 * , I 1 * , V 1 * ) ,其中,

S 1 * = μ 1 k 1 β 11 k 1 + β 12 r 1 ,  I 1 * = Λ 1 ( β 11 k 1 + β 12 r 1 ) μ 1 2 k 1 u 1 ( β 11 k 1 + β 12 r 1 ) ,  V 1 * = r 1 [ Λ 1 ( β 11 k 1 + β 12 r 1 ) μ 1 2 k 1 ] μ 1 k 1 ( β 11 k 1 + β 12 r 1 ) .

显然, S 1 * >0, I 1 * >0, V 1 * >0 当且仅当 R 10 >1 ,即地方病平衡点 E 1 * 存在当且仅当疾病传播的基本再生数大于1。

3.3. 平衡点的稳定性分析

对于平衡点的稳定性分析,包括无病平衡点的稳定性分析和地方病平衡点的稳定性分析。由文献[6]可以得到下面的定理是成立的。

定理1:当 R 10 <1 时,无病平衡点 E 0 是局部稳定的;当 R 10 >1 时,无病平衡点 E 0 是不稳定的。

定理2 R 10 >1 时,地方性平衡点 E 1 * 是局部稳定的。

借鉴文献[3]和文献[4]中的研究思路,下面利用李亚普诺夫函数法,验证无病平衡点 E 0 和地方病平衡点 E 1 * 的全局渐进稳定性。

定理3 R 10 <1 时,无病平衡点 E 0 是全局渐进稳定的。

证明:构造李雅普诺夫函数

L E 0 ( t )= L S 0 ( t )+ L I 0 ( t )+ L V 0 ( t ) L S 0 ( t )= S 0 g( S 1 S 0 ),  L I 0 ( t )= I 1 ( t ),  L V 0 ( t )= β 12 S 0 k 1 V 1 ( t ) (5)

其中 g( x )=x1lnx0 对所有 x>0 成立,则

d L S 0 dt = d S 1 dt ( 1 S 0 S 1 ),  d L I 0 dt = d I 1 ( t ) dt ,  d L V 0 dt = β 12 S 0 k 1 d V 1 ( t ) dt . (6)

所以当 R 10 <1 时有下列不等式成立

d L E 0 dt = d L S 0 dt + d L I 0 dt + d L V 0 dt = Λ 1 μ 1 S 1 μ 1 I 1 + S 0 ( β 12 r+ β 11 k 1 k 1 I 1 Λ 1 S 1 + μ 1 ) =2 Λ 1 μ 1 S 1 Λ 1 2 μ 1 S 1 + Λ 1 ( β 12 r+ β 11 k 1 ) μ 1 2 k 1 μ 1 k 1 I 1 =2 Λ 1 μ 1 S 1 Λ 1 2 μ 1 S 1 + μ 1 ( R 0 1 ) I 1 0 (7)

f 1 ( S 1 )= μ 1 S 1 + Λ 1 2 μ 1 S 1 ,可以得到当且仅当 S 1 = Λ 1 μ 1 时, f( S 1 ) 的最小值为 2 Λ 1 。从而当且仅当 S 1 = S 0 , I 1 =0, V 1 =0 时, d L E 0 dt =0 成立,故有李雅普诺夫函数法和定理1知无病平衡点 E 0 是全局渐进稳定的。

定理4:当 R 10 >1 时,地方病平衡点 E 1 * 是全局渐进稳定的。

证明:构造李雅普诺夫函数

L E ( t )= L S ( t )+ L I ( t )+ L V ( t ) L S ( t )= S 1 * g( S 1 S 1 * ), L I ( t )= I 1 * g( I 1 I 1 * ), L V ( t )= β 12 S 1 * V 1 * k 1 g( V 1 V 1 * ) (8)

d L S dt = d S 1 dt ( 1 S 1 * S 1 )=( Λ 1 μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1 )( 1 S 1 * S 1 ), d L I dt = d I 1 dt ( 1 I 1 * I 1 )=( μ 1 I 1 + β 11 S 1 I 1 + β 12 S 1 V 1 )( 1 I 1 * I 1 ), d L V dt = β 12 S 1 * k 1 d V 1 dt ( 1 V 1 * V 1 )= β 12 S 1 * k 1 ( r 1 I 1 k 1 V 1 )( 1 V 1 * V 1 ). (9)

由(4)得到

Λ 1 = μ 1 S 1 * + β 11 S 1 * I 1 * + β 12 S 1 * V 1 * μ 1 = β 11 S 1 * + β 12 S 1 * V 1 * I 1 * r 1 k 1 = V 1 * I 1 * (10)

分别带入(11)的 d L S dt , d L I dt , d L V dt 中,得到

d L S dt = d S 1 dt ( 1 S 1 * S 1 )=( Λ 1 μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1 )( 1 S 1 * S 1 ) =( μ 1 S 1 * + β 11 S 1 * I 1 * + β 12 S 1 * V 1 * μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1   )( 1 S 1 * S 1 ) =( 1 S 1 * S 1 )[ μ 1 ( S 1 S 1 * )+ β 11 S 1 * I 1 * β 11 S 1 I 1 + β 12 S 1 * V 1 * β 12 S 1 V 1 ] = μ 1 S 1 ( S 1 S 1 * ) 2 + β 11 S 1 * I 1 * ( 1 S 1 * S 1 )( 1 S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 S 1 * S 1 )( 1 S 1 V 1 S 1 * V 1 * ) = μ 1 S 1 ( S 1 S 1 * ) 2 + β 11 S 1 * I 1 * ( 1 S 1 * S 1 + I 1 I 1 * S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 S 1 * S 1 + V 1 V 1 * S 1 V 1 S 1 * V 1 * ).

d L I dt = d I 1 dt ( 1 I 1 * I 1 )=( μ 1 I 1 + β 11 S 1 I 1 + β 12 S 1 V 1 )( 1 I I 1 ) =( β 11 S 1 I 1 + β 12 S 1 V 1 β 11 S 1 * I 1 β 12 S 1 * V 1 * I 1 I 1 * )( 1 I I 1 ) = β 11 S 1 * I 1 * ( 1 I I 1 )( S 1 I 1 S 1 * I 1 * I 1 I 1 * )+ β 12 S 1 * V 1 * ( 1 I I 1 )( S 1 V 1 S 1 * V 1 * I 1 I 1 * ) = β 11 S 1 * I 1 * ( 1 S 1 S 1 * I 1 I 1 * S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 I 1 I 1 * S 1 V 1 I 1 * S 1 * V 1 * I 1 + S 1 V 1 S 1 * V 1 * ).

d L V dt = β 12 S 1 * k 1 d V 1 dt ( 1 V 1 * V 1 )= β 12 S 1 * k 1 ( r 1 I 1 k 1 V 1 )( 1 V 1 * V 1 ) = β 12 S 1 * ( 1 V 1 * V 1 )( r 1 I 1 k 1 V 1 )= β 12 S 1 * ( 1 V 1 * V 1 )( V 1 * I 1 I 1 * V 1 ) = β 12 S 1 * V 1 * ( I 1 I 1 * V 1 V 1 * V 1 * I 1 V 1 I 1 * +1 ).

所以当 R 10 >1 时有下列不等式成立

d L E dt = d L S dt + d L I dt + d L V dt = μ 1 S 1 ( S 1 S 1 * ) 2 + β 11 S 1 * I 1 * ( 1 S 1 * S 1 + I 1 I 1 * S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 S 1 * S 1 + V 1 V 1 * S 1 V 1 S 1 * V 1 * ) + β 11 S 1 * I 1 * ( 1 S 1 * S 1 I 1 I 1 * S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 I 1 I 1 * S 1 V 1 I 1 * S 1 * V 1 * I 1 + S 1 V 1 S 1 * V 1 * )+ β 12 S 1 * V 1 * ( I 1 I 1 * V 1 V 1 * V 1 * I 1 V 1 I 1 * +1 ) = μ 1 S 1 ( S 1 S 1 * ) 2 + β 11 S 1 * I 1 * ( 2 S 1 * S 1 S 1 S 1 * )+ β 12 S 1 * V 1 * ( 1 S 1 * S 1 +1 S 1 V 1 I 1 * S 1 * V 1 * I 1 +1 V 1 * I 1 V 1 I 1 * ) = μ 1 S 1 ( S 1 S 1 * ) 2 + β 11 S 1 * I 1 * ( 2 S 1 * S 1 S 1 S 1 * )+ β 12 S 1 * V 1 * [ g( S 1 * S 1 )g( S 1 V 1 I 1 * S 1 * V 1 * I 1 )g( V 1 * I 1 V 1 I 1 * ) ] 0

f 2 ( S 1 )= S 1 * S 1 + S 1 S 1 * ,则当且仅当 S 1 = S 1 * 时, f 2 ( S 1 ) 的最小值为2。由此可以看出,当且仅当 S 1 = S 1 * , I 1 = I 1 * , V 1 = V 1 * 时, d L E dt =0 成立,所以由定理2知地方病平衡点 E 是全局渐进稳定的。

4. 两个斑块中口蹄疫病毒的动力学行为

本节考虑两个斑块口蹄疫模型(1)的动力学行为。显然 E=( Λ 1 μ 1 ,0,0, Λ 2 μ 2 ,0,0 ) 为模型的无病平衡点。另外,在假设两个斑块对称,且 R 10 >1, R 20 >1 的情况下, E * =( S 1 * , I 1 * , V 1 * , S 2 * , I 2 * , V 2 * ) 为模型的一个地方病平衡点。将模型(1)在无病平衡点 E 处线性化,并写出变量 I 1 , V 1 , I 2 , V 2 对应的系数矩阵,得

J ¯ =( μ 1 + β 11 Λ 1 μ 1 β 12 Λ 1 μ 1 0 0 r 1 k 1 d 21 0 d 12 0 0 μ 2 + β 21 Λ 2 μ 2 β 22 Λ 2 μ 2 0 d 21 r 2 k 2 d 12 )

根据下一代矩阵法[7],得到计算两个斑块模型的基本再生数的矩阵为

F ¯ =( β 11 Λ 1 μ 1 β 12 Λ 1 μ 1 0 0 0 0 0 0 0 0 β 21 Λ 2 μ 2 β 22 Λ 2 μ 2 0 0 0 0 ),  V ¯ =( μ 1 0 0 0 r 1 k 1 + d 21 0 d 12 0 0 μ 2 0 0 d 21 r 2 k 2 + d 12 )

F ¯ ( V ¯ ) 1 =( β 11 Λ 1 μ 1 2 + β 12 Λ 1 r 1 μ 1 2 ( k 1 + d 21 ) β 12 Λ 1 μ 1 ( k 1 + d 21 ) β 12 Λ 1 2 d 12 ( k 1 + d 21 )( k 2 + d 12 ) 0 0 0 0 0 β 22 Λ 2 r 1 d 21 μ 1 μ 2 ( k 1 + d 21 )( k 2 + d 12 ) β 22 Λ 2 d 21 μ 2 ( k 1 + d 21 )( k 2 + d 12 ) β 21 Λ 2 μ 2 2 + β 22 Λ 2 r 2 μ 2 2 ( k 2 + d 12 ) β 22 Λ 2 k 2 + d 12 0 0 0 0 )

计算得到两个斑块模型的基本再生数为

R 0 =max{ β 11 Λ 1 μ 1 2 + β 12 Λ 1 r 1 μ 1 2 ( k 1 + d 21 ) , β 21 Λ 2 μ 2 2 + β 22 Λ 2 r 2 μ 2 2 ( k 2 + d 12 ) }

分别记为

R 1 = β 11 Λ 1 μ 1 2 + β 12 Λ 1 r 1 μ 1 2 ( k 1 + d 21 ) R 2 = β 21 Λ 2 μ 2 2 + β 22 Λ 2 r 2 μ 2 2 ( k 2 + d 12 ) .

平衡点的稳定性分析

定理5:当 R 10 <1, R 20 <1 ,且模型参数满 Λ 1 = Λ 2 , μ 1 = μ 2 , β 12 = β 22 , d 12 = d 21 , k 1 = k 2 时,模型(1)的无病平衡点 E 是全局渐进稳定的。

证明:构造李雅普诺夫函数

L E ( t )= L S 0 ( t )+ L I 0 ( t )+ L V 0 ( t )+ L S 0 ( t )+ L I 0 ( t )+ L V 0 ( t ) L S 0 ( t )= S 0 g( S 1 S 0 ), L I 0 ( t )= I 1 ( t ), L V 0 ( t )= β 12 S 0 k 1 V 1 ( t ), L S 0 ( t )= S 0 g( S 2 S 0 ), L I 0 ( t )= I 2 ( t ), L V 0 ( t )= β 22 S 0 k 2 V 2 ( t ). (11)

d L S 0 dt = d S 1 dt ( 1 S 0 S 1 ), d L I 0 dt = d I 1 ( t ) dt , d L V 0 dt = β 12 S 0 k 1 d V 1 ( t ) dt , d L S 0 dt = d S 2 dt ( 1 S 0 S 2 ), d L I 0 dt = d I 2 ( t ) dt , d L V 0 dt = β 22 S 0 k 2 d V 2 ( t ) dt . (12)

所以,

f 2 ( S 1 )= μ 1 S 1 + Λ 1 2 μ 1 S 1 ,则当且仅当 S 1 = Λ 1 μ 1 时, f 2 ( S 1 ) 取最小值为 2 Λ 1 。同理,记 f 3 ( S 2 )= μ 2 S 2 + Λ 2 2 μ 2 S 2 ,则当且仅当 S 2 = Λ 2 μ 2 时, f 3 ( S 2 ) 的最小值为 2 Λ 2 。故仅当在无病平衡点处, d L E dt =0 。所以,无病平衡点 E 是全局渐进稳定的。

定理6:当 R 10 >1, R 20 >1 ,且两斑块上对应参数都相等时,地方病平衡点 E * 是全局渐进稳定的。

证明:构造李雅普诺夫函数

L E ( t )= L S 1 * ( t )+ L I 1 * ( t )+ L V 1 * ( t )+ L S 2 * ( t )+ L I 2 * ( t )+ L V 2 * ( t ) L S 1 * ( t )= S 1 * g( S 1 S 1 * ), L I 1 * ( t )= I 1 * g( I 1 I 1 * ), L V 1 * ( t )= β 12 S 1 * V 1 * r 1 g( V 1 V 1 * ), L S 2 * ( t )= S 2 * g( S 2 S 2 * ), L I 2 * ( t )= I 2 * g( I 2 I 2 * ), L V 2 * ( t )= β 22 S 2 * V 2 * r 2 g( V 2 V 2 * ).

d L S 1 * dt = d S 1 dt ( 1 S 1 * S 1 )=( Λ 1 μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1 )( 1 S 1 * S 1 ) =( μ 1 S 1 * + β 11 S 1 * I 1 * + β 12 S 1 * V 1 * μ 1 S 1 β 11 S 1 I 1 β 12 S 1 V 1   )( 1 S 1 * S 1 ) = μ 1 S 1 ( S 1 S 1 * ) 2 + β 11 S 1 * I 1 * ( 1 S 1 * S 1 + I 1 I 1 * S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 S 1 * S 1 + V 1 V 1 * S 1 V 1 S 1 * V 1 * )

d L I 1 * dt = d I 1 dt ( 1 I 1 * I 1 )=( μ 1 I 1 + β 11 S 1 I 1 + β 12 S 1 V 1 )( 1 I I 1 ) =( β 11 S 1 I 1 + β 12 S 1 V 1 β 11 S 1 * I 1 β 12 S 1 * V 1 * I 1 I 1 * )( 1 I I 1 ) = β 11 S 1 * I 1 * ( 1 S 1 * S 1 I 1 I 1 * S 1 I 1 S 1 * I 1 * )+ β 12 S 1 * V 1 * ( 1 I 1 I 1 * S 1 V 1 I 1 * S 1 * V 1 * I 1 + S 1 V 1 S 1 * V 1 * )

d L V 1 * dt = β 12 S 1 * k 1 d V 1 dt ( 1 V 1 V 1 )= β 12 S 1 * k 1 ( r 1 I 1 k 1 V 1 + d 12 V 2 d 21 V 1 )( 1 V 1 * V 1 ) = β 12 S 1 * V 1 * ( I 1 I 1 * V 1 V 1 * V 1 * I 1 V 1 I 1 * +1 )+ β 12 S 1 * V 1 * ( d 12 V 2 d 21 V 1 k 1 d 12 V 2 d 21 V 1 k 1 V 1 )

同理可得

d L S 2 * dt = μ 2 S 2 ( S 2 S 2 * ) 2 + β 21 S 2 * I 2 * ( 1 S 2 * S 2 + I 2 I 2 * S 2 I 2 S 2 * I 2 * )+ β 22 S 2 * V 2 * ( 1 S 2 * S 2 + V 2 V 2 * S 2 V 2 S 2 * V 2 * ) d L I 2 * dt = β 21 S 2 * I 2 * ( 1 S 2 * S 2 I 2 I 2 * S 2 I 2 S 2 * I 2 * )+ β 22 S 2 * V 2 * ( 1 I 2 I 2 * S 2 V 2 I 2 * S 2 * V 2 * I 2 + S 2 V 2 S 2 * V 2 * ) d L V 2 * dt = β 22 S 2 * V 2 * ( I 2 I 2 * V 2 V 2 * V 2 * I 2 V 2 I 2 * +1 )+ β 22 S 2 * V 2 * ( d 21 V 1 d 12 V 2 k 2 d 21 V 1 d 12 V 2 k 2 V 2 )

所以,

f 4 ( S 1 )= S 1 * S 1 + S 1 S 1 * ,则当且仅当 S 1 = S 1 * 时, f 4 ( S 1 ) 的最小值为2。类似的,记 f 5 ( S 2 )= S 2 * S 2 + S 2 S 2 * ,则当且仅当 S 2 = S 2 * 时, f 5 ( S 2 ) 的最小值为2。故仅当在地方病平衡点 E * 处, d L E ( t ) dt =0 。所以,地方病平衡点 E * 是全局渐进稳定的。

5. 数值模拟

首先,对单斑块的环境下进行数值模拟。取初值为

S 1 =100, I 1 =5, V 1 =1,

其他参数分别为 Λ 1 =100, μ 1 =5, β 11 =0.3, β 12 =0.4, r 1 =2, k 1 =10 ,此时基本再生数 R 10 =1.52>1 。当基本再生数 R 10 >1 时,在单板块的情况下,口蹄疫病毒会在斑块中流行,患病的动物也会随之增加,且患病动物会一直存在,如图1

Figure 1. When R 10 >1 , the temporal variations of each compartment

1. R 10 >1 时,各仓室随时间的变化

当初值不变时,其他参数分别 Λ 1 =100, μ 1 =5, β 11 =0.15, β 12 =0.2, r 1 =2, k 1 =10 ,此时基本再生数 R 10 =0.76<1 。当基本再生数 R 10 <1 时,在单板块的情况下,口蹄疫病毒会在斑块中受到抑制,患病的动物的数量也会随之先增加,但随着病毒传播被抑制,患病动物的数量也会随之减少直至消失,见图2

两个斑块要考虑基本再生数的具体情况。取初值分别为

S 1 =90, I 1 =10, V 1 =5, S 2 =95, I 2 =5, V 2 =3,

以及其他参数分别为

Figure 2. When R 10 <1 , the temporal variations of each compartment

2. R 10 <1 时,各仓室随时间的变化

Λ 1 =100, Λ 2 =20, μ 1 = μ 2 =1, r 1 = r 2 =1, k 1 =0.1, k 2 =1 β 11 =0.1, β 12 =0.2, β 21 =0.01, β 22 =0.02, d 21 =0.1, d 12 =1.

基本再生数 R 1 = 110 >1, R 2 = 10 5 <1.

R 1 >1 R 2 <1 时,此时两个斑块的基本再生数 R 0 = R 1 ,可以得到斑块1中的病毒传播得比斑块2中更加剧烈,同时斑块1中,患病动物由于病毒传播数量迅速增加,又因自然死亡而数量慢慢减少。斑块2中由于病毒受到抑制,传播速度小于斑块1中的传播速度,患病动物的数量也要小于斑块1中的患病动物,如图3。反之,当 R 1 <1 R 2 >1 时,此时两个斑块的基本再生数 R 0 = R 2 ,会存在类似的情况。

当初值不变,其他参数分别为

Λ 1 =100, Λ 2 =80, μ 1 = μ 2 =1, r 1 = r 2 =1, k 1 =0.1, k 2 =1 β 11 =0.1, β 12 =0.2, β 21 =0.01, β 22 =0.02, d 21 =0.1, d 12 =1.

基本再生数分别为 R 1 = 110 >1, R 2 = 2 10 5 >1.

R 1 > R 2 >1 ,此时两个斑块的基本再生数 R 0 = R 1 ,可以得到两个斑块中的病毒都在传播,但由于 R 1 > R 2 ,所以斑块1中的病毒传播速度要大于斑块2中病毒传播速度。同时,斑块1中患病动物数量也是要大于斑块2中的患病动物数量,如图4。同理,当 R 2 > R 1 >1 时,此时两个斑块的基本再生数 R 0 = R 2 ,会存在类似的情况。

当初值不变,其他参数分别为

Figure 3. When R 1 >1 and R 2 <1 , the temporal dynamics of each compartment

3. R 1 >1 R 2 <1 时,各仓室随时间的变化

Figure 4. When R 1 > R 2 >1 , the temporal variations of each compartment

4. R 1 > R 2 >1 时,各仓室随时间的变化

Λ 1 =40, Λ 2 =20, μ 1 = μ 2 =1, r 1 = r 2 =1, k 1 = k 2 =1 β 11 = β 21 =0.01, β 12 = β 22 =0.02, d 21 = d 12 =1.

基本再生数分别为 R 1 = 2 5 5 <1, R 2 = 10 5 <1

R 2 < R 1 <1 时,此时两个斑块的基本再生数 R 0 = R 1 ,可以得到病毒的传播都受到了抑制,但由于 R 1 > R 2 ,所以斑块1中的患病动物的数量是要略多于斑块2中的患病动物的数量,如图5。同理,当 R 1 < R 2 <1 时,此时两个斑块的基本再生数 R 0 = R 2 ,会存在类似的情况。

Figure 5. When R 2 < R 1 <1 , the temporal variations of each compartment

5. R 2 < R 1 <1 时,各仓室随时间的变化

运用PRCC敏感性分析,可以更好地看出各个参数和基本再生数的相关性。对于单板块,如图6,可以看出参数 Λ 1 , β 11 , β 12 , k 1 , r 1 R 10 呈正相关,参数 μ 1 R 10 呈负相关。

对两个斑块进行相关性分析要考虑具体的情况。当 R 1 > R 2 时,此时两个斑块基本再生数 R 0 = R 1 。为了更直观地观察基本再生数与由易感动物转变为患病动物和通过病毒转变为患病动物的关系,绘制了三视图以便观察。可以看出如图7,当动物不管是因为受到同类的传染还是直接接触到了病毒而变成患病动物,对基本再生数都会产生较大的影响,这时会加剧病毒的传播。为了研究基本再生数与病毒的自然降解率和病毒从第一个斑块到第二个斑块的转移率,同样通过三视图观察它们的相关性。可以看出如图8,当病毒迅速向其他斑块传播以及快速降解时,基本再生数是逐渐变小的,此时病毒在斑块内传播的速度要减缓。反之,当病毒向其他斑块传播减缓以及病毒不能快速分解,此时病毒在斑块内传播加剧。同理,当 R 2 > R 1 时,此时两个斑块基本再生数 R 0 = R 2 ,会存在类似情况。

Figure 6. The relationship between R 10 and parameters

6. R 10 与参数的相关性

Figure 7. The relationship between   R 1 and parameters   β 11 , β 12

7.   R 1 与参数 β 11 , β 12 相关性

另外,动物的自然死亡率也对病毒的传播产生影响。如图9,当自然死亡率增大时,病毒的传播速度会逐渐减小。斑块环境中动物的数量也会影响病毒的传播。如图10,当动物的数量得到控制时,病毒的传播也会得到抑制。

Figure 8. The relationship between   R 1 and parameters   k 1 , d 21

8.   R 1 与参数 k 1 , d 21 的相关性

Figure 9. The relationship between   R 1 and parameter   μ 1

9.   R 1 与参数 μ 1 的相关性

Figure 10. The relationship between   R 1 and parameter   Λ 1

10.   R 1 与参数 Λ 1 的相关性

6. 结论

本文建立了两个斑块环境下的口蹄疫传染病模型,定义和计算了模型的基本再生数。对单斑块情形,验证了模型地方病平衡点存在的充分必要条件为基本再生数 R 10 >1 ;通过利用雅克比矩阵特征值的符号和李雅普诺夫法证明了模型平衡点的局部和全局稳定性。对两斑块情形,在假设斑块对称的情况下,通过构造合适的李雅普诺夫函数证明了系统平衡点的全局渐进稳定性。基于数值模拟对理论结果进行了验证,在单斑块环境中,当基本再生数 R 10 >1 时,病毒以及患病动物的数量相比于当 R 10 <1 时,传播速度更快,数量更多。在两个斑块的环境中,验证了病毒的传播速度也是会从基本再生数更大的斑块中,迅速传播的其他斑块。从两个斑块基本再生数和参数的相关性中,可以得到需要减少患病动物和其他动物以及动物和病毒源头的接触,或者对斑块环境中动物的数量进行控制,这样才能降低病毒传播的速度。

基金项目

天津职业技术师范大学科研启动项目(KYQD202324)。

参考文献

[1] Yang, J., Wang, X. and Li, K. (2022) Temporal-spatial Analysis of a Foot-and-Mouth Disease Model with Spatial Diffusion and Vaccination. Frontiers in Veterinary Science, 9, Article ID: 952382.
https://doi.org/10.3389/fvets.2022.952382
[2] Wang, J., Wu, S., Huang, M. and Zhao, H. (2024) Spatial Spread for a Delayed and Nonlocal Foot-and-Mouth Disease Model. Nonlinear Analysis: Real World Applications, 76, Article 104006.
https://doi.org/10.1016/j.nonrwa.2023.104006
[3] Wang, J. and Zhang, R. (2023) A Note on the Global Dynamics for a Diffusive Foot-and-Mouth Disease Model. Applied Mathematics Letters, 145, Article 108737.
https://doi.org/10.1016/j.aml.2023.108737
[4] Yang, J., Gong, M. and Sun, G. (2023) Asymptotical Profiles of an Age-Structured Foot-and-Mouth Disease with Nonlocal Diffusion on a Spatially Heterogeneous Environment. Journal of Differential Equations, 377, 71-112.
https://doi.org/10.1016/j.jde.2023.09.001
[5] Shuai, Z. and van den Driessche, P. (2014) Modelling and Control of Cholera on Networks with a Common Water Source. Journal of Biological Dynamics, 9, 90-103.
https://doi.org/10.1080/17513758.2014.944226
[6] van den Driessche, P. and Watmough, J. (2002) Reproduction Numbers and Sub-Threshold Endemic Equilibria for Compartmental Models of Disease Transmission. Mathematical Biosciences, 180, 29-48.
https://doi.org/10.1016/s0025-5564(02)00108-6
[7] Diekmann, O., Heesterbeek, J.A.P. and Roberts, M.G. (2009) The Construction of Next-Generation Matrices for Compartmental Epidemic Models. Journal of The Royal Society Interface, 7, 873-885.
https://doi.org/10.1098/rsif.2009.0386