基于脉冲接种扰动下随机传染病模型的稳定性分析

Stability Analysis Based on Random Infectious Disease Model under Pulse Vaccination Disturbance

摘要: 脉冲接种是一种有效控制疾病传播的方式，对传染病研究有重要意义。本文建立了一个具有标准发生率和垂直传播的噪声干扰脉冲接种随机SIR传染病模型，并研究了其在理论分析和数值模拟两个方面的动力学性质。首先构造辅助函数证明系统等价于一个不含脉冲的随机模型，并证明其正解的存在唯一性，其次利用伊藤公式给出疾病灭绝的充分条件，然后通过随机比较原理证明边界周期解的全局稳定性，最后用Matlab数值模拟来验证理论结果的正确性。

Abstract: Pulse vaccination is an effective way to control the spread of diseases and is of great significance for infectious disease research. In this paper, a random SIR infectious disease model of noise interfer-ence pulse inoculation with standard incidence and vertical propagation is established, and its ki-netic properties in both theoretical analysis and numerical simulation are studied. First, the auxil-iary function is constructed to prove that the system is equivalent to a random model without puls-es, and proves the existence uniqueness of its positive solution, and then uses Ito’s formula to give sufficient conditions for disease extinction, then proves the global stability of the boundary-periodic solution through the principle of random comparison, and finally uses Matlab numerical simulation to verify the correctness of the theoretical results.

1. 引言

在现实世界中，短暂的干扰会使事物的发展发生变化，传染病在传播的过程中也不可避免地会受到诸多脉冲现象的影响。靳桢，马知恩 [1] 首次证明了脉冲预防接种下SIR传染病模型无病周期解的存在性和全局渐近稳定性。Gao [2] 与Chen研究了一类具有潜伏期的脉冲接种SEIRS模型。大多数经典的传染病模型研究中发生率一般采用双线性发生率以及饱和发生率，Anderson与May [3] 的研究则证实了在某些时候采用标准发生率比线性发生率会更能符合一些疾病的传播机理。在传染病的传播方式中垂直传播指的是子女因遗传被已经感染了疾病的父母亲传染。Stevens [4] 研究分析了几种带有垂直感染的传染病模型，Ainseba [5] 等人分别对具有水平传播和垂直传播的SIS模型进行了对比。同时，环境噪声在不同程度下也影响着疾病的传播和发展，近年来有不少学者对随机传染病模型进行了研究。Ji [6] 等人研究了一个接触率扰动的随机SIR模型，讨论了疾病灭绝性和持久性。Dieu [7] 等人研究了一个具有退化扩散项的随机SIR模型，他们证明了其系统的不变分布的遍历性且估计了分布的收敛速度。此外还有很多学者研究了不同噪声下的随机传染病模型，但目前具有脉冲效应的随机传染病模型的研究较少，也很少考虑白噪声对系统动力学行为的影响。

考虑到疾病传播路径的多元化，本文综合的考察了具有标准发生率以及垂直感染的脉冲接种随机扰动SIR传染病模型，分析随机脉冲环境对疾病传播系统的稳定性及其相关动力学性质的影响。借助伊藤公式和随机比较原理等方法，通过构造Lyapunov函数，证明疾病传播系统正解存在唯一性，给出疾病灭绝的充分条件，并讨论边界周期解的全局稳定性，最后运用Matlab软件进行数值模拟，验证相关理论结果的有效性。

2. 模型描述

靳 [1] 讨论了脉冲SIR模型

$S^{\prime }\left(t\right)=\mu K-\beta S\left(t\right)I\left(t\right)-\mu S$,

$I^{\prime }\left(t\right)=\beta S\left(t\right)I\left(t\right)-\left(\mu +\alpha +\lambda \right)I\left(t\right)$,

$R^{\prime }\left(t\right)=\lambda I\left(t\right)-\mu R\left(t\right)$, $t\ne k$, $k\in N$, (1)

$S\left(k^{+}\right)=\left(1-p\right)S\left(k\right)$,

$I\left(k^{+}\right)=I\left(k\right)$,

$R\left(k^{+}\right)=R\left(k\right)+pS\left(k\right)$, $t=k$, $k\in N$,

其中 $S\left(t\right)$ 为易感者在时刻t的人数， $I\left(t\right)$ 为感染者在时刻t的人数， $R\left(t\right)$ 为康复者在时刻t的人数，K表示环境最大承载总人数， $\mu$ 表示出生率以及自然死亡率， $\beta$ 表示接触率， $\alpha$ 表示因病死亡率， $\lambda$ 表示恢复率，p表示预防接种率，以上参数都为正数。

本文在此基础上，考虑具有垂直感染 $q\mu I$ [2] [3]，标准发生率 $\beta SI/N$ [4] [5]，其中q为感染者的新生儿为感染者的比率( $0\le q\le 1$ )， $N=S\left(t\right)+I\left(t\right)+R\left(t\right)\le K$，得如下系统

$S^{\prime }\left(t\right)=\mu K-\frac{\beta S\left(t\right)I\left(t\right)}{N}-\mu S-q\mu I\left(t\right)$,

$I^{\prime }\left(t\right)=\frac{\beta S\left(t\right)I\left(t\right)}{N}-\left(\mu +\alpha +\lambda \right)I\left(t\right)+q\mu I\left(t\right)$,

$R^{\prime }\left(t\right)=\lambda I\left(t\right)-\mu R\left(t\right)$, $t\ne k$, $k\in N$, (2)

$S\left(k^{+}\right)=\left(1-p\right)S\left(k\right)$,

$I\left(k^{+}\right)=I\left(k\right)$,

$R\left(k^{+}\right)=R\left(k\right)+pS\left(k\right)$, $t=k$, $k\in N$.

但在现实生活中系统不可避免地会受到噪声干扰，假设噪声与系统变量成正比，其中 $R\left(t\right)$ 独立，考虑如下系统 [6] [7]

$\text{d}S\left(t\right)=\left[\mu K-\frac{\beta S\left(t\right)I\left(t\right)}{N}-\mu S\left(t\right)-q\mu I\left(t\right)\right]\text{d}t+{\sigma }_{1}S\left(t\right)\text{d}{B}_1(\; t\; )$

$\text{d}I\left(t\right)=\left[\frac{\beta S\left(t\right)I\left(t\right)}{N}-\left(\mu +\alpha +\lambda \right)I\left(t\right)+q\mu I\left(t\right)\right]\text{d}t+{\sigma }_{2}I\left(t\right)\text{d}{B}_2\left(t\right)$, $t\ne k$, $k\in N$, (3)

$S\left(k^{+}\right)=\left(1-p\right)S\left(k\right)$, $t=k$, $k\in N$,

其中 $B_i\left(t\right)\left(i=1,2\right)$ 为相互独立的标准布朗运动，其初值为 $B_i\left(0\right)=0$， ${\sigma }_1$， ${\sigma }_2$ 为相应噪声强度。

3. 全局正解的存在唯一性

定理1对于任意给定的初值 $S\left(0\right)>0$， $I\left(0\right)>0$，当 $t\ge 0$ 时，系统(3)有唯一正解 $\left(S\left(t\right),I\left(t\right)\right)$。

证明首先考虑如下无脉冲的系统

$\begin{array}{l}\text{d}x\left(t\right)=\left[\mu K{A}^{-1}\left(t\right)-\frac{\beta x\left(t\right)y\left(t\right)}{N}-\mu x\left(t\right)-q\mu y\left(t\right)\right]\text{d}t+{\sigma }_{1}x\left(t\right)\text{d}{B}_1\left(t\right)\\ \text{d}y\left(t\right)=\left[\frac{\beta A\left(y\right)x\left(t\right)y\left(t\right)}{N}-\left(\mu +\alpha +\lambda \right)y\left(t\right)+q\mu y\left(t\right)\right]\text{d}t+{\sigma }_{2}y\left(t\right)\text{d}{B}_2\left(t\right)\end{array}$ (4)

初值 $\left(x\left(0\right),y\left(0\right)\right)=\left(S\left(0\right),I\left(0\right)\right)$，其中 $A\left(t\right)=\{\begin{array}{l}{\left(1-p\right)}^{\left[t\right]-t},t\ne k\in N\\ {\left(1-p\right)}^{-1},t=k\in N\end{array}$， $A\left(t\right)$ 左连续且是1-周期的。

由随机微分方程解的存在定理 [8] [9] 知，系统(4)关于初值 $\left(x\left(0\right),y\left(0\right)\right)=\left(S\left(0\right),I\left(0\right)\right)$ 存在唯一的局部解 $\left(x\left(t\right),y\left(t\right)\right)$， $t\in \left[0,{\tau }_e\right)$，其中 ${\tau }_e$ 为爆破时间。这里定义停时 ${\tau }_k=\inf \left\{t\in \left[0,{\tau }_e\right]:x\left(t\right)\le 0或y\left(t\right)\le 0\right\}$，假设 $\inf \varphi =\infty$ ( $\varphi$ 表空集)，则由停时定义可知，如果可以证明 ${\tau }_k=\infty$，那么当 $t>0$ 时就有 ${\tau }_e=\infty$，且正解 $\left(x\left(t\right),y\left(t\right)\right)$ 必存在；否则存在常数 $T>0$ 和 $\epsilon \in \left(0,1\right)$ 使得 $P\left\{{\tau }_k\le T\right\}>\epsilon$。

下面证明 ${\tau }_k=\infty$，考虑如下Lyapunov函数 [10] $V\left(x\left(t\right),y\left(t\right)\right)=-\ln \left(\frac{\sigma x}{b}\right)-\ln \left(\frac{\sigma y}{b}\right)$，其中b为正常

数，由Itô公式得

$\begin{array}{c}\text{d}V\left(x,y\right)=-\frac{1}{x}\text{d}x-\frac{1}{y}\text{d}y+\frac{1}{2x^2}{\left(\text{d}x\right)}^2+\frac{1}{z{y}^2}{\left(\text{d}y\right)}^2\\ =\left[-\frac{\mu k{A}^{-1}}{x}+\frac{\beta y}{N}+\mu +q\mu \frac{y}{x}+\frac{1}{2}{\sigma }_1^2{A}^2\right]\text{d}t\\ +\left[-\frac{\beta Ax}{N}+\left(\mu +\alpha +\lambda \right)-q\mu +\frac{1}{2}{\sigma }_2^2\right]\text{d}t-{\sigma }_{1}A\text{d}{B}_1-{\sigma }_2\text{d}{B}_2\end{array}$,

则 $\text{d}V\left(x,y\right)\le \left[\frac{\beta y}{N}+\mu +q\mu \frac{y}{x}+\frac{1}{2}{\sigma }_1^2{A}^2\right]\text{d}t+\left[\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2\right]\text{d}t=Q\left(x,y\right)$，对其左右两边同时从0到t积分得 $V\left(x,y\right)\le V\left(x\left(0\right),y\left(0\right)\right)+{\displaystyle {\int }_0^{t}Q\left(x\left(s\right),y\left(s\right)\right)\text{d}s}=P\left(x,y\right)$，又因为 $V\left(x\left(t\right),y\left(t\right)\right)=-\ln \left(\frac{\sigma x}{b}\right)-\ln \left(\frac{\sigma y}{b}\right)$，且由停时的定义 [11] 可知 $\underset{t\to {\tau }_k}{\lim }V\left(x,y\right)=\infty$，可得出矛盾 $\infty \le P\left(x,y\right)<\infty$，故必有 ${\tau }_k=\infty$，得证。

下证系统(4)与(3)等价。

令 $\left(S\left(t\right),I\left(t\right)\right)=\left(A\left(t\right)x\left(t\right),y\left(t\right)\right)$，可知 $\left(x\left(t\right),y\left(t\right)\right)$ 在 $\left(k,k+1\right)\subset \left[0,+\infty \right)$ 是连续的， $k\in N$。

当 $t\ne k\in N$ 时，有

$\begin{array}{c}\text{d}S\left(t\right)=A^{-1}\left(t\right)x\left(t\right)\text{d}t+A\left(t\right)\text{d}x\left(t\right)\\ =A\left(t\right)\left[\mu K{A}^{-1}\left(t\right)-\frac{\beta x\left(t\right)y\left(t\right)}{N}-\mu x\left(t\right)-q\mu y\left(t\right)\right]\text{d}t+{\sigma }_{1}A\left(t\right)x\left(t\right)\text{d}{B}_1\left(t\right)\\ =\left[\mu K-\frac{\beta S\left(t\right)I\left(t\right)}{N}-\mu S\left(t\right)-q\mu I\left(t\right)\right]\text{d}t+{\sigma }_{1}S\left(t\right)\text{d}{B}_1\left(t\right)\end{array}$.

当 $t=k\in N$ 时，

$S\left(k^{-}\right)=\underset{t\to k^{-}}{\lim }A\left(t\right)x\left(t\right)={\left(1-p\right)}^{\left(k-1\right)-k}x\left(k\right)=S\left(k\right)$,

$S\left(k^{+}\right)=\underset{t\to k^{+}}{\lim }A\left(t\right)x\left(t\right)={\left(1-p\right)}^{k-k}x\left(k\right)=x\left(k\right)=\left(1-p\right)S\left(k\right)$.

同理可证，当 $t\ne k\in N$ 时，

$\text{d}I\left(t\right)=\left[\frac{\beta A\left(y\right)x\left(t\right)y\left(t\right)}{N}-\left(\mu +\alpha +\lambda \right)y\left(t\right)+q\mu y\left(t\right)\right]\text{d}t+{\sigma }_{2}y\left(t\right)\text{d}{B}_2\left(t\right)$.

当 $t=k\in N$ 时有 $I\left(k^{+}\right)=\underset{t\to k^{+}}{\lim }y\left(t\right)=\underset{t\to k^{-}}{\lim }y\left(t\right)=I\left(k\right)$。

综上，系统(4)与(3)等价。对于任意初值 $S\left(0\right)>0$， $I\left(0\right)>0$，我们由系统(4)正解的存在唯一性就能得到系统(3)的正解存在唯一性。

4. 系统的灭绝性

定理2假设 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 成立，令 $\left(S\left(t\right),I\left(t\right)\right)$ 为满足系统(3)任意初值的解，若 $R_0<1$，其中 $R_0=\frac{\beta {\langle S^{\ast }\rangle }_1}{\left(1-q\right)\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2}$，则 $\underset{t\to \infty }{\lim }\sup \frac{1}{t}\ln I\left(t\right)\le \left[\left(1-q\right)\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2\right]\left(R_0-1\right)<0$ a.s.，即疾病会灭绝。

证明首先构造一个1-周期的辅助函数 $h\left(t\right)$，它满足方程 $h^{\prime }\left(t\right)-h\left(t\right)\left(\mu -\ln \left(1-p\right)\right)=-A\left(t\right)$ 的正1-周期解。经过计算得 $h\left(t\right)=\frac{{\text{e}}^{\mu t}{\left(1-p\right)}^{-t}{\displaystyle {\int }_t^{t+1}{\text{e}}^{-\mu s}{\left(1-p\right)}^{\left[s\right]}\text{d}s}}{1-{\text{e}}^{-\mu }\left(1-p\right)},t\ge 0$，当 $t\in \left[0,1\right]$ 时有 $h\left(t\right)=\frac{1}{\mu }{\left(1-p\right)}^{-t}\left(1-\frac{p{\text{e}}^{\mu t}}{{\text{e}}^{\mu }-1+p}\right),t\in \left[0,1\right]$，又 $A\left(t\right)$ 与 $h\left(t\right)$ 都是1-周期函数，计算得到： ${\displaystyle {\int }_0^t{A}^{-1}\left(s\right)h\left(s\right)\text{d}s}={\displaystyle {\int }_0^1{A}^{-1}\left(s\right)h\left(s\right)\text{d}s}=\frac{1}{\mu }{\displaystyle {\int }_0^1\left(1-\frac{p{\text{e}}^{\mu s}}{{\text{e}}^{\mu }-1+p}\right)\text{d}s}=\frac{1}{\mu }\left[1-\frac{p\left({\text{e}}^{\mu }-1\right)}{\mu \left({\text{e}}^{\mu }-1+p\right)}\right]=\frac{1}{\mu K}{\langle S^{*}\rangle }_1$。

因为系统(2)与系统(3)无病周期解相同，用 $\left(S^{*}\left(t\right),0,R^{*}\left(t\right)\right)$ 表示， $S^{*}\left(t\right)=K\left[1-\frac{p{\text{e}}^{-\mu \left(t-k-1\right)}}{{\text{e}}^{\mu }-1+p}\right]$， $R^{*}\left(t\right)=K-S^{*}\left(t\right)$， $k<t\le k+1$，经过计算得到：

$\frac{1}{n}{\displaystyle {\int }_0^n{S}^{*}\left(t\right)\text{d}t}={\displaystyle {\int }_0^1{S}^{*}\left(t\right)\text{d}t}=K-\frac{Kp\left({\text{e}}^{\mu }-1\right)}{\mu \left({\text{e}}^{\mu }-1+p\right)}\triangleq {\langle S^{*}\rangle }_1$, $n\in {{\rm N}}^{+}$.

从而可以得到基本再生数 $R_0=\frac{\beta {\langle S^{\ast }\rangle }_1}{\left(1-q\right)\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2}$。

当 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$，则有

$\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}x\left(s\right)h\left(s\right)\text{d}{B}_1\left(s\right)}=0$, a.s.

$\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}y\left(s\right)h\left(s\right)\text{d}{B}_2\left(s\right)}=0$, $\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}y\left(s\right)\text{d}{B}_2\left(s\right)}=0$, a.s. (5)

$\underset{t\to \infty }{\lim }\frac{x\left(t\right)}{t}=0$, $\underset{t\to \infty }{\lim }\frac{y\left(t\right)}{t}=0$, a.s.

由Itô公式 [11] 得，

$\begin{array}{l}\text{d}\left[h\left(t\right)x\left(t\right)\right]=h^{\prime }\left(t\right)x\left(t\right)\text{d}t+h\left(t\right)\text{d}x\left(t\right)\\ =\left(\mu K{A}^{-1}\left(t\right)h\left(t\right)-A\left(t\right)\left(x\left(t\right)+qy\left(t\right)\right)-\frac{\beta h\left(t\right)x\left(t\right)y\left(t\right)}{N}\right)\text{d}t+{\sigma }_{1}h\left(t\right)x\left(t\right)\text{d}{B}_1\left(t\right)\\ =\left\{\mu K{A}^{-1}\left(t\right)h\left(t\right)+h^{\prime }\left(t\right)x\left(t\right)-h\left(t\right)\left[\mu -\ln \left(1-p\right)\right]\left(x\left(t\right)+qy\left(t\right)\right)-\frac{\beta h\left(t\right)x\left(t\right)y\left(t\right)}{N}\right\}\text{d}t\\ +{\sigma }_{1}h\left(t\right)x\left(t\right)\text{d}{B}_1\left(t\right)\end{array}$.

结合(5)式有

$\begin{array}{c}0=\underset{t\to \infty }{\lim }\frac{h\left(t\right)x\left(t\right)}{t}\\ =\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^t\mu K{A}^{-1}\left(s\right)x\left(s\right)\text{d}s}-\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}A\left(s\right)x\left(s\right)\text{d}s}\\ -\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}qA\left(s\right)y\left(s\right)\text{d}s}-\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^t\frac{\beta h\left(s\right)x\left(s\right)y\left(s\right)}{N}\text{d}s}\\ ={\langle S^{*}\rangle }_1-\underset{t\to \infty }{\lim }{\langle Ax\rangle }_t-q\underset{t\to \infty }{\lim }{\langle Ay\rangle }_t-\underset{t\to \infty }{\lim }{\langle \frac{\beta hxy}{N}\rangle }_t\end{array}$. (6)

类似地，可以得到

$0=\underset{t\to \infty }{\lim }\frac{y\left(t\right)}{t}=\underset{t\to \infty }{\lim }{\langle \frac{\beta Axy}{N}\rangle }_t-\left[\left(1-q\right)\mu +\alpha +\lambda \right]\underset{t\to \infty }{\lim }{\langle y\rangle }_t$. (7)

由Itô公式得

$\text{d}\ln y\left(t\right)=\left\{\frac{\beta A\left(t\right)x\left(t\right)}{N}-\left[\left(1-q\right)\mu +\alpha +\lambda \right]-\frac{1}{2}{\sigma }_2^2\right\}\text{d}t+{\sigma }_2\text{d}{B}_2\left(t\right)$,

将上式两端同时从0到t积分并乘 $\frac{1}{t}$ 有

$\begin{array}{c}\frac{1}{t}\ln y\left(t\right)={\langle \frac{\beta Ax}{N}\rangle }_t-\left[\left(1-q\right)\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2\right]+\frac{1}{t}\left({\sigma }_2{B}_2\left(t\right)+\ln y\left(0\right)\right)\\ \le \beta {\langle S^{*}\rangle }_1+\left[\left(1-q\right)\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2\right]+\frac{1}{t}\Phi \left(t\right)\end{array}$,(8)

其中

$\Phi \left(t\right)=\left\{\left[\left(1-q\right)\mu +\alpha +\lambda \right]{\langle y\rangle }_t-{\langle \frac{\beta Axy}{N}\rangle }_t\right\}+\beta t\left({\langle Ax\rangle }_t+q{\langle Ay\rangle }_t+{\langle \frac{\beta hxy}{N}\rangle }_t-{\langle S^{*}\rangle }_1\right)+{\sigma }_2{B}_2\left(t\right)+\ln y\left(0\right)$,

从而由布朗运动性质、式子(5)及(7)知 $\underset{t\to \infty }{\lim }\frac{1}{t}\Phi \left(t\right)=0$，a.s.。

若有 $R_0<1$ 成立，对(8)式两边取极限有

$\underset{t\to \infty }{\lim }\sup \frac{1}{t}\ln I\left(t\right)\le \left(\mu +\alpha +\lambda +\frac{1}{2}{\sigma }_2^2\right)\left(R_0-1\right)<0$, a.s. (9)

又因为 $y\left(t\right)=I\left(t\right)$，也就是说

$\underset{t\to \infty }{\lim }I\left(t\right)=0$, a.s. (10)

综上，疾病将以指数趋势灭绝。

5. 边界周期解的全局稳定性

定理3 假设 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 成立，若 $R_0<1$，则存在系统(3)的一个边界周期解 $\left(S_p\left(t\right),0\right)$，且

这个边界周期解具有全局稳定性。

证明由于 $\left(S\left(t\right),I\left(t\right)\right)=\left(A\left(t\right)x\left(t\right),y\left(t\right)\right)$，只需证明系统(4)有一边界周期解 $\left(x_p\left(t\right),0\right)$ 且其是全局吸引

的即可。

若满足 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 和 $R_0<1$，那么由定理2可得 $\underset{t\to \infty }{\lim }y\left(t\right)=0$，a.s.。对于充分小的 $\tau >0$，存在一正常数 $t_0=t_0\left(\omega \right)$ 和集合 ${\Omega }_{\tau }\subset \Omega$ 满足 $P\left({\Omega }_{\tau }\right)>1-\tau$ 且对 $t>t_0$， $\omega \in {\Omega }_{\tau }$ 有 $y\left(t\right)<{\text{e}}^{-ct}$ 与 $x\left(t\right)y\left(t\right)<\tau$。

对于 $t>t_0$， $\omega \in {\Omega }_{\tau }$ 有

$\begin{array}{c}\text{d}x\left(t\right)=\left(\mu K{A}^{-1}\left(t\right)-\mu x\left(t\right)+\mu qy\left(t\right)-\frac{\beta x\left(t\right)y\left(t\right)}{N}\right)\text{d}t+{\sigma }_{1}x\left(t\right)\text{d}{B}_1\left(t\right)\\ \ge \left(\mu K{A}^{-1}\left(t\right)-\mu x\left(t\right)+\mu qy\left(t\right)-\beta \tau \right)\text{d}t+{\sigma }_{1}x\left(t\right)\text{d}{B}_1\left(t\right)\end{array}$.

令 $\text{d}\underset{\_}{X}=\left(\mu K{A}^{-1}\left(t\right)-\mu \left(t\right)\underset{\_}{X}\left(t\right)+\mu qy\left(t\right)-\beta \tau \right)\text{d}t+{\sigma }_1\underset{\_}{X}\left(t\right)\text{d}{B}_1\left(t\right)$，其具有初值 $\underset{\_}{X}\left(0\right)=x\left(0\right)$ 的解 $\bar{X}\left(t\right)$。

由随机微分方程的比较定理 [12] 可知，对 $\omega \in {\Omega }_{\tau }$ 有

$\underset{\_}{X}\left(t\right)\le x\left(t\right)\le \bar{X}\left(t\right)$, $t>t_0$,

令 $\tau$ 趋于0有 $\underset{t\to \infty }{\lim }\left|\underset{\_}{X}\left(t\right)-\bar{X}\left(t\right)\right|=0$，那么 $\underset{t\to \infty }{\lim }\left|x\left(t\right)-\bar{X}\left(t\right)\right|=0$，则 $\underset{t\to \infty }{\lim }\left|x\left(t\right)-x_p\left(t\right)\right|=0$。

因此 $\left(x_p\left(t\right),0\right)$ 是全局吸引的，由于 $\left(S\left(t\right),I\left(t\right)\right)=\left(A\left(t\right)x\left(t\right),y\left(t\right)\right)$，那么系统(3)的边界周期解 $\left(S_p\left(t\right),0\right)=\left(A\left(t\right)x_p\left(t\right),0\right)$ 也是全局吸引的。

6. 数值模拟

通过研究具有垂直感染和脉冲接种的SIR模型的动力学性质，使用Matlab来验证所得结论。假定人口单位为百万，其中系统初值 $\left(S\left(0\right),I\left(0\right)\right)=\left(0.7,0.2\right)$，为演示白噪声的影响以及疾病的阈值，选择参数： $K=1$ ； $\mu =0.2$ ； $\beta =0.8$ ； $\alpha =0.05$ ； $\lambda =0.2$。

有以下四种情况：

1) 当随机系统(3)满足定理2的假设条件 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 时，令 $q=0.\text{3}$， ${\sigma }_1=0.03$， ${\sigma }_2=0.25$， $p=0.1$，此时 $R_0=1.2450>1$，由定理2，疾病 $I\left(t\right)$ 将持续，见图1中的(b)图。

2) 当随机系统(3)满足定理2的假设条件 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 时，令 $q=0.\text{3}$， ${\sigma }_1=0.03$， ${\sigma }_2=0.25$， $p=0.3$，此时 $R_0=0.6895<1$，由定理2，疾病 $I\left(t\right)$ 将灭绝，见图1中的(c)图。

3) 当随机系统(3)不满足定理2的假设条件 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 时， $q=0.4$， ${\sigma }_1=0.03$， ${\sigma }_2=0.5$， $p=0.3$，此时 $R_0=0.5868<1$，由图1中的(e)图看出疾病将灭绝。

4) 当随机系统(3)不满足定理2的假设条件 $\left(1-q\right)\mu >\frac{1}{2}\left({\sigma }_1^2\vee {\sigma }_2^2\right)$ 时， $q=0.4$， ${\sigma }_1=0.03$， ${\sigma }_2=0.5$， $p=0.1$，此时 $R_0=1.0595>1$，由图1中的(f)图可以看出疾病也将灭绝。

图1分别考察了易感人群( $S\left(t\right)$ )以及感染人群( $I\left(t\right)$ )随接种率变化(p分别等于0 (a, d)，0.1 (b, e)，0.3 (c, f))系统的动态轨线。从图中可以看出，当随机系统(3)满足定理2时，疾病以概率1灭绝，这与结论是一致的。当满足假设条件但 $R_0>1$ 时，疾病持续，如图2(b)所示。图2(b)与图2(c)除了参数p不同，其余参数保持不变，说明提高易感者预防接种率有利于疾病的灭绝。而当随机系统(3)不满足定理2假设条件时，即 $\left(1-q\right)\mu <\left({\sigma }_1^2\vee {\sigma }_2^2\right)/2$，此时白噪声较大，可以从图2(e)和图2(f)看出，无论 $R_0$ 是否大于1，疾病都会以指数速率灭绝，即系统白噪声对于抑制疾病的传播是有积极作用的。

图2考察了垂直感染率(q)以及噪声强度( $\sigma$ )对基本再生数( $R_0$ )的定量影响。显然，当外界噪声强度保持不变时，垂直感染率的增加将会有效的减少基本再生数的阈值效应；反之，外界环境噪声强度增大则将增加基本再生数，使得疾病传播性增强。

Figure 1. Analysis of persistence and extinction of diseases in system

图1. 系统中疾病持续与灭绝性分析

Figure 2. Diagram of the relationship between vertical infection rate and noise intensity of basic reproduction number R₀

图2. 基本再生数R₀关于垂直感染率和噪声强度关系图

7. 结束语

本文综合的考察了具有标准发生率以及垂直感染的脉冲接种随机扰动SIR传染病模型的动力学行为，主要借助伊藤公式和随机比较原理等方法，通过构造Lyapunov函数，证明疾病传播系统正解存在唯一性，分析系统的灭绝性、持久性以及讨论边界周期解的全局稳定性。理论研究发现白噪声在一定程度上对模型产生了影响，这对疾病传播的研究有着应用意义。但模型和研究方法都比较基础单一，对于复杂的混合脉冲以及有彩色噪声等的随机模型都值得关注和研究。

基金项目

海南省重点研发项目(ZDYF2021SHFZ231)，海南省自然科学基金(120RC451)，国家自然科学基金(11761025, 11961018)，海南省研究生创新创业项目(Qhys2021-208)。

NOTES

^*通讯作者。