一維格點鏈上食餌—捕食者模型的蒙特卡洛模擬

一維格點鏈上食餌—捕食者模型的蒙特卡洛模擬

胡玉彪，喬文華

(包頭師范學(xué)院物理科學(xué)與技術(shù)學(xué)院，內(nèi)蒙古包頭014030)

摘要：在一維格點鏈上定義了改進(jìn)的兩種群Lotka-Volterra模型，在只考慮近鄰相互作用的情況下，分析了種群中每個生命個體的空間分布，及環(huán)境噪聲對食餌—捕食者生命系統(tǒng)隨時間進(jìn)化的影響，利用蒙特卡洛方法，模擬計算了其動力學(xué)演化過程，給出了食餌—捕食者系統(tǒng)中兩生命體能夠長期共存，種群生態(tài)系統(tǒng)穩(wěn)定的基本條件。

關(guān)鍵詞：一維格點鏈；食餌—捕食者模型；環(huán)境噪聲；蒙特卡洛模擬

收稿日期：2014-11-04

作者簡介：胡玉彪(1988-)，吉林省德惠市人，碩士研究生，研究方向：計算物理。

中圖分類號：O186.13文獻(xiàn)標(biāo)識碼：A

0引言

眾所周知，生命現(xiàn)象的多樣性決定了生態(tài)模型和生物數(shù)學(xué)的復(fù)雜性。當(dāng)研究生命問題時，有許多因素會影響到生態(tài)模型的動力學(xué)性質(zhì)，包括物種間的相互競爭，生存環(huán)境的變遷，人工的干預(yù)，突發(fā)性自然災(zāi)害的發(fā)生，都會對生物種群的生存與發(fā)展產(chǎn)生影響。正是由于生態(tài)系統(tǒng)的復(fù)雜性，需要借助于各種不同的數(shù)學(xué)模型來描述處在同一生態(tài)環(huán)境下不同生物種群間通過相互競爭，共生共存的自然現(xiàn)象。在描述由食餌-捕食者構(gòu)成的兩種群模型中，就有各種各樣的數(shù)學(xué)模型[1-9]，但是最基本的Lotka-Volterra[2,3]模型仍然是討論這類問題的出發(fā)點，為此，我們首先對由Lotka-Volterra模型建立的數(shù)學(xué)方程作一個簡要說明。由食餌-捕食者構(gòu)成的兩種群Lotka-Volterra非線性方程如(1)式所示

(1)

其中，a(t)是A(食餌)物種在時刻t的濃度，b(t)是B(捕食者)物種在時刻t的濃度，σ是物種A獨立生存時的增長率，μ是物種B獨自存在時的自然死亡率，比例系數(shù)λ1表示捕食者掠取食餌的強度(能力)，比例系數(shù)λ2表示食餌的豐度(或食餌對捕食者的供養(yǎng)能力)。

對于標(biāo)準(zhǔn)的Lotka-Volterra模型，因為σ、μ、λ1和λ2都是常數(shù)，很容易得到食餌—捕食者種群濃度，具有持續(xù)震蕩的準(zhǔn)周期模式的數(shù)值解，但是，這樣的結(jié)果往往被認(rèn)為和實際的生態(tài)系統(tǒng)存在較大的差距[1,7,8]。認(rèn)真分析其原因，我們認(rèn)為由于λ1和λ2代表捕食者與食餌的相互作用，在自然情況下，發(fā)生捕食行為的過程中，捕食者掠取食餌的概率應(yīng)該與食餌和捕食者的空間分布和相對位置有關(guān)。例如，某一時刻，一個食餌個體若處在捕食者的包圍之中，此時它被捕食的可能性就高，反之，若食餌處在本物種的中央位置，它沒有暴露在與捕食者的直接接觸中，此時這個生命個體被捕食的可能性就很??；同理，某一時刻，若一個捕食者個體其四周近鄰都有食餌存在，此時它捕獲食餌的可能性就高，反之，若其四周近鄰無食餌存在，此時這個捕食者捕獲食餌的可能性就小很多。由于生活在同一空間區(qū)域里兩個種群中的每個生命個體都在運動，隨著時間的增長，它們會不斷的移動到新的位置；更由于新生命體的產(chǎn)生、原有的某些生命體的死亡、以及被掠殺等現(xiàn)象的發(fā)生，必然使每個生命體自身的近鄰不斷變化，即其生態(tài)環(huán)境發(fā)生改變，所以，代表捕食者與食餌相互作用的λ1和λ2也應(yīng)該不斷隨時間而變，任一時刻表示捕食者掠取食餌的強度λ1和表示食餌豐度的λ2應(yīng)該根據(jù)該時刻的近鄰條件計算和更新，而不該是一個恒定的數(shù)值。另一方面，自然環(huán)境、生態(tài)環(huán)境噪聲的影響也是不可忽略的因素，因此，根據(jù)以上分析，我們將對兩種群的Lotka-Volterra模型做出新的詮釋，本文中，我們考慮在一維空間中，由食餌(A)和捕食者(B)構(gòu)成的兩種群系統(tǒng)中，每個生命個體所處空間位置對種群生存的影響，為使模型更接近于自然界的實際情況，我們將隨機的給出生命個體的空間位置，然后在一定的演化機制下，用蒙特卡洛方法模擬其演化行為。

1一維格點鏈模型

(2)

考慮食餌種群(A)不受資源限制的一維鏈，初始時刻隨機配置兩種群個體的空間分布，考慮任意時刻tn時的兩種群中每個生命體個體的空間分布，對第i個生命個體，只考慮最近鄰相互作用的情況下：

(i) 若s(i)=0，且s(i-1)=0與s(i+1)=0，表示在第i格點位上有一個A物種的生命個體，且它的兩側(cè)近鄰都是A物種，可以認(rèn)為此刻它沒有被捕食的風(fēng)險。

(ii) 若s(i)=0，且s(i-1)=0而s(i+1)=1，或s(i-1)=1而s(i+1)=0，表示在第i格點位上有一個A物種的生命個體，而它的一側(cè)近鄰是A物種，另一側(cè)近鄰有B物種出現(xiàn)，它被捕食的可能性將比(i)增大，存在一定的風(fēng)險。

(iii) 若s(i)=0，且s(i-1)=1與s(i+1)=1，表示在第i格點位上有一個A物種的生命個體，而它的兩側(cè)近鄰都有B物種存在，它被捕食的風(fēng)險最高。

(iv) 若s(i)=1，且s(i-1)=0與s(i+1)=0，表示在第i格點位上有一個B物種的生命個體，而它的兩側(cè)近鄰都是A物種，提供給此捕食者的食餌豐度最大。

(v) 若s(i)=1，且s(i-1)=0而s(i+1)=1，或s(i-1)=1而s(i+1)=0，表示在第i格點位上有一個B物種的生命個體，它的一側(cè)近鄰是B物種，而另一側(cè)近鄰有A物種出現(xiàn)，提供給此捕食者的食餌豐度將比(iv)降低。

(vi) 若s(i)=1，且s(i-1)=1與s(i+1)=1，表示在第i格點位上有一個B物種的生命個體，它的兩側(cè)近鄰都是B物種，提供給此捕食者的食餌豐度最低。

到下一時刻tn+1，由于捕食行為和種群間的相互作用，種群A和種群B中生命個體數(shù)量都發(fā)生改變，將要根據(jù)(2)式計算每個種群的個體數(shù)量(或濃度)，并依照上述方法重新配置每個種群個體的空間分布。依此，可以計算系統(tǒng)在任意時刻的演化結(jié)果。

2蒙特卡洛模擬

2.1數(shù)據(jù)結(jié)構(gòu)和初始條件

定義在時刻t，種群A生命體個數(shù)和種群B生命體個數(shù)占生物體總數(shù)的比率為pA(t)=nA(t)/N(t) 和 pB(t)=nB(t)/N(t)=1-pA(t)，其中nA(t)代表A種群個體數(shù)，nB(t)代表B種群個體數(shù)，兩種生物的總個數(shù)為N(t)=nA(t)+nB(t)。給定常數(shù)σ、μ的值，設(shè)初始時刻兩種生物的總個數(shù)為N0，以及pA(0)=pB(0)=0.5，即有nA(0)=0.5N0個A種群個體，有nB(0)=0.5N0個B物種個體，它們隨機的占據(jù)N0個格點位，設(shè)兩種群各自的初始濃度a(0)=0.5、b(0)=0.5。

設(shè)在t時刻，已求得a(t)和b(t)及種群A和B中每個生命個體的空間分布，即A物種中的個體占據(jù)了nA(t)個格點位， B物種中的個體占據(jù)了nB(t)個格點位，兩種群生命體的總個數(shù)為N(t)=nA(t)+nB(t)， 共占據(jù)了N(t)格點位。

2.2模擬步驟

(1)逐次考察每一個格點位上生命體的自旋態(tài)，對第i個格點位上的生命體。

(i)若s(i)=0，且s(i-1)=0與s(i+1)=0，則令k1(t)=0×λ1。其中λ1是方程(1)中的比例系數(shù)，表示捕食者掠取食餌的強度(能力)。

(ii) 若s(i)=0，且s(i-1)=0而s(i+1)=1，或s(i-1)=1而s(i+1)=0，則令k1(t)=(1/2)×λ1×(0~1)。其中(0~1)表示取一個0到1之間的隨機數(shù)。(下同)

(iii) 若s(i)=0，且s(i-1)=1與s(i+1)=1，則令k1(t)=(0~1)×λ1。

(iv) 若s(i)=1，且s(i-1)=0與s(i+1)=0，則令k2(t)=(0~1)×λ2。其中λ2是方程(1)中的比例系數(shù)，表示食餌的豐度(或食餌對捕食者的供養(yǎng)能力)。

(v) 若s(i)=1，且s(i-1)=0而s(i+1)=1，或s(i-1)=1而s(i+1)=0，則令k2(t)=(1/2)×(0~1)×λ2。

(vi) 若s(i)=1，且s(i-1)=1與s(i+1)=1，則令k2(t)=0×λ2。

(2) 考察遍歷一維鏈上所有的N(t)個格點后，分別求出nA(t)個k1(t)，及nB(t)個k2(t)的平均值

重復(fù)以上步驟(1)～(4)。

3模擬結(jié)果

圖1、圖2及圖3是根據(jù)方程(2)的計算，給出了不同噪聲、不同的反應(yīng)速率(保持σ、μ不變，λ1、λ2取不同的值)情況下的結(jié)果，可以看出食餌濃度a(t)和捕食者濃度b(t)作為時間的函數(shù)在起始階段有一個振蕩過程，但是振蕩的振幅隨時間逐漸減弱，慢慢趨于平穩(wěn)恒定的過程；從圖中可以看出噪聲越強，不但振蕩衰減的越快，而且兩生物種群最終保留的生命體數(shù)量也會受其影響，當(dāng)噪聲超過一定的強度(D≥0.001)時，對種群的生存產(chǎn)生不可逆轉(zhuǎn)的影響，兩種群都將隨時間不斷衰退，只至消亡。但是在相同的噪聲強度下，較大的掠食系數(shù)λ1及食餌豐度系數(shù)λ2會使得這種衰退在一定程度上減弱，甚至使種群得以能夠延續(xù)。這意味著在相互競爭和掠殺行為不斷發(fā)生的過程中，兩種群生命個體數(shù)量經(jīng)歷增長和死亡的振蕩起伏，逐步趨向于增長和衰減互相平衡，每個種群中生命個體數(shù)趨于基本穩(wěn)定，而較小的外界的各種干擾(噪聲)有利于生態(tài)種群快速趨于某種分布“穩(wěn)定”，較大噪聲則會帶來毀滅性的災(zāi)難；另一方面，較大的掠食系數(shù)λ1及食餌豐度系數(shù)λ2，代表食物的豐富和食餌種群數(shù)量的龐大，對于外界影響的抵御較強，所以遲緩了(或減弱了)噪聲的作用，這是一種在外界環(huán)境因素干擾下，兩種群系統(tǒng)自適應(yīng)和自調(diào)整的過程。

圖1

圖2

圖3

不同噪聲強度下，食餌、捕食者濃度隨時間的演化結(jié)果，紅色線表示捕食者濃度b(t)，藍(lán)色線表示食餌濃度a(t)，D為噪聲強度(D=0表示無噪聲)。

圖4、圖5及圖6是食餌濃度和捕食者濃度在a-b平面內(nèi)的軌跡，每幅圖都描述了同一反應(yīng)速率在不同的噪聲下的相圖(在同一幅圖中的計算中我們使用了相同的λ1、λ2)，從上述圖中可以看出，當(dāng)沒有噪聲或噪聲較小時(D≤0.00001)，食餌濃度和捕食者濃度在最初的振蕩和漲落后，最終趨向于一個中心的固定點(數(shù)學(xué)上的焦點)，這是一種穩(wěn)定的結(jié)構(gòu)，說明在此條件下食餌和捕食者兩物種可以長期共存。但是當(dāng)噪聲較大時(D≥0.001)，長時間的演化行為表明，不再能形成長時間的準(zhǔn)周期振蕩過程，在一個短暫的振蕩漲落后(螺旋圈數(shù)明顯減少)，a(t)、b(t)很快就趨向于一個緩慢衰減區(qū)域，其在相平面內(nèi)的軌跡不再穩(wěn)定于數(shù)學(xué)上的焦點。

不同噪聲強度下， 在a—b平面內(nèi)的運行軌跡(相圖)，軌跡前進(jìn)的正方向為逆時針方向， D為噪聲強度(D=0表示無噪聲)。

圖4

圖5

圖6

圖7

圖8

圖9

4結(jié)論

參考文獻(xiàn)〔〕

[1]Washerberger.M.J, Mobilia.M, T.uber.U.C. Influence of Local Carrying Capacity Restrictions on Stochastic Predator-Prey Models[J]. J. Phys. Condens. Matter,2007,19 (6):1088-1105.

[2]Lotka.A.J. Analtical. Note on Certain Rhythmic Relations in Organic Systems[J]. Biology,1920,6:410-415.

[3]Volterra.V. Lecons sur la. Theorie Mathematique De La Lutte Pour La Vie Redigees Par Marcel Brelot[J]. Book Review,1931,30:157-159.

[4]Javier.S, Tania.T. Stochastic Lattice Gas Model for a Predator-Prey System[J]. Phys. Rev. E,1993,49(10):5073-5079.

[5]Boccara.N, Roblin.O, Roger.M. Automata Network Predator-Prey Model with Pursuit and Evation[J]. Phys.Rev.E, 1994,50(6):4531-4541.

[6]Rozefeld.A.F, Albano.E.V. Study of a Lattice-Gas Model for a Prey-Predator System[J].Physica.A,1999,266(4):322-329.

[7]Lipowski.A. Oscillatory Behavior in a Lattice Prey-Predator System[J]. Phys. Rev. E.1999,60(5):5179-5184.

[8]Monetti.R, Rozenfeld.A, Albano.E. Study of Interacting Particle Systems: the Transition to the Oscillatory Behavior of A Prey-Predator Model[J].Physica.A.2000,283(2):52-58.

[9]Antal.T, Droz.M, Lipowski.A. Critical Behavior of a Lattice Prey-Predator Model[J]. Phys. Rev.E.2001, 64(3):1-6.

Monte Carlo Simulation of One-Dimensional Lattice

Chains Model for a Predator-Prey System

HU Yu-biao,QIAO Wen-hua

(Faculty of Physics and Technology,Baotuo Teachers College,Baotuo,014030)

Abstract:In one dimensional lattice chain we defines the improved Lotka-Volterr model of two populations,in the case of only considering the interaction neighbor,we study the spatial distribution of each individual life in population and the environmental noise of life how to effect on evolution of predator-prey system with time,by using Monte Carlo method,the dynamic evolution process are simulated,the predator-prey system is given in two organisms can co-exist for a long time and the population ecological basic conditions for the robust stability of the system.

Key words:One dimensional lattice chain；The predator-prey system model；Environmental noise；Monte Carlo simulation