受扰动的捕食者–被捕食者模型的随机分岔分析

Stochastic Bifurcation Analysis of a Predator-Prey Model with Perturbations

摘要: 本研究针对一类具有随机增长率的捕食者–被捕食者系统，深入探讨了两种捕食者竞争单一被捕食者资源的动态行为及其随机演化规律。通过在不变曲面集上进行多阶段合理变形，推导出Fokker-Plank方程的封闭形式稳态解。在确立系统持久性的前提下，创新性地以环境噪声强度为分岔控制参数，揭示了所有种群的现象学分岔特性。本研究建立的随机分岔判据与竞争量化模型，为生物多样性保护与入侵物种控制策略的优化设计提供了动力学依据。

Abstract: In this paper, we focused on a predator-prey system with a random growth rate and deeply explored the dynamic behavior and stochastic evolution of two types of predators competing for a single prey resource. By performing multi-stage reasonable deformation on the invariant surface set, the closed-form steady-state solution of the Fokker-Plank equation is derived. Under the premise of establishing the persistence of the system, the environmental noise intensity was innovatively used as the bifurcation control parameter to reveal the phenomenological bifurcation characteristics of all populations. The stochastic bifurcation criterion and competition quantitative model established in this study provide a dynamic basis for the optimal design of biodiversity conservation and invasive species control strategies.

1. 引言

关于研究背景与意义，捕食者–被捕食者关系作为生态学中最基础的相互作用形式，其动力学研究可追溯至Lotka [1] and Volterra [2]的开创性工作，二者提出的经典模型为后续研究奠定了理论基础。近几十年来，随着生态动力学理论的不断发展，研究者们建立了包括确定性模型与随机模型在内的多种理论框架，用于揭示生物系统中的复杂动态行为(如竞争排斥原理、多物种共存、周期性振荡等[3]-[12])。这些研究成果不仅深化了人们对生态系统稳定机制的理解，更为生物防治策略的制定提供了理论支持。

针对两种捕食者竞争单一被捕食者的三维系统，学界提出了典型的确定性动力学模型[13]：

$\{\begin{array}{l}\text{d}N=\left[N\left(r-N/K\right)-{\alpha }_{1}Nx-{\alpha }_{2}Ny\right]\text{d}t,\\ \text{d}x=\left({\beta }_{1}N-{\mu }_1\right)x\text{d}t,\\ \text{d}y=\left({\beta }_{2}N-{\mu }_2\right)y\text{d}t.\end{array}$ (1)

其中，$N\left(t\right)$ 表示被捕食者种群密度，与其有关的参数有：固有增长率$r$ (被捕食者理想状态条件下的最大增长速率)、环境容纳量$K$ (环境中被捕食者能维持的最大种群密度)、捕食压力${\alpha }_1,{\alpha }_2$ (单位时间内捕食者$x\left(t\right)$ 、$y\left(t\right)$ 对被捕食者的捕食效率)；$x\left(t\right)$ 与$y\left(t\right)$ 为两种捕食者密度，与其有关的参数有：能量转化率${\beta }_1,{\beta }_2$ (捕食者通过捕食获得的能量转化为自身繁殖的效率)、死亡率${\mu }_1,{\mu }_2$ (捕食者在无被捕食者时的自然死亡率)。Llibre与Xiao [13]假设被捕食者资源是有限的，证明了竞争排斥原理的充分和必要条件，并揭示了周期振荡与稳态共存这两种三个物种共存的状态。后续研究进一步拓展至非线性函数情形，相关成果见Hsu [14] [15]、Groll [9]等学者的工作。

随机噪声因其理论和实践意义，在建模种群动态中被视为一个重要因素。通过将布朗运动引入确定性模型参数作为白噪声，可构建更贴近实际的随机微分方程。研究表明[10] [16]-[19]，随机扰动不仅能抑制种群数量的爆炸性增长，还可诱导系统产生新的稳定态。例如，Mao等人[16]发现微弱的环境噪声即可改变系统的全局稳定性；Zhu与Yin [11]通过长期平均分析，提出了随机竞争系统的部分普适性规律。特别地，Zhang与Jiang [12]在Lotka-Volterra模型中引入白噪声扰动，利用马尔可夫半群理论证明了共存原理成立的充分条件。这些进展凸显了将随机因素引入生态系统模型在生态系统建模中的理论价值与应用潜力。

尽管随机捕食者–被捕食者系统的研究已取得显著进展，但现有文献在持久性条件下系统稳态分布的分岔特性上依然缺乏定量分析。本文旨在引入环境噪声到确定性模型(1)中，以构建随机Lotka-Volterra模型，并分析噪声的影响，解决以上问题。

本文在第二节中展示了模型的构建以及目标；第三节证明了闭式稳态分布以便于进行后续分析；第四节通过选择噪声强度作为参数研究随机分岔。

2. 建立模型

2.1. 被捕食者增长率的随机项

模型(1)是一个确定性种群模型，通过理论分析可以将其改写为以下式子：

$\{\begin{array}{l}x\left(t\right)=x_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}},\\ y\left(t\right)=y_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}},\end{array}$ (2)

另外，被捕食者$N$ 满足以下方程，该方程与捕食者$x$ 和$y$ 无关：

$dN=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_1{x}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}}-{\alpha }_2{y}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}\right]\text{d}t,N\left(0\right)=N_0.$ (3)

也就是说，如果对(3)中定义的被捕食者的特性有充分的认识，那么就可以对捕食者$x$ 和$y$ 进行充分的分析。另外，对于某些位于食物链下游的被捕食者来说，其更容易受到温度、湿度、日照时间等自然环境的影响，因此无论是从理论上还是实践上，研究随机环境中被捕食者种群的演化规律都具有重要意义。以草履虫为例，众所周知它是位于食物链下游的动物，其平均生长率为0.4666%，数据来自实验[20] (见图1)。

Figure 1. Natural growth rate of large Paramecium

图1. 大型草履虫的自然增长速率

受到以上启发，在本文中，我们主要假设被捕食者的增长率$r$ 受到随机因素的影响，即$r\text{d}t\to r\text{d}t+\sigma \text{d}B\left(t\right)$ ，其中$B\left(t\right)$ 是标准布朗运动，$\sigma >0$ 表示噪声强度，从而将模型(3)转化为随机形式：

$\text{d}N=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_1{x}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}}-{\alpha }_2{y}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}\right]\text{d}t+\sigma N\text{d}B\left(t\right),N\left(0\right)=N_0.$ (4)

综上我们得到了由(2)和(4)定义的随机Lotka-Volterra模型，接下来我们将对其展开讨论。为了叙述方便，我们将它们合并为：

$\{\begin{array}{l}x\left(t\right)=x_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}},\\ y\left(t\right)=y_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}},\\ \text{d}N=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_1{x}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}}-{\alpha }_2{y}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}\right]\text{d}t+\sigma N\text{d}B\left(t\right),N\left(0\right)=N_0\end{array}$ (5)

理论上，Llibre和Xiao [13]证明了模型(1)中三个物种的两种共存状态：周期振荡和稳态。然而，与之前确定性的讨论不同，由于噪声的扰动，系统状态不会收敛到确定性的周期振荡或稳态，尤其是在持续性的情况下。因此，在本文中，发现环境噪声对捕食者与被捕食者之间的动态影响是我们的主要目标。具体可以分为以下几个方面：在噪声扰动下，求出封闭形式的稳态分布；在持续的情况下，通过有限的密度分布来显示由噪声引起的分岔；研究初始值对物种生存水平的影响。

定理2.1 对于任意初始值$\left(N_0,x_0,y_0\right)\in R_{+}^3$ ，方程(5)的解$\left(N\left(t\right),x\left(t\right),y\left(t\right)\right)$ 对于任意$t>0$ 都以概率1保持在$R_{+}^3$ 中。

证明：验证解的存在唯一性：

系统(5)的随机微分方程(SDE)为：

$\text{d}N=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_1{x}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}}-{\alpha }_2{y}_{0}N{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}\right]\text{d}t+\sigma N\text{d}B\left(t\right)$

初始条件为$N\left(0\right)=N_0>0$ 。

局部Lipschitz条件：漂移项和扩散项对$N$ 的依赖性为多项式与指数函数的组合。对于任意紧集$N\in \left[0,M\right]$ ，漂移项和扩散项关于$N$ 是局部Lipschitz连续的。线性增长条件：漂移项的增长速率被多项式$N^2/K$ 和指数项控制，扩散项为$\sigma N$ ，均满足线性增长条件。因此，根据随机微分方程的存在唯一性定理，系统在局部时间内存在唯一解。

全局解的延拓：

为了证明解在全局时间$t>0$ 内存在且保持正性，需验证解不会在有限时间内爆炸或到达零。

构造Lyapunov函数$V\left(N\right)=\ln N$ ，应用伊藤公式：

$\text{d}\left(\ln N\right)=\left[\frac{1}{N}\cdot 漂移项-\frac{{\sigma }^2}{2}\right]\text{d}t+\sigma \text{d}B\left(t\right).$

漂移项中的负项和捕食者项会抑制$N$ 的过度增大。因此，$\ln N$ 的增长率被控制，解不会在有限的时间内趋向无穷大。

当$N\left(t\right)$ 接近零时，方程的主导项为：

$\text{d}N\approx \left[rN-{\alpha }_{1}x\left(t\right)N-{\alpha }_{2}y\left(t\right)N\right]\text{d}t+\sigma N\text{d}B\left(t\right).$

由于$x\left(t\right)$ 和$y\left(t\right)$ 的表达式为指数积分，当$N$ 很小时，指数的积分可能趋向负值，导致$x\left(t\right)$ 和$y\left(t\right)$ 指数衰减。此时捕食者的影响趋近于零，漂移项由$rN$ 主导(正项)，推动$N\left(t\right)$ 远离零。同时，扩散项$\sigma N$ 在$N\to 0$ 时趋近于零，减少随机扰动的影响。

因此，当$N\left(t\right)$ 接近零时，漂移项使其无法达到零，即$\underset{t\to \tau }{\lim \inf }N\left(t\right)>0$ a.s.，其中$\tau$ 为首次触及零的时间，这说明$\tau =\infty$ 。

捕食者的动态由以下表达式给出：

$x\left(t\right)=x_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}},y\left(t\right)=y_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}.$

由于指数函数的非负性，且初始值$x_0,y_0>0$ ，根据之前步骤中$N\left(t\right)$ 的全局存在性和正性，积分路径有限，因此$x\left(t\right)$ 和$y\left(t\right)$ 保持严格正。

接下来，我们始终假设系统(5)具有全局唯一正解。为方便起见，我们记

${\lambda }_1=\frac{{\mu }_1}{{\beta }_1}$ , ${\lambda }_2=\frac{{\mu }_2}{{\beta }_2}$ , $\langle x\left(t\right)\rangle =\frac{1}{t}{\displaystyle {\int }_0^{t}x\left(s\right)\text{d}s}$ .

2.2. 多参数随机性扩展

在现有模型中，仅考虑了被捕食者增长率$r$ 的随机扰动。然而，实际生态系统中，捕食率${\alpha }_1$ 、${\alpha }_2$ ，捕食者死亡率${\mu }_1$ 、${\mu }_2$ 等参数同样可能受到环境随机性的影响。例如，捕食效率可能因猎物分布不均或气候波动而随机变化，而捕食者的死亡率可能因疾病暴发或资源短缺而呈现随机性。为此，本节提出一种将模型(5)扩展为多参数随机驱动系统的思路，以更全面地刻画生态系统的动态行为。

这样的思路面临的理论挑战与分析方法为：

1) 高维随机系统：扩展后的模型为多维随机微分方程，解析求解稳态分布或分岔条件将更为复杂。需借助降维技术(如不变流形理论)或近似方法(如小噪声摄动展开)。

2) 噪声耦合效应：不同参数的噪声可能通过非线性项产生耦合效应。例如，捕食率噪声与种群密度相乘，其影响随被捕食者数量动态变化，需分析噪声协方差矩阵对系统稳定性的调制作用。

3) 持久性与灭绝概率：多参数噪声可能改变系统的持久性阈值。例如，捕食者死亡率噪声的增加可能降低其生存概率，进而影响被捕食者的灭绝风险。需重新定义综合噪声强度阈值，并分析其与单参数噪声阈值的差异。

3. 封闭形式的稳定分布

在本文的讨论和分析中，平稳分布是一个重要的工具，因此首先证明平稳分布的存在性。在本节中，我们假设${\lambda }_1={\lambda }_2$ 且$r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}>0$ 。对于任意$c>0$ ，我们将研究系统(5)在表面集$F\left(N,x,y\right)=c$ 上的动力学问题，其中$F\left(N,x,y\right)=x^{{\beta }_2}{y}^{-{\beta }_1}$ 。

引理3.1 如果${\lambda }_1={\lambda }_2$ ，则存在常数$c$ 使得$x^{{\beta }_2}\left(t\right)y^{-{\beta }_1}\left(t\right)=c$ 对于任意的$t\ge 0$ 。

证明：由方程(5)可得：

$\begin{array}{c}{x}^{{\beta }_2}\left(t\right)y^{-{\beta }_1}\left(t\right)={\left(x_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u}}\right)}^{{\beta }_2}\cdot {\left(y_0{\text{e}}^{{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}\right)}^{-{\beta }_1}\\ =x_0^{{\beta }_2}{y}_0^{-{\beta }_1}{\text{e}}^{{\beta }_2{\displaystyle {\int }_0^t\left[{\beta }_{1}N\left(u\right)-{\mu }_1\right]\text{d}u-{\beta }_1{\displaystyle {\int }_0^t\left[{\beta }_{2}N\left(u\right)-{\mu }_2\right]\text{d}u}}}\end{array}$

其中指数项为：

${\beta }_2{\beta }_1{\displaystyle {\int }_0^{t}N\left(u\right)\text{d}u}-{\beta }_2{\mu }_{1}t-{\beta }_1{\beta }_2{\displaystyle {\int }_0^{t}N\left(u\right)\text{d}u}+{\beta }_1{\mu }_{2}t=\left(-{\beta }_2{\mu }_1+{\beta }_1{\mu }_2\right)t=0$ .

令$c=x_0^{{\beta }_2}{y}_0^{-{\beta }_1}$ ，则有$x^{{\beta }_2}\left(t\right)y^{-{\beta }_1}\left(t\right)=c,\forall t\ge 0$ 。

当${\lambda }_1={\lambda }_2$ 时，捕食者种群$x\left(t\right)$ 和$y\left(t\right)$ 的组合$x^{{\beta }_2}\left(t\right)y^{-{\beta }_1}\left(t\right)$ 保持恒定，表明两种捕食者在资源利用上达到平衡。

综上对于任意$c>0$ ，$F\left(N,x,y\right)=c$ 是(5)的不变集。在表面集合$F\left(N,x,y\right)=c$ 中，模型(5)可以简化为模型

$\{\begin{array}{l}\text{d}N=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_{1}x-{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]\text{d}t+\sigma N\text{d}B\left(t\right),\\ \text{d}x={\beta }_1\left(N-{\lambda }_1\right)x\text{d}t.\end{array}$ (6)

如果(6)存在平稳分布密度函数，则将其定义为$\psi \left(N,x\right)$ ，进一步我们定义其平稳分布的边缘密度函数：$\psi \left(N\right)={\displaystyle {\int }_0^{\infty }\psi \left(N,x\right)\text{d}x}$ 和$\psi \left(x\right)={\displaystyle {\int }_0^{\infty }\psi \left(N,x\right)\text{d}N}$ 。同样地，我们可以定义函数$\psi \left(N,y\right)$ 和$\psi \left(y\right)$ 。

注1. 在非退化情形下，随机Lotka-Volterra食物链模型的研究已由Hening与Nguyen在文献[8] [21]中完成，他们提出了系统随机持续性及物种灭绝的条件。近期，BenaIm、Bourquin与Nguyen在文献[22]中针对退化情形下的随机Lotka-Volterra食物链模型展开研究，证明了当顶端或底端物种受随机扰动时系统存在唯一不变概率测度。值得指出的是，文献[22]中的方法具有普适性，可推广至本文模型(5)以证明其平稳分布的存在性。因此，本文的研究重点将聚焦于平稳分布显式解的推导，而不再赘述其存在性及遍历性的证明过程。

定理3.2 在不变曲面集合$F\left(N,x,y\right)=c$ 中，如果${\lambda }_1={\lambda }_2$ 且$r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}＞0$ ，则(5)的边缘平稳分布是唯一的，并且有$\psi \left(N,x\right)=\psi \left(N\right)\psi \left(x\right)$ 和$\psi \left(N,y\right)=\psi \left(N\right)\psi \left(y\right)$ ，其中：

$\psi \left(N\right)={\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(N-1)}{N}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(u-1)}{u}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}}\text{d}u\right]}^{-1}$

$\psi \left(x\right)=\frac{{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)-1}}{{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}u+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{u}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{u}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)-1}}\text{d}u}$

$\psi \left(y\right)=\frac{{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}{c}^{-\frac{1}{{\beta }_2}}{y}^{\frac{{\beta }_1}{{\beta }_2}}+\frac{{\alpha }_2}{{\beta }_2}y\right]}{y}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{2}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)-1}}{{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}{c}^{-\frac{1}{{\beta }_2}}{u}^{\frac{{\beta }_1}{{\beta }_2}}+\frac{{\alpha }_2}{{\beta }_2}u\right]}{u}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{2}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)-1}}\text{d}u}$

证明：令${\theta }_1=\ln N$ ，${\theta }_2=\ln x$ ，则有$\text{d}{\theta }_1=\frac{1}{N}\text{d}N-\frac{1}{2N^2}{\left(\text{d}N\right)}^2$ ，$\text{d}{\theta }_2=\frac{1}{x}\text{d}x$ 。

代入(6)式可得：

$\{\begin{array}{l}\text{d}{\theta }_1=\left[\left(r-\frac{1}{K}{\text{e}}^{{\theta }_1}\right)-\frac{{\sigma }^2}{2}-{\alpha }_1{\text{e}}^{{\theta }_2}-{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{\text{e}}^{\frac{{\beta }_2}{{\beta }_1}{\theta }_2}\right]\text{d}t+\sigma \text{d}B\left(t\right),\\ \text{d}x={\beta }_1\left({\text{e}}^{{\theta }_1}-{\lambda }_1\right)\text{d}t.\end{array}$

设$\psi \left({\theta }_1\left(t\right),{\theta }_2\left(t\right)|{\theta }_1\left(s\right)={\theta }_1,{\theta }_2\left(s\right)={\theta }_2\right)$ 为(6)式从时间$s$ 到时间$t$ 的解过程的转移概率密度函数，那么$\psi$ 满足以下Fokker-Planck方程：

$\frac{\partial \psi }{\partial t}+\frac{\partial }{\partial {\theta }_1}\left[\left(r-\frac{1}{K}{\text{e}}^{{\theta }_1}\right)-\frac{{\sigma }^2}{2}-{\alpha }_1{\text{e}}^{{\theta }_2}-{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{\text{e}}^{\frac{{\beta }_2}{{\beta }_1}{\theta }_2}\right]\psi +\frac{\partial }{\partial {\theta }_2}\left[{\beta }_1\left({\text{e}}^{{\theta }_1}-{\lambda }_1\right)\psi \right]-\frac{{\sigma }^2}{2}\frac{{\partial }^2\psi }{\partial {\theta }_1^2}=0$

由注1我们可以注意到$\psi \left({\theta }_1,{\theta }_2\right)=\underset{t\to \infty }{\lim }\psi \left({\theta }_1\left(t\right),{\theta }_2\left(t\right)|{\theta }_1\left(s\right)={\theta }_1,{\theta }_2\left(s\right)={\theta }_2\right)$ 。因此，在稳态下，$\psi$ 满足平稳的Fokker-Planck方程[23]：

$\frac{\partial }{\partial {\theta }_1}\left[\left(r-\frac{1}{K}{\text{e}}^{{\theta }_1}\right)-\frac{{\sigma }^2}{2}-{\alpha }_1{\text{e}}^{{\theta }_2}-{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{\text{e}}^{\frac{{\beta }_2}{{\beta }_1}{\theta }_2}\right]\psi +\frac{\partial }{\partial {\theta }_2}\left[{\beta }_1\left({\text{e}}^{{\theta }_1}-{\lambda }_1\right)\psi \right]-\frac{{\sigma }^2}{2}\frac{{\partial }^2\psi }{\partial {\theta }_1^2}=0$ (7)

上式可以继续化简为解耦形式：

$L_1\left({\theta }_1,{\theta }_2\right)\left[G_1\left({\theta }_1\right)\psi +H_1\left({\theta }_1\right)\frac{\partial \psi }{\partial {\theta }_1}\right]+L_2\left({\theta }_1,{\theta }_2\right)\left[G_2\left({\theta }_2\right)\psi +H_2\left({\theta }_2\right)\frac{\partial \psi }{\partial {\theta }_2}\right]=0$

其中，$L_1\left({\theta }_1,{\theta }_2\right)=\frac{\partial }{\partial {\theta }_1}-{\beta }_{1}K\frac{\partial }{\partial {\theta }_2}$ ，$L_2\left({\theta }_1,{\theta }_2\right)=\frac{\partial }{\partial {\theta }_1}$ 是一阶偏微分算子，另外有：

$G_1\left({\theta }_1\right)=\frac{1}{K}\left({\lambda }_1-{\text{e}}^{{\theta }_1}\right)$ , $H_1\left({\theta }_1\right)=-\frac{{\sigma }^2}{2}$ , $G_2\left({\theta }_2\right)=r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}-{\alpha }_1{\text{e}}^{{\theta }_2}-{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{\text{e}}^{\frac{{\beta }_2}{{\beta }_1}{\theta }_2}$ , $H_2\left({\theta }_2\right)=-{\beta }_{1}K\frac{{\sigma }^2}{2}$ .

根据Liu的工作[24]，我们可以通过求解以下方程得到平稳概率密度函数$\psi$ ：

$G_1\left({\theta }_1\right)\psi +H_1\left({\theta }_1\right)\frac{\partial \psi }{\partial {\theta }_1}=0$ 和$G_2\left({\theta }_2\right)\psi +H_2\left({\theta }_2\right)\frac{\partial \psi }{\partial {\theta }_2}=0$

其中第一个式子可化简为：$\frac{\partial {\psi }_1}{\partial {\theta }_1}=\frac{2}{{\sigma }^2}\frac{1}{K}\left({\lambda }_1-{\text{e}}^{{\theta }_1}\right){\psi }_1$ 。

解得：${\psi }_1\left({\theta }_1\right)\propto \exp \left[\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}{\theta }_1-\frac{2}{{\sigma }^2}\frac{1}{K}\left({\text{e}}^{{\theta }_1}-1\right)\right]$ 。

同理第二个式子可化简为：$\frac{\partial {\psi }_2}{\partial {\theta }_2}=\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}-{\alpha }_1{\text{e}}^{{\theta }_2}-{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{\text{e}}^{\frac{{\beta }_2}{{\beta }_1}{\theta }_2}\right){\psi }_2$ 。

解得：${\psi }_2\left({\theta }_2\right)\propto \exp \left[\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right){\theta }_2-{\alpha }_1\left({\text{e}}^{{\theta }_2}-1\right)-\frac{{\beta }_1}{{\beta }_2}{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}\left({\text{e}}^{\frac{{\beta }_2}{{\beta }_1}{\theta }_2}-1\right)\right)\right]$ 。

从而$\psi \left({\theta }_1,{\theta }_2\right)=C{\psi }_1\left({\theta }_1\right){\psi }_2\left({\theta }_2\right)$ ，其中$C$ 是一个归一化常数。注意到$\underset{{\theta }_1\to +\infty }{\lim }\psi \left({\theta }_1,{\theta }_2\right)=0$ ，$\underset{{\theta }_1\to -\infty }{\lim }\psi \left({\theta }_1,{\theta }_2\right)=0$ ，$\underset{{\theta }_2\to +\infty }{\lim }\psi \left({\theta }_1,{\theta }_2\right)=0$ ，$\underset{{\theta }_2\to -\infty }{\lim }\psi \left({\theta }_1,{\theta }_2\right)=0$ 。

则(7)的边界条件是明确的。此外，注意到对于任意常数$A>0$ ，有：

$\exp \left|Az-{\text{e}}^z\right|\le \{\begin{array}{l}\exp \left[-z+A\ln \left(1+A\right)-\left(1+A\right)\right],z\ge 0\\ \exp \left[Az\right],z<0\end{array}$

那么${\displaystyle {\int }_{-\infty }^{+\infty }{\displaystyle {\int }_{-\infty }^{+\infty }\psi \left({\theta }_1,{\theta }_2\right)\text{d}{\theta }_1\text{d}{\theta }_2}}$ 是收敛的，并且对于某个常数$C$ 其值为1。根据Gray的唯一性定理[25]，概率密度$\psi \left(N,x\right)$ 是唯一的。由于$N={\text{e}}^{{\theta }_1}$ 和$x={\text{e}}^{{\theta }_2}$ ，则其雅可比行列式为$\frac{1}{Nx}$ ，则$\left(N,x\right)$ 的联合稳态概率函数具有如下形式：

$\begin{array}{c}\psi \left(N,x\right)=\frac{C}{Nx}{\mathrm{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left(N-1\right)}{N}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}}{\mathrm{e}}^{-\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}\left(x-1\right)+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}\left(x^{\frac{{\beta }_2}{{\beta }_1}}-1\right)\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)}\\ =C{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left(N-1\right)}{N}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}\left(x-1\right)+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}\left(x^{\frac{{\beta }_2}{{\beta }_1}}-1\right)\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)-1}\end{array}$ (8)

另外，通过公式$\psi \left(N\right)={\displaystyle {\int }_0^{\infty }\psi \left(N,x\right)\text{d}x}$ 和$\psi \left(x\right)={\displaystyle {\int }_0^{\infty }\psi \left(N,x\right)\text{d}N}$ ，我们可以得到$\psi \left(N\right)$ 和$\psi \left(x\right)$ 。对于$\psi \left(y\right)$ ，可以从$\psi \left(x\right)$ 和$F\left(N,x,y\right)=c$ 中推导出来，证明完成。

最后，我们证明了(5)的遍历性。

定理3.3 如果${\lambda }_1={\lambda }_2$ 且$r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}＞0$ ，那么(5)具有唯一的平稳分布，并且它具有遍历性质：

1) $\underset{t\to \infty }{\lim }EN\left(t\right)=\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}N\left(s\right)\text{d}s}={\displaystyle {\int }_0^{+\infty }N\psi \left(N\right)\text{d}N}={\lambda }_1$ ；

2) $\underset{t\to \infty }{\lim }Ex\left(t\right)=\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}x\left(s\right)\text{d}s}={\displaystyle {\int }_0^{+\infty }x\psi \left(x\right)\text{d}N}=x^{*}$ ；

3) $\underset{t\to \infty }{\lim }Ey\left(t\right)=\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}y\left(s\right)\text{d}s}={\displaystyle {\int }_0^{+\infty }y\psi \left(y\right)\text{d}N}=y^{*}$ ，

其中$x^{*}$ 和$y^{*}$ 满足${\alpha }_1{x}^{*}+{\alpha }_2{y}^{*}=r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}$ 。

证明：由于注1和(8)的唯一性保证了(6)具有唯一的平稳分布，因此根据定理3.2.6 [26]，(N(t), x(t))是遍历的。由此可知：

$\underset{t\to \infty }{\lim }EN\left(t\right)=\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}N\left(s\right)\text{d}s}={\displaystyle {\int }_0^{+\infty }N\psi \left(N\right)\text{d}N}$ ;

$\underset{t\to \infty }{\lim }Ex\left(t\right)=\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}x\left(s\right)\text{d}s}={\displaystyle {\int }_0^{+\infty }x\psi \left(x\right)\text{d}N}$ ;

$\underset{t\to \infty }{\lim }Ey\left(t\right)=\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}y\left(s\right)\text{d}s}={\displaystyle {\int }_0^{+\infty }y\psi \left(y\right)\text{d}N}$ .

根据定理3.2，我们可以得到：

$\begin{array}{c}\underset{t\to \infty }{\text{lim}}E\left(N\left(t\right)\right)={\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(u-1)}{u}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}}\text{d}u}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(u-1)}{u}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}}\text{d}u\right]}^{-1}\\ =\frac{{\sigma }^2}{2}\frac{K}{r}{\displaystyle {\int }_0^{+\infty }{\text{e}}^{-v}{v}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}}}\text{d}v\left[{\displaystyle {\int }_0^{+\infty }{\text{e}}^{-v}{u}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}}\text{d}u\right]\\ =\frac{{\sigma }^2}{2}K\frac{\Gamma \left(\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}+1\right)}{\Gamma \left(\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}\right)}={\lambda }_1\end{array}$

其中，$\Gamma \left(a\right)={\displaystyle {\int }_0^{\infty }{x}^{a-1}{\text{e}}^{-x}}\text{d}x$ 为伽马函数，满足$\Gamma \left(a+1\right)=a\Gamma \left(a\right)$ 。

$\underset{t\to \infty }{\lim }Ex\left(t\right)=C_0{\displaystyle {\int }_0^{\infty }x{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)-1}}\text{d}x\equiv x^{*}$ (9)

其中，$C_0={\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}u+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{u}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{u}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)-1}}\text{d}u\right]}^{-1}$ ，再通过$F\left(N,x,y\right)=c$ 我们可以得到：

$\begin{array}{c}\underset{t\to \infty }{\lim }Ey\left(t\right)=C_0{c}^{-\frac{1}{{\beta }_1}}\underset{t\to \infty }{\lim }E{x}^{\frac{{\beta }_2}{{\beta }_1}}\left(t\right)\\ =C_0{c}^{-\frac{1}{{\beta }_1}}{\displaystyle {\int }_0^{\infty }{x}^{\frac{{\beta }_2}{{\beta }_1}}{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)-1}}\text{d}x\equiv y^{*}\end{array}$ (10)

联立(9)、(10)则可得：

$\begin{array}{c}{\alpha }_1{x}^{*}+{\alpha }_2{y}^{*}=\underset{t\to \infty }{\lim }E\left({\alpha }_{1}x\left(t\right)+{\alpha }_{2}y\left(t\right)\right)=\underset{t\to \infty }{\lim }E\left({\alpha }_{1}x\left(t\right)+{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\left(t\right)\right)\\ =C_0{\displaystyle {\int }_0^{\infty }\left({\alpha }_{1}x+{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right){\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)-1}}\text{d}x\\ =C_0{\displaystyle {\int }_0^{\infty }\left({\alpha }_{1}x+{\alpha }_2{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}-1}\right){\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left[{\alpha }_{1}x+\frac{{\alpha }_2{\beta }_1}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)}}\text{d}x\end{array}$

$\begin{array}{c}=-C_0\frac{{\sigma }^2{\beta }_{1}K}{2}{\displaystyle {\int }_0^{\infty }{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)}}{\text{de}}^{-\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left[{\alpha }_{1}x+\frac{{\alpha }_2{\beta }_1}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}\\ =-C_0\frac{{\sigma }^2{\beta }_{1}K}{2}{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left[{\alpha }_{1}x+\frac{{\alpha }_2{\beta }_1}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)}|}_0^{+\infty }\\ +C_0\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right){\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left[{\alpha }_{1}x+\frac{{\alpha }_2{\beta }_1}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{x}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)-1}}\text{d}x\\ =r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}.\end{array}$

以上过程运用了分部积分法，证明完成。

引理3.4 对于任意初始值$\left(N_0,x_0,y_0\right)\in R_{+}^3$ ，方程(5)的解记为$\left(N\left(t\right),x\left(t\right),y\left(t\right)\right)$ ，对任意$p>1$ ，有：$\underset{t\to \infty }{\lim }\frac{\ln N\left(t\right)}{t^q}=0$ 和$\underset{t\to \infty }{\lim }\frac{N^p\left(t\right)+N\left(t\right)+x\left(t\right)+y\left(t\right)}{t^q}=0$ a.s.对于某个正常数$q\le 1$ 。

证明：在$R_{+}^3$ 上定义Lyapunov函数：

$V\left(N,x,y\right)=\frac{1}{p}{N}^p+N+\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}y-g\ln N$ ,

其中，$g\in \left(0,1\right)$ 是一个待定常数。在$V\left(N,x,y\right)$ 上应用伊藤公式可得：

$\text{d}V\left(N,x,y\right)=LV\left(N,x,y\right)\text{d}t+\sigma \left(N^p+N-g\right)\text{d}B\left(t\right)$ ,

其中

$\begin{array}{c}LV=N^p\left[\left(r-\frac{N}{K}\right)-{\alpha }_{1}x-{\alpha }_{2}y+\frac{p-1}{2}{\sigma }^2\right]+N\left(r-\frac{N}{K}\right)-{\alpha }_{1}Nx-{\alpha }_{2}Ny\\ +\frac{{\alpha }_1}{{\beta }_1}\left({\beta }_{1}N-{\mu }_1\right)x+\frac{{\alpha }_2}{{\beta }_2}\left({\beta }_{2}N-{\mu }_2\right)y+g\left(\frac{N}{K}+{\alpha }_{1}x+{\alpha }_{2}y\right)-g\left(r-\frac{{\sigma }^2}{2}\right)\\ \le N^p\left[\left(r-\frac{N}{K}\right)+\frac{p-1}{2}{\sigma }^2\right]+\left(r+\frac{g}{K}\right)N-\frac{{\alpha }_1{\mu }_1}{{\beta }_1}x-\frac{{\alpha }_2{\mu }_2}{{\beta }_2}y+g\left({\alpha }_{1}x+{\alpha }_{2}y\right)-g\left(r-\frac{{\sigma }^2}{2}\right)\end{array}$

接下来由于对任意$N>0$ 都有$\ln N\le N-1$ ，通过选择$g$ 使得$g\left({\alpha }_1+{\alpha }_2\right)\le \frac{c}{2}\min \left\{\frac{{\alpha }_1}{{\beta }_1},\frac{{\alpha }_2}{{\beta }_2}\right\}$ ，我们得到：

$LV\le -\frac{c}{2}\left(\frac{1}{p}{N}^p+N+\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}y-g\ln N\right)+Q\left(N\right)$ ,

其中，$Q\left(N\right)$ 满足以下约束条件：

$\begin{array}{c}Q\left(N\right)=N^p\left[\left(r-\frac{N}{K}\right)+\frac{p-1}{2}{\sigma }^2\right]+\left(r+\frac{g}{K}\right)N+\frac{c}{2}\left(\frac{1}{p}{N}^p+N\right)+\frac{cg}{2}\ln N-g\left(r-\frac{{\sigma }^2}{2}\right)\\ \le \underset{x＞0}{\sup }\left\{-\frac{x^{p+1}}{K}+\left(r+\frac{c}{2p}+\frac{p-1}{2}{\sigma }^2\right)x^p+\left(r+\frac{g}{K}+\frac{c\left(1+g\right)}{2}\right)x-g\left(r-\frac{{\sigma }^2}{2}\right)\right\}\equiv Q^{*}.\end{array}$

对于某个常数$\gamma >1$ 使得$\frac{\gamma -1}{2}{\left(\sigma p\right)}^2\le \frac{c}{4}$ ，通过定义$V_1\left(N,x,y\right)=V^{\gamma }\left(N,x,y\right)$ ，则有：

$\begin{array}{c}L{V}_1\le -\frac{c}{2}\gamma V^{\gamma }+Q^{*}\gamma V^{\gamma -1}+\frac{\gamma \left(\gamma -1\right)}{2}{\sigma }^2{V}^{\gamma -2}{\left(N^p+N-g\right)}^2\\ \le -\frac{c}{2}\gamma V^{\gamma }+Q^{*}\gamma V^{\gamma -1}+\gamma \left(\gamma -1\right){\sigma }^2{V}^{\gamma -2}{\left(N^p+N\right)}^2+\gamma \left(\gamma -1\right){\sigma }^2{g}^2{V}^{\gamma -2}\left(N,x,y\right)\\ \le -\frac{c}{2}\gamma V^{\gamma }+Q^{*}\gamma V^{\gamma -1}+\gamma \left(\gamma -1\right){\sigma }^2{V}^{\gamma -2}{\left(N^p+N+\left(1+\frac{1}{p}\right)\left(N-g\ln N\right)\right)}^2+\gamma \left(\gamma -1\right){\sigma }^2{g}^2{V}^{\gamma -2}\end{array}$

由于对任意$N>0$ 都有$N\le \frac{1}{p}{N}^p+\frac{p-1}{p}$ ，通过选择$\gamma$ 使得$\left(\gamma -1\right){\sigma }^2{\left(1+\frac{1}{p}\right)}^2\le \frac{c}{4}$ ，则有：

$\begin{array}{c}L{V}_1\le -\frac{c}{2}\gamma V^{\gamma }+Q^{*}\gamma V^{\gamma -1}+\gamma \left(\gamma -1\right){\sigma }^2{\left(1+\frac{1}{p}\right)}^2{V}^{\gamma -2}\left(N,x,y\right){\left(N^p+N-g\ln N\right)}^2+\gamma \left(\gamma -1\right){\sigma }^2\left(g^2+\frac{p-1}{p}\right)V^{\gamma -2}\\ \le -\frac{c}{4}\gamma V^{\gamma }+Q^{*}\gamma V^{\gamma -1}+\gamma \left(\gamma -1\right){\sigma }^2\left(g^2+\frac{p-1}{p}\right)V^{\gamma -2}.\end{array}$

由于$N-g-g\ln N\ge g\left(N-1-\ln N\right)\ge 0$ ，那么有：

$V\left(N,x,y\right)=\frac{1}{p}{N}^p+N+\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_2}y-g\ln N\ge g$ ,

因此：

$L{V}_1\le -\frac{c}{8}\gamma V^{\gamma }\left(N,x,y\right)+Q^{**}$ ,

其中$Q^{**}=\underset{v＞g}{\sup }\left\{-\frac{c}{8}\gamma v^{\gamma }+Q^{*}\gamma v^{\gamma -1}+\gamma \left(\gamma -1\right){\sigma }^2\left(g^2+\frac{p-1}{p}\right)V^{\gamma -2}\right\}<\infty$ 。通过重新定义$V_1\left(N,x,y,t\right)=V^{\gamma }\left(N,x,y\right){\text{e}}^{\frac{c\gamma }{8}t}$ 并对两边取期望，我们得到有界性性质，即：

$\begin{array}{c}E\left(V^{\gamma }\left(N\left(t\right),x\left(t\right),y\left(t\right)\right)\right)\le \left(V^{\gamma }\left(N\left(0\right),x\left(0\right),y\left(0\right)\right)\right){\text{e}}^{-\frac{c\gamma }{8}t}+Q^{**}{\displaystyle {\int }_0^t{\text{e}}^{-\frac{c\gamma }{8}(t-s)}}\text{d}s\\ \le V^{\gamma }\left(N\left(0\right),x\left(0\right),y\left(0\right)\right)+\frac{8}{c\gamma }{Q}^{**}.\end{array}$

注意到$\gamma >1$ ，通过使用经典Burkholder-Davis-Gundy不等式和Borel-Cantelli引理，并遵循[27]中的引理2.1证明的方案[28]，我们可以证明：

$\underset{t\to \infty }{\lim }\frac{V^{\sqrt{\gamma }}\left(N\left(t\right),x\left(t\right),y\left(t\right)\right)}{t}=0$ .

令$q=\frac{1}{\sqrt{\gamma }}$ ，则$\underset{t\to \infty }{\lim }\frac{\frac{1}{p}{N}^p\left(t\right)+N\left(t\right)+\frac{{\alpha }_1}{{\beta }_1}x\left(t\right)+\frac{{\alpha }_2}{{\beta }_2}y\left(t\right)-g\ln N\left(t\right)}{t^q}=0$ a.s.，从而得到所要的结论。

引理3.5 对于由模型(5)描述的两个具有正初始值的捕食者，在${\lambda }_1\ne {\lambda }_2$ 条件下，至少有一个捕食者会灭绝。具体有：

如果${\lambda }_1<{\lambda }_2$ ，则当$t$ 趋于无穷时，$y\left(t\right)$ 呈指数级趋于零。

如果${\lambda }_1>{\lambda }_2$ ，则当$t$ 趋于无穷时，$x\left(t\right)$ 呈指数级趋于零。

证明：根据系统(5)中的前两个方程，我们可以得到$\text{d}\left(\frac{\ln x}{{\beta }_1}-\frac{\ln y}{{\beta }_2}\right)=-\left({\lambda }_1-{\lambda }_2\right)\text{d}t$ ，其解的形式为：

$\frac{1}{{\beta }_2}\ln y\left(t\right)=\frac{1}{{\beta }_1}\ln x\left(t\right)+\ln x_0^{-\frac{1}{{\beta }_1}}{y}_0^{\frac{1}{{\beta }_2}}+\left({\lambda }_1-{\lambda }_2\right)t$ (11)

如果${\lambda }_1<{\lambda }_2$ ，那么根据引理3.4可得：$\underset{t\to \infty }{\lim \sup }\frac{1}{t}\ln y\left(t\right)\le {\beta }_2\left({\lambda }_1-{\lambda }_2\right)$ ，因此(1)成立。(2)的证明类似。

从定理3.3和引理3.5中可以很容易推导出以下结果：

定理3.6 令$\left(N\left(t\right),x\left(t\right),y\left(t\right)\right)$ 是方程(5)在$r-\frac{{\sigma }^2}{2}-\frac{\min \left\{{\lambda }_1,{\lambda }_2\right\}}{K}>0$ 条件下的解，则竞争排斥原理成立，即：

1) 如果${\lambda }_1<{\lambda }_2$ ，则有$\underset{t\to \infty }{\lim }y\left(t\right)=0$ 和$\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}x\left(s\right)\text{d}s}=x^{*}:=\frac{1}{{\alpha }_1}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)$ ；

2) 如果${\lambda }_1>{\lambda }_2$ ，则有$\underset{t\to \infty }{\lim }x\left(t\right)=0$ 和$\underset{t\to \infty }{\lim }\frac{1}{t}{\displaystyle {\int }_0^{t}y\left(s\right)\text{d}s}=y^{*}:=\frac{1}{{\alpha }_2}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_2}{K}\right)$ 。

定理3.7 设$\left(N\left(t\right),x\left(t\right),y\left(t\right)\right)$ 是方程(5)在$\left[0,\infty \right)$ 上的解，且初始值$\left(N_0,x_0,y_0\right)\in R_{+}^3$ ，则两个捕食者共存当且仅当${\lambda }_1={\lambda }_2$ 且$r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}>0$ ，此外$x^{*}>0$ ，$y^{*}>0$ 。

4. 随机分岔分析

注意到在满足$r-\frac{\lambda }{K}>0$ 的确定性模型中没有出现分岔行为，所以接下来我们主要讨论噪声强度对(5)式的影响。关于随机分岔，我们遵循[29]-[31]中的定义，其中P-分岔与随机动力系统中参数变化时稳态概率密度的变化有关，这是一种静态行为。

当${\lambda }_1\ne {\lambda }_2$ 时，其中一种捕食者将灭绝，系统简化为二维系统：

$\{\begin{array}{l}\text{d}N=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_{1}Nx\right]\text{d}t+\sigma N\text{d}B\left(t\right),\\ \text{d}x=\left({\beta }_{1}N-{\mu }_1\right)x\text{d}t.\end{array}$ (12)

$\{\begin{array}{l}\text{d}N=\left[N\left(r-\frac{N}{K}\right)-{\alpha }_{2}Ny\right]\text{d}t+\sigma N\text{d}B\left(t\right),\\ \text{d}y=\left({\beta }_{2}N-{\mu }_2\right)y\text{d}t.\end{array}$ (13)

根据定理3.2和定理3.3的结果，我们得出以下推论。

推论4.1 对于模型(12)，如果$r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}>0$ ，则其边际平稳分布为：

$\psi \left(N\right)={\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(N-1)}{N}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(u-1)}{u}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_1}{K}-1}}\text{d}u\right]}^{-1}$

$\psi \left(x\right)={\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\frac{{\alpha }_1}{{\beta }_1}x}{x}^{\frac{1}{{\beta }_1}\left[\frac{2}{{\sigma }^2}\frac{1}{K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)\right]-1}{\left\{{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\frac{{\alpha }_1}{{\beta }_1}u}{u}^{\frac{1}{{\beta }_1}\left[\frac{2}{{\sigma }^2}\frac{1}{K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)\right]-1}}\text{d}u\right\}}^{-1}$ .

此外，有$\underset{t\to \infty }{\lim }E\left(N\left(t\right)\right)={\lambda }_1$ 以及$\underset{t\to \infty }{\lim }E\left(x\left(t\right)\right)=\frac{1}{{\alpha }_1}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_1}{K}\right)$ 。

推论4.2 对于模型(12)，如果$r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_2}{K}>0$ ，则其边际平稳分布为：

$\psi \left(N\right)={\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(N-1)}{N}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_2}{K}-1}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(u-1)}{u}^{\frac{2}{{\sigma }^2}\frac{{\lambda }_2}{K}-1}}\text{d}u\right]}^{-1}$

$\psi \left(y\right)={\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\frac{{\alpha }_2}{{\beta }_2}y}{y}^{\frac{1}{{\beta }_2}\left[\frac{2}{{\sigma }^2}\frac{1}{K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_2\right)\right]-1}{\left\{{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\frac{{\alpha }_2}{{\beta }_2}u}{u}^{\frac{1}{{\beta }_2}\left[\frac{2}{{\sigma }^2}\frac{1}{K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_2\right)\right]-1}}\text{d}u\right\}}^{-1}$ .

此外，有$\underset{t\to \infty }{\lim }E\left(N\left(t\right)\right)={\lambda }_2$ 以及$\underset{t\to \infty }{\lim }E\left(y\left(t\right)\right)=\frac{1}{{\alpha }_2}\left(r-\frac{{\sigma }^2}{2}-\frac{{\lambda }_2}{K}\right)$ 。

由于$\left(N,x\right)$ 和$\left(N,y\right)$ 的平稳分布在$r-\frac{{\sigma }^2}{2}-\frac{\min \left\{{\lambda }_1,{\lambda }_2\right\}}{K}>0$ 下是不相关的，因此我们可以根据它们的边缘平稳分布$\psi \left(N\right)$ 、$\psi \left(x\right)$ 和$\psi \left(y\right)$ 分别研究$N$ 、$x$ 和$y$ 的分岔行为。

定义4.3 (定义2，[31]) 随机现象学分岔(随机P-分岔)是关于参数变化时，一个随机动力系统族的密度(平稳概率密度)形状的变化。如果存在一个常数${\alpha }_0$ 满足在${\alpha }_0$ 的任意邻域内，存在另外两个常数${\alpha }_1$ 和${\alpha }_2$ ，使得${\alpha }_1<{\alpha }_0<{\alpha }_2$ ，且它们对应的不变测度$P_{{\alpha }_1}$ 和$P_{{\alpha }_2}$ 是不相同的，则常数${\alpha }_0$ 是随机P-分岔的一个点。

4.1. 被捕食者中噪声引起的分岔

通过定义$\lambda =\min \left\{{\lambda }_1,{\lambda }_2\right\}$ ，定理3.2和推论4.1、4.2表明$\psi \left(N\right)$ 具有如下形式：

$\psi \left(N\right)={\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(N-1)}{N}^{\frac{2}{{\sigma }^2}\frac{\lambda }{K}-1}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}(u-1)}{u}^{\frac{2}{{\sigma }^2}\frac{\lambda }{K}-1}}\text{d}u\right]}^{-1}$ .

定理4.4 对于模型(5)，如果$r-\frac{{\sigma }^2}{2}-\frac{\lambda }{K}>0$ 且$r-\frac{2\lambda }{K}>0$ ，则在由${\sigma }_C=\sqrt{\frac{2\lambda }{K}}$ 定义的临界点处发生随机P-分岔。

证明：为方便起见，我们将其表示为：

$\psi \left(N\right)={\text{e}}^{\frac{2}{{\sigma }^2}\frac{1}{K}h(N)}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}u}{u}^{\frac{2}{{\sigma }^2}\frac{\lambda }{K}-1}}\text{d}u\right]}^{-1}.$ 其中$h\left(N\right)=-N+\left(\lambda -\frac{{\sigma }^{2}K}{2}\right)\ln N$ .

对N求导，我们可以得到：

$h^{\prime }\left(N\right)=-1+\left(\lambda -\frac{{\sigma }^{2}K}{2}\right)\frac{1}{N}=-\frac{1}{N}\left(N-\left(\lambda -\frac{{\sigma }^{2}K}{2}\right)\right)$ .

如果$\sigma >{\sigma }_C$ ，则有$h^{\prime }\left(N\right)<0$ 对所有$N>0$ 成立。由此可知，$h\left(N\right)$ 单调递减且有：

$\underset{N\to 0}{\lim }h\left(N\right)=\underset{N\to 0}{\lim }\left(-N+\left(\lambda -\frac{{\sigma }^{2}K}{2}\right)\ln N\right)=+\infty$ .

如果$\sigma ={\sigma }_C$ ，则$h\left(N\right)=-N$ ，并且在$N=0$ 时取得最大值0。

如果$\sigma <{\sigma }_C$ ，那么则有$\lambda -\frac{{\sigma }^{2}K}{2}>0$ ，从而$h^{\prime }\left(N\right)=0$ 有唯一解$N_0=\lambda -\frac{{\sigma }^{2}K}{2}$ ，在这种情况下，$h\left(N\right)$ 在$N_0$ 处达到最大值。

显然，当$\sigma$ 通过临界值${\sigma }_C$ 时，稳态概率密度函数的形状会发生突变，即$\sigma ={\sigma }_C$ 是模型(5)的一个随机P-分岔点。

注2 显然$\psi \left(N\right)$ 不包含参数$r$ ，但是其对$\psi \left(N\right)$ 的形状有显著影响。如果$r-\frac{2\lambda }{K}\le 0$ 并且$r-\frac{{\sigma }^2}{2}-\frac{\lambda }{K}>0$ ，则$\frac{{\sigma }^2}{2}-\frac{\lambda }{K}<r-\frac{2\lambda }{K}\le 0$ ，此时有$\lambda -\frac{{\sigma }^{2}K}{2}>0$ ，这意味着$h^{\prime }\left(N\right)=0$ 拥有一个唯一的正解，记为$\overline{N}$ ，此时$h\left(N\right)$ 在$\overline{N}$ 处达到峰值。从另一个角度来看，在这种情况下，$\sigma$ 不能诱导分岔。

例1 将模型(5)的参数设置如下：$r=0.5$ ，$K=10$ ，${\alpha }_1=0.2$ ，${\alpha }_2=0.3$ ，${\beta }_1=0.05$ ，${\beta }_2=0.02$ ，${\mu }_1=0.09$ ，${\mu }_2=0.1$ ，$\lambda =1.8$ ，${\sigma }_C=0.6$ 。为了验证噪声对被捕食者分岔的影响，我们分别选择$\sigma =0.55$ 、$\sigma =0.6$ 和$\sigma =0.65$ 。采用欧拉丸山法(Euler-Maruyama)对随机微分方程(5)进行离散化，步长固定为0.01。当$\sigma$ 超过临界值${\sigma }_C=0.6$ 时，稳态概率密度$\psi \left(N\right)$ 发生急剧变化，即在${\sigma }_C=0.6$ 处发生P-分岔，这证实了理论结果(见图2)。

Figure 2. Numerical simulation at $\sigma =0.55$ , $\sigma =0.6$ , and $\sigma =0.65$

图2. $\sigma =0.55$ 、$\sigma =0.6$ 、$\sigma =0.65$ 时的数值模拟

4.2. 捕食者中噪声引起的分岔

定理4.5 对于模型(5)，在${\lambda }_1={\lambda }_2$ 和$r-\frac{{\sigma }^2}{2}-\frac{\lambda }{K}>0$ 条件下有：

1) 如果$\sigma$ 的值自下而上经过${\sigma }_C^x=\sqrt{\frac{2}{{\beta }_{1}K+1}\left(r-\frac{\lambda }{K}\right)}$ ，那么捕食者$x$ 将发生随机P-分岔；

2) 如果$\sigma$ 的值自下而上经过${\sigma }_C^y=\sqrt{\frac{2}{{\beta }_{2}K+1}\left(r-\frac{\lambda }{K}\right)}$ ，那么捕食者$y$ 将发生随机P-分岔。

证明：根据定理3.2，我们有

$\psi \left(x\right)={\text{e}}^{\frac{2}{{\sigma }^2}\frac{1}{K}h(x)}{\left[{\displaystyle {\int }_0^{\infty }{\text{e}}^{-\frac{2}{{\sigma }^2}\frac{1}{K}\left[\frac{{\alpha }_1}{{\beta }_1}u+\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{u}^{\frac{{\beta }_2}{{\beta }_1}}\right]}{u}^{\frac{2}{{\sigma }^2}\frac{1}{{\beta }_{1}K}\left(r-\frac{{\sigma }^2}{2}-\frac{1}{K}{\lambda }_1\right)-1}}\text{d}u\right]}^{-1}$

其中，$h\left(x\right)=-\frac{{\alpha }_1}{{\beta }_1}x-\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}+B\ln x$ ，$B=\frac{1}{{\beta }_1}\left(r-\frac{{\sigma }^2}{2}-\frac{\lambda }{K}\right)-\frac{{\sigma }^{2}K}{2}$ 。

如果$\sigma >{\sigma }_C^x$ 成立，那么$B<0$ ，因此对于所有$x>0$ 有$h^{\prime }\left(x\right)<0$ ，故$h\left(x\right)$ 单调递减。进一步有：

$\underset{x\to 0^{+}}{\lim }h\left(x\right)=\underset{x\to 0^{+}}{\lim }\left(-\frac{{\alpha }_1}{{\beta }_1}x-\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}+B\ln x\right)=+\infty$

如果$\sigma ={\sigma }_C^x$ 成立，那么$h\left(x\right)=-\frac{{\alpha }_1}{{\beta }_1}x-\frac{{\alpha }_2}{{\beta }_2}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}$ 且在$x=0$ 时取得最大值0。

如果$\sigma <{\sigma }_C^x$ 成立，那么$B>0$ ，由于$h^{\prime }\left(x\right)=-\frac{{\alpha }_1}{{\beta }_1}-\frac{{\alpha }_2}{{\beta }_1}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}-1}+\frac{B}{x}$ ，令$g\left(x\right)=-x{h}^{\prime }\left(x\right)$ ，则可得$g\left(x\right)=\frac{{\alpha }_1}{{\beta }_1}x+\frac{{\alpha }_2}{{\beta }_1}{c}^{-\frac{1}{{\beta }_1}}{x}^{\frac{{\beta }_2}{{\beta }_1}}-B$ 是单调函数，因此$h^{\prime }\left(x\right)=0$ 有唯一解，记为$\overline{x}$ 。此时$h\left(x\right)$ 在$\overline{x}$ 处取得最大值。

由于当$\sigma$ 穿过${\sigma }_C^x$ 时，稳态概率密度函数会发生突变，因此$\sigma ={\sigma }_C^x$ 是捕食者$x$ 的随机P-分岔点，同理可证捕食者$y$ 的情况。

由此定理可知，对于共存的捕食者$x$ 和$y$ ，分岔点是不同的，在证明的过程中可以得知这种差异主要取决于捕食者在存在被捕食者时的效率和传播率。当传播率相同时，二者的分岔点相同，从而可得到以下结论。

推论4.6 若${\lambda }_1<{\lambda }_2$ 且$r-\frac{{\sigma }^2}{2}-\frac{\lambda }{K}>0$ ，则在系统(12)中将在临界值${\sigma }_C^x=\sqrt{\frac{2}{{\beta }_{1}K+1}\left(r-\frac{\lambda }{K}\right)}$ 处发生随机P-分岔。

推论4.7 若${\lambda }_1>{\lambda }_2$ 且$r-\frac{{\sigma }^2}{2}-\frac{\lambda }{K}>0$ ，则在系统(13)中将在临界值${\sigma }_C^y=\sqrt{\frac{2}{{\beta }_{2}K+1}\left(r-\frac{\lambda }{K}\right)}$ 处发生随机P-分岔。

例2 在模型(5)中，我们选择$c=1$ ，$r=1.26$ ，$K=10$ ，${\alpha }_1=0.2$ ，${\alpha }_2=0.3$ ，${\beta }_1=0.05$ ，${\beta }_2=0.02$ ，${\mu }_1=0.09$ ，${\mu }_2=0.36$ ，${\lambda }_1={\lambda }_2=1.8$ ，${\sigma }_C^x=1.2$ ，${\sigma }_C^y=1.342$ 。为了验证噪声对捕食者分岔的影响，我们分别将$\psi \left(x\right)$ 中的$\sigma$ 设置为1.1、1.2和1.25，将$\psi \left(y\right)$ 中的$\sigma$ 设置为1.3、1.342和1.4。采用欧拉丸山法(Euler-Maruyama)对随机微分方程(5)进行离散化，步长固定为0.01，验证了定理4.5给出的分岔结果(见图3、图4)。

Figure 3. Numerical simulations of $\psi \left(x\right)$ at $\sigma =1.1$ , $\sigma =1.2$ and $\sigma =1.25$ were used, respectively

图3. 分别使用$\sigma =1.1$ 、$\sigma =1.2$ 、$\sigma =1.25$ 时的$\psi \left(x\right)$ 数值模拟

Figure 4. Numerical simulations of $\psi \left(x\right)$ at $\sigma =1.3$ , $\sigma =1.342$ and $\sigma =1.4$ were used, respectively

图4. 分别使用$\sigma =1.3$ 、$\sigma =1.342$ 、$\sigma =1.4$ 时的$\psi \left(x\right)$ 数值模拟

注3 本文的分岔分析主要基于噪声强度，称之为噪声诱发分岔。从所得结果可以看出，每一种分岔行为都是由一个唯一的条件引起的。这意味着我们可以将模型中的另一个参数视为分岔参数来分析系统的动力学，从而可以直接获得类似的分岔结果。

4.3. 被捕食者灭绝概率分析

在确定性模型(1)中，若满足共存条件，被捕食者和捕食者可长期共存。然而在随机扰动下，被捕食者可能有概率灭绝。本节将基于模型(5)推导被捕食者的灭绝概率，并分析噪声强度对其的影响。

定理4.8 (被捕食者灭绝概率) 假设系统(5)满足持久性条件(即$\sigma <{\sigma }_C$ ，其中${\sigma }_C$ 为噪声强度阈值)，则被捕食者$N\left(t\right)$ 的灭绝概率可表示为：

$P_{ext}={\displaystyle {\int }_0^{+\infty }{\displaystyle {\int }_0^{+\infty }{p}^{*}\left(0,x,y\right)\text{d}x\text{d}y}}$

其中，$P_{ext}=p^{*}\left(N,x,y\right)$ 为系统(5)的联合稳态概率密度函数。进一步地：

1) 当噪声强度$\sigma \ge {\sigma }_C$ 时，$P_{ext}\to 1$ ，即捕食者以概率1灭绝；

2) 当$\sigma <{\sigma }_C$ 时，$P_{ext}=0$ ，系统具有持久性。

对于捕食者$N\left(t\right)$ ，其动态方程为：

$\text{d}N=N\left(r\left(1-\frac{N}{K}\right)-{\alpha }_{1}x-{\alpha }_{2}y\right)\text{d}t+\sigma N\text{d}B\left(t\right)$

当$\sigma \ge {\sigma }_C$ 时，漂移项主导方向为负，扩散项加剧$N\left(t\right)$ 的波动，导致其以概率1灭绝。而当$\sigma <{\sigma }_C$ 时，系统持久性保证$N\left(t\right)$ 远离边界，故$P_{ext}=0$ 。

4.4. 其他参数作为分岔参数的扩展分析

尽管噪声强度是本文的核心分岔参数，但在实际生态系统中，捕食率、捕食者死亡率以及环境容纳量等参数的变化同样可能诱导系统的分岔行为。比如，捕食率反映了捕食者对被捕食者的捕食效率，其变化直接影响捕食者能量获取能力，进而改变种群的竞争平衡；捕食者死亡率的增加会削弱其种群维持能力，当其超过临界值时，捕食者可能灭绝，从而改变被捕食者的动态。本节将其作为概念上的拓展，不作展开描述。

5. 总结

我们的工作是研究了一个随机捕食者–被捕食者模型，其中两个捕食者竞争一个被捕食者。引入噪声的做法对系统动力学有显著影响，当噪声强度增大时，捕食者难以生存，导致系统内所有物种灭绝，而这种现象在确定性系统中并不存在。在假设被捕食者持续存在的前提下，通过降低噪声强度来分析捕食者的生存水平。在这种情况下，资源需求最低的捕食者将持续存在并排除需要更高资源需求的其他捕食者，所有资源需求最低的捕食者将在系统中共存。

正如Huang等人[31]指出的，平稳概率密度的极值包含了随机系统最重要的本质。如果平稳概率密度具有单峰结构，则峰值在概率意义上是稳定的，种群大小最终将以非常高的概率分布在该点附近。为了在模型(5)上证明这点，我们从稳态下的相应Fokker-Plank方程推导出封闭形式的稳态概率密度。结果表明，噪声可以诱导所有物种的分岔。当噪声较小时，稳态概率密度具有单调结构，在0处取无穷大值。当噪声强度增加并超过阈值时，稳态概率密度变为单峰结构。稳态概率密度的定性变化表明了P-分岔的发生，这种变化有助于我们从统计学角度更好地理解噪声对物种的动态影响。