一類蚜蟲種群突變模型的建立與實(shí)證分析

那 穎,趙立純,劉敬娜

(1.遼寧師范大學(xué) 數(shù)學(xué)學(xué)院,遼寧 大連 116029;2.鞍山師范學(xué)院 數(shù)學(xué)與信息科學(xué)學(xué)院,遼寧 鞍山 114007)

突變理論在生態(tài)學(xué)中的應(yīng)用已覆蓋7種突變模型中的5種,其中,尖角突變模型的應(yīng)用最為廣泛.如趙惠燕研究了以蚜蟲為主要控制對(duì)象的麥田生態(tài)系統(tǒng),發(fā)現(xiàn)麥蚜突變形成危害的過程符合尖角突變模型及其性質(zhì),從而用突變理論對(duì)麥蚜造成危害的現(xiàn)象給予科學(xué)的解釋與分析[1];她們還提出棉花苗蚜是棉花苗期的主要害蟲,每年都有發(fā)生,且受氣候、天敵等多因素影響,由于用傳統(tǒng)方法很難完全描述和預(yù)測(cè),所以利用尖角突變模型解釋和預(yù)測(cè)上述現(xiàn)象[2].這為突變理論在蚜蟲種群爆發(fā)的深入研究奠定了基礎(chǔ).

回歸分析方法在生態(tài)系統(tǒng)中的應(yīng)用也十分廣泛.王秀美等人針對(duì)小麥蚜蟲預(yù)測(cè)預(yù)警準(zhǔn)確率不高的問題,使用局部支持向量回歸分析構(gòu)建預(yù)測(cè)模型,并對(duì)未知樣本進(jìn)行預(yù)測(cè)[3].王向陽(yáng)等根據(jù)小麥蚜蟲及其天敵的調(diào)查數(shù)據(jù),采用逐步回歸分析,篩選出相關(guān)性強(qiáng)的預(yù)報(bào)因子,建立中長(zhǎng)期預(yù)測(cè)模型并給出預(yù)測(cè)結(jié)論[4].田潔基于實(shí)驗(yàn)數(shù)據(jù),利用回歸分析方法,建立擬合度較高的尖角突變模型,刻畫蚜蟲種群的生長(zhǎng)情況[5].這些研究的應(yīng)用均為多元線性回歸,而多元非線性回歸在突變理論中的成果較少見.

本文結(jié)合突變理論和回歸分析方法,分析數(shù)據(jù)的突變特征并判別控制變量類型,建立尖角突變模型.再利用R2值和F′統(tǒng)計(jì)量,以線性模型和Logistic模型為參考,分析所建模型的擬合度.最后,根據(jù)突變理論對(duì)模型的平衡點(diǎn)和穩(wěn)定性進(jìn)行分析,結(jié)果可為蚜蟲種群的預(yù)測(cè)預(yù)報(bào)提供參考.

1 模型建立

1.1 數(shù)據(jù)處理

為使所描述的模型更符合實(shí)際,以西北農(nóng)林科技大學(xué)某年蚜蟲種群生態(tài)系統(tǒng)調(diào)查資料為基礎(chǔ)數(shù)據(jù),其中，包括百株蚜量、天敵量和以月、旬或候?yàn)橐粋€(gè)時(shí)段單位的平均氣溫、平均相對(duì)濕度、平均降水、平均日照時(shí)數(shù).

根據(jù)蚜蟲種群生態(tài)系統(tǒng)的調(diào)查資料,百株蚜量的數(shù)據(jù)很全面,所以把蚜蟲種群作為本文研究的狀態(tài)變量,記為x;由于蜘蛛量的數(shù)據(jù)也很全面,所以把蜘蛛量作為與天敵種群相關(guān)的控制變量,記為a;而另一個(gè)控制變量是復(fù)合平均氣溫、平均相對(duì)濕度、平均降水、平均日照時(shí)數(shù)的一個(gè)復(fù)合的氣象因素,記為b.利用主成分分析法,得到主成分得分和各個(gè)主成分貢獻(xiàn)率的結(jié)果如下：

其中，矩陣的第一列代表b1的主成分得分,第二列代表b2的主成分得分,第三列代表b3的主成分得分,第四列代表b4的主成分得分.

b1的貢獻(xiàn)率為42.935 5%,b2的貢獻(xiàn)率為40.653 9%,b3的貢獻(xiàn)率為12.144 4%,b4的貢獻(xiàn)率為4.266 2%.

根據(jù)b1,b2,b3和b4的主成分得分,并用其貢獻(xiàn)率加權(quán),得到基于氣象因素的復(fù)合控制變量b的表達(dá)式為:

b=42.935 5%b1+40.653 9%b2+12.144 4%b3+4.266 2%b4,

其中,b1代表平均氣溫,b2代表平均相對(duì)濕度,b3代表平均降水,b4代表平均日照時(shí)數(shù).

注:由于控制變量b能夠較好地反映氣象因素的變化信息,因此本文中氣象因素可以看成蚜蟲種群發(fā)展變化的主導(dǎo)因素之一.

1.2 突變特征分析

考慮到以往研究從理論上說明突變特征問題的主觀盲目性[6-7],本文基于實(shí)際數(shù)據(jù),繪制狀態(tài)變量分別隨兩個(gè)控制變量變化的散點(diǎn)圖(見圖1),來判斷蚜蟲種群的突變特征.

圖1 蚜蟲種群生態(tài)系統(tǒng)散點(diǎn)圖

由文獻(xiàn)[8]得到,當(dāng)控制變量連續(xù)變化時(shí),如果狀態(tài)變量有不連續(xù)的跳躍行為,則說明蚜蟲種群具有突變特征.

圖1中存在5個(gè)不連續(xù)行為,其中1是由于控制變量a突然變化引起的不連續(xù)行為,所以1不是突變現(xiàn)象;2是由于控制變量a和控制變量b突然變化引起的不連續(xù)行為,所以2不是突變現(xiàn)象;3是由于控制變量a和控制變量b漸變引起的不連續(xù)行為,所以3是突變現(xiàn)象;4是由于控制變量a和控制變量b突然變化引起的不連續(xù)行為,所以4不是突變現(xiàn)象;5是由于控制變量a和控制變量b漸變引起的不連續(xù)行為,所以5是突變現(xiàn)象.

綜上可知蚜蟲種群具有突變特征,由于這個(gè)突變特征是由兩個(gè)控制變量a和b連續(xù)變化引起的,所以可以考慮用尖角突變模型來刻畫蚜蟲種群的突變現(xiàn)象.但哪個(gè)變量對(duì)模型突變起主要作用,需要進(jìn)一步思考.

1.3 控制變量類型判別

文獻(xiàn)[9]提出,尖角突變模型的一般形式為:

x3+ux+v=0,

(1)

其中,x代表尖角突變模型的狀態(tài)變量,u代表剖分變量,v代表正則變量.

為了準(zhǔn)確判別正則變量和剖分變量,假設(shè)u=a和v=b,得到如下模型:

x3+ax+b=0,

(2)

或者相反,假設(shè)u=b和v=a,得到如下模型:

x3+bx+a=0.

(3)

再借助SPSS軟件,分別對(duì)尖角突變模型(2)和(3)進(jìn)行多元非線性回歸分析,得到結(jié)果如下:

SSResidualcusp1=7.397,R2cusp1=0.171,

SSResidualcusp2=138.364,R2cusp2=0.166.

由于R2cusp1>R2cusp2且SSResidualcusp1

(4)

2 擬合度分析

Malthusian線性模型[11]和Logistic模型[12]在刻畫種群增長(zhǎng)方面應(yīng)用極其廣泛,用本文提出的尖角突變模型與線性模型和Logistic模型對(duì)比,來分析擬合度更具有說服力.下面基于實(shí)際數(shù)據(jù)分析3個(gè)模型的擬合度,具體分析如下:

根據(jù)實(shí)測(cè)資料,許多生物和環(huán)境研究量沒有負(fù)值,例如天敵和氣象因素,因此,需要將原來的坐標(biāo)系進(jìn)行坐標(biāo)變換,將原來的a,b軸平移和旋轉(zhuǎn),使之具有生物意義.這樣將模型(4)平移旋轉(zhuǎn)后得:

p1xi+q13+p2ai+q2xi+q1+bi+q3=0.

(5)

根據(jù)Malthusian模型,在其它條件最優(yōu)時(shí),天敵種群數(shù)量a的變化直接影響蚜蟲種群的變化率,而氣象因素對(duì)蚜蟲種群的影響具有一個(gè)最優(yōu)值取為ε,基于此,得到如下模型:

(6)

將模型(6)平移旋轉(zhuǎn)后得:

p1xi+q1ai+q2-|b+q3-ε|=0，

(7)

對(duì)模型(7)進(jìn)行回歸分析,得到結(jié)果如下:

SSResidualplane=8.920,R2plane=0.082.

根據(jù)模型Logistic,一般狀態(tài)下,蜘蛛量a是與內(nèi)稟增長(zhǎng)率r相關(guān)的變量,氣象因素b是與環(huán)境容納量K相關(guān)的變量,得到如下模型:

，

(8)

將模型(8)平移旋轉(zhuǎn)后得:

(9)

對(duì)模型(9)進(jìn)行回歸分析,得到結(jié)果如下:

SSResidualLogistic=603.783,R2Logistic=0.161.

上述ai,bi,xi(i=1,2,…,16)是實(shí)測(cè)數(shù)據(jù)，qi(i=1,2,3)為平移的距離,pj(j=1,2)為旋轉(zhuǎn)的角度.

從得到的R2值和SSR可知,尖角突變模型的擬合度最優(yōu).為進(jìn)一步驗(yàn)證,根據(jù)顯著性檢驗(yàn)[8,13],設(shè)F′統(tǒng)計(jì)量如下:

(10)

其中，n表示樣本總數(shù),在本文中n=16.若F′>0,說明SSResidualplane>SSResidualcusp,意味著尖角突變模型擬合度優(yōu)于線性模型;同理若F′<0,說明SSResidualplane

(11)

顯然F1′=1.132>0,說明尖角突變模型擬合度優(yōu)于線性模型.同理構(gòu)造Logistic模型和尖角突變模型的統(tǒng)計(jì)量F2′,通過計(jì)算F2′=0.287>0,說明尖角突變模型擬合度優(yōu)于Logistic模型.

3 模型分析

考慮尖角突變模型(4)，把由x′(t)=0確定的曲面稱為平衡曲面，即

x3+ax+b=0.

(12)

把x′(t)=0與x″(t)=3x2+a=0聯(lián)立消去x，可得到控制空間中的方程

4a3+27b2=0,

(13)

式(13)所確定的曲線稱為分歧點(diǎn)集.

方程(12)是一個(gè)3次方程，它或是有1個(gè)實(shí)根，或者有3個(gè)實(shí)根，實(shí)根數(shù)目由判別式Δ=4a3+27b2決定.如果Δ≤0，則有3個(gè)實(shí)根；否則只有1個(gè)實(shí)根.

Δ=0時(shí)，或者兩個(gè)實(shí)根相同a,b≠0，或者3個(gè)根全都相同a=b=0.

這說明當(dāng)Δ<0時(shí)，x(t)有兩個(gè)極小值，一個(gè)極大值，即兩個(gè)平衡點(diǎn).Δ>0時(shí)，x(t)只有一個(gè)平衡點(diǎn).

利用卡爾丹諾公式求解x3+ax+b=0得：

圖2 尖角突變分區(qū)仿真圖

根據(jù)各平衡點(diǎn)的Jacobi矩陣,可得Dxfxi,a,b=3xi2+ai=1,2,3，根據(jù)Dxfxi,a,b的符號(hào)判斷，當(dāng)3xi2+a<0時(shí)，平衡點(diǎn)穩(wěn)定；當(dāng)3xi2+a>0時(shí)，平衡點(diǎn)不穩(wěn)定.

利用Matlab繪制尖角突變分區(qū)仿真圖(見圖2).

由圖2可以看出:當(dāng)Δ<0時(shí),滿足條件的控制變量a,b值位于尖角區(qū)域內(nèi),處于潛在危害區(qū),這是一不穩(wěn)定區(qū)域,隨時(shí)可向兩極發(fā)展;當(dāng)Δ>0時(shí),a,b值所決定的點(diǎn)落在尖角區(qū)域外,分別為穩(wěn)定安全區(qū)和穩(wěn)定危害區(qū),其中,凡是控制變量(a,b)控制的點(diǎn)落在穩(wěn)定安全區(qū)內(nèi),它們都處于尖角區(qū)域的左側(cè);當(dāng)Δ=0時(shí),a=b,則x=0,這時(shí)對(duì)應(yīng)的正是尖角的地方.如果平衡點(diǎn)位于潛在危害區(qū),這表明蚜蟲種群可能發(fā)生突變,此時(shí)可采用生物防治的辦法抑制蚜蟲種群爆發(fā);如果平衡點(diǎn)位于穩(wěn)定安全區(qū),這表明該時(shí)間段蚜蟲種群沒有威脅到小麥生長(zhǎng),此時(shí)可做常規(guī)防治工作;如果平衡點(diǎn)位于穩(wěn)定危害區(qū),這表明蚜蟲種群處于大爆發(fā)期,已經(jīng)嚴(yán)重危害小麥生長(zhǎng),此時(shí)應(yīng)及時(shí)進(jìn)行化學(xué)防治.