基于自由擴散下被捕食者-捕食者模型的圖靈不穩定分析

賴 婷，袁 權，楊喜艷，王浩華，4

（1.海南大學 理學院，海南 ?？?，570228；2.廣東金融學院 金融數學與統計學院，廣東 廣州，510521；3.海南省工程建模與統計計算重點實驗室，海南 ?？?，570228；4.海南大學 熱帶特色林木花卉遺傳與種質創新教育部重點實驗室，海南 ?？?570228）

被捕食者-捕食者模型通過定義所謂的功能性函數來描述捕食者在單位時間內消耗被捕食者的數量，以此來刻畫生態系統中不同物種之間存在的相互作用［1-4］.作為該類模型的一個關鍵性因素，Holling［5］提出了3種不同類型的功能反應函數來模擬被捕食者-捕食者模型，大量結果表明，Beddington-DeAngelis功能反應更能貼近真實的種群行為特征［6-9］.此外，生態系統中物種間相互作用在空間上的行為特征引起了學者的廣泛關注.自Turing［10］提出反應-擴散系統以來，反應-擴散系統已被廣泛研究.交叉擴散問題首先由Kerner［11］提出，Shigesada 等［12］將交叉擴散問題應用于競爭種群系統中，刻畫種群在相位空間中的相對聯系.在自然界中，考慮到被捕食者之間存在自由競爭而導致自相殘殺或同類相食的情況，大多數雄性動物為了種群繁衍，擴大種群數量，會找多個雌性動物繁衍后代，因此會導致同類爭奪配偶的情況發生，例如，眼鏡蛇.由于生物種群之間的自擴散和交叉擴散的普遍存在性，不同于文獻［13］，筆者考慮了具有Beddington-DeAngelis 功能反應且種群之間具有自擴散和交叉擴散模式下的被捕食者-捕食者模型.假定被捕食者具有Logistic增長及種群內部自由競爭，且該種群在空間中自由擴散.首先給出了系統滿足耗散性和一致持久性的充分條件，并對系統的非負平衡點進行穩定性分析和圖靈不穩定性分析，然后確定產生圖靈斑圖的分岔參數滿足的條件，并由此識別出圖靈區域，最后通過Matlab軟件進行數值模擬，選取不同的分岔參數可以觀察到，分岔參數會影響種群的動力學形態.通過數值模擬結論表明，系統產生的圖靈模式主要分為2 種：1）斑點模式遍布整個空間；2）條紋模式和斑點模式共存.但后者經過時間的推移，條紋分裂成斑點，最終以斑點模式遍布于整個空間為穩定的動力學形態.

1 模型建立

在自然界中，物種的生長與其自身的數量以及相關物種的數量有關.前人對形如以下模型的動力學行為的研究已做了大量的工作，主要考慮Beddington-DeAngelis 功能反應，被捕食者具有Logistic增長［14-15］：

基于系統（1），考慮了被捕食者存在種群內部自由競爭的情形，即被捕食者除了在Logistic 增長條件下受到環境承載率約束而導致食物、空間等資源的不足而產生的種內競爭外，還存在種群內部的自由競爭，得到如下模型

其中，U(t)，V(t)分別代表在t時刻被捕食者和捕食者的密度，r，K，h，e，c，d 為正常數，分別代表物種的內在增長率，被捕食者的環境承載率，被捕食者之間的種內競爭，被捕食者轉化為捕食者的轉化率，被捕食者的死亡率以及捕食者的死亡率稱為Beddington-DeAngelis 功能反應，參數a為被捕食者的最大消耗率，用參數b來衡量被捕食者干擾的影響.

對于模型（2），令

則模型（2）簡化成以下無量綱形式的模型（3）

通過改變自變量dτ→(β+u+v)dt，得到與模型（3）對應的模型（4）

分別記模型（4）的第一個方程為f(u，v)，第二個方程為g(u，v)，為了研究模型（2）的空間動力學，將模型（2）

轉化為研究模型（4），基于模型（4），考察如下自擴散和交叉擴散模型

其中，?2=?2/?x2+?2/?y2是二維空間中的拉普拉斯算子，描述了物種間的隨機運動，其中非負常數d11和d22分別為被捕食者和捕食者的自擴散系數，d12和d21分別為被捕食者和捕食者的交叉擴散系數，可能取正值或零，取正值表示一個物種向另一個較低密度的物種方向擴散［16］，主要考慮交叉擴散和自擴散同時存在的被捕食者-捕食者系統.

對于模型（5），在給定的初始條件和零通量條件［17］下進行分析.

其中，LX和LY分別表示系統在X，Y方向上的大小，n表示?Ω的向外單位法向量，零通量意味著沒有其他種群通過邊界.

2 系統的耗散性

從生物學的角度來看，耗散性意味著所有的種群都是有界的.將證明系統（3）的耗散性.在證明耗散性之前先給出2個引理，引理的證明前人已經給出［18］.

定理1當δ-μ＞0時，系統（3）是耗散的.

證明對于系統（3）的第一個方程，可得

設δ-μ＞0，然后使用引理2，并假設存在一個正實數M1，可得

因此，對于ε1＞0，存在一個T1＞0，對于任意的t＞T1有，u(t) ≤M1+ε1.

現在，對于任意的t＞T1有，

其中，G=(δ-μ+ω)(M1+ε1).

因此，使用引理1，可得

現在假設存在正實數M2，可得

因此，對于ε2＞0，存在一個正實數T2＞T1，對于任意t＞T2可得v(t) ≤M2+ε2.

綜上可知，當δ-μ＞0時，系統（3）是耗散的.

3 系統的一致持久性

一致持久性在生物學意義上保證了物種的長期生存.從分析學的角度來看，一致持久性的定義如下：

定義1如果系統（3）的每個初始條件為(u(0)，v(0)) ∈int()的解(u(t)，v(t))都滿足以下條件，則說系統（3）是一致持久的.

1）u(t) ≥0，v(t) ≥0對于?t≥0.

2）存在ε＞0，使得成立.

定理2當δ-μ-1 ＞0和(γ-ω)(m1-ε3) ＞ωβ時，系統（3）是一致持久的.

證明從系統（3）的第一個方程可得，對于任意的t＞T2，

4 非負平衡點的穩定性分析

為了研究圖靈不穩定性，考慮模型（4），通過計算，得出3個非負平衡點：

定理3當δ-μ＞0時系統（4）在平衡點E1處不穩定，反之，系統（4）在E1處穩定.

證明在穩定點E1處的雅可比矩陣為

由此可知，J(0，0)的特征多項式為f(λ)=(λ-δβ+μβ)(λ+ωβ)，特征值分別為λ1=(δ-μ)β，λ2=-ωβ，當δ-μ＞0 時，系統（4）有一個正根(δ-μ)β和一個負根-ωβ，所以系統（4）在E1點處不穩定.從生物學的意義上來說，在系統（4）中，所有物種都不可能全部滅絕.反之，若δ-μ＜0，則系統（4）的特征多項式有2個負根，故系統（4）在E1點處穩定.

定理4當β＞時，系統（4）在平衡點E2處局部漸近穩定，反之，系統（4）在平衡點E2處不穩定.

證明通過計算系統（4）在平衡點E2處的雅可比矩陣可知

根據前人的研究，在不帶有擴散項的常微分方程中，穩態是穩定的.在帶有交叉擴散的偏微分方程中穩態是經過非均勻擾動從穩定狀態到不穩定狀態，從而產生圖靈模式，給出系統（4）穩態穩定的條件［19］

定理5當tr(JE*) ＜0，且det(JE*) ＞0時，系統（4）在正平衡點E*=(u*，v*)處局部漸近穩定.

證明求出系統（4）在穩定點E*處的雅可比矩陣，雅可比矩陣在上文中已經給出，當

可以得到，矩陣A的2 個特征根都為負值，由此可知，系統（4）在穩定點E*處局部漸近穩定，通過計算，可以得到

綜上可知，當a11+a22＜0，a11a22-a12a21＞0時系統（4）在正平衡點E*處局部漸近穩定.

5 圖靈不穩定性分析

考慮系統（5）在正穩態E*處的圖靈不穩定性，將系統（5）在E*周圍線性化，對于依賴于空間和時間較小的擾動［19］

圖1 波色數限制產生圖靈不穩定性的范圍

為了觀察圖靈區域，設置參 數α=0.1，δ=0.5，γ=0.8，ω=0.3，μ=0.03，d11=0.01，d21=0.001，d22=1.

現在討論由參數β以及d12張成的空間，如圖2所示.

圖2 系統（5）的分岔圖

圖2 中綠色的直線代表Hopf 分岔，對應的Hopf 分岔參數值為βT=0.027 4，紅色的曲線代表Turing 分岔.區域Ⅳ表示圖靈區域，此區域中找到的穩定點在非均勻擾動下是不穩定的，因此可以觀察到圖靈斑圖.區域Ⅰ表示具有齊次均衡的系統是無條件穩定的.區域Ⅱ表示只存在Hopf 不穩定性.區域Ⅲ表示穩定點有可能發生圖靈不穩定性也有可能發生Hopf不穩定性.

保持其他參數不變，改變β值，得到不同分叉參數β對應的色散關系，如圖3所示.

圖3 不同的分叉參數對應的色散關系

圖3 中從下到上的第三條曲線對應的臨界參數值βT=0.110 311 177 473 585，當β＜βT時，圖靈不穩定性發生，反之，當β＞βT時，不會發生圖靈不穩定性，即系統的穩態是穩定的.β的值從圖3中的曲線從下到上分別取β=0.15，0.13，0.110 311 177 473 585，0.08，0.05，從圖3 中也可以看出，波色數在(0.055 9，2.342 1)區間內，會發生圖靈不穩定性，這與圖1相對應.

6 數值模擬

通過繪制圖4 觀察系統（4）的穩定性，設置參數為α=0.1，δ=0.5，γ=0.8，ω=0.3，μ=0.03，β=0.05，藍色和紅色曲線的初始值分別為(0.05，0.08)，(0.15，0.20).通過計算得出a11a22-a12a21=0.001 9，a11+a22=-0.005 4，從圖4 可以觀察到，無論被捕食者與捕食者種群密度的初始值取何值，系統最終都會趨于穩定，定理5得到驗證.

圖4 系統（4）趨于穩定

通過繪制圖5 觀察2 個種群數量關系的變化，圖5a、b 和c 為系統（3）的種群密度關系變化圖，圖5d 為系統（1）的種群密度關系圖.圖5a結果表明，控制被捕食者種內競爭的變量α=0.1時，被捕食者與捕食者密度變化出現同增同減，且以時間t為周期的周期震蕩現象.圖5c結果表明，當α=20，由于被捕食者種內競爭強度增強導致種群數量急劇下降，從而使得捕食者由于缺乏食物導致種群數量減少，當捕食者數量減少到一定程度時，被捕食者數量開始少量增長，并趨于平衡，捕食者由于食物嚴重匱乏而出現衰亡現象.對比圖5a與圖5c發現，被捕食者種內自由競爭的強度增強，捕食者由于食物匱乏導致衰亡，種群密度以時間t為周期的周期震蕩現象消失.

圖5 2個種群的數量關系

圖5b 結果表明，隨著時間的推移，種群密度變化趨于平穩，系統逐漸穩定，被捕食者與捕食者共生存.圖5d結果表明，被捕食者與捕食者密度變化也會出現同增同減，以時間t為周期的周期變化現象.對比圖5b與d發現，圖5d中被捕食者與捕食者的種群密度都比圖5b中2個種群的種群密度大，說明系統由于存在被捕食者的種內自由競爭從而導致被捕食者與捕食者的種群數量減少.

在圖靈斑圖的數值模擬中，需要將系統的空間和時間離散化，即從無限維變換為有限維的形式.在實際應用中，反應擴散系統在二維空間中的連續問題在M×N網格點的離散域內求解.網格的長度設置為常數h，時間步長設置為常數τ，描述擴散的Laplacian 算子采用有限差分法進行求解，對于模型（5），在(xi，yi)位置上的進行迭代，有如下形式

其中，Laplacian算子為

同樣地，也可以寫出?2Φ2.設置網格步長h=0.15，時間步長τ=0.001 以及網格大小M=N=220，并且設置固定的參數值為α=0.1，δ=0.5，γ=0.8，ω=0.3，μ=0.03，d11=0.01，d12=0.58，d21=0.001，d22=1，初始值為在穩定點E*附近均勻分布的隨機擾動，記為

其中，η1(x)，η2(x) ∈[-5×10-5，5×10-5]，通過改變β的值進行數值模擬，直到其行為特性不再改變時停止，觀察被捕食者種群的斑圖.數值模擬的結果如圖6所示.

圖6 β=0.015時被捕食者的圖靈斑圖

從圖6 可以看出，圖6a 由于迭代時間較少，在隨機擾動的作用下，形成不規則的圖案，隨著迭代時間的逐漸增大，系統逐漸形成斑點模式且最終斑點模式不再改變.

圖7為β=0.02，迭代時間分別為t=2×104、t=1×106、t=2×106、t=5×106所呈現的圖靈斑圖.

圖7 β=0.02時被捕食者的圖靈斑圖

同樣地，在圖7 a 中，由于迭代時間較少，在隨機擾動的作用下形成不規則的圖案，隨著時間的推移，系統開始形成長度較短的條紋模式，條紋隨著時間的增加長度逐漸增長，并趨于分裂，且空間中出現斑點，最后形成形態穩定的斑點圖.對比圖6 和圖7，可以看出，改變分岔參數β，會導致系統呈現出不同的圖靈模式.

7 小 結

研究了具有Beddington-DeAngelis 功能反應和自擴散以及交叉擴散的被捕食者-捕食者模型的空間動力學模式，在經典的被捕食者只具有Logistic 增長模型上，考慮了被捕食者存在種群內部自由競爭的情形，即被捕食者除了在Logistic 增長條件下受到環境承載力約束而導致的食物、空間等資源不足外，種群內部還存在自由的競爭形態.首先分析了系統滿足耗散性以及一致持久性的條件，其次對系統的3 個非負平衡點進行穩定性分析以及圖靈不穩定性分析，給出了產生圖靈不穩定的條件，通過圖靈不穩定性分析發現，波色數、分叉參數、擴散系數都會影響系統圖靈不穩定性的發生，從而影響圖靈斑圖的空間形態，在數值模擬中也驗證了這一點.同時，在數值模擬中可以清晰地看到，2個種群的數量關系呈現出以時間為周期的周期震蕩現象，隨著時間的推移，2 個種群數量最終趨于穩定，被捕食者與捕食者共生存.除此之外，通過圖5 的結果及分析表明，被捕食者存在種內自由競爭會導致被捕食者與捕食者種群數量減少，且當被捕食者種內競爭強度增大時，會導致捕食者因食物嚴重匱乏而衰亡，且種群數量以時間為周期的周期震蕩現象消失.