包虫病传播模型的最优控制分析

doi:10.12677/aam.2025.147358

期刊菜单

包虫病传播模型的最优控制分析
Optimal Control Analysis of an Echinococcosis Transmission Model

DOI: 10.12677/aam.2025.147358, PDF, HTML, XML, 科研立项经费支持
作者: 俞语涵, 李金辉^*：阜阳师范大学数学与统计学院，安徽阜阳
关键词: 包虫病模型；平衡点；稳定性；持续性；最优控制；Echinococcosis Model； Equilibrium； Stability； Persistence； Optimal Control

摘要: 包虫病是一种人畜共患病，主要在最终宿主、中间宿主和最终宿主之间进行循环。例如羊作为中间宿主，由于摄入被寄生虫卵污染过的食物和水而被感染，而狗作为最终宿主可能在吃了有寄生虫幼虫的羊的内脏后被感染。本文提出了一个具有标准发病率的包虫病传播模型，得到了模型的基本再生数，研究了疾病的灭绝性与持续性。考虑屠宰场卫生检测、给狗驱虫、给绵羊接种疫苗和环境消毒等各种控制，建立控制系统并证明了最优解的存在性、刻画。

Abstract: Echinococcosis is a zoonotic disease that mainly circulates among definitive hosts, intermediate hosts and definitive hosts. For example, sheep, as intermediate hosts, are infected by ingesting food and water contaminated with parasite eggs, while dogs, as definitive hosts, may be infected after consuming the internal organs of sheep with parasitic larvae. In this paper, we propose an echinococcosis transmission model with a standard incidence rate. We obtain the basic reproduction number of the model, and study the extinction and persistence of the disease. Considering various control measures such as hygiene inspection in slaughterhouses, deworming of dogs, vaccination of sheep and environmental disinfection, a control system is established, and the existence and characterization of the optimal solution are shown.

文章引用：俞语涵, 李金辉. 包虫病传播模型的最优控制分析[J]. 应用数学进展, 2025, 14(7): 214-225. https://doi.org/10.12677/aam.2025.147358

1. 引言

包虫病是由棘球绦虫引起的一种常见人畜共患病[1]。据报道，包虫病有四种类型：由细粒棘球绦虫的幼虫引起的囊型包虫病，由多房棘球绦虫的幼虫引起的泡型包虫病。另外，有两种热带棘球蚴病类型：由伏氏棘球绦虫的幼虫引起的多囊性棘球蚴病和由少节棘球绦虫的幼虫引起的单囊型棘球蚴病。其中，囊型包虫病与人类关系最为密切而且分布广泛，存在于除南极洲以外的各大洲，泡型包虫病主要分布在北半球，特别是中国、俄罗斯、欧洲大陆和北美洲的一些国家。

包虫病的生命周期分为卵、幼虫和成虫三个发展阶段。绵羊和山羊等草食性和杂食性动物作为中间宿主。它们通过食用含有寄生虫卵的食物和水而被感染，之后卵就会在它们的内脏中逐渐发育为幼虫。狗等食肉动物(最终宿主)通过食用含有寄生虫幼虫的中间宿主的内脏而被感染，然后幼虫将在它们的肠道中会含有成虫，成虫在最终宿主体内产卵，随后虫卵会随最终宿主的粪便一起排出。通过风、水或者苍蝇，这些虫卵会传播到不同的地方，之后可能会被中间宿主摄取。人类会作为偶然的中间宿主，被感染的方式与其他中间宿主相同，但是不会传染给最终宿主。

传染病模型是在有限资源情况下寻找根除感染的控制措施的有效工具。现在已经有很多关于使用数学模型研究包虫病传播的工作成果[2]-[5]。特别地，在文献[5]中，Wang等人提出了一个具有双线性发病率的包虫病模型，并定义了基本再生数 $R_{0}$ ，并证明了如果 $R_{0} < 1$ ，包虫病的无病平衡点是全局渐近稳定的，如果 $R_{0} > 1$ ，无病平衡点是不稳定的但是地方病平衡点是稳定的。此外他们还发现人的参数并不会影响疾病的行为。

Rong等[6] [7]在分析高维包虫病动力学模型后将流浪犬认定为预防和控制包虫病的一个关键因素，并且得出了能够更好地控制包虫病的传播——减少流浪犬数量及提高羊疫苗的免疫。更多的相关研究成果参见文献[8]-[10]。在关于彭阳县的囊型包虫病的传播机理的研究之后，崔倩倩[11]建立了包虫病在人、羊、犬及环境中传播的数学模型，同时预测了彭阳县的包虫病流行趋势评估了防控措施对疾病传播的影响。在2021年，Zhao等人[12]提出了一个具有最优控制策略的包虫病传播的动力学模型，研究出了该模型的基本再生数，并利用其研究了无病平衡点和地方病平衡点的全局稳定性。接着进行了成本效益分析，得出采用家庭屠宰卫生检查和驱虫药物治疗相结合的方法，能够提供最具成本效益的控制包虫病传播的策略。本文研究具有标准发生率的包虫病传播模型，分析疾病的灭绝性与持续性及最优控制的存在性、刻画。

2. 模型建立

我们将最终宿主(通常是狗)分为三类：易感者、潜伏者和感染者，分别用 $S_{D}$ 、 $E_{D}$ 和 $I_{D}$ 表示。我们用相同的方式将中间宿主(通常是家畜)分为三类：易感者 $S_{L}$ 、潜伏者 $E_{L}$ 和感染者 $I_{L}$ 。因此，狗和家畜的总数为 $N_{D} = S_{D} + E_{D} + I_{D}$ 和 $N_{L} = S_{L} + E_{L} + I_{L}$ 。然后，建立具有标准发生率的包虫病动力学的模型：

${\begin{cases} \frac{d S_{D}}{d t} = A_{1} - \frac{β_{1} S_{D}}{N_{D}} I_{L} - d_{1} S_{D} + σ I_{D}, \\ \frac{d E_{D}}{d t} = \frac{β_{1} S_{D}}{N_{D}} I_{L} - (d_{1} + γ_{1}) E_{D}, \\ \frac{d I_{D}}{d t} = γ_{1} E_{D} - (d_{1} + α_{1} + σ) I_{D}, \\ \frac{d S_{L}}{d t} = A_{2} - \frac{β_{2} S_{L}}{N_{L}} x - d_{2} S_{L}, \\ \frac{d E_{L}}{d t} = \frac{β_{2} S_{L}}{N_{L}} x - (d_{2} + γ_{2}) E_{L}, \\ \frac{d I_{L}}{d t} = γ_{2} E_{L} - (d_{2} + α_{2}) I_{L}, \\ \frac{d x}{d t} = a I_{D} - d x, \end{cases}$ (2.1)

Figure 1. Flow chart of the model

图1. 模型示意图

其中参数均为正数，描述如表1所示。包虫病传播过程示意图如图1所示。

Table 1. Parameters and meanings

表1. 参数及含义

参数	参数含义	参数	参数含义
$A_{1}$ 、 $A_{2}$	狗、家畜的年出生量	$β_{1}$ 、 $β_{2}$	狗、牲畜的感染率
$d_{1}$ 、 $d_{2}$ 、 $d$	狗、家畜、卵的自然死亡率	$γ_{2}$	由牲畜暴露者到染病者的发展率
$γ_{1}$	由狗暴露者到染病者的发展率	$α_{1}$ 、 $α_{2}$	狗、家畜因病的死亡率
$σ$	狗的免疫失去率	$a$	感染狗的卵释放率

3. 模型分析

3.1. 解的正性和有界性

由于模型(2.1)的生物学背景，我们需要证明模型(2.1)中解的正性和有界性。

引理3.1 模型(2.1)的具有正初始条件的解 $(S_{D}, E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x)$ 对于所有 $t \geq 0$ 都是正的，且满足

$\underset{t \to \infty}{\lim \sup} N_{D} \leq A_{1} / d_{1}, \underset{t \to \infty}{\lim \sup} N_{L} \leq A_{2} / d_{2}, \underset{t \to \infty}{\lim \sup} x \leq a A_{1} / d d_{1} .$

证：设 $m (t) = \min {S_{D} (t), E_{D} (t), I_{D} (t), S_{L} (t), E_{L} (t), I_{L} (t), x (t)}$ ，则 $m (0) > 0$ 。设 $\bar{t} > 0$ ，对于所有的 $t \in [0, \bar{t}]$ ， $m (t) > 0$ 且 $m (\bar{t}) = 0$ 。如果 $m (\bar{t}) = S_{D} (\bar{t})$ ，则有 $d S_{D} (\bar{t}) / d t \leq 0$ 。

但是由模型(2.1)的第一个方程可以知道 $d S_{D} (\bar{t}) / d t = A_{1} + σ I_{D} (\bar{t}) > 0$ ，这就导致了一个矛盾。同理，当 $m (\bar{t}) = E_{D} (\bar{t}), I_{D} (\bar{t}), S_{L} (\bar{t}), E_{L} (\bar{t}), I_{L} (\bar{t})$ 或 $x (\bar{t})$ 时，仍然有矛盾。因此，对于所有的 $t \in [0, \bar{t}]$ 都有 $m (t) > 0$ 。

设 $N_{D} = S_{D} + E_{D} + I_{D}$ ， $N_{L} = S_{L} + E_{L} + I_{L}$ ，则由模型(2.1)可知

$d N_{D} / d t = A_{1} - d_{1} N_{D} - α_{1} I_{D} \leq A_{1} - d_{1} N_{D},$

这意味着

$\underset{t \to \infty}{\lim \sup} N_{D} \leq A_{1} / d_{1}$ 。

同理，我们也有

$\underset{t \to \infty}{\lim \sup} N_{L} \leq A_{2} / d_{2}, \underset{t \to \infty}{\lim \sup} x \leq a A_{1} / d d_{1}$ 。

因此，对于所有的 $t \in [0, \bar{t}]$ ， $m (t)$ 都是有界的。此外，我们有 $\bar{t} = + \infty$ 。这就完成了证明。

3.2. 平衡点与稳定性分析

显然可见模型(2.1)具有无病平衡点 $P_{0} (A_{1} / d_{1}, 0, 0, A_{2} / d_{2}, 0, 0, 0)$ 。使用文献[13]中的下一代矩阵方法，我们能够得出模型(2.1)的基本再生数

$R_{0} = ρ (F V^{- 1}) = \sqrt{\frac{a β_{1} β_{2} γ_{1} γ_{2}}{d (d_{1} + γ_{1}) (d_{1} + α_{1} + σ) (d_{2} + γ_{2}) (d_{2} + α_{2})}}$

其中

$F = [\begin{matrix} 0 & 0 & 0 & β_{1} & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & β_{2} \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \end{matrix}], V = [\begin{matrix} d_{1} + γ_{1} & 0 & 0 & β_{1} & 0 \\ γ_{1} & d_{1} + α_{1} + σ & 0 & 0 & 0 \\ 0 & 0 & d_{2} + γ_{2} & 0 & 0 \\ 0 & 0 & - γ_{2} & d_{2} + α_{2} & 0 \\ 0 & - a & 0 & 0 & d \end{matrix}] .$

根据文献[13]，我们得到以下结论。

定理3.3 如果 $R_{0} < 1$ ，那么模型(2.1)的无病平衡点 $P_{0}$ 是局部渐近稳定的；如果 $R_{0} > 1$ ，那么无病平衡点是不稳定的。

此外，我们还得到了无病平衡点的全局渐近稳定性的结果。证明见附录B。

由于模型(2.1)的复杂性，其地方病平衡点的存在性是很难分析的。因此，假设 $α_{1}$ 为零。从模型(2.1)的前三个方程中可以得到 $d N_{D} / d t = A_{1} - d_{1} N_{D}$ 。不失一般性，令 $N_{D} = A_{1} / d_{1}$ 。在模型(2.1)中用 $A_{1} / d_{1} - E_{D} - I_{D}$ 来表示 $S_{D}$ 。那么，模型(2.1)将变成模型(3.2)，进而有下述结论，证明过程见附录A。

定理3.4 当 $R_{0} \leq 1$ 且 $α_{1} = 0$ 时，模型(2.1)的无病平衡点是全局渐近稳定的。

${\begin{cases} \frac{d E_{D}}{d t} = \frac{β_{1} d_{1}}{A_{1}} (\frac{A_{1}}{d_{1}} - E_{D} - I_{D}) I_{L} - (d_{1} + γ_{1}) E_{D}, \\ \frac{d I_{D}}{d t} = γ_{1} E_{D} - (d_{1} + σ) I_{D}, \\ \frac{d S_{L}}{d t} = A_{2} - \frac{β_{2} S_{L}}{N_{L}} x - d_{2} S_{L}, \\ \frac{d E_{L}}{d t} = \frac{β_{2} S_{L}}{N_{L}} x - (d_{2} + γ_{2}) E_{L}, \\ \frac{d I_{L}}{d t} = γ_{2} E_{L} - (d_{2} + α_{2}) I_{L}, \\ \frac{d x}{d t} = a I_{D} - d x . \end{cases}$ (3.1)

设 $P^{*} = (E_{D}^{*}, I_{D}^{*}, S_{L}^{*}, E_{L}^{*}, I_{L}^{*}, x^{*})$ 是模型(3.2)的一个地方病平衡点，则可以得到以下表达式：

$\begin{array}{l} E_{D}^{*} = \frac{A_{1} β_{1} (d_{1} + σ) I_{L}^{*}}{A_{1} q_{1} (d_{1} + σ) + β_{1} d_{1} (d_{1} + γ_{1} + σ) I_{L}^{*}}, I_{D}^{*} = \frac{γ_{1}}{d_{1} + σ} E_{D}^{*}, x^{*} = \frac{a}{d} I_{D}^{*}, \\ S_{L}^{*} = \frac{A_{2}}{λ_{L}^{*} + d_{2}}, E_{L}^{*} = \frac{A_{2} λ_{L}^{*}}{(λ_{L}^{*} + d_{2}) q_{3}}, I_{L}^{*} = \frac{A_{2} γ_{2} λ_{L}^{*}}{(λ_{L}^{*} + d_{2}) q_{3} q_{4}}, \end{array}$

其中 $λ_{L}^{*} = \frac{β_{2} x^{*}}{N_{L}^{*}}$ ，满足方程：

$a_{1} λ_{L}^{*^{2}} + b_{1} λ_{L}^{*} + c_{1} = 0$ (3.2)

其中

$\begin{array}{l} a_{1} = (q_{4} + γ_{2}) (β_{1} d_{1} A_{2} γ_{2} (q_{2} + γ_{1}) + A_{1} q_{1} q_{2} q_{3} q_{4}), c_{1} = A_{1} d_{2} q_{1} q_{2} q_{3}^{2} q_{4}^{2} (1 - R_{0}^{2}), \\ b_{1} = q_{3} q_{4} (A_{1} d_{2} q_{1} q_{2} (q_{4} + γ_{2}) + A_{2} β_{1} γ_{2} d_{1} (q_{2} + γ_{1}) + A_{1} q_{1} q_{2} q_{3} q_{4} (1 - R_{0}^{2})) . \end{array}$

由等式(3.3)可得如果 $R_{0} \geq 1$ ，则 $c_{1} \geq 0$ 和 $b_{1} \geq 1$ 。因此，当 $R_{0} \leq 1$ 时，模型(3.2)没有地方病平衡点，当 $R_{0} \leq 1$ 时，模型(3.2)具有唯一的地方病平衡点 $P^{*}$ 。

定理3.5 当 $R_{0} > 1$ 且 $α_{1} = 0$ 时，模型(2.1)具有唯一的地方病平衡点 $P^{*}$ 。

此外，应用动力学系统的持续性理论[14]，我们得到以下结果。证明见附录C。

定理3.6 当 $R_{0} > 1$ 且 $α_{1} = 0$ 时，模型(2.1)是一致连续的，即存在一个正常数 $p$ ，使得(2.1)的所有具有正初始条件的解都满足

$\underset{t \to \infty}{\lim \inf} (S_{D} (t), E_{D} (t), I_{D} (t), S_{L} (t), E_{L} (t), I_{L} (t), x (t)) \geq (p, p, p, p, p, p, p) .$

4. 最优控制问题

4.1. 最优控制的存在性

在本小节当中，考虑 $u_{1}$ 屠宰场卫生监测， $u_{2}$ 家犬的驱虫控制， $u_{3}$ 家畜疫苗接种， $u_{4}$ 环境消毒这四个干预措施，我们将包虫病模型(2.1)扩展为如下形式控制模型，

${\begin{cases} \frac{d S_{D}}{d t} = A_{1} - \frac{β_{1} (1 - u_{1} (t)) S_{D}}{N_{D}} I_{L} - d_{1} S_{D} + σ I_{D} + c δ_{1} u_{2} (t) I_{D}, \\ \frac{d E_{D}}{d t} = \frac{β_{1} (1 - u_{1} (t)) S_{D}}{N_{D}} I_{L} - (d_{1} + γ_{1}) E_{D}, \\ \frac{d I_{D}}{d t} = γ_{1} E_{D} - (d_{1} + α_{1} + σ) I_{D} - δ_{1} u_{2} (t) I_{D}, \\ \frac{d S_{L}}{d t} = (1 - δ_{2} u_{3} (t)) A_{2} - \frac{β_{2} S_{L}}{N_{L}} x - d_{2} S_{L}, \\ \frac{d E_{L}}{d t} = \frac{β_{2} S_{L}}{N_{L}} x - (d_{2} + γ_{2}) E_{L}, \\ \frac{d I_{L}}{d t} = γ_{2} E_{L} - (d_{2} + α_{2}) I_{L}, \\ \frac{d x}{d t} = a I_{D} - d x - δ_{3} u_{4} (t) x . \end{cases}$ (4.1)

定义控制集合

$U = {(u_{1} (t), u_{2} (t), u_{3} (t), u_{4} (t)) : 0 \leq u_{i} (t) \leq 1, 勒 � 格可 �, t \in [0, T], i = 1, 2, 3, 4}$ 。记 $u = (u_{1} (t), u_{2} (t), u_{3} (t), u_{4} (t))$ ，设 $φ = (S_{D}, E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x)$ 是模型(4.1)对于控制 $u$ 的状态解。我们的目的是使家畜中的潜伏者 $E_{L}$ 和家畜中的感染者 $I_{L}$ 的数量尽可能的少，同时花费尽可能少。为此，我们考虑目标函数

$\begin{matrix} J (u) : = \int_{0}^{T} (ξ_{1} E_{D} + ξ_{2} I_{D} + ξ_{3} E_{L} + ξ_{4} I_{L} + ξ_{5} x + B_{1} u_{1} + C_{1} u_{1}^{2} \\ + B_{2} I_{D} u_{2} + C_{2} u_{2}^{2} + B_{3} S_{L} u_{3} + C_{3} u_{3}^{2} + B_{4} x u_{4} + C_{4} u_{4}^{2}) d t . \end{matrix}$

其中 $ξ_{i} > 0 (i = 1, 2, 3, 4, 5), B_{i}, C_{i} > 0 (i = 1, 2, 3, 4)$ 是平衡系数，将积分转化为在 $T$ 天内的花费。其中控制的二次表达式，表明在控制策略实施过程中可能产生的非线性成本。

定理4.1 给定 $J (u)$ 满足具有正初值的系统(4.1)，则存在一个最优控制 $u^{*}$ 和相应的 $φ$ ，使 $J (u)$ 最小。

为了证明最优控制的存在性，需要检验假设文献[15]中的假设(具体证明过程见附录C)

1) 控制集 $U$ 是凸闭的。

2) 存在常数 $C_{1}$ 和 $C_{2}$ 使得

$f (t, φ, u) \leq C_{1} (1 + | φ | + | u |), | f (t, φ_{1}, u) - f (t, φ, u) | \leq C_{2} (1 + | u |) | φ_{1} - φ |$ (4.2)

3) 目标泛函的被积函数 $g (t, φ, u)$ 是凸的。

4) 存在常数 $D_{1}, D_{2} > 0$ 且 $β^{*} > 1$ 使得目标泛函的被积函数满足

$g (t, φ, u) \geq D_{1} {({| u_{1} |}^{2} + {| u_{2} |}^{2} + {| u_{3} |}^{2} + {| u_{4} |}^{2} + {| u_{5} |}^{2})}^{\frac{β^{*}}{2}} - D_{2}$

4.2. 最优控制的刻画

为了得到最优控制的必要条件，即最优控制的表达式，应用庞特里亚金最大值原理。因此，需要将最优控制问题转变为关于 $u$ 的逐点哈密顿量最小化问题，其中哈密顿量定义为

$H (t, φ, u, λ) = g (t, φ, u) + λ_{1} \frac{d S_{D}}{d t} + λ_{2} \frac{d E_{D}}{d t} + λ_{3} \frac{d I_{D}}{d t} + λ_{4} \frac{d S_{L}}{d t} + λ_{5} \frac{d E_{L}}{d t} + λ_{6} \frac{d I_{L}}{d t} + λ_{7} \frac{d x}{d t},$

这里 $λ (t) = (λ_{1} (t), \dots, λ_{7} (t))$ 是伴随变量，满足方程组

$\begin{array}{l} \frac{d λ_{1}}{d t} = - \frac{\partial H}{\partial S_{D}} = d λ_{1} + \frac{β_{1} (1 - u_{1}) (E_{D} + I_{D}) I_{L}}{N_{D}^{2}} (λ_{1} - λ_{2}) \\ \frac{d λ_{2}}{d t} = - \frac{\partial H}{\partial E_{D}} = - ξ_{1} + (d_{1} + γ_{1}) λ_{2} - γ_{1} λ_{3} + \frac{β_{1} (1 - u_{1}) S_{D} I_{L}}{N_{D}^{2}} (λ_{2} - λ_{1}) \\ \frac{d λ_{3}}{d t} = - \frac{\partial H}{\partial I_{D}} = - ξ_{2} - B_{2} u_{2} - (σ + c δ_{1} u_{2}) λ_{1} + \frac{β_{1} (1 - u_{1}) S_{D} I_{L}}{N_{D}^{2}} (λ_{2} - λ_{1}) + χ λ_{3} - a λ_{7} \\ \frac{d λ_{4}}{d t} = - \frac{\partial H}{\partial S_{L}} = - B_{3} u_{3} + \frac{β_{2} x (E_{L} + I_{L})}{N_{L}^{2}} (λ_{4} - λ_{5}) + d_{2} λ_{4} \\ \frac{d λ_{5}}{d t} = - \frac{\partial H}{\partial E_{L}} = - ξ_{3} + \frac{β_{2} S_{L} x}{N_{L}^{2}} (λ_{5} - λ_{4}) + (d_{2} + γ_{2}) λ_{5} - γ_{2} λ_{6} \\ \frac{d λ_{6}}{d t} = - \frac{\partial H}{\partial I_{L}} = - ξ_{4} + \frac{β_{1} (1 - u_{1}) S_{D}}{N_{D}} (λ_{1} - λ_{2}) + \frac{β_{2} S_{L} x}{N_{L}^{2}} (λ_{5} - λ_{4}) + (d_{2} + α_{2}) λ_{6} \\ \frac{d λ_{7}}{d t} = - \frac{\partial H}{\partial x} = - ξ_{5} - B_{4} u_{4} + \frac{β_{2} S_{L}}{N_{L}} (λ_{5} - λ_{4}) + (d + δ_{4} u_{5}) λ_{7} \\ λ (T) = 0 \end{array}$

其中 $χ = d_{1} + α_{1} + σ + δ_{1} u_{2}$ 。此外，利用庞特里亚金原理的必要条件，我们得到 $u^{*}$ 是 $H (t, φ, u, λ)$ 的一个临界点，即 $u^{*}$ 满足 $\partial H / \partial u_{i} = 0, i = 1, 2, 3, 4$ 。因此，我们有

$\begin{array}{l} B_{1} + 2 C_{1} u_{1} + \frac{β_{1} S_{D} I_{L}}{N_{D}} λ_{1} - \frac{β_{1} S_{D} I_{L}}{N_{D}} λ_{2} = 0, \\ B_{2} I_{D} + 2 C_{2} u_{2} + λ_{1} δ_{1} I_{D} - δ_{1} I_{D} λ_{3} = 0, \\ B_{3} S_{L} + 2 C_{3} u_{3} - λ_{4} δ_{2} A_{2} = 0, \\ B_{4} x + 2 C_{4} u_{4} - λ_{7} δ_{3} x = 0. \end{array}$

计算可得最优控制的表达式：

$\begin{array}{l} u_{1}^{*} = \min {1, \max {0, \frac{- B_{1} N_{D} + β_{1} S_{D} I_{L} (λ_{2} - λ_{1})}{2 C_{1} N_{D}}}}, u_{2}^{*} = \min {1, \max {0, \frac{- B_{2} I_{D} + δ_{1} I_{D} (λ_{3} - λ_{1})}{2 C_{2}}}}, \\ u_{3}^{*} = \min {1, \max {0, \frac{- B_{3} S_{L} + δ_{2} A_{2} λ_{4}}{2 C_{3}}}}, u_{4}^{*} = \min {1, \max {0, \frac{- B_{4} x + δ_{3} x λ_{7}}{2 C_{4}}}} . \end{array}$

(a) (b)

Figure 2. Solutions of the model with $E_{D} (0)$ values 500, 1000, 1500, represented by dashed lines, solid lines, and dotted lines

图2. 模型解曲线图，虚线、实线、点线的 $E_{D} (0)$ 初值分别为500，1000，1500

5. 数值模拟

取参数 $A_{1} = 2.1 \times 10^{5}, β_{1} = 2.16 \times 10^{- 2}, d_{1} = 0.08, σ = 2, γ_{1} = 0.2, A_{2} = 5.4 \times 10^{5}, β_{2} = 7.4, d_{2} = 0.9, α_{1} = 0.001,$ $α_{2} = 0.001, a = 970$ 计算可得基本再生数 $R_{0} = 1.4428$ 。取初始值 $(3 \times 10^{5}, 500, 1000, 5.386 \times 10^{5}, 2.5 \times 10^{5}, 2.5 \times 10^{4}, 2.5 \times 10^{5})$ 值，再改变 $E_{D} (0)$ 为取值为1000，1500，使用MATLAB可以看出解最终都趋于正平衡点，见图2。

6. 结论

本文一方面建立了一个具有标准发生率的包虫病传播模型，分析了模型的灭绝性与持续性。另一方面考虑一些控制策略，构造控制系统，应用最优控制理论证明了最优控制的存在性并给出了刻画。

基金项目

安徽省高等学校自然科学研究重点项目(2022AH051320, 2023AH050415)，阜阳师范大学博士启动基金项目(2021KYQD0001, 2021KYQD0002)，国家级大学生创新创业项目“基于流行病学模型的包虫病传播与最优控制研究”(202410371033)。

附录A 定理3.4的证明

对于模型(3.2)，我们考虑李雅普诺夫函数

$V (E_{D}, I_{D}, E_{L}, I_{L}, x) = ξ_{1} E_{D} + ξ_{2} I_{D} + ξ_{3} E_{L} + I_{L} + ξ_{4} x$

其中 $ξ_{1} = \frac{α β_{2} γ_{1} γ_{2}}{d q_{1} q_{2} q_{3}}, ξ_{2} = \frac{α β_{2} γ_{2}}{d q_{2} q_{3}}, ξ_{3} = \frac{γ_{2}}{q_{3}}, ξ_{4} = \frac{β_{2} γ_{2}}{d q_{3}}$ 。沿着(3.2)的解计算

$V (E_{D}, I_{D}, E_{L}, I_{L}, x)$ 的全导数，我们得到

$\begin{matrix} \frac{d V}{d t} = ξ_{1} (\frac{β_{1} d_{1}}{A_{1}} (\frac{A_{1}}{d_{1}} - E_{D} - I_{D}) I_{L} - (d_{1} + γ_{1}) E_{D}) + ξ_{2} (γ_{1} E_{D} - (d_{1} + σ) I_{D}) \\ + ξ_{3} (\frac{β_{2} S_{L}}{N_{L}} x - (d_{2} + γ_{2}) E_{L}) + (γ_{2} E_{L} - (d_{2} + α_{2}) I_{L}) + ξ_{4} (a I_{D} - d x) \\ = - \frac{β_{2} γ_{2}}{q_{3}} (1 - \frac{S_{L}}{N_{L}}) x + (ξ_{1} \frac{β_{1} d_{1}}{A_{1}} (\frac{A_{1}}{d_{1}} - E_{D} - I_{D}) I_{L} - (d_{2} + α_{2})) I_{L} \\ \leq - \frac{β_{2} γ_{2}}{q_{3}} (1 - \frac{S_{L}}{N_{L}}) x + (ξ_{1} β_{2} - q_{4}) I_{L} = - \frac{β_{2} γ_{2}}{q_{3}} (1 - \frac{S_{L}}{N_{L}}) x + q_{4} (R_{0}^{2} - 1) I_{L} . \end{matrix}$

因此，当 $R_{0} \leq 1$ 且 $α_{1} = 0$ 时，对于所有的 $(E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x) \in ℝ_{+}^{6}$ 都有 $d V / d t \leq 0$ 。当 $d V / d t = 0$ 时， $x = 0$ 或 $S_{L} = N_{L}$ 。如果 $x = 0$ ，则由模型(3.2)，我们能够得到 $E_{D} = I_{D} = E_{L} = I_{L} = 0$ 和 $S_{L} = A_{2} / d_{2}$ ，如果 $S_{L} = N_{L}$ ，则可以得到 $E_{L} = I_{L} = 0$ 和 $S_{L} = A_{2} / d_{2}$ 。因此，由模型(2.1)可以得到 $E_{D} = I_{D} = x = 0$ 。因此，在 ${(E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x) \in ℝ_{+}^{6} : d V / d t = 0}$ 中的最大正不变集合是单点集 $(0, 0, A_{2} / d_{2}, 0, 0, 0)$ 。由拉萨尔不变性原理可知，无病平衡点是全局逐渐稳定的。

附录B 定理3.6的证明

定义集合：

$\begin{array}{l} {(S_{D}, E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x) : S_{D} > 0, E_{D} \geq 0, I_{D} \geq 0, S_{L} > 0, E_{L} \geq 0, I_{L} \geq 0, x \geq 0}, \\ X_{0} = {(S_{D}, E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x) \in X : S_{D}, E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x > 0}, \\ \partial X_{0} = {(S_{D}, E_{D}, I_{D}, S_{L}, E_{L}, I_{L}, x) \in X : E_{D} I_{D} E_{L} I_{L} x = 0} . \end{array}$

设 $(S_{D} (t), E_{D} (t), I_{D} (t), S_{L} (t), E_{L} (t), I_{L} (t), x (t)) = φ (t)$ 是具有初始条件为

$(S_{D} (0), E_{D} (0), I_{D} (0), S_{L} (0), E_{L} (0), I_{L} (0), x (0)) = φ (0)$ 的模型(2.1)的解。定义如下的一个集合 $M_{\partial} = {φ (0) \in \partial X_{0} : φ (t) \in \partial X_{0}, t \geq 0}$ 。显然，具有初始条件为 $(S_{D} (0), 0, 0, S_{L} (0), 0, 0, 0)$ 的模型(2.1)的解形如 $(S_{D} (t), 0, 0, S_{L} (t), 0, 0, 0)$ 。因此，我们有

$A : = {(S_{D} (0), 0, 0, S_{L} (0), 0, 0, 0) : S_{D} (0), S_{L} (0) > 0} \subset M_{\partial}$ 。

设存在 $φ (0) \in M_{\partial}$ 使得 $φ (0) \notin A$ ，则在 $E_{D} (0), I_{D} (0), E_{L} (0), I_{L} (0)$ 和 $x (0)$ 中至少有一个数大于零。不失一般性，假设 $E_{D} (0) > 0$ ，则对于所有 $t > 0$ 时，我们有 $E_{D} > E_{D} (0) e^{- (d_{1} + γ_{1}) t} > 0$ 。由模型(2.1)的第二个方程，对于所有的 $t > 0$ 可得 $I_{D} > I_{D} (0) e^{- (d_{1} + α_{1} + σ) t} > 0$ 。由模型(2.1)的第六个方程，对于所有的 $t > 0$ 有 $x > x (0) e^{- d t} > 0$ 。由模型(2.1)的第四个和第五个方程，对于所有的 $t > 0$ ，同理得到 $E_{L} \geq E_{L} (0) e^{- (d_{2} + γ_{2}) t} > 0, I_{L} > I_{L} (0) e^{- d_{2} t} \geq 0$ 。因此，我们得到 $φ (0) \notin M_{\partial}$ ，这就导致了一个矛盾。所以我们也得到 $M_{\partial} \subset A$ 。因此，我们最终能够得到 $M_{\partial} = A$ 。

很显然，对于 $M_{\partial}$ 的初值即模型(2.1)的无病平衡 $P_{0}$ 是具有全局吸引的。这表明 $M_{\partial}$ 中的 ${P_{0}}$ 是孤立不变的且不具有周期性的。由定理3.3可知，如果 $R_{0} > 1$ ，则 $P_{0}$ 是不稳定的。因此，当 $t \to \infty$ 时，对于任意的 $φ (0) \in X_{0}$ ， $φ (t)$ 的解不趋向 $P_{0}$ 。因此，我们有 $W^{s} (P_{0}) = {φ (0) : \lim_{t \to \infty} φ (t) = P_{0}}$ 。

由文献[14]中的定理4.2可知存在一个数 $p > 0$ ，使得

$\underset{t \to \infty}{\lim \inf} d (φ (t), \partial X_{0}) \geq p, \forall φ (0) \in X_{0} .$

这意味着对于初值为正值的解，每个分量都满足

$\begin{array}{l} \underset{t \to \infty}{\lim \inf} S_{D} (t) \geq p, \underset{t \to \infty}{\lim \inf} E_{D} (t) \geq p, \underset{t \to \infty}{\lim \inf} I_{D} (t) \geq p, \\ \underset{t \to \infty}{\lim \inf} S_{L} (t) \geq p, \underset{t \to \infty}{\lim \inf} E_{L} (t) \geq p, \underset{t \to \infty}{\lim \inf} I_{L} (t) \geq p, \underset{t \to \infty}{\lim \inf} x (t) \geq p . \end{array}$

附录C 定理4.1的证明

1) 对于任意 $y = (y_{1}, y_{2}, y_{3}, y_{4}), z = (z_{1}, z_{2}, z_{3}, z_{4}) \in U$ 。因为 $0 \leq v y_{i} + (1 - v) z_{i} \leq 1$ ，对于任意的 $v \in [0, 1]$ 可以得到 $v y + (1 - v) z \in U, \forall v \in [0, 1]$ 。因此， $U$ 是凸的。此外，对于当 $n \to \infty$ 时满足 $u^{n} \to u$ 的一个数列 ${u^{n}} \subseteq U$ ，若能够得到 $u \in U$ ，我们称 $U$ 是有界的。事实上，由于当 $n \to \infty$ 时 $u_{i}^{n} \to u_{i}$ ，即当 $n \to \infty$ 时， $u_{i}^{n} \to u_{i}, i = 1, 2, 3, 4$ ，从而对于任意 $ε > 0$ ，存在 $N > 0$ ，使得对于任意的 $n > N, t \in [0, T]$ 时，有 $| u_{i}^{n} - u_{i} | < ε$ ，即 $ε \leq ε - u_{i}^{n} \leq u_{i} \leq u_{i}^{n} + ε \leq 1 + ε$ 。令 $ε \to 0$ ，则对于 $i = 1, 2, 3, 4, 5$ 时，有 $0 \leq u_{i} \leq 1$ ，故 $u \in U$ 。

2) 显然系统(4.1)可以改写为 $d φ / d t = F_{1} (t, φ) + F_{2} (t, φ) u = : f (t, φ, u)$ ，其中

$F_{1} (t, φ) = [\begin{matrix} A_{1} - \frac{β_{1} S_{D}}{N_{D}} I_{L} - d_{1} S_{D} + σ I_{D} \\ \frac{β_{1} S_{D}}{N_{D}} I_{L} - (d_{1} + γ_{1}) E_{D} \\ γ_{1} E_{D} - (d_{1} + α_{1} + σ) I_{D} \\ A_{2} - \frac{β_{2} S_{L}}{N_{L}} x - d_{2} S_{L} \\ \frac{β_{2} S_{L}}{N_{L}} x - (d_{2} + γ_{2}) E_{L} \\ γ_{2} E_{L} - (d_{2} + α_{2}) I_{L} \\ σ I_{D} - d x \end{matrix}], F_{2} (t, φ) = [\begin{matrix} \frac{β_{1} S_{D}}{N_{D}} I_{L} & c δ_{1} I_{D} & 0 & 0 \\ - \frac{β_{1} S_{D}}{N_{D}} I_{L} & 0 & 0 & 0 \\ 0 & - c δ_{1} I_{D} & 0 & 0 \\ 0 & 0 & - A_{2} & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & v \end{matrix}] .$

因为

$‖ u ‖ = \sqrt{\sum_{i = 1}^{4} {| u_{i} |}^{2}} \leq \sqrt{5},$

即 $u$ 是有界的，所以为了证明不等式(4.2)，只需证明

$‖ f (t, φ, u) ‖ \leq C_{1} (1 + ‖ φ ‖), ‖ f (t, φ_{1}, u) - f (t, φ, u) ‖ \leq C_{2} ‖ φ_{1} - φ ‖ .$

为了证明上述不等式，只需证明存在一个常数C，使得

$‖ f (t, 0, 0) ‖ \leq C, ‖ f_{φ} (t, φ, u) ‖ \leq C, ‖ f_{u} (t, φ, u) ‖ \leq C .$

在 $f (t, φ, u) - f (t, φ, 0), f (t, φ, 0) - f (t, 0, 0), f (t, φ_{1}, u_{1}) - f (t, φ, u)$ 中应用中值定理，则(4.2)立即成立。

因为 $f (t, 0, 0) = {[A_{1}, 0, 0, A_{2}, 0, 0, 0]}^{T}$ ，所以我们有 $‖ f (t, 0, 0) ‖ \leq \sqrt{A_{1}^{2} + A_{2}^{2}} : = C^{'}$ 此外，我们有 $f_{u} (t, φ, u) = F_{2} (t, φ)$ ，因此，

$‖ f_{u} (t, φ, u) ‖ \leq \sqrt{2 β_{1}^{2} I_{L}^{2} + 2 c^{2} δ_{1}^{2} I_{D}^{2} + A_{2}^{2} + v^{2}} \leq \sqrt{(\frac{2 β_{2}^{2}}{d_{2}} + 1) A_{2}^{2} + \frac{2 c^{2} δ_{1}^{2}}{d_{1}^{2}} + A_{1}^{2} + δ_{3}^{2}} : = C^{″}$

此外，我们还有

$f_{φ} (t, φ, u) = [\begin{matrix} D_{1} & 0 \\ 0 & D_{2} \end{matrix}]$

其中

$D_{1} = [\begin{matrix} - ϑ_{4} I_{L} - d_{1} & ϑ_{1} S_{D} I_{L} & ϑ_{1} S_{D} I_{L} + σ + δ_{1} u_{2} \\ ϑ_{4} I_{L} & - ϑ_{1} S_{D} I_{L} - q_{1} & ϑ_{1} S_{D} I_{L} \\ 0 & γ_{1} & ϑ_{0} \end{matrix}]$ ， $D_{2} = [\begin{matrix} - ϑ_{2} x & ϑ_{3} S_{L} x & ϑ_{3} S_{L} x & - \frac{β_{2} S_{L}}{N_{L}} \\ ϑ_{2} x & - ϑ_{3} S_{L} x - q_{3} & - ϑ_{3} S_{L} x & \frac{β_{2} S_{L}}{N_{L}} \\ 0 & γ_{2} & - q_{4} & 0 \\ 0 & 0 & 0 & ϑ_{5} \end{matrix}]$

并且

$\begin{array}{l} ϑ_{0} = d_{1} + α_{1} + σ + c δ_{1} u_{2}, ϑ_{1} = \frac{β_{1} (1 - u_{1})}{N_{D}^{2}}, ϑ_{4} = \frac{β_{1} (E_{D} + I_{D}) (1 - u_{1})}{N_{D}^{2}}, \\ ϑ_{2} = \frac{β_{2} (E_{L} + I_{L})}{N_{D}^{2}}, ϑ_{3} = \frac{β_{2}}{N_{L}^{2}}, ϑ_{5} = - d - δ_{3} u_{4} . \end{array}$

注意矩阵 $f_{φ} (t, φ, u)$ 的每一个元素都是有界的，即 $‖ f_{φ} (t, φ, u) ‖ \leq C$ 。因此，只要取 $C = \max {C^{'}, C^{″}, C^{‴}}$ ，就能使得(4.2)成立。

设

$h_{1} (u_{1}) = B_{1} u_{1} + C_{1} u_{1}^{2}, h_{2} (u_{2}) = B_{2} I_{D} u_{2} + C_{2} u_{2}^{2},$

$h_{3} (u_{3}) = B_{3} S_{L} u_{3} + C_{3} u_{3}^{2}, h_{4} (u_{4}) = B_{4} x u_{4} + C_{4} u_{4}^{2} .$

要证明 $g (t, φ, u)$ 关于 $u$ 是凸的，我们需要证明 $h_{i} (i = 1, 2, 3, 4)$ 是凸的。不失一般性，我们只给出 $h_{1} (u_{1})$ 关于 $u_{1}$ 是凸性。设 $0 \leq y, z \leq 1$ 和 $v \in [0, 1]$ ，则

$h_{1} (v y + (1 - v) z) \leq v h_{1} (y) + (1 - v) h_{1} (z) = - C_{1} v (1 - v) {(y - z)}^{2} \leq 0$ 。

因此， $h_{1} (v y + (1 - v) z) \leq v h_{1} (y) + (1 - v) h_{1} (z)$ ，也就说明了 $h_{1} (u_{1})$ 是凸的。

4) 对于被积函数 $g (t, φ, u)$ ，我们有

$\begin{matrix} g (t, φ, u) = a_{1} E_{L} + a_{2} E_{L} + B_{1} u_{1} + B_{2} I_{D} u_{2} + B_{3} u_{3} + B_{4} x u_{4} + \sum_{i = 1}^{4} c_{i} u_{i}^{2} \\ \geq \sum_{i = 1}^{4} c_{i} u_{i}^{2} \\ \geq D_{1} {(\sum_{i = 1}^{4} u_{i}^{2})}^{β^{*} / 2} - D_{2} . \end{matrix}$

这样就完成了证明。

NOTES

^*通讯作者。

参考文献

[1]	Azlaf, R. and Dakkak, A. (2006) Epidemiological Study of the Cystic Echinococcosis in Morocco. Veterinary Parasitology, 137, 83-93.
[2]	Gemmell, M.A., Lawson, J.R. and Roberts, M.G. (1986) Population Dynamics in Echinococcosis and Cysticercosis: Biological Parameters of Echinococcosis granulosus in Dogs and Sheep. Parasitology, 92, 599-620.
[3]	Huang, L., Huang, Y., Wang, Q., et al. (2011) An Agent-Based Model for Control Strategies of Echinococcosis granulosus. Veterinary Parasitology, 179, 84-91.
[4]	Torgerson, P.R., Ziadinov, I., Aknazarov, D., Nurgaziev, R. and Deplazes, P. (2009) Modelling the Age Variation of Larval Protoscoleces of Echinococcosis granulosus in Sheep. International Journal for Parasitology, 39, 1031-1035.
[5]	Wang, K., Zhang, X. and Jin, Z. (2013) Modeling and Analysis of Transmission of Echinococcosis with Application to Xingjiang Uygur Autonomous Region of China. Journal of Theoretical Biology, 333, 78-90.
[6]	Rong, X., Fan, M., Sun, X., et al. (2018) Impact of Disposing Stray Dogs on Risk Assessment and Control of Echinococcosis in Inner Mongolia. Mathematical Biosciences, 299, 85-96.
[7]	Rong, X., Fan, M., Zhu, H. and Zheng, Y. (2021) Dynamic Modeling and Optimal Control of Cystic Echinococcosis. Infectious Diseases of Poverty, 10, Article No. 38.
[8]	Roberts, M.G., Lawson, J.R. and Gemmell, M.A. (1986) Population Dynamics in Echinococcosis and Cysticercosis: Mathematical Model of the Life-Cycle of Echinococcus granulosus. Parasitology, 92, 621-641.
[9]	赵瑜. 一类包虫病传播动力学模型的研究[J]. 生物数学学报, 2011, 26(3): 441-450.
[10]	Zhu, G., Chen, S., Shi, B., et al. (2019) Dynamics of Echinococcosis Transmission among Multiple Species and Acase Study in Xinjiang, China. Chaos, Solitons & Fractals, 127, 103-109.
[11]	崔倩倩. 宁夏彭阳县包虫病传播动力学建模及分析[J]. 应用数学, 2023, 36(2): 277-289.
[12]	Zhao, J. and Run, Y. (2021) A Dynamical Model of Echinococcosis with Optimal Control and Cost-Effectiveness. Nonlinear Analysis: Real World Applications, 62, Article 103388.
[13]	van den Dressche, P. and Watmough, J. (2002) Reproduction Numbers and Sub-Throshold Endemic Equilibria for Compartmental Models of Disease Transmission. Mathematical Biosciences, 180, 29-48.
[14]	Hale, J.K. and Waltman, P. (1989) Persistence in Infinite-Dimensional Systems. SIAM Journal on Mathematical Analysis, 20, 388-395.
[15]	Lenhart, S.M. and Workman, J.T. (2007) Optimal Control Applied to Biological Models. Chapman et Hall/CRC Books.

为你推荐

友情链接