馬永剛,張啟敏,劉俊梅
(1.寧夏大學數學統計學院,寧夏銀川 750021)
(2.榆林學院數學與統計學院,陜西榆林 719000)
隨機微分方程是近幾年熱門的數學學科,它是微分方程、動力系統及隨機分析相互交叉的學科.在實際應用中,隨機微分方程的參數一般是未知的,因此參數估計問題成為一個重要的研究方向.關于此類方程參數估計問題已有很多結果,見文獻[1–3,9]等.在文獻[1]中,Bishwal詳細介紹了參數估計的理論方法和技巧,并運用到各種隨機模型中.在文獻[3]中,作者基于采樣的時間序列討論了非線性隨機微分方程的參數估計問題.
1926年意大利數學家Volterra提出著名的Lotka-Volterra微分方程模型,該模型引起眾多生物學家、數學家極大的興趣.Lotka-Volterra模型所表現的生態學現象,在自然界中處處可見.例如,病蟲害的周期爆發,農作物的豐收年與低產年的周期循環,動物之間的捕食與被捕食等,這些現象都是Lotka-Volterra模型變化規律的具體展示.Lotka-Volterra模型的一般形式為[4]

其中αi,βi,γi(i=1,2)均為常數,α1與α2分別表示兩種群的內稟增長率,β1與γ2分別表示種內作用系數,γ1與β2分別表示種間作用系數.就其生態意義可分為三類:競爭模型(αi> 0,βi< 0,γi< 0,i=1,2)、互惠模型 (αi> 0,βi< 0,γi> 0,i=1,2)、捕食 -被捕食模型(αi> 0,βi< 0,γ1< 0,γ2> 0,i=1,2).基于此模型,各種形式的Lotka-Volterra生態模型被提出,其中一類是受外界隨機因素影響的Lotka-Volterra生態模型.下面研究白噪聲擾動的Lotka-Volterra競爭模型[5]

其中r1,r2,a11,a12,a21,a22均為正常數,ω1(t),ω2(t)為白噪聲且相互獨立.許多學者對此模型進行了研究,見文獻[4,5,6]等.其中最重要的結果之一,種群隨機持久生存,需要滿足下面條件

其中

當上述條件變為

兩種群中至少有一個種群隨機滅絕.
本文主要研究隨機種群Lotka-Volterra競爭模型的參數估計,即對方程組(1.2)中參數r1,a11,a12,r2,a21,a22,σ1,σ2應用最小二乘法理論進行估計,得到較好的估計結果.主要結果是參數的點估計值和區間估計,同時獲得影響估計區間長度的主要因素.最后給出了數值模擬,結果表明了該方法的可行性與有效性.

其中方程組(2.1)可轉化為

其中

進一步方程組(2.2)表示為

方程組(2.3)是簡單的線性回歸模型,應用線性回歸理論的方法估計參數,估計過程基于最小二乘法,使下面兩個式子的值越小越好

為計算方便,方程(2.3)可寫成分塊矩陣的形式

其中

在塊矩陣中,子矩陣分別為

由多元線性回歸理論的公式估計參數β和η:

記

故(3.1)式可進一步化簡為

為了獲得參數的置信區間,需要估計參數β和η的方差,方差估計公式為


估計得到,p為被估計參數個數.由(3.1),(3.3)式可得

等式(3.4),(3.5)化簡得

定理3.1方程組(3.6)中分別為方程組(2.1)中參數的漸進無偏估計,即當n→∞,

證類似文獻[7]證明.利用分別估計(3.2)式中參數得到被估參數的方差

如果σ2已知,由最小二乘法回歸理論[8],參數估計值?β,?η的各分量都滿足正態分布.當觀察值的數量n足夠大時,可由?σ2代替σ2,因此(1?α)的置信區間(CIs)分別為

下面給出具體數據,進行參數的估計.
例1在模型(1.2)式中,我們選擇滿足種群持久生存的數據進行模擬,即兩種群隨機持久存在的情形,其它狀態的數據也可獲得類似結論.假設參數真實值分別為r1=1,a11=0.3,a12=0.2,σ1=0.04,r2=1.2,a21=0.3,a22=0.4,σ2=0.05, 初始值為x1(0)=1.5,x2(0)=2.應用這些參數的值,通過Euler-Maruyama方法離散,我們獲得三組模擬數值xk1(t),xk2(t),數量分別為(A)n=5000;(B)n=10000;(C)n=50000;保存這些數據為模擬數據集.對每一組數據,給步長?t=0.02,0.05,0.1進行模擬,給出一些模擬結果,以便比較真實值與估計值.在表1–4中,樣本的大小以符號“Size n”表示,且在表中第一列給出.其中表1,2樣本大小從5000增加到50000,在每一組模擬數據中,分別以三種不同步長進行模擬.模擬數據表明參數的估計值與真實值沒有明顯區別,絕對誤差與樣本大小有關,與步長?t無關,且隨著樣本數量的增大,絕對誤差越來越小.

表1:(r1,a11,a12)的估計模擬結果(真實值:r1=1,a11=0.3,a12=0.2)

表2:(r2,a21,a22)的估計模擬結果(真實值:r2=1,a21=0.3,a22=0.2)
其次給出各個參數置信水平為0.95的置信區間的一些模擬結果.表3,4中,給出參數的估計區間及區間長度,模擬數據表明隨著樣本量的增大,置信水平長度越來越小.

表3:置信水平為0.95的參數r1,a11,a12置信區間模擬結果(真實值:r1=1,a11=0.3,a12=0.2)

表4:置信水平為0.95的參數r2,a21,a22置信區間模擬結果(真實值:r2=1,a21=0.3,a22=0.2)
本文對隨機兩種群Lotka-Volterra競爭模型應用最小二乘法理論進行參數估計,獲得參數r1,a11,a12,r2,a21,a22,σ1,σ2的估計值及估計區間.通過觀察結果,得到影響置信區間長度的主要因素.結果表明,參數估計區間的長度隨著樣本數量n的增加而減小,不依賴步長?t的長度.例1給出了具體的數值模擬.極大似然估計與貝葉斯估計方法是另外兩種參數估計的方法,在將來的工作中,將應用這些方法到隨機兩種群Lotka-Volterra模型.
[1]Bishwal JPN.Parameter estimation in stochastic differential equations[M].Berlin:Springer,2008.
[2]Kristensen NR,Madsen H,Young PC.Parameter estimation in stochastic grey-box model[J].Automatica,2004,40(2):225–237.
[3]Timmer J.Parameter estimation in nonlinear stochastic differential equatuons[J].Chaos Sol.Fract.,2000,11(15):2571–2578.
[4]Zhu C,Yin G.On competitive Lotka-Volterra model in random environments[J].J.Math.Anal.Appl.,2009,357(1):154–170.
[5]Jiang D Q,Ji C Y,Li X,O’Regan D.Analysis of autonomous Lotka-Volterra competition systems with random perturbation[J].J.Math.Anal.Appl.,2012,390:582–595.
[6]王佳,丁潔麗.Logistic回歸模型中參數極大似然估計的二次下界算法及其應用[J].數學雜志,2015,35(6):1521–1532
[7]Pan J F,Alison G,David G,Mao X R.Parameter estimation for the stochastic SIS epidemic model[J].Stat.Infer.Stoch Proc.,2014,17(1):75–98.
[8]Rawlings JO.Applied regression analysis:a research tool[M].Belmont,CA:Wadsworth,1988.
[9]洪志敏,閆在在.Volterra積分方程的蒙特卡羅數值求解方法[J].數學雜志,2016,36(2):425–436.