雙獵物種群合作對抗捕食者種群的建模與分析

摘要： 基于三維捕食-獵物模型討論獵物種群間的合作及其對捕食過程的影響. 首先， 利用動力系統理論分析系統平衡點的存在性和穩定性， 并證明無獵物種群平衡點的周期震蕩特征. 其次， 通過數值模擬驗證理論結果. 數值模擬結果表明， 無獵物種群y的平衡點E1和無獵物種群z的平衡點E2都是周期震蕩的， 并且隨著合作效應的增大震蕩幅度也逐漸增大.

關鍵詞： 捕食-獵物模型； 合作； 平衡點； 數值模擬

中圖分類號： O29文獻標志碼： A文章編號： 1671-5489（2025）02-0417-08

Modeling and Analysis of Cooperation between Double Prey Populations against Predator Population

FENG Shanshan， ZHANG Yongxin

（School of Mathematics， North University of China， Taiyuan 030051， China）

Abstract： Based on a three-dimensional predator-prey model， we discussed the cooperation between prey populations and their impact on the predation process. Firstly，

we analyzed the existence and stability of equilibrium points of the system by using dynamical system theory， and proved the periodic oscillation characteristics of equilibrium point

s in prey free populations. Secondly， the theoretical results were verified through numerical simulation. The numerical simulation results show that the equilibrium point E1 of the prey free populations y

and the equilibrium point E2 of the prey free populations z are both periodic oscillations， and the amplitude of the oscillations gradually increases with the increase of the cooperative effect.

Keywords： predator-prey model; cooperation; equilibrium point; numerical simulation

收稿日期： 2024-08-05.

第一作者簡介： 馮姍姍（1992—）， 女， 漢族， 博士， 副教授， 從事網絡傳播動力學的研究， E-mail： fengshanshan@nuc.edu.cn.

基金項目： 國家自然科學基金青年科學基金（批準號： 12201583）和中國博士后科學基金（批準號： 2024M752989）.

0 引 言

隨著自然環境的變化和人類活動的干擾， 許多物種的生存面臨嚴峻挑戰. 其中， 獵物種群面臨捕食者的威脅， 而獵物種群之間合作對抗捕食者已成為一種普遍的生存策略. 在這種策略下， 不同種群的獵物相互合作， 共同抵抗捕食者的攻擊， 從而提高自身的生存能力和繁衍成功率^［1］.雙獵物種群合作對抗捕食者是一種特殊的合作現象， 它通常發生在兩個或多個獵物種群之間.在這種情況下， 獵物種群之間通過一系列的信號和行為相互協作， 共同抵抗捕食者的威脅^［2-3］.例如， 一些昆蟲種群會與其他昆蟲種群合作， 以對抗同一種捕食者. 這種合作不僅可以提高獵物種群的生存率， 而且還可以促進它們之間的交流和合作^［4-5］.對于兩個獵物種群的合作， 可以有多種形式. 例如， 通過共享資源、 互相幫助等方式提高自身的安全性^［6］.研究兩個獵物種群合作對抗捕食者的現象具有重要的生態學意義.

目前， 關于種群之間合作的研究已有很多成果. 例如： Alves等^［7］在兩種群Lotka-Volterra捕食模型的基礎上， 考慮在捕食者的攻擊率上增加一個合作項， 提出了一類具有合作狩獵的功能性反應模型； Capone等^［8］研究具有擴散的合作狩獵捕食系統， 討論了平衡點的穩定性和Turing不穩定性條件； Singh等^［9］通過具有合作狩獵的空間捕食系統， 研究了不同的合作率、 捕食者的基本再生數和擴散系數對捕食系統的影響； Enatsu等^［10］分析了合作狩獵和具有成長時滯的捕食-被捕食模型的動力學性態； Shivam等^［11］研究了合作狩獵和隨機移動對捕食系統的影響； Banerjee等^［12］研究了一個由兩種獵物物種和一種捕食者種群組成的群落模塊， 其中獵物種群考慮了相互合作影響出生率， 而捕食者表現出捕獵合作. 但上述研究成果都更關注捕食者之間的合作捕食， 而忽略了不同獵物種群之間的合作也會對抗捕食者的捕食行為. 基于此， 本文從兩個獵物種群之間相互合作影響捕食者的捕食率出發， 建立相應的數學模型并分析其性態.設獵物種群之間的合作函數為φM=ωHM1+M， 表示獵物種群H在合作效應中獵物種群M對它的幫助， 其中合作效應隨合作獵物種群的豐度增加而飽和， ωgt;0表示獵物從合作中獲得的最大利益^［13-14］.本文建立如下模型：

dxdt=-dx+cp1x-ωz1+zy+bp2x-ωy1+yz，dydt

=r1y-p1x-ωz1+zy，dzdt=r2z-p2x-ωy1+yz，（1）

其中x表示捕食者種群， y，z表示兩個合作的獵物種群， r1，r2分別表示兩個獵物種群的內稟增長率， d表示捕食者種群的自然死亡率， b，c分別表示捕食者種群對獵物種群z和y的轉換率， p1，p2分別表示捕食者種群對獵物種群y和z的捕食率. 首先， 討論系統（1）邊界平衡點和正平衡點的存在性； 其次， 對平衡點的性態進行理論分析； 最后， 通過數值模擬驗證理論結果.

1 平衡點的存在性

顯然， 系統（1）一直存在滅絕平衡點E0=（0，0，0）. 下面分析其他情形. 當捕食者種群為零， 其他種群不為零， 即x=0， y≠0， z≠0時， 系統（1）平衡點滿足如下方程：

cωyz1+z+bωyz1+y=0，r1y+ωyz1+z=0，r2z+ωyz1+y=0.（2）

由式（2）中第二和第三個方程得y=-r2r2ωlt;0， z=-r1r1+ωlt;0， 因此在第一象限內不存在無捕食者種群的平衡點.

當獵物種群y=0， 其他種群不為零， 即x≠0， y=0， z≠0時， 系統（1）平衡點滿足

-dx+bp2xz=0，r2z-p2xz=0.（3）

解得系統存在無獵物種群y的平衡點為E1=r2p2，0，dbp2. 當獵物種群z=0， 其他種群不為零， 即x≠0， y≠0， z=0時， 系統（1）平衡點滿足

-dx+cp1xy=0，r1y-p1xy=0.（4）

解得系統（1）存在無獵物種群z的平衡點為E2=r1p1，dcp1，0. 當所有種群都不為零， 即x≠0， y≠0， z≠0時， 系統（1）平衡點滿足

-dx+cp1xy-cωyz1+z+bp2xz-bωyz1+y=0，r1-p1x+ωz1+z=0，r2-p2x+ωy1+y=0.（5）

由式（5）中第二和第三個方程可得

y=p2x-r2r2-p2x+ω，z=p1x-r1r1-p1x+ω.

若ygt;0， 則r2p2lt;xlt;r2+ωp2； 若zgt;0，則r1p1lt;xlt;r1+ωp1. 從而正平衡點E3存在的定義域為

maxr1p1，r2p2lt;xlt;minr1+ωp1，r2+ωp2.

最后將y=p2x-r2r2-p2x+ω和z=p1x-r1r1-p1x+ω代入式（5）中第一個方程并化簡， 得

1（r1-r2+ω）（r1-p1x+ω）（r2-p2x+ω）（Ax3+Bx2+Cx+D）=0，（6）

其中

A=（cp1+d）（r2p1p2+r1）+dω，B=p1（bp2-cp1-d）r22+{[cp21+[（-b-c）p2+d]p1-dp2]r1-[（b+2c）p1+d]p2ω}r2+（ω+r

1）[p2（cp1+d）r1+dω（p1+p2）]，C=-b（p1+p2）r32+{[（b+c）p1+bp2+d]r1+ω（bp2+2cp1+d）}r22+（ω+r1）[（-cp1-d）r1+（bp1+cp

2）ω]r2-dω（ω+r1）2，D=-r22[-br22+br1r2+ω（ω+r1）（b+c）].（7）

由方程（6）可知， x的值由三次方程

Ax3+Bx2+Cx+D=0（8）

在區間maxr1p1，r2p2lt;xlt;minr1+ωp1，r2+ωp2

的正實根決定. 下面求解三次方程的根. 將方程（8）除以A， 并令x=ξ-B3A， 則可化為以下形式：

ξ3+pξ+q=0，（9）

其中

p=3AC-B23A2，q=27A2D-9ABC+2B327A3.（10）

令

u=3-q2+q24+p327，v=3-q2-q24+p327.（11）

則方程（8）的3個根可表示為

ξk=u+v，k=1，-1+3i2u+-1-3i2v，k=2，

-1-3i2u+-1+3i2v，k=3.（12）

特別地， 方程Δ=q22+p33是三次方程（9）的判別式. 當Δlt;0時， 3個根均為不同的實根. 因此， 在定義域

maxr1p1，r2p2lt;xlt;minr1+ωp1，r2+ωp2

內模型存在共存平衡點E3=ξk-B3A，p2x-r2r2-p2x+ω，p1x-r1r1-p1x+ω.

綜上， 可得如下定理.

定理1 系統（1）的平衡點特征如下：

1） 滅絕平衡點E0=（0，0，0）；

2） 無獵物種群y的平衡點E1=r2p2，0，dbp2；

3） 無獵物種群z的平衡點E2=r1p1，dcp1，0；

4） 共存平衡點E3=ξk-B3A，p2x-r2r2-p2x+ω，p1x-r1r1-p1x+ω.

2 平衡點的穩定性

2.1 滅絕平衡點E0

定理2 系統（1）的滅絕平衡點E0是一個鞍點.

證明： 系統（1）在E0處的Jacobi矩陣為

JE0=-d000r1000r2，

JE0的3個特征值分別為λ1=-d， λ2=r1， λ3=r2. 顯然， λ1lt;0， λ2gt;0， λ3gt;0， 故E0是鞍點.

2.2 無獵物種群y的平衡點E1

定理3 無獵物種群y的平衡點E1是不穩定的.

證明： 系統（1）在E1處的Jacobi矩陣為

JE1=0cp1r2p2-cωdbp2+d-ωdp2br2

0r1-p1r2p2+ωdbp2+d0-dbωdbp20，

JE1對應的特征方程為

（λ2+br2）λ-r1-p1r2p2+ωdbp2+d=0，（13）

其特征值分別為

λ1=λ2=±br2i，λ3=r1-p1r2p2+ωdbp2+d，

故平衡點E1是不穩定的. 證畢.

進一步， 集合{（x，y，z）y=0}顯然是系統（1）的一個不變集. 為此， 令y=0， 則此時可得系統（1）的一個等價系統為

dxdt=-dx+bp2xz，dzdt=r2z-p2xz.（14）

易得系統（14）的奇點M1r2p2，dbp2. 此時M1對應于系統（1）的平衡點E1. 對M1做平移變換

x=x-r2p2，z=z-dbp2，

則系統（14）變為

dxdt=br2z+bp2xz，dzdt=-dbx-p2x

z.（15）

奇點M1變為新坐標下的原點1（0，0）. 由軌線滿足的方程dxdz=br2+bp

2xz-db+p2zx

求得首次積分為

bz+x-dlnp2z+db-r2lnp2x+r2=k.（16）

易知， 在1的足夠小鄰域內式（16）左端連續. 根據文獻［15］可知， M1為中心點. 其外圍被一閉軌線族{Γk}環繞. 說明在M1點附近， x與z呈周期變化. 固定參數取值為r1=0.5， r2=0.8， b=0.2

， c=0.5， d=0.3， p1=0.5， p2=0.5， ω=0.3， 則系統（1）在平面（x，z）內的軌線為一簇周期解， 如圖1所示.

2.3 無獵物種群z的平衡點E2

定理4 無獵物種群z的平衡點E2是不穩定的.

證明： 系統（1）在E2處的Jacobi矩陣為

JE2=0cr1ωdp1+br1p2p1-bωdcp1+d-

dc0ωdcp100r2-p2r1p1+ωdcp1+d，

其對應的特征方程為

（λ2+r1d）λ-r2-p2r1p1+ωdcp1+d=0，（17）

特征值分別為λ1=λ2=±dr1i， λ3=r2-p2r1p1+ωdcp1+d，

故平衡點E2是不穩定的. 證畢.

進一步， 集合{（x，y，z）z=0}是系統（1）的一個不變集. 為此， 令z=0， 則可得系統（1）的一個等價系統為

dxdt=-dx+cp1xy，dydt=r1y-p1xy.（18）

易得系統（18）的奇點為M2r1p1，dcp1. 此時M2對應于系統（1）的平衡點E2. 對M2做平移變換

x=x-r1p1，y=y-dcp1，

則系統（18）變為

dxdt=cr1y+cp1xy，dydt=-dcx-p1x

y.（19）

奇點M2變為新坐標下的原點2（0，0）. 由軌線滿足的方程

dxdy=cr1+cp1xy-dc+p1yx

求得首次積分為

cy+x-dlnp1y+dc-r1lnp1x+r1=m.（20）

易知在2的足夠小鄰域內， 式（20）左端連續. 根據文獻［15］可知， M2為中心點. 其外圍被一閉軌線族{Γm}環繞. 說明在M2點附近， x與y呈周期變化. 固定參數取值為r1=0.5， r2=0.8， b=0.2

， c=0.5， d=0.3， p1=0.5， p2=0.5， ω=0.3， 則系統（1）在平面（x，y）內的軌線為一簇周期解， 如圖2所示.

圖1 （x，z）平面內的周期解

Fig.1 Periodic solutions in plane （x，z）

圖2 （x，y）平面內的周期解

Fig.2 Periodic solutions in plane （x，y）

3 數值模擬

下面利用MATLAB軟件模擬系統（1）的動力學性態， 考察獵物合作對整個系統的影響， 其中的數據參考文獻［16-17］.

首先， 模擬無獵物種群y的平衡點E1的性態. 由定理2可得系統（1）在無獵物種群y的平衡點E1處的Jacobi矩陣為一對純虛根和一個實根， 因此無獵物種群y的平衡點E1必是不穩定的. 取一組參數r1=0.6， r2=0.8， b=0.6， c=0.5， d=0.1， p1=0.5， p2=0.5， ω=0.3， 初始點為（x，y，z）=（2，2，2）， 模擬結果如圖3（A）所示. 由圖3（A）可見， 獵物種群y趨于滅絕， 捕食者種群和獵物種群z一直處于周期震蕩水平. 為考察獵物合作對系統性態的影響， 改變參數ω， 分別取ω=0.3， ω=0.6， ω=0.8， 模擬結果如圖3所示.

由圖3可見， 隨著參數ω的增大， 捕食者種群和獵物種群z的震蕩幅度逐漸增大（以物種x為例， 其值由5逐漸增大為5.9）. 因此， 獵物從合作中獲得的最大利益ω可以提高捕食者種群和獵物種群z的最大種群規模.

其次， 模擬無獵物種群z的平衡點E2的性態. 由定理3可得系統（1）在無獵物種群z的平衡點E2處的Jacobi矩陣為一對純虛根和一個實根， 因此無獵物種群z的平衡點E2必是不穩定的. 取一組參數r1=0.7， r2=0.8， b=0.6， c=0.5， d=0.1， p1=0.5， p2=0.7， ω=0.3， 初始點為（x，y，z）=（2，2，2）， 模擬結果圖4（A）所示. 由圖4（A）可見， 獵物種群z趨于滅絕， 捕食者種群和獵物種群y一直處于周期震蕩水平. 為考察獵物合作對系統性態的影響， 改變參數ω， 分別取ω=03， ω=06， ω=0.8

， 模擬結果如圖4所示. 由圖4可見， 隨著參數ω的增大， 捕食者種群和獵物種群y的震蕩幅度逐漸增大（以物種x為例， 其值由4.8逐漸增大為6.5）. 因此， 獵物從合作中獲得的最大利益ω可以提高捕食者種群和獵物種群y的最大種群規模. 綜上可見， 無論是無獵物種群y的平衡點E1還是無獵物種群z的平衡點E2， 增加獵物之間合作的影響都會提高捕食者種群的最大規模.

最后， 對系統（1）正平衡點的穩定性進行模擬. 固定參數取值為r1=05， r2=08， b=02， c=0.5， d=0.3， p1=0.5， p2=0.5， ω=0.3， 模擬結果圖5所示. 由圖5可見， 系統（1）存在唯一的正平衡點（123，108，062）， 且該正平衡點是不穩定的.

綜上所述， 本文基于兩個獵物種群之間相互合作影響捕食者的捕食過程， 建立了雙獵物種群合作對抗捕食者種群的三維模型并進行了分析. 結果表明， 系統存在滅絕平衡點E0、 無獵物種群y的平衡點E1、 無獵物種群z的平衡點E2和共存平衡點E3， 且滅絕平衡點E0是一個鞍點， 無獵物種群y的平衡點E1和無獵物種群z的平衡點E2是不穩定的. 特別地， 系統（1）在平面（x，z）內的軌線為一簇周期解， 即當物種y的取值為0時， 物種x和z呈現周期性變化; 系統（1）在平面（x，y）內的軌線為一簇周期解， 即當物種z的取值為0時， 物種x和y呈現周期性變化. 數值模擬結果表明， 無獵物種群y的平衡點E1和無獵物種群z的平衡點E2都是周期震蕩的， 并且隨著合作效應的增大震蕩的幅度也逐漸增大.

參考文獻

［1］ABERNETHY G M， MULLAN R， GLASS D H， et al. A Multiple Phenotype Predator-Prey

Model with Mutation［J］.Physica A： Statistical Mechanics and Its Applications， 2017， 465： 762-774.

［2］BULL J J， RICE W R. Distinguishing Mechanisms for the Evolution of Co-operation［J］.Journal of Theoretical Biology， 1991， 149（1）： 63-74.

［3］SACHS J L， MUELLER U G， WILCOX T P， et al. The Evolution of Cooperation［J］.The Quarterly Review of Biology， 2004， 79（2）： 135-160.

［4］HUGHES D P， PIERCE N E， BOOMSMA J J. Social Insect Symbionts： Evolution in Hom

eostatic Fortresses［J］.Trends in Ecology and Evolution， 2008， 23（12）： 672-677.

［5］TUTTLE L J， LAMB R W， STRINGER A L. Differential Learning by Native versus Inv

asive Predators to Avoid Distasteful Cleaning Mutualists［J］.Functional Ecology， 2021， 35（7）： 1481-1490.

［6］Ellison A M. Ecological Neetworks： Linking Structure and Dynamics in Food Webs［J］.Ecology， 2007， 88（1）： 265-266.

［7］ALVES M T， HILKER F M. Hunting Cooperation and Allee Effects in Predators［J］.Journal of Theoretical Biology， 2017， 419： 13-22.

［8］CAPONE F， CARFORA M F， LUCA R D， et al. Turing Patterns in a Reaction-Diffusi

on System Modeling Hunting Cooperation［J］.Mathematics and Computers in Simulation， 2019， 165（C）： 172-180.

［9］SINGH T， DUBEY R， MISHRA V N. Spatial Dynamics of Predator-Prey System with H

unting Cooperation in Predators and Type Ⅰ Functional Response［J］.AIMS Mathematics， 2020， 5（1）： 673-684.

［10］ENATSU Y， ROY J， BANERJEE M. Hunting Cooperation in a Prey-Predator Model wit

h Maturation Delay［J］.Journal of Biological Dynamics， 2024， 18（1）： 2332279-1-2332279-26.

［11］SHIVAM， SINGH T， KUMAR M. Effect of Random Movement and Cooperative Hunting in

the Prey-Predator System： A Dynamical Approach［J］.International Journal of Biomathematics， 2024， 17（3）： 2350029-1-2350029-30.

［12］BANERJEE S， SHA A， CHATTOPADHYAY J. Cooperative Pr

edation on Mutualistic Prey Communities［J］.Journal of Theoretical Biology， 2020， 490： 110156-1-110156-12.

［13］HOLLAND J N， DEANGELIS D L， BRONSTEIN J L. Population Dynamics and Mutualism： F

unctional Responses of Benefits and Costs［J］.The American Naturalist， 2002， 159（3）： 231-244.

［14］MOUGI A， KONDOH M. Instability of a Hybrid Module of Antagonistic and Mutualistic Interactions［J］.Population Ecology， 2014， 56（2）： 257-263.

［15］馬知恩， 周義倉， 李承治. 常微分方程定性與穩定性方法［M］.2版. 北京： 科學出版

社， 2015： 144-147. （MA Z E， ZHOU Y C， LI C Z. Qualitative and Stability Methods

of Ordinary Differential Equations［M］.2nd ed. Beijing： Science and Technology Press， 2015： 144-147.）

［16］KLEBANOFF A， HASTINGS A. Chaos in One-Predator， Two-Prey Models： General Res

ults from Bifurcation Theory［J］.Mathematical Biosciences， 1994， 122（2）： 221-233.

［17］KRIVAN V， EISNER J. Optimal Foraging and Predator-Prey Dynamics Ⅲ［J］.Theoretical Population Biology， 2003， 63（4）： 269-279.

（責任編輯： 李 琦）