具有两类靶细胞和潜伏感染的HIV感染模型

期刊菜单

具有两类靶细胞和潜伏感染的HIV感染模型
A Model of HIV Infection with Two Types of Target Cells and Latent Infection

DOI: 10.12677/aam.2024.136257, PDF, HTML, XML, 下载: 23 浏览: 45
作者: 张佳豪, 郭黄涛, 王艳^*：中国石油大学(华东)理学院，山东青岛
关键词: HIV病毒；CD4⁺T细胞；巨噬细胞；基本再生数；稳定性分析；HIV； CD4⁺T Cells； Macrophages； Basic Reproduction Number； Stability Analysis

摘要: 本文研究具有潜伏感染阶段的HIV在两种不同类型靶细胞(活化的CD4⁺T细胞和巨噬细胞)中的传播动力学，考虑HIV感染的多个方面，包括病毒复制、细胞感染、药物治疗干预，以及潜伏感染细胞的作用，建立具有两类靶细胞的数学模型。该模型存在无病平衡点和感染平衡点，通过对模型特征方程根的分析得到各平衡点的局部渐近稳定性。我们关注抗逆转录病毒药物(RTIs和PIs)对HIV传播的抑制效果，并探讨这些药物在不同靶细胞中的作用。数值模拟结果表明，模型的动态行为受到药物抑制效果以及其它关键生物参数的显著影响。特别是，我们展示在特定的药物效能和病毒感染参数下，系统能够达到无病状态或维持低病毒负荷的长期稳定状态。此外，本研究还通过数值模拟揭示不同药物组合和剂量下模型的响应，强调优化HIV治疗策略的潜在重要性。

Abstract: This paper studies the dynamics of HIV transmission in two different types of target cells (activated CD4⁺T cells and macrophages) with latent infection. It considers various aspects of HIV infection, including viral replication, cell infection, antiretroviral intervention, and the role of latently infected cells, establishing a mathematical model with two types of target cells. The model exhibits both disease-free and infected equilibrium points, and the local asymptotic stability of each equilibrium point is determined by analyzing the roots of the model’s characteristic equation. We focus on the inhibitory effects of antiretroviral drugs (RTIs and PIs) on HIV transmission and explore the effects of these drugs in different target cells. Numerical simulations indicate that the model’s dynamics are significantly influenced by the drug inhibition efficacy and other key biological parameters. In particular, we demonstrate that under certain drug efficacy and viral infection parameters, the system can achieve a disease-free state or maintain a long-term stable state with a low viral load. Moreover, this study reveals the model's response under different drug combinations and dosages through numerical simulations, emphasizing the potential importance of optimizing HIV treatment strategies.

文章引用：张佳豪, 郭黄涛, 王艳. 具有两类靶细胞和潜伏感染的HIV感染模型[J]. 应用数学进展, 2024, 13(6): 2687-2700. https://doi.org/10.12677/aam.2024.136257

1. 引言

感染人类免疫缺陷病毒(HIV-1)的细胞大多是活化的CD4⁺T细胞。一旦感染，病毒就会在这些细胞上开始复制，从而延长感染时间。只要检测到了这种感染，免疫系统就会相应地做出反应，在有限的范围内控制病毒数量。

抗逆转录病毒药物可以进一步协助控制病毒数量，通常有两类抗逆转录病毒药物用于降低病毒载量和限制受感染的T细胞数量。其中一类被称为逆转录酶抑制剂(RTIs)，可以降低活化的CD4⁺T细胞的感染率。另一类是蛋白酶抑制剂(PIs)，可阻止HIV-1蛋白酶将HIV多聚蛋白裂解为功能单元，从而使受感染细胞产生不具感染性的未成熟病毒颗粒，减少新产生的传染性病毒的数量。如果长期以足够的频率服用这些药物，病毒的数量就会受到很大的抑制。

但是抗逆转录病毒治疗并不能完全根除病毒，因为一旦中断治疗，病毒数量就会反弹。研究学者们提出了许多思路来解释这种反弹。部分研究学者认为，HIV可以作为前病毒在静息记忆CD4⁺T细胞中持续存在。这些潜伏感染的细胞数量比正常的易感T细胞少得多，在受相关抗原激活前不会产生新的病毒[1] [2]。但是它们的寿命周期相对较长，并且可以为病毒提供一个藏身之处以逃避免疫系统或抗逆转录病毒药物的控制，因此成为了一个可以长期产生病毒的贮存库。人们已经建立了包含潜伏感染的数学模型来研究病毒和CD4⁺T细胞在感染过程中的动态变化，分析抗逆转录病毒疗法的影响、耐药性的出现等问题[3] [4]。

${\begin{cases} \dot{T} = λ - d_{T} T - k (1 - ε_{R T}) T V_{I}, \\ \dot{L} = p k (1 - ε_{R T}) T V_{I} - α L - d_{L} I, \\ \dot{I} = (1 - p) k (1 - ε_{R T}) T V_{I} + α L - d_{I} I, \\ {\dot{V}}_{I} = N (1 - ε_{P I}) d_{I} I - d_{V} V_{I}, \\ {\dot{V}}_{N I} = N ε_{P I} d_{I} I - d_{V^{'}} V_{N I}, \end{cases}$

其中T代表初始易感染的靶细胞，L代表潜伏感染细胞，I代表已感染靶细胞，V_I和V_NI分别代表具有感染性和不具感染性的病毒颗粒。

然而，HIV病毒除了入侵CD4⁺T细胞，还可能攻击巨噬细胞以及树突状细胞等。因此，考虑病毒粒子与多类靶细胞的互相作用是合理并且必要的。从生物学角度来看，HIV感染人体的靶细胞可以分为两类，分别为整合速度快的和整合速度慢的。其中整合速度快的以CD4⁺T细胞为代表，而整合速度慢的则以巨噬细胞为代表。人们已经根据这一病毒感染机制提出了考虑两种靶细胞的HIV感染模型[5] [6]。

${\begin{cases} {\dot{T}}_{1} = λ_{1} - μ_{1} T_{1} - k (1 - ε_{1 R T}) T_{1} V_{I}, \\ {\dot{I}}_{1} = k (1 - ε_{1 R T}) T_{1} V_{1} - δ_{1} I_{1}, \\ {\dot{T}}_{2} = λ_{2} - μ_{2} T_{2} - k_{2} (1 - ε_{2 R T}) T_{2} V_{I}, \\ {\dot{I}}_{2} = k_{2} (1 - ε_{2 R T}) T_{2} V_{2} - δ_{2} I_{2}, \\ {\dot{V}}_{I} = N_{1} (1 - ε_{1 P I}) δ_{1} I_{1} + N_{2} (1 - ε_{2 P I}) δ_{2} I_{2} - c_{1 V} V_{I}, \\ {\dot{V}}_{N I} = N_{1} ε_{1 P I} δ_{1} I_{1} + N_{2} ε_{2 P I} δ_{2} I_{2} - c_{2 V} V_{N I}, \end{cases}$

其中，T₁和T₂分别代表了初始易感染的两种靶细胞：CD4⁺T细胞和巨噬细胞，I₁和I₂代表已感染靶细胞，V_I和V_NI分别代表具有感染性和不具感染性的病毒颗粒。

2. 包含潜伏感染阶段的具有两种靶细胞的HIV感染模型

在前人研究的基础上，我们提出一种具有两种靶向细胞和潜伏感染阶段的HIV感染模型，并考虑两类抑制剂(逆转录酶抑制剂和蛋白酶抑制剂)对HIV病毒的治疗效果。建立的模型如下：

${\begin{cases} {\dot{T}}_{1} = λ_{1} - d_{T_{1}} T_{1} - β_{1} (1 - ε_{1 R T}) T_{1} V_{I}, \\ {\dot{L}}_{1} = f_{1} β_{1} (1 - ε_{1 R T}) T_{1} V_{I} + α_{1} L_{1} - δ_{L_{1}} L_{1}, \\ {\dot{I}}_{1} = (1 - f_{1}) β_{1} (1 - ε_{1 R T}) T_{1} V_{I} + α_{1} L_{1} - δ_{I_{1}} I_{1}, \\ {\dot{T}}_{2} = λ_{2} - d_{T_{2}} T_{2} - β_{2} (1 - ε_{2 R T}) T_{2} V_{I}, \\ {\dot{L}}_{2} = f_{2} β_{2} (1 - ε_{2 R T}) T_{2} V_{I} - α_{2} L_{2} - δ_{L_{2}} L_{2}, \\ {\dot{I}}_{2} = (1 - f_{2}) β_{2} (1 - ε_{2 R T}) T_{2} V_{I} + α_{2} L_{2} - δ_{I_{2}} I_{2}, \\ {\dot{V}}_{I} = N_{1} (1 - ε_{1 P I}) δ_{I_{1}} I_{1} + N_{2} (1 - ε_{2 P I}) δ_{I_{2}} I_{2} - c_{1 V} V_{I}, \\ {\dot{V}}_{N I} = N_{1} ε_{1 P I} δ_{I_{1}} I_{1} + N_{2} ε_{2 P I} δ_{I_{2}} I_{2} - c_{2 V} V_{N I}, \end{cases}$ (2-1)

其中，T₁和T₂分别代表了初始易感染的两种靶细胞：CD4⁺T细胞和巨噬细胞，L₁和L₂代表潜伏感染细胞，I₁和I₂代表已感染靶细胞，V_I和V_NI分别代表具有感染性和不具感染性的病毒颗粒。其余参数含义见表1。

Table 1. Meaning of parameters

表1. 参数含义

参数	含义	参数	含义
$λ_{1}, λ_{2}$	未感染靶细胞的产生速率	$α_{1}, α_{2}$	潜伏感染细胞的激活率
$d_{T_{1}}, d_{T_{2}}$	未感染靶细胞死亡率	$δ_{L_{1}}, δ_{L_{2}}$	潜伏感染细胞的死亡率
$β_{1}, β_{2}$	感染率	$δ_{I_{1}}, δ_{I_{2}}$	感染细胞的死亡率
$ε_{1 R T}, ε_{2 R T}$	逆转录酶抑制剂疗效	$N_{1}, N_{2}$	两种靶细胞内病毒的裂解量
$f_{1}, f_{2}$	转为潜伏感染细胞比例	$c_{1 V}, c_{2 V}$	病毒清除速率
$ε_{1 P I}, ε_{2 P I}$	蛋白酶抑制剂疗效

为了方便书写，我们令 ${\bar{β}}_{1} = β_{1} (1 - ε_{1 R T})$ ， ${\bar{β}}_{2} = β_{2} (1 - ε_{2 R T})$ ， ${\bar{N}}_{1} = N_{1} (1 - ε_{1 P I})$ ， ${\bar{N}}_{2} = N_{2} (1 - ε_{2 P I})$ 。

正解的存在唯一性及有界性

定理1: 当模型(2-1)的初始条件满足 $T_{1} (0) > 0$ ， $L_{1} (0) > 0$ ， $I_{1} (0) > 0$ ， $T_{2} (0) > 0$ ， $L_{2} (0) > 0$ ， $I_{2} (0) > 0$ ， $V_{I} (0) > 0$ ， $V_{N I} (0) > 0$ 时，模型(2-1)存在唯一有界的正解。

证明根据非线性微分方程组解的存在唯一性定理[7]，在 $T_{1} (0) > 0$ ， $L_{1} (0) > 0$ ， $I_{1} (0) > 0$ ， $T_{2} (0) > 0$ ， $L_{2} (0) > 0$ ， $I_{2} (0) > 0$ ， $V_{I} (0) > 0$ ， $V_{N I} (0) > 0$ 的初始条件下，模型(2-1)存在唯一解。

接下来证明解是非负的。

假设 $\exists t_{1}, t_{2}$ ，并且满足 $t_{1} = \inf {t | T_{1} (t) = 0, t > 0}$ ， $t_{2} = \inf {t | T_{2} (t) = 0, t > 0}$ ，这就使得 ${{\dot{T}}_{1} |}_{t = t_{1}} = λ_{1} < 0$ ， ${{\dot{T}}_{2} |}_{t = t_{2}} = λ_{2} < 0$ ，这与 $λ_{1} > 0$ ， $λ_{2} > 0$ 矛盾，故 $T_{1} (t) \geq 0$ ， $T_{2} (t) \geq 0$ 。

假设 $\exists t_{3}, t_{4}$ ，并且满足 $t_{3} = \inf {t | L_{1} (t) = 0, t > 0}$ ， $t_{4} = \inf {t | L_{2} (t) = 0, t > 0}$ ，这就使得 ${{\dot{L}}_{1} |}_{t = t_{3}} = f_{1} {\bar{β}}_{1} T_{1} (t_{3}) V_{I} (t_{3}) < 0$ ， ${{\dot{L}}_{2} |}_{t = t_{4}} = f_{2} {\bar{β}}_{2} T_{2} (t_{4}) V_{I} (t_{4}) < 0$ ，故 $V_{I} (t_{3}) < 0$ ， $V_{I} (t_{4}) < 0$ 。因为 $V_{I} (0) > 0$ ，所以 $\exists t_{5}$ 使得 $t_{5} = \inf {t | V_{I} (t) = 0, t \in (0, t_{\min}), t_{\min} = \min {t_{3}, t_{4}}}$ ，故 ${{\dot{V}}_{I} |}_{t = t_{5}} = {\bar{N}}_{1} δ_{I_{1}} I_{1} (t_{5}) + {\bar{N}}_{2} δ_{I_{2}} I_{2} (t_{5}) < 0$ 。

不失一般性，假设 $I_{1} (t_{5}) < 0$ 。因为 $I_{1} (0) > 0$ ，所以 $\exists t_{6} = \inf {t | I_{1} (t) = 0, t \in (0, t_{5})}$ ，使得 ${{\dot{I}}_{1} |}_{t = t_{6}} = (1 - f_{1}) {\bar{β}}_{1} T_{1} (t_{6}) V_{I} (t_{6}) < 0$ ，故 $V_{I} (t_{6}) < 0$ 。而 $t_{6} < t_{5}$ ，故 $V_{I} (t_{6}) > 0$ ，出现矛盾，故 $I_{1} (t) \geq 0$ ( $t \in (0, t_{5}]$ )，同理 $I_{2} (t) \geq 0$ ( $t \in (0, t_{5}]$ )，这与 ${{\dot{V}}_{I} |}_{t = t_{5}} < 0$ 矛盾，故 $V_{I} (t) \geq 0$ ( $t \in (0, t_{\min}]$ ， $t_{\min} = \min {t_{3}, t_{4}}$ )，这又与 ${{\dot{L}}_{1} |}_{t = t_{3}} < 0$ ， ${{\dot{L}}_{2} |}_{t = t_{4}} < 0$ 矛盾，故不存在这样的 $t_{3}$ ， $t_{4}$ 。因此，在 $t > 0$ 时 $L_{1} (t)$ ， $L_{2} (t)$ ， $I_{1} (t)$ ， $I_{2} (t)$ ， $V_{I} (t) \geq 0$ 。

假设 $\exists t_{7}$ ，并且满足 $t_{7} = \inf {t | V_{N I} (t) = 0, t > 0}$ ，故 ${{\dot{V}}_{N I} |}_{t = t_{7}} = N_{1} ε_{1 P I} δ_{I_{1}} I_{1} (t_{7}) + N_{2} ε_{2 P I} δ_{I_{2}} I_{2} (t_{7}) < 0$ ，这与 $I_{1} (t)$ ， $I_{2} (t) \geq 0$ 矛盾，故 $V_{N I} (t) \geq 0$ 。

综上，解的非负性得证。下面证明解的有界性。我们定义

$Z (t) = T_{1} (t) + L_{1} (t) + I_{1} (t) + T_{2} (t) + L_{2} (t) + I_{2} (t) + V_{I} (t) + V_{N I} (t),$

$\dot{Z} = λ_{1} + λ_{2} - d_{T_{1}} T_{1} - δ_{L_{1}} L_{1} - d_{T_{2}} T_{2} - δ_{L_{2}} L_{2} - (1 - N_{1}) δ_{I_{1}} I_{1} - (1 - N_{2}) δ_{I_{2}} I_{2} - c_{1 V} V_{I} - c_{2 V} V_{N I},$

不妨记 $γ = \min {d_{T_{1}}, d_{T_{2}}, δ_{L_{1}}, δ_{L_{2}}, (1 - N_{1}) δ_{I_{1}}, (1 - N_{2}) δ_{I_{2}}, c_{1 V}, c_{2 V}}$ ,

则可得

$\dot{Z} \leq λ_{1} + λ_{2} - γ (T_{1} + T_{2} + L_{1} + L_{2} + I_{1} + I_{2} + V_{I} + V_{N I}) = λ_{1} + λ_{2} - γ Z,$

$\lim_{t \to \infty} \sup Z \leq \frac{λ_{1} + λ_{2}}{γ},$

因此 $Z (t)$ 是有界的，故 $T_{1} (t)$ ， $L_{1} (t)$ ， $I_{1} (t)$ ， $T_{2} (t)$ ， $L_{2} (t)$ ， $I_{2} (t)$ ， $V_{I} (t)$ ， $V_{N I} (t)$ 是有界的。

综上，当模型(2-1)的初始条件满足 $T_{1} (0) > 0$ ， $L_{1} (0) > 0$ ， $I_{1} (0) > 0$ ， $T_{2} (0) > 0$ ， $L_{2} (0) > 0$ ， $I_{2} (0) > 0$ ， $V_{I} (0) > 0$ ， $V_{N I} (0) > 0$ 时，模型存在唯一有界的正解。

3. 模型平衡点

下面我们来求解模型(2-1)的无病平衡点和感染平衡点，由于 $V_{N I}$ 代表不具感染性的病毒颗粒，所以不会对模型的稳定性产生影响，所以在本节稳定性分析的过程中我们暂不考虑 $V_{N I}$ 。先考虑无病平衡点，我们令模型(2-1)的右边都等于0，并假设 $V_{I} = 0$ ，得到以下方程：

${\begin{cases} λ_{1} - d_{T_{1}} T_{1} - {\bar{β}}_{1} T_{1} V_{I} = 0, \\ f_{1} {\bar{β}}_{1} T_{1} V_{I} - α_{1} L_{1} - δ_{L_{1}} L_{1} = 0, \\ {(1 - f}_{1}) {\bar{β}}_{1} T_{1} V_{I} + α_{1} L_{1} - δ_{I_{1}} I_{1} = 0, \\ λ_{2} - d_{T_{2}} T_{2} - {\bar{β}}_{2} T_{2} V_{I} = 0, \\ f_{2} {\bar{β}}_{2} T_{2} V_{I} - α_{2} L_{2} - δ_{L_{2}} L_{2} = 0, \\ {(1 - f}_{2}) {\bar{β}}_{2} T_{2} V_{I} + α_{2} L_{2} - δ_{I_{2}} I_{2} = 0, \\ {\bar{N}}_{1} δ_{I_{1}} I_{1} + {\bar{N}}_{2} δ_{I_{2}} I_{2} - c_{1 V} V_{I} = 0, \end{cases}$ (3-1)

由此得到 $I_{1} = I_{2} = 0$ ， $L_{1} = L_{2} = 0$ ， $T_{1} = \frac{λ_{1}}{d_{T_{1}}}$ ， $T_{2} = \frac{λ_{2}}{d_{T_{2}}}$ ，所以无病平衡点 $E^{0} = (\frac{λ_{1}}{d_{T_{1}}}, 0, 0, \frac{λ_{2}}{d_{T_{2}}}, 0, 0, 0)$ 。

接着我们定义基本再生数 $R_{0}$ ，它表示每个HIV病毒的存活期间所能感染健康T细胞和健康巨噬细胞的平均数量：

$R_{0} = \frac{{\bar{N}}_{1} {\bar{β}}_{1} [α_{1} + (1 - f_{1}) δ_{L_{1}}] T_{1}^{0}}{c_{1 V} (δ_{L_{1}} + α_{1})} + \frac{{\bar{N}}_{2} {\bar{β}}_{2} [α_{2} + (1 - f_{2}) δ_{L_{2}}] T_{2}^{0}}{c_{2 V} (δ_{L_{2}} + α_{2})},$ (3-2)

下面求解感染平衡点，此时 $V_{I} \neq 0$ ，表示为 $V_{I}^{*}$ ，通过化简求解方程组(3-1)得到：

$T_{1}^{*} = \frac{λ_{1}}{d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}}$ ； $T_{2}^{*} = \frac{λ_{2}}{d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}}$ ；

$L_{1}^{*} = \frac{f_{1} {\bar{β}}_{1} T_{1}^{*} V_{I}^{*}}{α_{1} + δ_{L_{1}}}$ ； $L_{2}^{*} = \frac{f_{2} {\bar{β}}_{2} T_{2}^{*} V_{I}^{*}}{α_{2} + δ_{L_{2}}}$ ；

$I_{1}^{*} = \frac{λ_{1} α_{1} {\bar{β}}_{1} V_{I}^{*} + λ_{1} δ_{L_{1}} (1 - f_{1}) {\bar{β}}_{1} V_{I}^{*}}{δ_{I_{1}} (α_{1} + δ_{L_{1}}) [d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}]}$ ； $I_{2}^{*} = \frac{λ_{2} α_{2} {\bar{β}}_{2} V_{I}^{*} + λ_{2} δ_{L_{2}} (1 - f_{2}) {\bar{β}}_{2} V_{I}^{*}}{δ_{I_{2}} (α_{2} + δ_{L_{2}}) [d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}]}$ ；

其中 $V_{I}^{*}$ 为方程：

${\bar{N}}_{1} \frac{λ_{1} α_{1} {\bar{β}}_{1} + λ_{1} δ_{L_{1}} (1 - f_{1}) {\bar{β}}_{1}}{(α_{1} + δ_{L_{1}}) [d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}]} + {\bar{N}}_{2} \frac{λ_{2} α_{2} {\bar{β}}_{2} + λ_{2} δ_{L_{2}} (1 - f_{2}) {\bar{β}}_{2}}{(α_{2} + δ_{L_{2}}) [d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}]} - c_{1 V} = 0,$ (3-3)

的解。

记

$G (V_{I}) = {\bar{N}}_{1} \frac{λ_{1} α_{1} {\bar{β}}_{1} + λ_{1} δ_{L_{1}} (1 - f_{1}) {\bar{β}}_{1}}{(α_{1} + δ_{L_{1}}) [d_{T_{1}} + {\bar{β}}_{1} V_{I}]} + {\bar{N}}_{2} \frac{λ_{2} α_{2} {\bar{β}}_{2} + λ_{2} δ_{L_{2}} (1 - f_{2}) {\bar{β}}_{2}}{(α_{2} + δ_{L_{2}}) [d_{T_{2}} + {\bar{β}}_{2} V_{I}]} - c_{1 V},$

下面证明当 $R_{0} > 1$ 时模型(2-1)在 $(0, + \infty)$ 仅有一个感染平衡点 $E^{*}$ ：

$\begin{matrix} G (0) = {\bar{N}}_{1} \frac{λ_{1} α_{1} {\bar{β}}_{1} + λ_{1} δ_{L_{1}} (1 - f_{1}) {\bar{β}}_{1}}{(α_{1} + δ_{L_{1}}) d_{T_{1}}} + {\bar{N}}_{2} \frac{λ_{2} α_{2} {\bar{β}}_{2} + λ_{2} δ_{L_{2}} (1 - f_{2}) {\bar{β}}_{2}}{(α_{2} + δ_{L_{2}}) d_{T_{2}}} - c_{1 V} \\ = \frac{{\bar{N}}_{1} {\bar{β}}_{1} [α_{1} + (1 - f_{1}) δ_{L_{1}}] T_{1}^{0}}{(δ_{L_{1}} + α_{1})} + \frac{{\bar{N}}_{2} {\bar{β}}_{2} [α_{2} + (1 - f_{2}) δ_{L_{2}}] T_{2}^{0}}{(δ_{L_{2}} + α_{2})} - c_{1 V} = c_{1 V} (R_{0} - 1) \end{matrix}$

所以当 $R_{0} > 1$ 时，有 $T (0) > 0$ ，又

$G^{'} (V_{I}) = - {\bar{N}}_{1} (α_{1} + δ_{L_{1}}) {\bar{β}}_{1}^{2} \frac{λ_{1} α_{1} + λ_{1} δ_{L_{1}} (1 - f_{1})}{{(α_{1} + δ_{L_{1}})}^{2} {[d_{T_{1}} + {\bar{β}}_{1} V_{I}]}^{2}} - {\bar{N}}_{2} (α_{2} + δ_{L_{2}}) {\bar{β}}_{2}^{2} \frac{λ_{2} α_{2} + λ_{2} δ_{L_{2}} (1 - f_{2})}{{(α_{2} + δ_{L_{2}})}^{2} {[d_{T_{2}} + {\bar{β}}_{2} V_{I}]}^{2}} < 0,$

所以 $R_{0} > 1$ 时，在 $(0, + \infty)$ 内 $V_{I}^{*}$ 有唯一解，即当 $R_{0} > 1$ 时模型(2-1)在 $(0, + \infty)$ 仅有一个感染平衡点 $E^{*}$ 。

所以当 $R_{0} > 1$ 时，模型(2-1)的感染平衡点为： $E^{*} = (T_{1}^{*}, L_{1}^{*}, I_{1}^{*}, T_{2}^{*}, L_{2}^{*}, I_{2}^{*}, V_{I}^{*})$ 。

4. 平衡点稳定性分析

4.1. 无病平衡点的局部稳定性

定理2: 当 $R_{0} < 1$ 时，无病平衡点 $E^{0}$ 是局部渐近稳定的。当 $R_{0} > 1$ 时无病平衡点 $E^{0}$ 是不稳定的。

证明：模型(2-1)在无病平衡点 $E^{0}$ 处的Jacobi矩阵为：

$J (E^{0}) = [\begin{matrix} - d_{T_{1}} & 0 & 0 & 0 & 0 & 0 & - {\bar{β}}_{1} T_{1}^{0} \\ 0 & - α_{1} - δ_{L_{1}} & 0 & 0 & 0 & 0 & f_{1} {\bar{β}}_{1} T_{1}^{0} \\ 0 & α_{1} & - δ_{I_{1}} & 0 & 0 & 0 & (1 - f_{1}) {\bar{β}}_{1} T_{1}^{0} \\ 0 & 0 & 0 & - d_{T_{2}} & 0 & 0 & - {\bar{β}}_{2} T_{2}^{0} \\ 0 & 0 & 0 & 0 & - α_{2} - δ_{L_{2}} & 0 & f_{2} {\bar{β}}_{2} T_{2}^{0} \\ 0 & 0 & 0 & 0 & α_{2} & - δ_{I_{2}} & (1 - f_{2}) {\bar{β}}_{2} T_{2}^{0} \\ 0 & 0 & {\bar{N}}_{1} δ_{I_{1}} & 0 & 0 & {\bar{N}}_{2} δ_{I_{2}} & - c_{1 V} \end{matrix}]$

令

$| ξ E - J (E^{0}) | = | \begin{matrix} ξ + d_{T_{1}} & 0 & 0 & 0 & 0 & 0 & {\bar{β}}_{1} T_{1}^{0} \\ 0 & ξ + α_{1} + δ_{L_{1}} & 0 & 0 & 0 & 0 & - f_{1} {\bar{β}}_{1} T_{1}^{0} \\ 0 & - α_{1} & ξ + δ_{I_{1}} & 0 & 0 & 0 & - (1 - f_{1}) {\bar{β}}_{1} T_{1}^{0} \\ 0 & 0 & 0 & ξ + d_{T_{2}} & 0 & 0 & {\bar{β}}_{2} T_{2}^{0} \\ 0 & 0 & 0 & 0 & ξ + α_{2} + δ_{L_{2}} & 0 & - f_{2} {\bar{β}}_{2} T_{2}^{0} \\ 0 & 0 & 0 & 0 & - α_{2} & ξ + δ_{I_{2}} & - (1 - f_{2}) {\bar{β}}_{2} T_{2}^{0} \\ 0 & 0 & - {\bar{N}}_{1} δ_{I_{1}} & 0 & 0 & - {\bar{N}}_{2} δ_{I_{2}} & ξ + c_{1 V} \end{matrix} | = 0$

得到特征方程：

$\begin{array}{l} (ξ + d_{T_{1}}) (ξ + d_{T_{2}}) (ξ + α_{1} + δ_{L_{1}}) (ξ + α_{2} + δ_{L_{2}}) [(ξ + δ_{I_{1}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V}) \\ - {\bar{N}}_{1} δ_{I_{1}} (ξ + δ_{I_{2}}) (1 - f_{1}) {\bar{β}}_{1} T_{1}^{0} - {\bar{N}}_{2} δ_{I_{2}} (ξ + δ_{I_{1}}) (1 - f_{2}) {\bar{β}}_{2} T_{2}^{0}] \\ - (ξ + d_{T_{1}}) (ξ + d_{T_{2}}) (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) {\bar{N}}_{1} f_{1} {\bar{β}}_{1} α_{1} δ_{I_{1}} T_{1}^{0} \\ - (ξ + d_{T_{1}}) (ξ + d_{T_{2}}) (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) {\bar{N}}_{2} f_{2} {\bar{β}}_{2} α_{2} δ_{I_{2}} T_{2}^{0} = 0 \end{array}$

整理得：

$(ξ + d_{T_{1}}) (ξ + d_{T_{2}}) (ξ + α_{1} + δ_{L_{1}}) (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{1}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V}) q_{1} (ξ) = 0,$

其中

$\begin{matrix} q_{1} (ξ) = 1 - \frac{{\bar{N}}_{1} (1 - f_{1}) δ_{I_{1}} {\bar{β}}_{1} T_{1}^{0}}{(ξ + δ_{I_{1}}) (ξ + c_{1 V})} - \frac{{\bar{N}}_{2} (1 - f_{2}) δ_{I_{2}} {\bar{β}}_{2} T_{2}^{0}}{(ξ + δ_{I_{2}}) (ξ + c_{1 V})} - \frac{{\bar{N}}_{1} α_{1} f_{1} δ_{I_{1}} {\bar{β}}_{1} T_{1}^{0}}{(ξ + δ_{I_{1}}) (ξ + α_{1} + δ_{L_{1}}) (ξ + c_{1 V})} \\ - \frac{{\bar{N}}_{2} α_{2} f_{2} δ_{I_{2}} {\bar{β}}_{2} T_{2}^{0}}{(ξ + δ_{I_{2}}) (ξ + α_{2} + δ_{L_{2}}) (ξ + c_{1 V})} \end{matrix}$

显然特征方程的根取决于 $q_{1} (ξ) = 0$ 的根，由此确定 $E^{0}$ 的局部稳定性。

1) 当 $R_{0} > 1$ 时

$\begin{matrix} q_{1} (0) = 1 - \frac{{\bar{N}}_{1} (1 - f_{1}) {\bar{β}}_{1} T_{1}^{0}}{c_{1 V}} - \frac{{\bar{N}}_{2} (1 - f_{2}) {\bar{β}}_{2} T_{2}^{0}}{c_{1 V}} - \frac{{\bar{N}}_{1} α_{1} f_{1} {\bar{β}}_{1} T_{1}^{0}}{(α_{1} + δ_{L_{1}}) c_{1 V}} - \frac{{\bar{N}}_{2} α_{2} f_{2} {\bar{β}}_{2} T_{2}^{0}}{(α_{2} + δ_{L_{2}}) c_{1 V}} \\ = 1 - \frac{{\bar{N}}_{1} {\bar{β}}_{1} T_{1}^{0} [α_{1} + (1 - f_{1}) δ_{L_{1}}]}{c_{1 V} (α_{1} + δ_{L_{1}})} - \frac{{\bar{N}}_{2} {\bar{β}}_{2} T_{2}^{0} [α_{2} + (1 - f_{2}) δ_{L_{2}}]}{c_{1 V} (α_{2} + δ_{L_{2}})} = 1 - R_{0} < 0 \end{matrix}$

又 $\lim_{ξ \to + \infty} q_{1} (ξ) = 1 > 0$ ，故 $q_{1} (ξ)$ 在 $(0, + \infty)$ 上至少有一个正根，从而特征方程存在正实部的根，所以 $E^{0}$ 在 $R_{0} > 1$ 时是不稳定的。

2) 当 $R_{0} < 1$ 时

假设存在 $R e (ξ) \geq 0$ ，且满足 $q_{1} (ξ) = 0$ ，则

$\begin{matrix} 1 = \frac{{\bar{N}}_{1} (1 - f_{1}) δ_{I_{1}} {\bar{β}}_{1} T_{1}^{0}}{(ξ + δ_{I_{1}}) (ξ + c_{1 V})} + \frac{{\bar{N}}_{2} (1 - f_{2}) δ_{I_{2}} {\bar{β}}_{2} T_{2}^{0}}{(ξ + δ_{I_{2}}) (ξ + c_{1 V})} + \frac{{\bar{N}}_{1} α_{1} f_{1} δ_{I_{1}} {\bar{β}}_{1} T_{1}^{0}}{(ξ + δ_{I_{1}}) (ξ + α_{1} + δ_{L_{1}}) (ξ + c_{1 V})} + \frac{{\bar{N}}_{2} α_{2} f_{2} δ_{I_{2}} {\bar{β}}_{2} T_{2}^{0}}{(ξ + δ_{I_{2}}) (ξ + α_{2} + δ_{L_{2}}) (ξ + c_{1 V})} \\ \leq | \frac{{\bar{N}}_{1} (1 - f_{1}) {\bar{β}}_{1} T_{1}^{0}}{c_{1 V}} + \frac{{\bar{N}}_{2} (1 - f_{2}) {\bar{β}}_{2} T_{2}^{0}}{c_{1 V}} + \frac{{\bar{N}}_{1} α_{1} f_{1} {\bar{β}}_{1} T_{1}^{0}}{(α_{1} + δ_{L_{1}}) c_{1 V}} + \frac{{\bar{N}}_{2} α_{2} f_{2} {\bar{β}}_{2} T_{2}^{0}}{(α_{2} + δ_{L_{2}}) c_{1 V}} | \\ = | \frac{{\bar{N}}_{1} {\bar{β}}_{1} T_{1}^{0} [α_{1} + (1 - f_{1}) δ_{L_{1}}]}{c_{1 V} (α_{1} + δ_{L_{1}})} + \frac{{\bar{N}}_{2} {\bar{β}}_{2} T_{2}^{0} [α_{2} + (1 - f_{2}) δ_{L_{2}}]}{c_{1 V} (α_{2} + δ_{L_{2}})} | \\ \leq | \frac{{\bar{N}}_{1} {\bar{β}}_{1} T_{1}^{0} [α_{1} + (1 - f_{1}) δ_{L_{1}}]}{c_{1 V} (α_{1} + δ_{L_{1}})} | + | \frac{{\bar{N}}_{2} {\bar{β}}_{2} T_{2}^{0} [α_{2} + (1 - f_{2}) δ_{L_{2}}]}{c_{1 V} (α_{2} + δ_{L_{2}})} | = R_{0} \end{matrix}$

与 $R_{0} < 1$ 矛盾，由此可以得到 $R e (ξ) < 0$ ，其中 $q_{1} (ξ) = 0$ 。所以当 $R_{0} < 1$ 时无病平衡点 $E^{0}$ 是局部渐近稳定的[23]。

综上得出 $R_{0} < 1$ 时，无病平衡点 $E^{0}$ 是局部渐近稳定的；当 $R_{0} > 1$ 时无病平衡点 $E^{0}$ 是不稳定的。

4.2. 感染平衡点的局部稳定性

定理3： $R_{0} > 1$ 时，感染平衡点 $E^{*}$ 局部渐近稳定。

证明：模型(2-1)在感染平衡点 $E^{*}$ 处的Jacobi矩阵为：

$J (E^{*}) = [\begin{matrix} - d_{T_{1}} - {\bar{β}}_{1} V_{I}^{*} & 0 & 0 & 0 & 0 & 0 & - {\bar{β}}_{1} T_{1}^{*} \\ f_{1} {\bar{β}}_{1} V_{I}^{*} & - α_{1} - δ_{L_{1}} & 0 & 0 & 0 & 0 & f_{1} {\bar{β}}_{1} T_{1}^{*} \\ (1 - f_{1}) {\bar{β}}_{1} V_{I}^{*} & α_{1} & - δ_{I_{1}} & 0 & 0 & 0 & (1 - f_{1}) {\bar{β}}_{1} T_{1}^{*} \\ 0 & 0 & 0 & - d_{T_{2}} - {\bar{β}}_{2} V_{I}^{*} & 0 & 0 & - {\bar{β}}_{2} T_{2}^{*} \\ 0 & 0 & 0 & f_{2} {\bar{β}}_{2} V_{I}^{*} & - α_{2} - δ_{L_{2}} & 0 & f_{2} {\bar{β}}_{2} T_{2}^{*} \\ 0 & 0 & 0 & (1 - f_{2}) {\bar{β}}_{2} V_{I}^{*} & α_{2} & - δ_{I_{2}} & (1 - f_{2}) {\bar{β}}_{2} T_{2}^{*} \\ 0 & 0 & {\bar{N}}_{1} δ_{I_{1}} & 0 & 0 & {\bar{N}}_{2} δ_{I_{2}} & - c_{1 V} \end{matrix}]$

令

$\begin{array}{l} | ξ E - J (E^{*}) | \\ = | \begin{matrix} ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*} & 0 & 0 & 0 & 0 & 0 & {\bar{β}}_{1} T_{1}^{*} \\ - f_{1} {\bar{β}}_{1} V_{I}^{*} & ξ + α_{1} + δ_{L_{1}} & 0 & 0 & 0 & 0 & - f_{1} {\bar{β}}_{1} T_{1}^{*} \\ - (1 - f_{1}) {\bar{β}}_{1} V_{I}^{*} & - α_{1} & ξ + δ_{I_{1}} & 0 & 0 & 0 & - (1 - f_{1}) {\bar{β}}_{1} T_{1}^{*} \\ 0 & 0 & 0 & ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*} & 0 & 0 & {\bar{β}}_{2} T_{2}^{*} \\ 0 & 0 & 0 & - f_{2} {\bar{β}}_{2} V_{I}^{*} & ξ + α_{2} + δ_{L_{2}} & 0 & - f_{2} {\bar{β}}_{2} T_{2}^{*} \\ 0 & 0 & 0 & - (1 - f_{2}) {\bar{β}}_{2} V_{I}^{*} & - α_{2} & ξ + δ_{I_{2}} & - (1 - f_{2}) {\bar{β}}_{2} T_{2}^{*} \\ 0 & 0 & - {\bar{N}}_{1} δ_{I_{1}} & 0 & 0 & - {\bar{N}}_{2} δ_{I_{2}} & ξ + c_{1 V} \end{matrix} | = 0 \end{array}$

得到特征方程：

$\begin{array}{l} [ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V}) \\ - [ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (1 - f_{2}) {\bar{β}}_{2} T_{2}^{*} {\bar{N}}_{2} δ_{I_{2}} \\ - [ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] f_{2} {\bar{β}}_{2} T_{2}^{*} α_{2} {\bar{N}}_{2} δ_{I_{2}} \\ - [ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) {\bar{β}}_{2} T_{2}^{*} {\bar{N}}_{2} δ_{I_{2}} [α_{2} f_{2} {\bar{β}}_{2} V_{I}^{*} + (ξ + α_{2} + δ_{L_{2}}) (1 - f_{2}) {\bar{β}}_{2} V_{I}^{*}] \\ - [ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] {\bar{β}}_{1} T_{1}^{*} [α_{1} f_{1} - (ξ + α_{1} + δ_{L_{1}}) (1 - f_{1})] {\bar{N}}_{1} δ_{I_{1}} [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) \\ - {\bar{β}}_{1}^{2} T_{1}^{*} V_{I}^{*} [α_{1} f_{1} - (ξ + α_{1} + δ_{L_{1}}) (1 - f_{1})] {\bar{N}}_{1} δ_{I_{1}} [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) = 0 \end{array}$

整理得：

$[ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V}) p (ξ) = 0,$

其中

$\begin{matrix} p (ξ) = 1 + \frac{{\bar{N}}_{1} {\bar{β}}_{1} δ_{I_{1}} (1 - f_{1}) T_{1}^{*}}{(ξ + δ_{I_{1}}) (ξ + c_{1 V})} - \frac{{\bar{N}}_{2} {\bar{β}}_{2} δ_{I_{2}} (1 - f_{2}) T_{2}^{*}}{(ξ + δ_{I_{2}}) (ξ + c_{1 V})} \\ + \frac{{\bar{N}}_{1} {\bar{β}}_{1}^{2} δ_{I_{1}} (1 - f_{1}) T_{1}^{*} V_{I}^{*}}{[ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + δ_{I_{1}}) (ξ + c_{1 V})} - \frac{{\bar{N}}_{2} {\bar{β}}_{2}^{2} δ_{I_{2}} (1 - f_{2}) T_{2}^{*} V_{I}^{*}}{[ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + δ_{I_{2}}) (ξ + c_{1 V})} \\ - \frac{{\bar{N}}_{1} α_{1} {\bar{β}}_{1}^{2} δ_{I_{1}} f_{1} T_{1}^{*} V_{I}^{*}}{[ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) (ξ + c_{1 V})} \\ - \frac{{\bar{N}}_{2} α_{2} {\bar{β}}_{2}^{2} δ_{I_{2}} f_{2} T_{2}^{*} V_{I}^{*}}{[ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V})} \\ - \frac{{\bar{N}}_{1} α_{1} {\bar{β}}_{1} δ_{I_{1}} f_{1} T_{1}^{*}}{(ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) (ξ + c_{1 V})} - \frac{{\bar{N}}_{2} α_{2} {\bar{β}}_{2} δ_{I_{2}} f_{2} T_{2}^{*}}{(ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V})} \end{matrix}$

下面特征方程的根取决于 $p (ξ) = 0$ 的根，假设存在 $R e (ξ) \geq 0$ ，令 $p (ξ) = 0$ ，由前面计算得到的感染平衡点 $T_{1}^{*} = \frac{λ_{1}}{d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}}$ ， $T_{2}^{*} = \frac{λ_{2}}{d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}}$ ；代入整理得：

$\begin{matrix} 1 = \frac{{\bar{N}}_{1} δ_{I_{1}} {\bar{β}}_{1} [λ_{1} α_{1} + λ_{1} δ_{L_{1}} (1 - f_{1})] [ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}]}{[ξ + d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) (ξ + c_{1 V}) [d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}]} \\ + \frac{{\bar{N}}_{2} δ_{I_{2}} {\bar{β}}_{2} [λ_{2} α_{2} + λ_{2} δ_{L_{2}} (1 - f_{2})] [ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}]}{[ξ + d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V}) [d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}]} \\ = \frac{{\bar{N}}_{1} δ_{I_{1}} {\bar{β}}_{1} [λ_{1} α_{1} + λ_{1} δ_{L_{1}} (1 - f_{1})]}{(ξ + α_{1} + δ_{L_{1}}) (ξ + δ_{I_{1}}) (ξ + c_{1 V}) [d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}]} \\ + \frac{{\bar{N}}_{2} δ_{I_{2}} {\bar{β}}_{2} [λ_{2} α_{2} + λ_{2} δ_{L_{2}} (1 - f_{2})]}{(ξ + α_{2} + δ_{L_{2}}) (ξ + δ_{I_{2}}) (ξ + c_{1 V}) [d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}]} \\ < | \frac{{\bar{N}}_{1} λ_{1} δ_{I_{1}} {\bar{β}}_{1} [α_{1} + δ_{L_{1}} (1 - f_{1})]}{[d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (α_{1} + δ_{L_{1}}) δ_{I_{1}} c_{1 V}} | + | \frac{{\bar{N}}_{2} λ_{2} δ_{I_{2}} {\bar{β}}_{2} [α_{2} + δ_{L_{2}} (1 - f_{2})]}{[d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (α_{2} + δ_{L_{2}}) δ_{I_{2}} c_{1 V}} | \end{matrix}$

再根据(3-1)的第六个式子 ${\bar{N}}_{1} δ_{I_{1}} I_{1}^{*} + {\bar{N}}_{2} δ_{I_{2}} I_{2}^{*} - c_{1 V} V_{I}^{*} = 0$ 得到：

$\frac{{\bar{N}}_{1} δ_{I_{1}} I_{1}^{*} + {\bar{N}}_{2} δ_{I_{2}} I_{2}^{*}}{c_{1 V} V_{I}^{*}} = 1,$ (4-1)

而感染平衡点 $I_{1}^{*} = \frac{λ_{1} α_{1} {\bar{β}}_{1} V_{I}^{*} + λ_{1} δ_{L_{1}} (1 - f_{1}) {\bar{β}}_{1} V_{I}^{*}}{δ_{I_{1}} (α_{1} + δ_{L_{1}}) [d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}]}$ ， $I_{2}^{*} = \frac{λ_{2} α_{2} {\bar{β}}_{2} V_{I}^{*} + λ_{2} δ_{L_{2}} (1 - f_{2}) {\bar{β}}_{2} V_{I}^{*}}{δ_{I_{2}} (α_{2} + δ_{L_{2}}) [d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}]}$ 带入(4-1)得：

$\frac{{\bar{N}}_{1} λ_{1} δ_{I_{1}} {\bar{β}}_{1} [α_{1} + δ_{L_{1}} (1 - f_{1})]}{[d_{T_{1}} + {\bar{β}}_{1} V_{I}^{*}] (α_{1} + δ_{L_{1}}) δ_{I_{1}} c_{1 V}} + \frac{{\bar{N}}_{2} λ_{2} δ_{I_{2}} {\bar{β}}_{2} [α_{2} + δ_{L_{2}} (1 - f_{2})]}{[d_{T_{2}} + {\bar{β}}_{2} V_{I}^{*}] (α_{2} + δ_{L_{2}}) δ_{I_{2}} c_{1 V}} = 1,$

即推出 $1 < 1$ ，产生矛盾，因此有 $R e (ξ) < 0$ ，从而定理3得证。

综上，当 $R_{0} > 1$ 时，感染平衡点 $E^{*}$ 局部渐近稳定。

5. 数值模拟

在进行了相关理论证明后，为了更直观地观察模型的动力学特征，本节对模型(2-1)进行数值模拟实验。模型中相关参数的取值如表2所示

Table 2. Parameter table

表2. 参数表

参数	含义	单位	取值	参考
$λ_{1}$	未感染T细胞产生速率	mm⁻³ day⁻¹	10	[8] -[10]
$λ_{2}$	未感染巨噬细胞产生速率	mm⁻³ day⁻¹	0.04	[11] [12]
$d_{T_{1}}$	未感染T细胞死亡率	day⁻¹	0.01	[10] [13] [14]
$d_{T_{2}}$	未感染巨噬细胞死亡率	day⁻¹	0.01	[11] [12]
$β_{1}$	T细胞病毒感染率	mm³ day⁻¹	10⁻⁵~0.5	[15]
$β_{2}$	巨噬细胞病毒感染率	mm³ day⁻¹	10⁻⁵~0.5	[15] [16]
$f_{1}$	感染的T细胞转为潜伏感染的比例	day⁻¹	0.1	[17]
$f_{2}$	感染的巨噬细胞转为潜伏感染的比例	day⁻¹	0.1	[17]
$α_{1}$	潜伏感染T细胞激活率	day⁻¹	0.2	[18]
$α_{2}$	潜伏感染巨噬细胞激活率	day⁻¹	0.05	[19]
$δ_{L_{1}}$	潜伏感染T细胞死亡率	day⁻¹	0.02	[9] [19]
$δ_{L_{2}}$	潜伏感染巨噬细胞死亡率	day⁻¹	0.1	[19]
$δ_{I_{1}}$	感染T细胞死亡率	day⁻¹	0.5	[20][22] -
$δ_{I_{2}}$	感染巨噬细胞死亡率	day⁻¹	0.1	[22]
$N_{1}$	T细胞内病毒裂解量	virions cells⁻¹	15	[19]
$N_{2}$	巨噬细胞内病毒裂解量	virions cells⁻¹	10	[19]
$c_{1 V}, c_{2 V}$	病毒清除速率	day⁻¹	3	[10]
$ε_{1 R T}, ε_{2 R T}$	逆转录酶抑制剂抑制率	—	(0, 1)	—
$ε_{1 P I}, ε_{2 P I}$	蛋白酶抑制剂抑制率	—	(0, 1)	—

本章所选取不同的初值如下所示：

第一组初值：

$(T_{1} (0), L_{1} (0), I_{1} (0), T_{2} (0), L_{2} (0), I_{2} (0), V_{I} (0)) = (900, 2, 7, 2, 0.15, 0.5, 18)$

第二组初值：

$(T_{1} (0), L_{1} (0), I_{1} (0), T_{2} (0), L_{2} (0), I_{2} (0), V_{I} (0)) = (700, 1.5, 5, 1.5, 0.0005, 0.0005, 15)$

第三组初值：

$(T_{1} (0), L_{1} (0), I_{1} (0), T_{2} (0), L_{2} (0), I_{2} (0), V_{I} (0)) = (500, 1, 3, 1, 0.001, 0.001, 12)$

第四组初值：

$(T_{1} (0), L_{1} (0), I_{1} (0), T_{2} (0), L_{2} (0), I_{2} (0), V_{I} (0)) = (300, 0.5, 1.5, 0.5, 0.1, 0.3, 8)$

在上述初值条件下，选择不同的 $β_{1}$ ， $β_{2}$ ， $ε_{1 R T}$ ， $ε_{2 R T}$ ， $ε_{1 P I}$ ， $ε_{2 P I}$ 会得到不同的数值结果。

5.1. 平衡点稳定性的数值模拟

为验证模型(2-1)的无病平衡点 $E^{0}$ 与感染平衡点 $E^{*}$ ，首先不考虑抑制剂的作用，即取 $ε_{1 R T} = ε_{2 R T} = ε_{1 P I} = ε_{2 P I} = 0$ 。

当取 $β_{1} = 0.0001$ ， $β_{2} = 0.015$ 时， $R_{0} = 0.6821 < 1$ ，如图1所示，四组不同初值开始的曲线轨迹最后均趋于无病平衡点 $E^{0} (1000, 0, 0, 4, 0, 0, 0)$ 。

当取 $β_{1} = 0.0006$ ， $β_{2} = 0.001$ 时， $R_{0} = 2.9852 > 1$ ，如图2所示，四组不同初值开始的曲线轨迹最后均趋于点 $(335.4188, 3.0208, 13.1704, 0.9297, 0.0205, 0.2866, 33.0290)$ 。

将 $V_{I} = 33.0290$ 代入方程(3-3)左侧，结果为−3.8943*10⁻⁴，接近为0，可以说明V_I是方程(3-3)的解，即 $V_{I}^{*} = 33.0290$ 。将 $V_{I}^{*}$ 进一步代入其余变量的感染平衡点计算式计算，并与数值稳定点结果对比，验证得四组不同初值开始的曲线轨迹最后趋于感染平衡点 $E^{*} (335.4188, 3.0208, 13.1704, 0.9297, 0.0205, 0.2866, 33.0290)$ 。

Figure 1. Local stability of the equilibrium point when β₁ = 0.0001, β₂ = 0.015

图1. β₁ = 0.0001，β₂ = 0.015时平衡点的局部稳定性

Figure 2. Local stability of the equilibrium point when β₁ = 0.0006, β₂ = 0.001

图2. β₁ = 0.0006，β₂ = 0.001时平衡点的局部稳定性

5.2. 不同抑制剂效能的数值模拟

由式(3-2)，抑制剂的抑制效能也会影响R₀的大小，故接下来考虑逆转录酶抑制剂和蛋白酶抑制剂的抑制效能对R₀的影响。在取 $β_{1} = 0.0006$ ， $β_{2} = 0.001$ 的情况下，R₀的取值随 $ε_{R T}$ 与 $ε_{P I}$ 的变化如图3所示。其中 $ε_{R T}$ 与 $ε_{P I}$ 不同取值所形成的红色等高线为使 $R_{0} = 1$ 的边界线。

Figure 3. The value of R₀ varies with the changes in ε_RT and ε_PI

图3. R₀的取值随ε_RT与ε_PI的变化

如图3，四种 $ε_{R T}$ 与 $ε_{P I}$ 的取值组合如表3所示：

Table 3. Nhibitor efficacy values

表3. 抑制剂效能取值

序号	逆转录酶抑制剂ε_RT	蛋白酶抑制剂ε_PI	R₀
A	$ε_{1 R T} = ε_{2 R T} = 0.13$	$ε_{1 P I} = ε_{2 P I} = 0.23$	1.9998
B	$ε_{1 R T} = ε_{2 R T} = 0.23$	$ε_{1 P I} = ε_{2 P I} = 0.33$	1.5401
C	$ε_{1 R T} = ε_{2 R T} = 0.33$	$ε_{1 P I} = ε_{2 P I} = 0.43$	1.1400
D	$ε_{1 R T} = ε_{2 R T} = 0.43$	$ε_{1 R T} = ε_{2 R T} = 0.53$	0.7997

在 $(T_{1} (0), L_{1} (0), I_{1} (0), T_{2} (0), L_{2} (0), I_{2} (0), V_{I} (0)) = (500, 1, 3, 1, 0.001, 0.001, 12)$ 的初值条件下，V_I，T₁和T₂在四种抑制剂效能取值的轨迹如图4所示。

从数值结果可以发现，在同样的初值条件下，当 $R_{0} > 1$ 时，随着R₀的减小，V_I的感染平衡点逐渐降低；当 $R_{0} < 1$ 时，V_I的轨迹最终稳定在0，即代表着在该抑制剂效能的组合下，病毒可以被完全清除。对于靶向细胞CD4⁺T细胞T₁和巨噬细胞T₂，随着R₀的减小，T₁和T₂轨迹的稳定值越高，即代表着健康的靶细胞数量越多，当R₀减小到小于1时，T₁和T₂的轨迹最终都回归到无病平衡点的状态。

数值结果表明，只要取在 $R_{0} = 1$ 边界线右侧的不同 $ε_{R T}$ 与 $ε_{P I}$ 组合，各细胞和病毒就能够回归到无病平衡点的状态。

Figure 4. Trajectory plot of the four types of inhibitor efficacy values

图4. 四种抑制剂效能取值的轨迹图

6. 结论

在本研究中，我们从HIV的传播动力学出发，构建了一个包含潜伏感染的多维度常微分方程组模型，覆盖了健康和感染的CD4⁺T细胞及巨噬细胞。第一部分详细介绍了构建HIV动力学模型所需的生物学背景和数学理论，并回顾了相关研究成果，为模型的理解和后续开发奠定了基础。在第二部分中，模型通过引入非线性感染率和细胞恢复概念，对HIV的传播及其在药物治疗影响下的行为进行了深入分析，计算了无病和感染平衡点，并引出了基本再生数这一关键参数。第三部分通过求解平衡点处的Jacobi矩阵和特征方程，对模型的稳定性进行了严格的数学证明，确认了在不同基本再生数条件下，无病和感染平衡点的稳定性。第四部分使用MATLAB进行数值模拟，验证了理论分析的准确性，并根据模拟结果探讨了不同药物治疗策略的影响，突出了在有感染细胞可能恢复为健康细胞的情况下，药物治疗对抑制病毒传播的重要性。

整体而言，本文的研究不仅增进了对HIV传播动力学的理解，还为未来的研究提供了新的视角和方法，特别是在改进模型以包括更多类型的治疗策略，以期更全面地反映实际治疗场景中的HIV感染与控制过程。通过这种综合性的研究方法，可以为制定更有效的HIV预防和治疗策略提供更为科学的依据。

在未来的研究中，应更全面地考虑模型的拓展与改进，特别是纳入多种药物治疗策略，如新型抗病毒药物的影响，以更精确地模拟HIV在不同治疗条件下的传播动力学。此外，为了使模型更贴合实际治疗场景，将病人的经济状况和身体承受能力整合入模型也显得尤为重要。长期的临床数据收集与分析将有助于调整和验证模型，确保其预测结果的准确性和实用性。这样的综合研究不仅能深化我们对HIV传播机制的理解，还能为制定更有效的预防和治疗策略提供科学依据，从而在全球范围内更有效地控制HIV的传播。

NOTES

^*通讯作者。

参考文献

[1]	Clark, I.C., Wheeler, M.A., Lee, H.-G., Li, Z., Sanmarco, L.M., Thaploo, S., et al. (2023) Identification of Astrocyte Regulators by Nucleic Acid Cytometry. Nature, 614, 326-333. https://doi.org/10.1038/s41586-022-05613-0
[2]	Bacchus-Souffan, C., Fitch, M., Symons, J., Abdel-Mohsen, M., Reeves, D.B., Hoh, R., et al. (2021) Relationship between CD4 T Cell Turnover, Cellular Differentiation and HIV Persistence during Art. PLOS Pathogens, 17, e1009214. https://doi.org/10.1371/journal.ppat.1009214
[3]	Perelson, A.S., Essunger, P., Cao, Y., Vesanen, M., Hurley, A., Saksela, K., et al. (1997) Decay Characteristics of HIV-1-Infected Compartments during Combination Therapy. Nature, 387, 188-191. https://doi.org/10.1038/387188a0
[4]	Pankavich, S. (2015) The Effects of Latent Infection on the Dynamics of HIV. Differential Equations and Dynamical Systems, 24, 281-303. https://doi.org/10.1007/s12591-014-0234-6
[5]	Shen, L., Peterson, S., Sedaghat, A.R., McMahon, M.A., Callender, M., Zhang, H., et al. (2008) Dose-Response Curve Slope Sets Class-Specific Limits on Inhibitory Potential of Anti-HIV Drugs. Nature Medicine, 14, 762-766. https://doi.org/10.1038/nm1777
[6]	Wang, X., Song, X., Tang, S. and Rong, L. (2016) Dynamics of an HIV Model with Multiple Infection Stages and Treatment with Different Drug Classes. Bulletin of Mathematical Biology, 78, 322-349. https://doi.org/10.1007/s11538-016-0145-5
[7]	Edwards, C.H. and Penney, D.E. (2015) Differential Equations with Boundary-Value Problems. 5th Edition, Machinery Industry Press.
[8]	Li, F. and Ma, W. (2018) Dynamics Analysis of an HTLV-1 Infection Model with Mitotic Division of Actively Infected Cells and Delayed CTL Immune Response. Mathematical Methods in the Applied Sciences, 41, 3000-3017. https://doi.org/10.1002/mma.4797
[9]	Perelson, A.S., Kirschner, D.E. and De Boer, R. (1993) Dynamics of HIV Infection of CD4⁺ T Cells. Mathematical Biosciences, 114, 81-125. https://doi.org/10.1016/0025-5564(93)90043-a
[10]	Elaiw, A.M. and AlShamrani, N.H. (2020) HTLV/HIV Dual Infection: Modeling and Analysis. Mathematics, 9, Article 51. https://doi.org/10.3390/math9010051
[11]	Adams, B.M., Banks, H.T., Kwon, H.-D. and Tran, H.T. (2004) Dynamic Multidrug Therapies for HIV: Optimal and STI Control Approaches. Mathematical Biosciences and Engineering, 1, 223-241. https://doi.org/10.3934/mbe.2004.1.223
[12]	Adams, B.M., Banks, H.T., Davidian, M., Kwon, H., Tran, H.T., Wynne, S.N., et al. (2005) HIV Dynamics: Modeling, Data Analysis, and Optimal Treatment Protocols. Journal of Computational and Applied Mathematics, 184, 10-49. https://doi.org/10.1016/j.cam.2005.02.004
[13]	Callaway, D.S. and Perelson, A.S. (2002) HIV-1 Infection and Low Steady State Viral Loads. Bulletin of Mathematical Biology, 64, 29-64. https://doi.org/10.1006/bulm.2001.0266
[14]	Mohri, H., Bonhoeffer, S., Monard, S., Perelson, A.S. and Ho, D.D. (1998) Rapid Turnover of T Lymphocytes in SIV-Infected Rhesus Macaques. Science, 279, 1223-1227. https://doi.org/10.1126/science.279.5354.1223
[15]	Wang, Y., Zhou, Y., Brauer, F. and Heffernan, J.M. (2012) Viral Dynamics Model with CTL Immune Response Incorporating Antiretroviral Therapy. Journal of Mathematical Biology, 67, 901-934. https://doi.org/10.1007/s00285-012-0580-3
[16]	Hill, A.L., Rosenbloom, D.I.S., Nowak, M.A. and Siliciano, R.F. (2018) Insight into Treatment of HIV Infection from Viral Dynamics Models. Immunological Reviews, 285, 9-25. https://doi.org/10.1111/imr.12698
[17]	Tuckwell, H.C. and Le Corfec, E. (1998) A Stochastic Model for Early HIV-1 Population Dynamics. Journal of Theoretical Biology, 195, 451-463. https://doi.org/10.1006/jtbi.1998.0806
[18]	Wang, Y., Liu, J. and Liu, L. (2016) Viral Dynamics of an HIV Model with Latent Infection Incorporating Antiretroviral Therapy. Advances in Difference Equations, 2016, Article No. 225. https://doi.org/10.1186/s13662-016-0952-x
[19]	Elaiw, A.M., Raezah, A.A. and Azoz, S.A. (2018) Stability of Delayed HIV Dynamics Models with Two Latent Reservoirs and Immune Impairment. Advances in Difference Equations, 2018, Article No. 414. https://doi.org/10.1186/s13662-018-1869-3
[20]	Perelson, A.S., Neumann, A.U., Markowitz, M., Leonard, J.M. and Ho, D.D. (1996) HIV-1 Dynamics in vivo: Virion Clearance Rate, Infected Cell Life-Span, and Viral Generation Time. Science, 271, 1582-1586. https://doi.org/10.1126/science.271.5255.1582
[21]	Nelson, P.W., Murray, J.D. and Perelson, A.S. (2000) A Model of HIV-1 Pathogenesis that Includes an Intracellular Delay. Mathematical Biosciences, 163, 201-215. https://doi.org/10.1016/s0025-5564(99)00055-3
[22]	Perelson, A.S. and Nelson, P.W. (1999) Mathematical Analysis of HIV-1 Dynamics in Vivo. SIAM Review, 41, 3-44. https://doi.org/10.1137/s0036144598335107
[23]	Wang, X., Ge, Q. and Chen, Y. (2020) Threshold Dynamics of an HIV Infection Model with Two Distinct Cell Subsets. Applied Mathematics Letters, 103, Article 106242. https://doi.org/10.1016/j.aml.2020.106242

为你推荐

友情链接