噪聲和捕撈對(duì)捕食生態(tài)系統(tǒng)穩(wěn)定性的影響*

董慶國(guó) 寧麗娟

（陜西師范大學(xué)，數(shù)學(xué)與信息科學(xué)學(xué)院，西安 710062）

噪聲和捕撈對(duì)捕食生態(tài)系統(tǒng)穩(wěn)定性的影響*

董慶國(guó) 寧麗娟?

（陜西師范大學(xué)，數(shù)學(xué)與信息科學(xué)學(xué)院，西安 710062）

討論了在色噪聲激勵(lì)下，具有獨(dú)立常數(shù)率捕撈和庇護(hù)所效應(yīng)的捕食生態(tài)系統(tǒng)的穩(wěn)定性問題.在弱擾動(dòng)假設(shè)下應(yīng)用Stratonovich-Khasminskii隨機(jī)平均原理分別得到了兩個(gè)物種的穩(wěn)態(tài)概率密度，并研究了捕撈強(qiáng)度E1，色噪聲強(qiáng)度Kii，譜寬和噪聲相關(guān)時(shí)間對(duì)兩個(gè)物種的穩(wěn)態(tài)概率密度的影響.Monte-Carlo模擬驗(yàn)證理論求解的合理性.研究表明：1）隨著捕撈活動(dòng)的增大，隨機(jī)因素對(duì)生態(tài)系統(tǒng)的影響逐漸減弱；2）噪聲強(qiáng)度越大，生態(tài)系統(tǒng)越不穩(wěn)定；3）隨機(jī)激勵(lì)的譜帶越寬，生態(tài)系統(tǒng)越穩(wěn)定；4）隨機(jī)激勵(lì)的相關(guān)時(shí)間越小，生態(tài)系統(tǒng)越穩(wěn)定.

色噪聲，常數(shù)率捕撈，隨機(jī)平均方法，穩(wěn)態(tài)概率密度

引言

經(jīng)典的Lotka-Volterra模型奠定了利用數(shù)學(xué)建模方法來(lái)分析種間關(guān)系的理論基礎(chǔ)［1，2］.學(xué)者們對(duì)此模型進(jìn)行了研究并提出了許多改進(jìn)［3-7］，例如，考慮了種群間的自競(jìng)爭(zhēng)，被捕食者依賴型Holling I-III功能反應(yīng)［5］，比率依賴型功能反應(yīng)［6］，庇護(hù)所效應(yīng)［7］等.以上模型可以描述不同類型的確定性被捕食—捕食生態(tài)系統(tǒng)模型.一般地，由于環(huán)境中存在眾多的不確定因素，隨機(jī)模型比確定性模型更能準(zhǔn)確表達(dá)生態(tài)系統(tǒng)的狀態(tài).可以對(duì)隨機(jī)生態(tài)系統(tǒng)模型分析的途徑和方法有信息熵［8］，廣義圖胞映射法［9］，本征值和本征矢理論（對(duì)單種群Logistic模型）［10］，路徑積分法等.2004年，Cai和Lin利用隨機(jī)平均法研究了在白噪聲擾動(dòng)下的捕食生態(tài)模型，得到了系統(tǒng)中兩物種的穩(wěn)態(tài)概率密度［11］.2013年，Qi等進(jìn)一步建立了受色噪聲擾動(dòng)的含食物過剩和捕食者飽和的統(tǒng)一捕食類型的生態(tài)系統(tǒng)模型，并利用隨機(jī)平均法得到了該統(tǒng)一模型中兩物種的穩(wěn)態(tài)概率密度［12］.

以上研究的模型中忽略了人類活動(dòng)對(duì)模型的影響.事實(shí)上，研究人類活動(dòng)對(duì)生態(tài)系統(tǒng)穩(wěn)定性的影響具有重要的生態(tài)和經(jīng)濟(jì)價(jià)值.2006年，Kar研究了具有庇護(hù)所效應(yīng)的食餌—捕食者捕撈模型，發(fā)現(xiàn)捕撈活動(dòng)不僅可以破壞系統(tǒng)的周期行為，還可以控制系統(tǒng)的最終性態(tài)［13］.2010年，Lv等建立了在外界捕獲下水生浮游植物和浮游動(dòng)物構(gòu)成的被捕食—捕食模型，并對(duì)其穩(wěn)定性和最優(yōu)捕獲策略進(jìn)行了研究，發(fā)現(xiàn)外界捕獲可以導(dǎo)致種群滅絕，合理的捕獲可以維持種群數(shù)量的穩(wěn)定［14］.然而，在人類捕撈活動(dòng)下的生態(tài)系統(tǒng)的隨機(jī)模型的研究甚少.

綜上所述，本文對(duì)色噪聲激勵(lì)下，具有獨(dú)立常數(shù)率捕撈和庇護(hù)所效應(yīng)的隨機(jī)生態(tài)系統(tǒng)模型進(jìn)行了研究.為了更好的描述捕撈活動(dòng)對(duì)生態(tài)系統(tǒng)的影響，討論了被捕食者的捕撈強(qiáng)度E1變化時(shí)，整個(gè)模型的性態(tài)變化.首先，研究了確定性生態(tài)系統(tǒng)的穩(wěn)定性態(tài)，并得到了捕撈強(qiáng)度E1（0＜E1＜z/q1）的兩個(gè)閾值e1和e2.其次，將隨機(jī)模型轉(zhuǎn)化為It?隨機(jī)微分方程，并在e1＜E1＜e2區(qū)間內(nèi)應(yīng)用隨機(jī)平均法，從理論上推導(dǎo)了隨機(jī)系統(tǒng)中兩物種的穩(wěn)態(tài)概率密度.然后，利用低通濾過過程來(lái)模擬環(huán)境中的噪聲，并對(duì)真實(shí)環(huán)境中的噪聲譜密度進(jìn)行分析.最后，通過數(shù)值模擬，具體分析了人類捕撈和隨機(jī)環(huán)境對(duì)兩物種的穩(wěn)態(tài)概率密度的影響，并利用Monte-Carlo模擬來(lái)驗(yàn)證理論分析的準(zhǔn)確性.

1　確定性模型研究

1.1經(jīng)典Lotka-Volterra模型

描述物種間相互作用的Lotka-Volterra模型［1-2］如下：

其中，x，y分別代表被捕食者和捕食者的種群數(shù)量；a，b，c和f為正常數(shù)：a代表被捕食者自然生長(zhǎng)率（出生率），c代表捕食者自然死亡率，b和f分別代表捕食者捕食率和食物轉(zhuǎn)化效率.（1）式是一保守系統(tǒng)，有一個(gè)不穩(wěn)定平衡點(diǎn)（0，0），一個(gè)穩(wěn)定的非漸近平衡點(diǎn)（c/f，a/b）.具有如下形式的首次積分［15］

其中γ（x，y）=0對(duì)應(yīng)正平衡點(diǎn)（c/f，a/b），且對(duì)于γ（x，y）=R（R≥0），正數(shù)x＞0，y＞0，即γ（x，y）=R描述了動(dòng)力系統(tǒng)中兩物種的周期軌道.相應(yīng)軌道的周期定義為

其環(huán)路積分定義在由（2）式確定的周期軌道上.

1.2具有獨(dú)立常數(shù)率捕撈的生態(tài)系統(tǒng)模型

包含人類對(duì)被捕食者和捕食者獨(dú)立常數(shù)率捕撈，被捕食者具有庇護(hù)所的捕食生態(tài)系統(tǒng)的確定性模型可以表示為：

由動(dòng)力系統(tǒng)局部穩(wěn)定性判定條件，得到關(guān)于對(duì)被捕食者捕撈強(qiáng)度E1的兩個(gè)閥值，記為e1，e2，其值為

根據(jù)極限環(huán)存在的條件［17］，當(dāng)0＜E1＜e1時(shí)，系統(tǒng)（4）存在一個(gè)極限環(huán)；當(dāng)e1＜E1＜e2時(shí)，P2是漸近穩(wěn)定點(diǎn)；當(dāng)e2＜E1＜z/q1時(shí)，P1是漸近穩(wěn)定點(diǎn).

表1　正平衡點(diǎn)P2（）的部分值（e1＜E1＜e2）Table 1 The partially value of P2（x～，y～）（e1＜E1＜e2）

表1　正平衡點(diǎn)P2（）的部分值（e1＜E1＜e2）Table 1 The partially value of P2（x～，y～）（e1＜E1＜e2）

E12468 P2（143.3，57.80）（143.3，45.15）（143.3，32.50）（143.3，19.84）

圖1　不同E1取值下系統(tǒng)（4）的軌道：（a）E1=4；（b）E1=6.Fig.1 Trajectories of the system（4）with two different E1values：（a）E1=4 and（b）E1=6.

2　隨機(jī)生態(tài)系統(tǒng)模型研究

2.1建立具有獨(dú)立常數(shù)率捕撈的隨機(jī)生態(tài)系統(tǒng)模型

環(huán)境中存在的各種不確定性因素，會(huì)對(duì)食餌出生率和捕食者死亡率產(chǎn)生影響.本文建立色噪聲激勵(lì)下的隨機(jī)生態(tài)系統(tǒng)模型：

其中，X，Y分別代表被捕食者和捕食者的種群數(shù)量；ξi（t）（i=1，2）是相互獨(dú)立的平穩(wěn)色噪聲，分別作用于被捕食者自然增長(zhǎng)率z（忽略色噪聲對(duì)被捕食者間競(jìng)爭(zhēng)的影響）和捕食者自然死亡率d.其均值為零，相關(guān)函數(shù)為：

對(duì)模型（6）進(jìn)行如下變換：

根據(jù)文獻(xiàn)［18］，對(duì)于一般的隨機(jī)模型：

其中，fi（X，t）和gij（X，t）是關(guān)于X和時(shí)間t的函數(shù)，ξj（t）是零均值的隨機(jī)激勵(lì).若ξj（t）的相關(guān)時(shí)間遠(yuǎn)小于隨機(jī)過程Xi的觀察時(shí)間間隔（數(shù)值積分步長(zhǎng)）Δt，合理的選取Δt，隨機(jī)過程Xi能夠近似為一個(gè)It?隨機(jī)微分方程，其漂移系數(shù)aj（X，t），擴(kuò)散系數(shù)bik（X，t）為：

Cjs（τ）為ξj（t）和ξs（t）的相關(guān)函數(shù).對(duì)于相互獨(dú)立的隨機(jī)激勵(lì)，即有Cjs（τ）=Csj（τ）=0，（j≠s）.

將（11）～（13）式應(yīng)用到（9）式，可以轉(zhuǎn)化為It?隨機(jī)微分方程：

其中，Kjj如（13）式定義，即近似為噪聲強(qiáng)度. B1（t）和B2（t）為兩個(gè)獨(dú)立單位Wiener（維納）過程.

2.2穩(wěn)態(tài)概率密度研究

為了獲得隨機(jī)捕食生態(tài)系統(tǒng)中兩物種的穩(wěn)態(tài)概率密度，結(jié)合實(shí)驗(yàn)所得數(shù)據(jù)及（8），（10）式，對(duì)（14）式作如下說明：1）被捕食者環(huán)境容納量參數(shù)K遠(yuǎn)大于被捕食者自然增長(zhǎng)率參數(shù)z，故s是小量；2）當(dāng)e1＜E1＜e2時(shí)，利用四階Runge-Kutta算法計(jì)算（14）式得到G1（X，Y），G2（X，Y）均相對(duì)參數(shù)a，b，c和f是小量；3）合理假設(shè)系統(tǒng)受到的隨機(jī)擾動(dòng)ξi（t）（i=1，2）的噪聲強(qiáng)度是小量.

定義相應(yīng)于（2）式的隨機(jī)過程R（t）：

利用It?公式及（14）式可得到關(guān)于R（t）的It?隨機(jī)微分方程：

由以上說明可知，（16）式右邊是小量，即R（t）是一個(gè)慢變過程，而系統(tǒng)（14）可以作為保守系統(tǒng)（1）受擾動(dòng)后產(chǎn)生的新系統(tǒng).

根據(jù)Stratonovich-Khasminskii隨機(jī)平均原理，可以得到隨機(jī)過程R（t）的平均It?隨機(jī)微分方程為：

其中，m（R）和σ（R）分別為漂移系數(shù)和擴(kuò)散系數(shù)，可以由下式解得：

〈·〉t表示一個(gè)周期內(nèi)的時(shí)間平均，T（R）如（3）式定義.結(jié)合（1），（3）式可以得到：

其中

將式（20）～（22）代入（18）式和（19）式得到簡(jiǎn)化后的漂移系數(shù)和擴(kuò)散系數(shù)：

It?隨機(jī)微分方程（17）式對(duì)應(yīng)的穩(wěn)態(tài)概率密度滿足的FPK方程為：

其中，r是隨機(jī)過程R（t）的狀態(tài)變量，m（r）和σ2（r）對(duì)應(yīng)m（R）和σ2（R），求得（23）式對(duì)應(yīng)的穩(wěn)態(tài)概率密度［11］：

其中，C為歸一化常數(shù).假設(shè)當(dāng)R=r時(shí)隨機(jī)過程X的條件概率密度為p（x|r），則有：

由上式可得R和X的聯(lián)合概率密度為：

將Y作為R和X的函數(shù)，條件概率公式和Jacobi變換得到R和X的聯(lián)合概率密度p（r，x）為：

其中，C為歸一化常數(shù)，從而可以得到兩物種的穩(wěn)態(tài)概率密度為：

p（x），p（y）的存在說明隨機(jī)系統(tǒng)不再穩(wěn)定到某個(gè)穩(wěn)定點(diǎn)而是穩(wěn)定到一定的分布.

3　數(shù)值模擬

本文采用常用的低通濾過過程產(chǎn)生的色噪聲模型ξi（t）來(lái)模擬隨機(jī)環(huán)境.低通濾過過程由下列一階線性微分方程產(chǎn)生：

其中，W1（t），W2（t）是獨(dú)立的高斯白噪聲，其強(qiáng)度分別為D1（t），D2（t），可以描述色噪聲ξi（t）的相關(guān)時(shí)間.ξi（t）的譜密度Φi（ω）和相關(guān)函數(shù)Cii（τ）為：

圖2　αi取不同值時(shí)低通濾過過程的譜密度曲線Figure.2 Power spectral densities of the low-pass fltered processes（24）for different values of band width parameter αi

圖2給出了白噪聲強(qiáng)度為0.005，αi取不同值時(shí)低通濾過過程的譜密度曲線.可以得到，隨著αi的增大，譜帶變的平緩.

圖3給出了參數(shù)2D1=2D2=0.01，α1=α2=2，其余參數(shù)取值如（5）式定義，E1取不同值時(shí)，兩種群的穩(wěn)態(tài)概率密度的解析解（Analytical）和模擬（Simulation）曲線.從圖中可以得到：與確定性生態(tài)系統(tǒng)數(shù)據(jù)表1對(duì)應(yīng)，隨機(jī)系統(tǒng)漸近穩(wěn)定時(shí)，被捕食者種群數(shù)量隨著E1的增大在143.3附近振動(dòng)（圖3（a）），捕食者種群逐漸減少（圖3（b））；兩種群x，y均以固定概率在系統(tǒng)的平衡點(diǎn)附近振動(dòng)，當(dāng)E1=2時(shí)，振動(dòng)范圍最大，隨著E1的變大，振動(dòng)明顯減少，系統(tǒng)振動(dòng)范圍趨于穩(wěn)定.因此，在捕撈強(qiáng)度較小時(shí)，隨機(jī)因素不能忽略，隨著捕撈強(qiáng)度的加大，隨機(jī)因素對(duì)系統(tǒng)的影響減小.解析解曲線與Monte-Carlo模擬兩種群密度基本一致，說明了隨機(jī)平均方法的正確性.

在其余參數(shù)取值不變的前提下，圖4（a）和（b）給出了E1=6，α1=α2=2，隨機(jī)激勵(lì)強(qiáng)度2Di分別取0.005，0.01，0.05時(shí)兩種群的穩(wěn)態(tài)概率密度（其他參數(shù)如上）.從圖中可以得到，當(dāng)2D1=2D2=0. 005時(shí)，系統(tǒng)振動(dòng)范圍較小，隨著2Di的增大，兩種群x，y振動(dòng)的范圍明顯較少，即說明外界不確定性因素越強(qiáng)，生態(tài)系統(tǒng)越不穩(wěn)定.

圖3　不同E1取值下兩物種的穩(wěn)態(tài)概率密度.Fig.3 The stationary probability density of two species with two different E1values..

圖4　不同2Di取值下兩物種的穩(wěn)態(tài)概率密度.Fig.4 The stationary probability density of two species with two different 2Divalues.

圖5　不同αi取值下兩物種的穩(wěn)態(tài)概率密度.Fig.5 The stationary probability density of two species with two different αivalues

圖5（a）和（b）給出了在其余參數(shù)取值不變，2D1=2D2=0.01，E1=6，隨機(jī)激勵(lì)譜帶寬度αi分別取1，2，3時(shí)兩種群的穩(wěn)態(tài)概率密度.從圖中可以得到，隨著αi的增大，系統(tǒng)振動(dòng)范圍減少，即可以得到，噪聲譜帶越寬，生態(tài)系統(tǒng)越穩(wěn)定；另即噪聲的相關(guān)時(shí)間越小，生態(tài)系統(tǒng)越穩(wěn)定.

4　結(jié) 論

本文研究了被捕食者和捕食者具有獨(dú)立常數(shù)率捕撈，被捕食者常數(shù)比例進(jìn)入庇護(hù)所的捕食生態(tài)系統(tǒng)模型，主要分析了被捕食者捕撈強(qiáng)度E1在區(qū)間（e1，e2）內(nèi)系統(tǒng)的動(dòng)力學(xué)行為；建立了在更加接近實(shí)際環(huán)境的色噪聲激勵(lì)下的隨機(jī)生態(tài)模型，并將隨機(jī)模型通過近似手段轉(zhuǎn)化為It?隨機(jī)微分方程，在合理假設(shè)下，利用隨機(jī)平均原理得到兩物種的穩(wěn)態(tài)概率密度的理論解；運(yùn)用Monte-Carlo模擬驗(yàn)證了整個(gè)理論的正確性.對(duì)兩物種穩(wěn)態(tài)概率密度的研究發(fā)現(xiàn)：捕撈活動(dòng)較小時(shí)，環(huán)境因素對(duì)系統(tǒng)影響較大.隨著捕撈活動(dòng)的加大，隨機(jī)因素影響逐漸減少，捕撈活動(dòng)主要影響系統(tǒng)的穩(wěn)態(tài).噪聲強(qiáng)度越大，隨機(jī)因素對(duì)系統(tǒng)影響越大.噪聲相關(guān)時(shí)間越小，隨機(jī)因素對(duì)系統(tǒng)影響越小.另外，本文的模型研究捕食者捕撈強(qiáng)度E2不變時(shí)，被捕食者捕撈強(qiáng)度E1對(duì)系統(tǒng)的影響，對(duì)于研究E1不變，E2對(duì)系統(tǒng)影響的方法以及結(jié)果類似.本文所有結(jié)果與自然規(guī)律一致，可以有效地解釋并指導(dǎo)生產(chǎn)實(shí)踐.

本文的研究建立在具有常數(shù)率捕撈和常數(shù)比例被捕食者進(jìn)入庇護(hù)所的確定性模型的基礎(chǔ)上，并且限制在了（e1，e2）之間，對(duì)于其他類型的捕撈活動(dòng)，庇護(hù)所效應(yīng)模型，以及在隨機(jī)激勵(lì)下（0，e1），（e2，z/q1）范圍內(nèi)系統(tǒng)的運(yùn)動(dòng)性態(tài)及其方法，還有待研究.

1 Lotka A J.Element of Physical Biology.Baltimore：Williams and Wilkins Press，1925

2 Volterra V.Varizaioni e fluttuazioni del numero d’individui in specie d’animani conviventi.Mem Acad Lincei，1926，2：31～113

3 王洪禮，許佳，郭龍等.海洋赤潮藻類的生態(tài)動(dòng)力學(xué)穩(wěn)定性研究.動(dòng)力學(xué)與控制學(xué)報(bào)，2005，3（2）：40～43（Wang H L，Xu J，Guo L，et al.Research of red tide in ocean on dynamical stability.Journal of Dynamics and Control，2005，3（2）：40～43（in Chinese））

4 Saha T，Bandyopadhyay M.Dynamical analysis of toxin producing Phytoplankton-Zooplankton interactions.Nonlinear Analysis：Real World Applications，2009，10：314～332

5 Holling C S.Some characteristics of simple types of predation and parasitism.The Canadian Entomologist，1959，91（7）：385～398

6 Ak?akaya H R，Arditi R，Ginzburg L R.Ratio-dependent predation：an abstraction that works.Ecology，1995，76（3）：995～1004

7 Kar T K.Stability analysis of a prey-predator model incorporating a prey refuge.Communications in Nonlinear Science and Numerical Simulation，2005，10（6）：681～691

8 謝文賢，蔡力，岳曉樂.兩種群隨機(jī)動(dòng)力系統(tǒng)的信息熵和動(dòng)力學(xué)研究.物理學(xué)報(bào)，2012，61：170509（Xie W X，Cai L，Yue X L，et al.Information entropies and dynamics in the stochastic ecosystem of two competing species.Acta Physica Sinica，2012，61：170509（in Chinese））

9 Wu Y，Zhu W Q.Stochastic analysis of a pulse-type preypredator mode.Physical Review E，2008，77（4）：041911

10 楊會(huì)會(huì)，寧麗娟.非線性漂移的Fokker-Planck方程的近似非定態(tài)解.物理學(xué)報(bào)，2013，62（18）：180501（Yang H H，Ning L J.Approximate time-dependent solution of Fokker-Planck equation with non-linear drift force.Acta Physica Sinica，2013，62（18）：180501（in Chinese））

11 Cai G Q，Lin Y K.Stochastic analysis of the Lotka-Volterra model for ecosystems.Physical Review E，2004，70：041910

12 Qi L Y，Xu W，Gao W T.Stationary response of Lotka-Volterra system with real noises.Communications in Theoretical Physics，2013，4：503～509

13 Kumar Kar T.Modelling and analysis of a harvested preypredator system incorporating a prey refuge.Journal of Computational and Applied Mathematics，2006，185（1）：19～33

14 Lv Y，Pei Y，Gao S，et al.Harvesting of a phytoplankton-zooplankton model.Nonlinear Analysis：Real World Applications，2010，11（5）：3608～3619

15 馬知恩.種群生態(tài)學(xué)的數(shù)學(xué)建模與研究.合肥：安徽教育出版社，1996，43（Ma Z N.Methods of mathematics and research in population ecology.Hefei：Anhui Education Publishing House，1996，43（in Chinese））.

16 Bairagi N，Jana D.On the stability and Hopf bifurcation of a delay-induced predator-prey system with habitat complexity.Applied Mathematical Modelling，2011，35（7）：3255～3267

17 Kuang Y，F(xiàn)reedman H I.Uniqueness of limit cycles in Gause-type models of predator-prey systems.Mathematical Biosciences，1988，88（1）：67～84

18 Lin Y K.Some observations on the stochastic averaging method.Probabilistic Engineering Mechanics，1986，1（1）：23～27

EFFECTS OF NOISES AND HARVESTING ON THE STABILITY IN THE PREY-PREDATOR ECOSYSTEM*

Dong Qingguo Ning Lijuan?
（College of Mathematics and Information Science，Shaanxi Normal University，Xi’an 710062，China）

This paper discusses stability problem of a new stochastic predation type ecosystem with corporating a prey-refuge and independent harvesting in either species.To explore the prey-harvesting and colored noises effects on the stability of the ecosystem，with the assumption of weak disturbances，the stationary probability density functions for both species were obtained by applying the Stratonovich-Khasminskii averaging principle.The accuracy of the results obtained from theoretical method was demonstrated by those obtained from Monte Carlo simulation.Results obtained show that：1）the ecosystem with smaller harvesting is less stable when the system is disturbed by noises；2）the stronger the noise intensities are，the less stable the ecosystem will be；3）the narrower the band width is，the less stable the ecosystem will be；a narrower band width leads to a less stable system；4）a smaller correlation time leads to a more stable system.

colored noise，harvesting，stochastic averaging，the stationary probability density

6 December 2013，revised 7 January 2014.

E-mail：ninglijuan@snnu.edu.cn

10.6052/1672-6553-2014-033

2013-12-06收到第1稿，2014-01-07收到修改稿.

*國(guó)家自然科學(xué)基金資助項(xiàng)目（11202120和61273311）

E-mail：ninglijuan@snnu.edu.cn

*The project supported by the National Natural Science Foundation of China（11202120 and 61273311）.