江福才 范慶波 馬全黨 張 帆 馬 勇
(武漢理工大學航運學院1) 武漢 430063) (武漢理工大學內河航運技術湖北省重點實驗室2) 武漢 430063)
隨著水上交通在智能化、信息化等方面的發展,實時、準確的交通流量預測可為水上智能交通系統進行交通誘導及控制提供有效的基礎數據,便于海事主管機關制定合理、有效的調度、交通組織方案.同時港口部門可根據船舶交通流量的變化趨勢,合理規劃碼頭、航道、錨地等水域的布局,提高港口資源的利用效率[1].
目前國內外的學者在水上交通流預測方面取得了許多研究成果,但大多是借鑒陸地、航空等領域的預測方法,比較常用的方法有支持向量機(support vector machine,SVM)模型、灰色模型、時間序列模型、神經網路模型、組合模型等[2-5],這些模型有各自的優勢,但很少考慮序列的非平穩性,因此,預測精度相對較小.在進行月、季、年度交通流量的預測時,往往受到趨勢性、季節性等因素的影響,季節性自回歸移動平均(seasonal auto regressive integrated moving average,SARIMA)模型可以很好的解決這些問題.在SARIMA模型研究方面,李曉磊等[6]運用SARIMA模型對荊州港月度船舶交通流量進行了預測,并將該方法與Holt-Winters加法模型、自回歸移動平均(auto regressive integrated moving average,ARIMA)模型進行了對比,驗證了SARIMA模型的可行性;郝軍章等[7]將SARIMA模型運用到全國鐵路客運量預測中,通過將預測值與實際值進行對比,得到模型的預測能力為“優良”;王瑩等[8]選擇SARIMA模型對北京地鐵進站客流量進行時間序列建模,從預測結果看,該模型可以較準確的描述進站客流量的變化趨勢.
馬爾科夫(Markov)模型具有無后效性,廣泛應用到交通流量、瓦斯流量、交通事故等預測方面.Yang等[9]運用Markov模型與小波神經網絡模型結合,深入洞察樣本流量模式,提高預測的準確性與魯棒性,對德克薩斯州休斯敦的高速公路流量進行了預測;Lyu等[10]運用BP神經網絡與Markov模型組合,對防城港船舶交通流量進行了預測;陳釗等[11]運用GM(1,1)模型與Markov模型結合對礦井瓦斯流量進行了預測等等.因傳統的Markov模型在求取灰色區間的白化系數時存在不足,本文選用具備較強全局搜索能力的粒子群算法求解最佳的白化系數,在粒子群算法優化方面,樊冬雪[12]利用無偏灰色馬爾科夫模型對鐵路客流量進行預測,并利用粒子群算法(particle swarm optimization,PSO)輔助求解馬爾科夫模型中最佳的白化系數;Eseye[13]運用PSO進行全局搜索,對支持向量機中的參數進行優化,以提高預測精度等.
鑒于以上分析,本文將SARIMA-Markov模型引入船舶交通流量預測領域,發揮了SARIMA模型與Markov模型所具有的優勢,并且運用PSO進行優化,求取最優的白化系數,以期更好的表征船舶交通流量的發展趨勢.
ARIMA模型(B-J方法),常用于解決隨機時序問題[14].SARIMA模型是由ARIMA模型衍生而來的,基本形式為SARIMA(p,d,q)(P,D,Q)s.
設原始時間序列是它的當期和前期的隨機誤差項及前期值的線性函數,即可表示為
Xi=φ1xi-1+φ2xi-2+…+φpxi-p+
μi-θ1μi-1-θ2μi-2-…-θqμi-q
(1)
式中:時間序列Xi為自回歸移動平均序列;ui為相互獨立的白噪聲序列.上述為(p,q)階的自回歸移動平均模型,記為ARMA(p,q);φ1,φ2,…,φp為自回歸系數;θ1,θ2,…,θq為移動平均系數.
引入滯后算子,式(1)可簡記為
φ(B)xi=θ(B)μi
(2)
若序列在季節性差分前還進行了d階的逐期差分才平穩,則最終形成的SARIMA模型表達式為
(3)

軟件常被用于解決時間序列建模及數據處理問題[15].
Markov模型具有無后效性,對于一個具有無后效性的隨機序列{Y(t),t∈T},T為離散的時間集合,其無后效性可以表示為
P{Yn+1=in+1|Y1=i1,Y2=i2,…,Yn=in}=
P{Yn+1=in+1|Yn=in}
(4)
Markov預測模型的構建可以分為以下幾個步驟.
步驟1求取殘差序列e(i) 對原始數據Xi,運用SARIMA模型進行計算得到預測值fi,同時可以得到殘差序列e(i)
e(i)=xi-fi,i=0,1,2,…
(5)
步驟2狀態區間劃分 將殘差序列按照大小劃分為間隔寬度相等的r個狀態,殘差序列e(i)的第i步的第j個狀態的上下邊界分別用Qij和Hij表示,計算方法為
(6)
步驟3求取k步狀態轉移概率 設狀態空間為I={i1,i2,…,ir},k步狀態轉移概率為
{ii→i1,ii→i2,…,ii→ir}
(7)
步驟4構建狀態轉移概率矩陣 構建能反映系統間各狀態概率轉移的一步轉移矩陣,其大小為r·r.該轉移矩陣描述為
(8)
步驟5求取預測值 Markov模型的預測值表示為
(9)

選用粒子群算法進行全局搜索來確定Markov算法中的白化系數λi(i=1,2,…,r)[16-17].
對于某個粒子i,其粒子群搜索公式為
vi(t+1)=wvi(t)+c1rand1(xpbesti-xi(t))+
c2rand2(xgbesti-xi(t))
(10)
xi(t+1)=xi(t)+vi(t+1)
式中:xi為在t時刻的位置;vi為粒子在t時刻的速度;xpbest為群體微粒搜索到的自身最佳位置;xgbest為群體微粒經歷過的全局最好位置;rand1,rand2為服從U(0,1)的隨機數;w為慣性系數;c1,c2為加速度系數.
用均方差fit來衡量每個粒子的適應度函數
(11)
選擇軟件編寫相關程序求得最佳的白化系數λi,最終得到狀態區間中可能的預測值zi為
zi=(1-λi)Qij+λiHijλi∈[0,1]
(12)
為了更好的表征船舶交通狀況,實現對不同水域相同時段或同一水域不同時段的船舶交通流量的可比性,提出了標準船舶的概念[18].我國在評估海上交通安全時提出的標準船舶換算系數見表1.

表1 標準船舶換算系數表
運用以上換算系數對原始數據的船舶交通流量進行加權,轉化為標準船舶,原始數據為從海事相關部門調研得到的赤壁長江公路大橋斷面2010年1月—2016年6月的船舶交通流量數據,經過標準船舶轉化后的樣本數據見表2.

表2 赤壁長江公路大橋斷面標準船舶交通流量統計表
為驗證SARIMA-Markov模型的適用性,運用、軟件,選取赤壁長江公路大橋斷面2010年1月—2015年12月(編號1~36)的樣本數據用于建立模型,2016年1—6月(編號37~42)的樣本數據用于模型驗證,并運用相關評價指標對模型進行評價.
2.2.1SARIMA模型
依據前面介紹的建模步驟,首先對序列進行平穩性處理,見圖1.各觀測值沒有完全圍繞其均值上下波動,曲線前后波動具有不一致性;表3為原始數據、一階差分后ACF-PACF分析圖,由圖3可知,原始數據的自相關系數(ACF)起初都落在顯著水平α=0.05的置信區間以外,綜上所述,得到序列xi為非平穩時間序列.
對非平穩序列通常采用逐期差分的方法,ACF很快的趨于0(滯后階數k>2或3時),落入置信區間內;對序列xi進行單位根檢驗,結果見表4.顯著水平為1%的臨界值大于t統計值,得到序列xi轉化為平穩序列,因此d=1.

圖1 時間序列圖(2010年1月—2015年12月)
由表3可知,ACF在12階,24階等以12為周期的階數位置存在較大的峰值,表明序列存在以12個月為周期的季節波動,為消除模型的季節性,本文應用軟件中的seasonal adjustment-census X12功能.參照原始序列平穩性的分析步驟,得到D=1.
ACF-PACF分析圖常被用來對SARIMA(p,d,q)(P,D,Q)s模型進行識別和定階.根據以上分析,得到d=1,D=1,見表5,由表5可知,ACF在k=12,24處顯著超出置信區間,存在以12為周期的波動性,得到P=1、Q=1.參照相關文獻、書籍對p、q值的取值方法[19-20],對比表3、表5可知,得到可供選擇的(p、q)組合有(3,1)、(4,0)、(2,1)、(3,0).
表3原始數據、一階差分后ACF-PACF分析圖

表4 一階差分后ADF檢驗圖

根據不同的(p、q)組合建立SARIMA(3,1,1)(1,1,1)12,SARIMA(4,1,0)(1,1,1)12,SARIMA(2,1,1)(1,1,1)12,SARIMA(3,1,0)(1,1,1)12四組模型,運用軟件中的Equation Estimation功能對模型進行檢驗,得到選取SARIMA(3,1,0)(1,1,1)12是比較合適的.
SARIMA(3,1,0)(1,1,1)12模型表達式為
(1+B12)(1+0.621B+0.376B2+0.138B3)
(1-B)(1-B12)lg(xi)=(1-0.009B12})μi
表5一階季節差分后ACF-PACF分析圖

運用SARIMA(3,1,0)(1,1,1)12模型對2013年1月—2015年12月的月度標準船舶交通流量進行預測,結果見表6.
2.2.2SARIMA-Markov模型
根據Markov模型的建模步驟,首先確定殘差序列的范圍e(i)∈[-1 455,1 442],將殘差序列劃分為4個狀態區間,區間的上下邊界值為



表6 SARIMA模型預測值及殘差序列統計表 艘·次
同理,E2、E3、E4狀態的上下邊界區間分別為[-730.75,-6.5],[-6.5,717.75],[717.75,1 442].殘差序列狀態分布見表7.

表7 狀態分布表
根據表7中狀態劃分情況,對從Ei經一步轉移到Ej的頻數統計見表8.

表8 一步狀態轉移頻數統計表
狀態轉移概率:P11=P(E1→E1)=P(E1|E1)=0.4
同理可得到:P12=0.3,P13=0.2,P14=0.1,P21=P22=0.364,P23=0.091,P24=0.182,P31=0.2,P32=0.2,P33=0.3,P34=0.1,P43=0.75,P44=0.25.
狀態轉移概率矩為
運用編制粒子群算法程序,進行全局搜索,以找到最佳的白化系數,計算結果為
λ1=0.015 4,λ2=0.142 0,
λ3=0.327 8,λ4=0.988 4
根據式(12),得到:z1=-741.90、z2=-109.34、z3=480.34、z4=726.15.根據式(9)計算得到SARIMA-Markov模型的預測值,見表9和圖2.

表9 SARIMA-Markov模型預測值及殘差序列統計表 艘·次

圖2 SARIMA-Markov模型擬合示意圖
由圖2可知,觀察兩種模型的擬合曲線,SARIMA-Markov模型的擬合效果要優于SARIMA模型.
2.2.3模型驗證
為對所建模型的預測精度進行驗證,擬運用構建好的模型對赤壁長江公路大橋斷面2016年1—6月的月度標準船舶交通流量進行預測.
以第37、38組數據為例,運用軟件,計算SARIMA(3,1,0)(1,1,1)12模型2016年1月、2月的預測值為6 015,4 213/(艘·次),然后根據式(8)~(9)計算得到第36組數據轉移到37、38組數據的1步、2步狀態轉移概率矩陣,結果為
第36組數據處于E4,得到第37組數據的預測值為:5 800/(艘次),38組數據的預測值為:6 009/艘次.
同理可得到后面幾個月的預測值,結果見表10.

表10 預測結果統計表 艘·次
本文擬采用平均絕對誤差、平均絕對百分誤差對模型擬合效果進行評價,見表11.優化后模型的兩項指標均小于優化前,SARIMA-Markov模型的擬合精度達到了92.0849%,優于SARIMA模型84.755 2%的擬合精度.SARIMA-Markov模型預測精度為95.786 1%,優于單純利用SARIMA模型87.445 4%的預測精度,SARIMA-Markov模型預測精度較高,滿足要求.

表11 模型檢驗表
1) 將SARIMA-Markov模型應用到月度船舶交通流量預測,樣本數據為經過加權后的標準船舶交通流量.該模型綜合考慮了SARIMA模型對月度船舶交通流量季節相關性的表達,Markov模型運用殘差序列、狀態劃分、狀態轉移矩陣等方法處理樣本近期狀況對預測值的影響,并運用具有全局搜索能力的粒子群算法求取模型中的最佳白化系數,為船舶交通流量預測提供了新方法.
2) 本文模型與單純利用SARIMA模型相比,能確保更高精度的預測結果,驗證了該模型的科學性、可行性.
3) 對于該模型的可推廣性及最佳狀態劃分方案等,在后續的研究工作仍有待進一步完善.
4) 預測結果可為港口水域的合理布局、水上智能交通系統、海事主管機關等提供基礎性數據,該方法對受季節性影響的港口、航道等水域船舶交通流量預測也有很好的借鑒意義.