具有Beddington-DeAngelis功能反应和分布时滞的病毒动力学模型的稳定性

周欣然; 郑涛; 张龙

doi:10.6054/j.jscnun.2020084

具有Beddington-DeAngelis功能反应和分布时滞的病毒动力学模型的稳定性

新疆大学数学与系统科学学院，乌鲁木齐 830046

基金项目:

国家自然科学基金项目 11861065

国家自然科学基金项目 11361059

国家自然科学基金项目 11771373

国家自然科学基金项目 11702237

新疆维吾尔自治区自然科学基金项目 2019D01C076

新疆维吾尔自治区自然科学基金项目 2017D01C082

新疆维吾尔自治区教育厅重点实验室项目 2016D03022

新疆维吾尔自治区普通高等学校科研计划-创新团队项目 XJEDU2017T001

新疆大学博士创新项目 XJUBSCX-2017005

详细信息

通讯作者:
张龙, 教授, Email:longzhang_xj@sohu.com

中图分类号: O175
计量
- 文章访问数: 871
- HTML全文浏览量: 293
- PDF下载量: 83
出版历程
- 收稿日期: 2020-01-20
- 网络出版日期: 2020-11-02
- 刊出日期: 2020-10-24

The Stability of Virus Dynamics Model with Beddington-DeAngelis Functional Response and Distributed Time Delay

College of Mathematics and System Sciences, Xinjiang University, Urumqi 830046, China

摘要

摘要: 为了更真实地描述细胞被病毒感染的现象, 在常见的病毒时滞模型基础上建立了一类具有Beddington-DeAngelis功能反应和分布时滞的病毒动力学模型.首先，讨论了该模型解的正性和有界性，得到了该模型的基本再生数R₀；然后，证明了该模型的无病平衡点的局部渐近稳定性和全局渐近稳定性；最后，通过数值模拟验证了结论的有效性.
- Beddington-DeAngelis /
- 分布时滞 /
- 基本再生数 /
- Lyapunov函数 /
- 稳定性
Abstract: A class of virus dynamics model with Beddington-DeAngelis functional response and distributed delay is established on the basis of the common virus time delay model to describe the phenomenon of virus infection in cells more accurately. Firstly, the positivity and boundedness of the solution of the model are discussed, and the basic reproduction number R₀ of the model is obtained. Then, the local and global stability of the infection-free equilibrium of the model are proved. Further, numerical simulations are demonstrated to verify the theoretical results.
- Beddington-DeAngelis /
- distributed delay /
- basic reproduction number /
- Lyapunov function /
- stability

HTML全文

自20世纪以来，国内外生物动力学的研究进展迅速，许多专家学者都是通过动力学模型来分析病毒感染^[1-6], 这为更好地理解病毒感染和疾病的控制及不同药物治疗策略效应提供了至关重要的作用.

ANDERSON和MAY^[7]介绍了一种具有自由生存阶段的寄生虫模型, NOWAK和BANGHAM^[8]将其应用于病毒动力学.考虑到病毒感染靶细胞过程中普遍存在的非线性感染率, HUANG等^[9]提出了一类具有Beddington-DeAngelis功能反应的病毒动力学模型:

$\left\{ {\begin{array}{*{20}{l}} {\dot x(t) = \lambda - dx(t) - \frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}},}\\ {\dot y(t) = \frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}} - py(t),}\\ {\dot v(t) = ky(t) - uv(t),} \end{array}} \right.$

(1)

其中，x(t)为健康靶细胞在t时刻的细胞密度，y(t)为感染细胞在t时刻的细胞密度, v(t)为游离病毒在t时刻的细胞密度, λ为健康靶细胞的生成速率, d为健康靶细胞的自然死亡率, β为病毒感染健康靶细胞的感染率, p为感染细胞的自然死亡率, 常数k为每个感染细胞在单位时间内生成的游离病毒数量, u为游离病毒的自然死亡率；并得到了模型(1)的无病平衡点和地方病平衡点的全局渐近稳定性.又考虑到病毒颗粒接触靶细胞的时间可能与接触细胞受到积极影响的时间存在一定的滞后性, 即接触病毒体进入细胞需要时间, HUANG等^[10]将细胞内延迟纳入模型, 得到了以下离散时滞病毒模型：

$\left\{ {\begin{array}{*{20}{l}} {\dot x(t) = \lambda - dx(t) - \frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}},}\\ {\dot y(t) = {{\rm{e}}^{ - p\tau }}\frac{{\beta x(t - \tau )v(t - \tau )}}{{1 + ax(t - \tau ) + bv(t - \tau )}} - py(t),}\\ {\dot v(t) = ky(t) - uv(t),} \end{array}} \right.$

(2)

并通过对特征方程的分析和构造合适的李雅普诺夫函数, 建立了模型(2)的无病平衡点和地方病平衡点的局部稳定性；且证明了当R₀ < 1时无病平衡点全局渐近稳定, R₀>1时地方病平衡点全局渐近稳定.

基于以上模型，考虑到接触细胞被感染变成感染细胞的滞后时间可能由靶细胞自身细胞内物质决定，即不同靶细胞被感染变成感染细胞的滞后时间是可变的, 则可以理想化假设t－τ到t时刻接触病毒细胞的靶细胞都有可能被感染形成感染细胞，这一现象可以用分布时滞来描述^[11].又考虑到目前在全世界范围内仍缺乏根治HIV感染的有效药物, 现阶段的治疗目标是最大限度和持久地抑制患者体内的病毒复制, 而抗逆转录酶药物则是通过抑制病毒感染健康细胞来抑制病毒的复制.因此，引入参数γ∈(0, 1)来代表药物对病毒感染的抑制作用.为了更加细致地描述病毒感染健康靶细胞后的滞后效应以及了解药物治疗参数γ对艾滋病治疗的影响，本文建立了如下数学模型：

$\left\{ {\begin{array}{*{20}{l}} {\dot x(t) = s - dx(t) - (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}},}\\ {\dot y(t) = (1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} d\sigma - py(t),}\\ {\dot v(t) = ky(t) - uv(t),} \end{array}} \right.$

(3)

其中, x(t)为健康靶细胞在t时刻的细胞密度, y(t)为感染细胞在t时刻的细胞密度, v(t)为游离病毒在t时刻的细胞密度, s为健康靶细胞的生成速率, d为健康靶细胞的自然死亡率, β为描述细胞感染的速率, γ∈(0, 1)为抗逆转录酶对病毒感染健康靶细胞的抑制率, $\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \tau)}}\beta x(\sigma)v(\sigma)}}{{1 + ax(\sigma) + bv(\sigma)}}} {\rm{d}}\sigma$ 为在t－τ到t时刻接触病毒的健康靶细胞被感染变成感染细胞的个数, p为感染细胞的自然死亡率, 常数k为每个感染细胞在单位时间内生成的游离病毒数量, u为游离病毒的自然死亡率.

本文主要讨论了模型(3)解的正性和有界性, 得到了该模型的基本再生数R₀; 讨论了该模型的无病平衡点局部渐近稳定性和全局渐近稳定性的判定准则(R₀≤1, R₀ < 1);最后, 通过数值模拟验证本文结论的有效性.

1. 解的非负性与有界性

设C=C([－τ, 0]; $\mathbb{R}$ ³)是[－τ, 0]到 $\mathbb{R}$ ³上的连续函数全体所形成的具有上确界范数的Banach空间, 模型(3)的初始条件为：

$\begin{array}{l} {\kern 1pt} {\kern 1pt} {\kern 1pt} x(\theta ) = {\varphi _1}(\theta ),y(\theta ) = {\varphi _2}(\theta ),v(\theta ) = {\varphi _3}(\theta ),\\ \theta \in [ - \tau ,0],{\varphi _i}(\theta ) \ge 0,{\varphi _i}(0) > 0\quad (i = 1,2,3), \end{array}$

(4)

其中，φ=(φ₁, φ₂, φ₃)^T∈C.

由文献[12]知:模型(3)满足初始条件(4)的任意解(x(t), y(t), v(t))^T在[－τ, ∞)上存在且唯一.在此条件下，可得:

定理1 设(x(t), y(t), v(t))^T为模型(3)在初始条件(4)下的任意解, 则任意解(x(t), y(t), v(t))^T非负且最终有界.

证明首先证明模型(3)在初始条件(4)下的任意解都为正.先证x(t)≥0. (反证)假设对任意的t≥0, 存在t₁>0, x(t₁)=0且x(t)>0.当t∈[0, t₁)时, 有

$\dot x(t) \ge - dx(t) - (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}},$

因此,

$x({t_1}) \ge x(0){\rm{exp}}\left( {\int_0^{{t_1}} {\left( { - d - (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}}} \right){\rm{d}}t} } \right) > 0,$

这与假设x(t₁)=0矛盾.

同理, 假设对任意的t≥0, 存在t₁>0, y(t₁)=0且y(t)>0.当t∈[0, t₁)时, 有

$\dot y(t) = (1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma - py(t).$

(5)

由式(5)，有

$y({t_1}) = y(0){\rm{exp}}\left( {\int_0^{{t_1}} {\frac{{1 - \gamma }}{{y(t)}}} \int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma - p} \right) > 0,$

这与假设y(t₁)=0矛盾.

同理, 由于y(t)>0, 则有

$\dot v(t) = ky(t) - uv(t) > - uv(t).$

(6)

对式(6)两边同时积分, 可得

$v({t_1}) = v(0){\rm{exp}}\left( {\int_0^{{t_1}} {( - u)} {\rm{d}}t} \right) = v(0){{\rm{e}}^{ - u{t_1}}} > 0.$

(7)

这与假设v(t₁)=0矛盾，因此，模型(3)的任意解(x(t), y(t), v(t))^T的非负性得证.

接下来, 证明任意解(x(t), y(t), v(t))^T最终有界.对于模型(3)的第1个等式, 有

$\dot x(t) = s - dx(t) - (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}} \le s - dx(t),$

(8)

由比较定理^[12]可知,

$0 \le x(t) \le \frac{s}{d} + x(0){{\rm{e}}^{ - dt}},$

(9)

则 $\mathop {\lim \sup }\limits_{t \to + \infty } x(t) \le s/d$ , 所以, x(t)在[0，+∞)上有界.

定义

$F(t) = x(t) + y(t) + (1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma ,$

(10)

且δ=min{d－1, p}.对式(10)两边同时求导, 有

$\begin{array}{l} \frac{{{\rm{d}}F(t)}}{{{\rm{d}}t}} = \frac{{{\rm{d}}x(t)}}{{{\rm{d}}t}} + \frac{{{\rm{d}}y(t)}}{{{\rm{d}}t}} + (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax (t) + bv(t)}} - \\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} (1 - \gamma )\frac{{{{\rm{e}}^{ - d\tau }}\beta x(t - \tau )v(t - \tau )}}{{1 + ax(t - \tau ) + bv(t - \tau )}} - \\ \begin{array}{*{20}{l}} {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} d(1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma = }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} s - dx(t) + (1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma - } \end{array}\\ \begin{array}{*{20}{l}} {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} py(t) - \frac{{{{\rm{e}}^{ - d\tau }}\beta x(t - \tau )v(t - \tau )}}{{1 + ax(t - \tau ) + bv(t - \tau )}} - }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} d(1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma < s - dx(t) - }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} py(t) - (d - 1)(1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta x(\sigma )v(\sigma )}}{{1 + ax(\sigma ) + bv(\sigma )}}} {\rm{d}}\sigma < }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} s - \delta F(t)}. \end{array} \end{array}$

(11)

由比较定理可知F(t)有界，从而y(t)有界；再由模型(3)的第3个等式可知v(t)有界.故模型(3)的任意解(x(t), y(t), v(t))^T最终有界.

2. 平衡点及无病平衡点的局部渐近稳定性

模型(3)总存在一个无病平衡点E₀(x₀, 0, 0)，其中x₀=s/d.下面将分析模型(3)是否存在正的地方病平衡点.

若模型(3)存在正的地方病平衡点，则平衡点需满足以下等式：

$\left\{ {\begin{array}{*{20}{l}} {s - d{x^*} - (1 - \gamma )\frac{{\beta {x^*}{v^*}}}{{1 + a{x^*} + b{v^*}}} = 0,}\\ {(1 - \gamma )\int_{t - \tau }^t {\frac{{{{\rm{e}}^{ - d(t - \sigma )}}\beta {x^*}{v^*}}}{{1 + a{x^*} + b{v^*}}}} {\rm{d}}\sigma - p{y^*} = 0,}\\ {k{y^*} - v{u^*} = 0.} \end{array}} \right.$

(12)

由式(12)可得

${{x^*} = \frac{{dpu + (1 - {e^{ - d\tau }})kbs}}{{(1 - {e^{ - d\tau }})[bkd + (1 - \gamma )\beta k] - adpu}},}$

${{v^*} = \frac{k}{{pu}}(s - d{x^*})(\frac{1}{d} - \frac{1}{d}{{\rm{e}}^{ - d\tau }}),{y^*} = \frac{u}{k}{v^*}.}$

若要满足地方病平衡点的正性，则要求x^*>0，即(1－e^－dτ)[kbd+(1－γ)βk]－adpu>0，即要满足

$\frac{{(1 - {{\rm{e}}^{ - d\tau }})[kbd + (1 - \gamma )\beta k]}}{{adpu}} > 1.$

(13)

由于

$1 < \frac{{k\beta s(1 - \gamma )(1 - {{\rm{e}}^{ - d\tau }})}}{{adpu(d + as)}} < \frac{{(1 - {{\rm{e}}^{ - d\tau }})[kbd + (1 - \gamma )\beta k]}}{{adpu}},$

(14)

则可定义模型(3)的基本再生数为：

${R_0} = \frac{{k\beta s(1 - \gamma )(1 - {{\rm{e}}^{ - d\tau }})}}{{dpu(d + as)}}.$

因此，以R₀为阈值可得以下结论：

定理2 若R₀ < 1，则模型(3)的无病平衡点E₀局部渐近稳定；反之，E₀不稳定.

证明系统(3)在E₀处的特征方程为

$\begin{array}{*{20}{c}} {(\lambda + d)[(\lambda + p)(\lambda + u) - \frac{{k\beta s}}{{d + as}}(\frac{1}{d} - \frac{1}{d}{{\rm{e}}^{ - d\tau }}) \times }\\ {(1 - \gamma ){{\rm{e}}^{ - \lambda \tau }}] = 0,} \end{array}$

(15)

由于λ=－d < 0，故E₀的局部稳定性取决于

$f(\lambda ) = (\lambda + p)(\lambda + u) - \frac{{k\beta s}}{{d(d + as)}}(1 - {{\rm{e}}^{ - d\tau }})(1 - \gamma ){{\rm{e}}^{ - \lambda \tau }} = 0.$

(16)

下面分2种情况讨论.

情况1:假设 ${R_0} = \frac{{k\beta s(1 - \gamma)\left({1 - {{\rm{e}}^{ - d\tau }}} \right)}}{{dpu(d + as)}} > 1$ ，则

$f(0) = pu - \frac{{k\beta s(1 - \gamma )}}{{d(d + as)}}(1 - {{\rm{e}}^{ - d\tau }}) < 0,$

(17)

且 $\mathop {\lim }\limits_{\lambda \to \infty } f(\lambda) = + \infty$ .从而，f(λ)=0在0, +∞)上至少有1个正根.因此，若R₀>1，E₀不稳定.

情况2:假设 ${R_0} = \frac{{k\beta s(1 - \gamma)\left({1 - {{\rm{e}}^{ - {d_\tau }}}} \right)}}{{dpu(d + as)}} < 1$ ，则

$f(\lambda ) = (\lambda + p)(\lambda + u) - \frac{{k\beta s(1 - \gamma )}}{{d(d + as)}}(1 - {{\rm{e}}^{ - d\tau }}){{\rm{e}}^{ - \lambda \tau }} = 0,$

(18)

等价于

$G(\lambda ) = \frac{{k\beta s(1 - \gamma )(1 - {{\rm{e}}^{ - d\tau }}){{\rm{e}}^{ - \lambda \tau }}/d(d + as)}}{{(\lambda + p)(\lambda + u)}} = 1.$

(19)

令λ=x+iy(x, y∈ $\mathbb{R}$ ). (反证)假设x≥0，则

$\begin{array}{*{20}{l}} {1 = |G(\lambda )| = \frac{{|k\beta s(1 - \gamma )(1 - {{\rm{e}}^{ - d\tau }}){{\rm{e}}^{ - \lambda \tau }}/d(d + as)|}}{{|(\lambda + p)(\lambda + u)|}} \le }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \frac{{|k\beta s(1 - \gamma )(1 - {{\rm{e}}^{ - d\tau }}){{\rm{e}}^{ - x\tau }}/d(d + as)|}}{{|(x + p)(x + u)|}} \le }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \frac{{|k\beta s(1 - \gamma )(1 - {{\rm{e}}^{ - d\tau }})/d(d + as)|}}{{|pu|}} = }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \frac{{k\beta s(1 - \gamma )}}{{pud(d + as)}}(1 - {{\rm{e}}^{ - d\tau }}) = {R_0} < 1,} \end{array}$

矛盾，所以x < 0，即此时f(λ)的根具有负实部.因此，当R₀ < 1时, E₀是局部渐近稳定的.

3. 无病平衡点E₀的全局渐近稳定性

定理3 对于模型(3), 若R₀≤1, 则无病平衡点E₀全局渐近稳定.

证明设x₀=s/d，则模型(3)等价于：

$\left\{ {\begin{array}{*{20}{l}} {\dot x(t) = s - dx(t) - (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}},}\\ {\dot y(t) = (1 - \gamma )\int_0^\tau {\frac{{{{\rm{e}}^{ - d\sigma }}\beta x(t - \sigma )v(t - \sigma )}}{{1 + ax(t - \sigma ) + bv(t - \sigma )}}} {\rm{d}}\sigma - py(t),}\\ {\dot v(t) = ky(t) - uv(t),} \end{array}} \right.$

(20)

定义李雅普诺夫函数

$\begin{array}{*{20}{l}} {V = \frac{1}{d}(1 - {{\rm{e}}^{ - d\tau }})\frac{{{x_0}}}{{1 + a{x_0}}}\left( {\frac{{x(t)}}{{{x_0}}} - 1 - {\rm{ln}}\frac{{x(t)}}{{{x_0}}}} \right) + y(t) + }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \frac{p}{k}v(t) + (1 - \gamma ){U^ - }(t),} \end{array}$

(21)

其中

${U^ - }(t) = \int_0^\tau {\int_{t - \sigma }^t {\frac{{{{\rm{e}}^{ - d\sigma }}\beta x(\theta )v(\theta )}}{{1 + ax(\theta ) + bv(\theta )}}} } {\rm{d}}\theta {\rm{d}}\sigma ,$

则有

$\begin{array}{l} \begin{array}{*{20}{c}} {\frac{{{\rm{d}}{U^ - }(t)}}{{{\rm{d}}t}} = \int_0^\tau {\frac{{{{\rm{e}}^{ - d\sigma }}\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}}} {\rm{d}}\sigma - \int_0^\tau {\frac{{{{\rm{e}}^{ - d\sigma }}\beta x(t - \sigma )v(t - \sigma )}}{{1 + ax(t - \sigma ) + bv(t - \sigma )}}} {\rm{d}}\sigma = }\\ {\frac{1}{d}(1 - {{\rm{e}}^{ - d\tau }})\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}} - \int_0^\tau {\frac{{{{\rm{e}}^{ - d\sigma }}\beta x(t - \sigma )v(t - \sigma )}}{{1 + ax(t - \sigma ) + bv(t - \sigma )}}} {\rm{d}}\sigma ,} \end{array}\\ \begin{array}{*{20}{l}} {\frac{{{\rm{d}}V}}{{{\rm{d}}t}} = \frac{1}{d}(1 - {{\rm{e}}^{ - d\tau }})\frac{1}{{1 + a{x_0}}}\left( {1 - \frac{{{x_0}}}{{x(t)}}} \right)\dot x(t) + \dot y(t) + \frac{p}{k}\dot v(t) + }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} (1 - \gamma )\frac{{{\rm{d}}{U^ - }(t)}}{{{\rm{d}}t}} = \frac{1}{d}(1 - {{\rm{e}}^{ - d\tau }})\frac{1}{{1 + a{x_0}}}\left( {1 - \frac{{{x_0}}}{{x(t)}}} \right)\left( {d{x_0} - dx(t) - } \right.}\\ {\left. {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} (1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}}} \right) - \frac{{pu}}{k}v(t) + \frac{1}{d}(1 - {{\rm{e}}^{ - d\tau }}) \times } \end{array}\\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \begin{array}{*{20}{l}} {(1 - \gamma )\frac{{\beta x(t)v(t)}}{{1 + ax(t) + bv(t)}} = - \frac{{(1 - {{\rm{e}}^{ - d\tau }}){{(x(t) - {x_0})}^2}}}{{x(t)(1 + a{x_0})}} + }\\ {\frac{{puv(t)}}{k}\frac{{1 + ax(t)}}{{1 + ax(t) + bv(t)}}({R_0} - 1) - }\\ {\frac{{pub}}{{k(1 + ax(t) + bv(t))}}{v^2}(t).} \end{array} \end{array}$

(22)

由式(22)可知:当R≤1时, dV/dt≤0；当dV/dt=0时，v=0当且仅当x=x₀.由拉萨尔不变原理可知E₀是全局渐近稳定的.

4. 数值模拟

本节通过MATLAB软件进行数值模拟, 以验证该模型的无病平衡点具有全局渐近稳定性.在此取定相应系数值与初值：a=0.6, b=0.8, d=0.1, s=0.8, β=0.4, p=0.5, k=0.3, u=0.2, γ=0.5.(x(0), y(0), v(0))分别取为(9.5, 5.5, 6.6)、(7.4, 4.8, 5.3)、(6.8, 5, 6)、(4.4, 3.1, 4.8)、(5, 3.2, 4)、(14.5, 3.2, 2.6)、(2, 3.2, 2.6).

取τ=0.2，可得R₀=0.163 9 < 1.由此可以分别模拟得到模型(3)在R₀ < 1时的时间序列图(图 1)与三维空间相图(图 2).

图 1 模型(3)的时间序列图(R₀ < 1)

Figure 1. The time series plot of model(3) when R₀ < 1

下载: 全尺寸图片幻灯片

图 2 模型(3)的三维空间相图(R₀ < 1)

注：颜色不同的曲线代表初值不同的解曲线.

Figure 2. The 3-dimensional phase diagram of model (3) when R₀ < 1

下载: 全尺寸图片幻灯片

由图 1可知：x(t)随着时间的增加达到一定值, y(t)、v(t)随时间的增加而最终趋于0.由图 1和图 2可知：对任意给定的初值, x(t)、y(t)、v(t)最终趋于x轴的一定值E₀=(8, 0, 0).由此可知：无病平衡点在R₀ < 1时全局渐近稳定.

取τ=1.287 8, 可得R₀=1.由模拟得到的模型(3)在R₀=1时的时间序列图(图 3)和三维空间相图(图 4)可知：R₀=1时，x(t)、y(t)、v(t)的走向与R₀ < 1时的一致, 则可知无病平衡点在R₀=1时也是全局渐近稳定的.

图 3 模型(3)的时间序列图(R₀=1)

Figure 3. The time series plot of model(3) when R₀=1

下载: 全尺寸图片幻灯片

图 4 模型(3)的三维空间相图(R₀=1)

注：颜色不同的曲线代表初值不同的解曲线.

Figure 4. The 3-dimensional phase diagram of model(3) when R₀=1

下载: 全尺寸图片幻灯片

取τ=6.6, 则可得R₀=3.998 5>1.从而模拟得到模型(3)在R₀>1时的时间序列图(图 5)与三维空间相图(图 6).

图 5 模型(3)的时间序列图(R₀>1)

Figure 5. The time series plot of model(3) when R₀>1

下载: 全尺寸图片幻灯片

图 6 模型(3)的三维空间相图(R₀>1)

注：颜色不同的曲线代表初值不同的解曲线.

Figure 6. The 3-dimensional phase diagram of model(3) when R₀>1

下载: 全尺寸图片幻灯片

由图 5可知：x(t)随着时间的变化达到一定值, 但y(t)、v(t)并未随时间的变化而消亡, 故此时无病平衡点不稳定.由图 5和图 6可知：对任意给定的初值, x(t)、y(t)、v(t)最终趋于一定值E^*=(2.971 4, 4.859 1, 7.288 7).由此可知：地方病平衡点在R₀>1是全局渐近稳定.

5. 小结

本文建立了具有Beddington-DeAngelis功能反应和分布延迟的病毒动力学模型，并讨论了该模型的稳定性.为了研究药物对病毒感染的抑制作用，引入了药物治疗参数γ∈(0, 1), 并求解该模型在药物治疗下的基本再生数R₀, 通过分析该模型的特征方程, 证明了该模型的无病平衡点E₀的局部渐近稳定性; 通过构造合适的李雅普诺夫泛函和利用李亚普诺夫-拉萨尔型定理, 证明了该模型的无感染平衡点E₀在R₀≤1时的全局渐近稳定性; 通过Matlab验证了本文结论的有效性.本文虽然没有讨论地方病平衡点的全局渐近稳定性, 但由数值模拟可知地方病平衡点E^*的全局渐近稳定性, 由三维空间相图可知地方病平衡点E^*在R₀>1时是全局渐近稳定的.

图 1 模型(3)的时间序列图(R₀ < 1)

Figure 1. The time series plot of model(3) when R₀ < 1

下载: 全尺寸图片幻灯片

图 2 模型(3)的三维空间相图(R₀ < 1)

注：颜色不同的曲线代表初值不同的解曲线.

Figure 2. The 3-dimensional phase diagram of model (3) when R₀ < 1

下载: 全尺寸图片幻灯片

图 3 模型(3)的时间序列图(R₀=1)

Figure 3. The time series plot of model(3) when R₀=1

下载: 全尺寸图片幻灯片

图 4 模型(3)的三维空间相图(R₀=1)

注：颜色不同的曲线代表初值不同的解曲线.

Figure 4. The 3-dimensional phase diagram of model(3) when R₀=1

下载: 全尺寸图片幻灯片

图 5 模型(3)的时间序列图(R₀>1)

Figure 5. The time series plot of model(3) when R₀>1

下载: 全尺寸图片幻灯片

图 6 模型(3)的三维空间相图(R₀>1)

注：颜色不同的曲线代表初值不同的解曲线.

Figure 6. The 3-dimensional phase diagram of model(3) when R₀>1

下载: 全尺寸图片幻灯片

参考文献(12)

[1]	李盈科, 滕志东, 努尔别克·艾孜玛洪.时变多时滞连续系统的有限时间稳定性分析[J].新疆大学学报(自然科学版), 2016, 33(2):161-164. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=xjdxxb201602006 LI Y K, TENG Z D, NERBEK A. Finite-time stability of continuous-time systems with time-varying delays[J]. Journal of Xinjiang University(Natural Science Edition), 2016, 33(2):161-164. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=xjdxxb201602006
[2]	许明星, 滕志东, 张龙.具有密度制约的非自治SIRS传染病模型的持久性与灭绝性[J].新疆大学学报(自然科学版), 2015, 32(1):33-39. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=xjdxxb201501007 XU M X, TENG Z D, ZHANG L. Permanence and extinction for nonautonomous SIRS epidemic model with density dependence[J]. Journal of Xinjiang University(Natural Science Edition), 2015, 32(1):33-39. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=xjdxxb201501007
[3]	黑力军, 陈斯养.具有时滞三斑块扩散捕食模型的全局稳定性[J].工程数学学报, 2003, 20(4):129-133. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=gcsxxb200304026 HEI L J, CHEN S Y. Global stability for a three-patch diffusion predator-prey model with time delays[J]. Journal of Engineering Mathematics, 2003, 20(4):129-133. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=gcsxxb200304026
[4]	MCCLUSKEY C C. Complete global stability for an SIR epidemic model with delay distributed or discrete[J]. Nonlinear Analysis:Real World Applications, 2010, 11(1):55-59. doi: 10.1016/j.nonrwa.2008.10.014
[5]	毛卫华.二阶非线性具时滞脉冲微分方程振动性[J].华南师范大学学报(自然科学版), 2002(1):47-53. http://journal-n.scnu.edu.cn/article/id/1092 MAO W H. Oscillations of second-order nonlinear impulse differential equaitons with delays[J]. Journal of South China Normal University (Natural Science Edition), 2002(1):47-53. http://journal-n.scnu.edu.cn/article/id/1092
[6]	黄燕革, 姚晓洁, 黄勇.时间测度上具有Beddington-DeAngelis类功能反应和扩散的捕食系统的周期解[J].华南师范大学学报(自然科学版), 2010(4):19-26. http://journal-n.scnu.edu.cn/article/id/445 HUANG Y G, YAO X J, HUANG Y. Periodic solutions for a class of predator-prey system on time scales with Beddington-DeAngelis functional response and diffusion[J]. Journal of South China Normal University (Natural Science Edition), 2010(4):19-26. http://journal-n.scnu.edu.cn/article/id/445
[7]	ANDERSON R M, MAY R M. The population dynamics of microparasites and their invertebrate hosts[J]. Philosophi-cal Transactions of the Royal Society:Series B, 1981, 291:451-524. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=1e8f990c1ce9ec8d24fa63ad4bf10950
[8]	NOWAK M A, BANGHAM C R M. Population dynamics of immune responses to persistent virus[J]. Science, 1996, 272:74-79. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=PubMed000000973303
[9]	HUANG G, MA W B, TAKEUCHI Y. Global properties for virus dynamics model with Beddington-DeAngelis functional response[J]. Applied Mathematics Letters, 2009, 22:1690-1693. doi: 10.1016/j.aml.2009.06.004
[10]	HUANG G, MA W B, TAKEUCHI Y. Global analysis for delay virus dynamics model with Beddington-DeAngelis functional response[J]. Applied Mathematics Letters, 2011, 24:1199-1203. doi: 10.1016/j.aml.2011.02.007
[11]	XIAMIXINUER A, TENG Z D. Global dynamics in a he-roin epidemic model with different conscious stages and two distributed delays[J]. International Journal of Biomathematics, 2019, 12:1-19. http://www.zhangqiaokeyan.com/academic-journal-cn_journal-biomathematics-english_thesis/0201275812541.html
[12]	JACK H. Theory of functional differential equations[M]. Berlin:Springer, 1977:36-42.

施引文献(2)

期刊类型引用(1)

宋冰，张育茹，桑媛，张龙. 具有饱和发生率和分布时滞的HIV免疫模型的稳定性. 广西师范大学学报(自然科学版). 2023(02): 106-117 .

百度学术

其他类型引用(1)

资源附件(0)

图(6)

计量

文章访问数: 871
HTML全文浏览量: 293
PDF下载量: 83
被引次数: 2

1. 解的非负性与有界性
2. 平衡点及无病平衡点的局部渐近稳定性
3. 无病平衡点E0的全局渐近稳定性
4. 数值模拟
5. 小结

具有Beddington-DeAngelis功能反应和分布时滞的病毒动力学模型的稳定性

通讯作者: 张龙, 教授, Email:longzhang_xj@sohu.com

计量

出版历程