下载:
-
开放科学（资源服务）标志码（OSID）：
-
蚊子是世界上最致命的动物之一，由其传播的疾病每年导致数百万人死亡[1]. 因此，持续的蚊虫控制工作对于预防这些疾病的爆发非常重要.
Wolbachia是一种革兰氏阴性菌，是世界上分布最为广泛的共生菌，并且可以人工植入埃及伊蚊体内[2]. 研究发现Wolbachia可以抑制登革热病毒在蚊子体内的复制并且可以降低病毒的传染性[3]，故使用Wolbachia控制蚊媒传染病的传播备受关注. Wolbachia对蚊子的繁殖影响主要有：垂直传播和细胞质分离(CI)，即携带Wolbachia的雌性的下一代也会携带此种菌，正常雌性与携带Wolbachia的雄性的下一代会因细胞质分离而死亡.
近年来，许多文献利用数学模型研究了Wolbachia在蚊虫种群中的传播. 文献[4]提出并分析了Wolbachia感染蚊子种群与未感染蚊子种群之间的基本竞争模型，然后，利用反馈控制技术设计了Wolbachia的引入方案. 文献[5-7]分别建立了具有脉冲一般出生和死亡率函数的脉冲模型、具有脉冲出生和投放的性别结构和状态依赖脉冲的综合控制模型，对使用Wolbachia控制蚊媒传染病的各种控制策略进行了研究. 文献[8]考虑了环境的异质性，建立了两种机制随机切换的数学模型. 研究发现：在均匀环境中维持的Wolbachia的初始状态在异质环境中会灭绝，频繁的环境转换有利于Wolbachia的传播. 此外，文献[9-11]等都讨论了宿主、媒介和疾病之间的相互作用.
本文将建立具有成熟时滞的蚊子种群模型，通过理论分析和数值模拟讨论时滞对Wolbachia传播的影响.
全文HTML
-
文献[12]建立了如下数学模型描述Wolbachia在蚊子种群中的传播：
其中：I(t)表示被Wolbachia感染的蚊子数量，U(t)表示未被Wolbachia感染的蚊子数量；β∈[0, 1]是垂直传播的可能性；q∈[0, 1]是CI影响的可能性；b＞0是出生率；d＞0是死亡率；D≥0是适合度损失. 模型(1)考虑了Wolbachia的垂直传播和CI两种影响，由文献[3]可知，对于某种特定的Wolbachia菌株，其垂直传播率和CI影响发生的可能性都近似接近1，也即，如果雌性蚊子感染了Wolbachia，其后代也会感染；正常雌性蚊子与感染雄性蚊子的后代会因为发生CI而死亡. 所以，我们假设β=1，q=1，并忽略掉适合度的影响(D=0).蚊子的一生有4个阶段(卵、幼虫、蛹、成虫)，从卵到成虫大约需要14~21 d，雌性蚊子的寿命约3个月，故蚊子的成熟时间不能被忽略，此因素可用成熟时滞描述. 而成熟时滞又受季节性影响，故考虑为周期时间依赖的时滞，记为τ(t). 根据文献[13]的方法，将τ(t)引入模型(1)中得到如下模型：
其中：δ为幼年蚊子的死亡率，e-δτ(t)是幼年蚊子成熟的可能性. τ(t)是[0，+∞)上的非负、有界且连续可微的ω-周期函数，满足0≤τ(t)≤τ，常数τ=
$\mathop {\sup }\limits_{t \ge 0}$ {τ(t)}，τ≥0；μ(t)和ν(t)是关于t∈[-τ，0]的正连续函数.
-
首先引入如下定义：令χ：=C([-τ，0]，
$\mathbb{R}^{2}$ )，且具有最大模. 如果函数x(·)∈C([-τ，∞]，$\mathbb{R}^{2}$ )，则xt∈χ可被定义为xt(θ)=x(t+θ)，∀θ∈[-τ，0]. 对任意ϕ∈χ，定义f(t，ϕ)=(f1(t，ϕ)，f2(t，ϕ))，其中因为τ(t)是ω-周期的，故有f(t+ω，ϕ)=f(t，ϕ). 所以，模型(2)是一个ω-周期泛函微分系统.
假设g(t)是一个连续的ω-周期函数，令
则对于系统(2)解的适定性有如下结果.
引理1 对任意ϕ=(ϕ1，ϕ2)∈χ+：=C([-τ，0]，
$\mathbb{R}_+^{2}$ )，系统(2)在[0，∞)上有唯一的非负有界解S(t，ϕ)(S0=ϕ).证 在χ+的每个紧子集上，f(t，ϕ)关于ϕ是连续且Lipschitz的. 所以，对任意ϕ∈χ+，系统(2)在其最大存在区间上有唯一的具有初值S0=ϕ的解S(t，ϕ).
定义bc=(1-τ′(t))e-δτ(t)b. 由系统(2)可知：
我们令
则有
易知
又令
并由
可得到
再令
对任意给定的ρ≥1，令
则[0，ρx*]χ是关于χ的有序区间.
容易证明：对任意ψ∈[0，ρx*]χ(ψi(0)=0~(ψi(0)=ρxi*)，i=1，2)和t∈
$\mathbb{R}$ ，有fi(t，ψ)≥0~(fi(t，ψ)≤0). 进一步地，[0，ρx*]χ对于系统(2)是正不变的. 通过选择任意大的ρ即可得到解关于χ+的正性和有界性. 证毕.定理1 系统(2)有一个种群灭绝平衡态(0，0)，且是不稳定的.
证 (0，0)的存在性容易得到，此略去. 定义Y(t)=(I(t)，U(t))为系统(2)的任意解，且|Y(t)|=
$\sqrt{I_{(t)}^{2}+U_{(t)}^{2}}$ . 对于t∈[t0-τ，t0]，Y(t)=Φ(t). 定义Φ(t)的范数为||Φ||=sup{|Φ(t)|：t0-τ≤t≤t0}.假设(0，0)是稳定的. 则对任意ε＞0，存在δ=δ(t0，ε)＞0使得|Y(t)|＜ε对于任意的t≥t0和Φ(t)(||Φ||≤δ)成立. 令Y(t)是系统(2)满足I0=inf{I(t)：t0-τ≤t≤t0}＞0的解，则I0≤δ. 下一步，我们将证明I(t)＞
$\frac{I_{0}}{2}$ 对所有的t≥t0成立.假设对任意t≥t0，I(t)＞
$\frac{I_{0}}{2}$ 不成立. 则存在一个t1＞t0使得对于t∈[t0-τ，t1)有I(t)＞$\frac{I_{0}}{2}$ 和I(t1)=$\frac{I_{0}}{2}$ . 所以，I′(t1)≤0. 由系统(2)的第一个方程知：由于|Y(t)|＜ε，则有I(t1)＜ε和U(t1)＜ε. 所以，
我们取
则有
故对任意的t≥t0，有I(t)＞
$\frac{I_{0}}{2}$ 成立，矛盾. 所以(0，0)是不稳定的. 证毕.
-
因为τ(t)∈[0，τ]，所以，我们首先分析常数时滞的情况. 令τ(t)=r，系统(2)存在如下两个非平凡平衡态(I*，U*)：
下分析E1和E2稳定性，将系统(2)在(I*，U*)处线性化后可写为
其中x(t)=(I(t)，U(t))T且
这里
且有
假设系统(2)具有形如x(t)=c eλt(c≠0)的解，则可得到系统(3)的特征方程如下：
展开得
其中：a1，a2分别为矩阵A的第一列和第二列；b1，b2分别为矩阵B的第一列和第二列.
可以计算得到
则在E1=(0，
$\frac{\mathrm{e}^{-\delta r} b}{d}$ )处，特征方程为由特征根的符号易知E1是不稳定的.
在E2=(
$\frac{\mathrm{e}^{-\delta r} b}{d}$ ，0)处，特征方程为易知E2是渐近稳定的.
综上可得如下结果：
定理2 在区域{(I，U)：I，U≥0}中，当τ(t)=r时，系统(2)有两个非平凡平衡态E1和E2，且E1是不稳定的，E2是局部渐近稳定的.
由定理2可知：如果垂直传播是完全的，Wolbachia最终将完全入侵蚊子种群. 图 1a将系统(2)在有时滞((μ(t)，ν(t))=(20，100))和无时滞((I0，U0)=(20，100))两种情况下的动力学行为做了比较，比较结果见表 1. 表 1说明：无时滞的ODE系统将高估蚊子总量，这对蚊子种群的控制是不利的，因此成熟时滞对于蚊子种群的动力学行为有重要影响.
图 1b将不同时滞对系统(2)动力学行为的影响做了对比. 可以看出：时滞越小，系统(2)解趋于稳定的时间越短，也即，如果环境有利于蚊子的成熟，Wolbachia完全入侵蚊子种群的速度会加快. 但从蚊子总量看，相比于大时滞，小时滞会使环境中有更多的蚊子. 在实际中，我们应该寻找平衡这两种互反影响的方法.
当成熟时滞是周期函数τ(t)时，我们将通过数值模拟讨论时滞对系统(2)动力学行为的影响. 根据广州市的登革热数据，拟合出成熟时滞表达式为：
其中
并且
图 2显示了不同历史值对系统(2)周期解的影响. 可以看出：(μ(t)，ν(t))=(20，100)时，正常蚊子种群将会灭绝；(μ(t)，ν(t))=(10，200)时，正常蚊子种群与携带Wolbachia的蚊子种群共存，并呈现周期波动. 从取值来看，当携带Wolbachia的蚊子种群在历史值中所占比例较大时，Wolbachia会完全入侵野生蚊子种群.
图 3在图 2历史值的基础上，研究了不同出生率对系统(2)周期解的影响. 图 3a中(μ(t)，ν(t))=(10，200)，当b=20时两种蚊子种群共存，当b=40时正常蚊子种群灭绝，图 3b中(μ(t)，ν(t))=(20，100)，当b=10时两种蚊子种群共存，当b=20时正常蚊子种群灭绝，说明出生率越大，越有利于Wolbachia的传播.
-
本文主要研究成熟时滞对Wolbachia传播的影响. 研究发现：①对于常数时滞，时滞越大越有利于疾病的控制；②携带Wolbachia的蚊子种群在历史值中所占比例越大，越有利于Wolbachia的传播；③蚊子种群出生率越大，越有利于Wolbachia的传播.