薛 晗，邵哲平，潘家財，張 鋒

(集美大學(xué) 航海學(xué)院， 福建 廈門 361021)

由于海上交通越來越繁忙,船舶交通流量快速増長,導(dǎo)致海上交通事故頻繁發(fā)生,造成巨大的經(jīng)濟(jì)損失.船舶交通流量是海上交通工程學(xué)中一個重要的基本量，也是衡量海上交通基礎(chǔ)設(shè)施建設(shè)的一個重要指標(biāo),其預(yù)測結(jié)果能為相關(guān)主管部門制定科學(xué)航道管理規(guī)劃和做好船舶通航管理等工作提供基礎(chǔ)性依據(jù).因此，保障船舶交通流量預(yù)測的準(zhǔn)確性和合理性對于改善航道基礎(chǔ)設(shè)施建設(shè)與制定科學(xué)的航道管理策略等具有重要的意義.

船舶交通流量預(yù)測研究采用多種先進(jìn)人工智能優(yōu)化算法的預(yù)測模型,如人工神經(jīng)網(wǎng)絡(luò)[1]、極限學(xué)習(xí)機(jī)[2]、模糊控制[3]、自回歸移動平均模型[4]、低秩稀疏分解[5]、支持向量機(jī)[6]、支持向量機(jī)回歸[7]、遺傳算法[8]及灰色理論[9]等人工智能算法.由于能處理復(fù)雜非線性問題，神經(jīng)網(wǎng)絡(luò)取得了一定的成果.李俊等[10]采用遺傳算法改善神經(jīng)網(wǎng)絡(luò)進(jìn)行了船舶交通流量的預(yù)測.鈕浩東等[11]采用果蠅優(yōu)化算法優(yōu)化廣義回歸神經(jīng)網(wǎng)絡(luò)(GRNN)進(jìn)行了船舶交通流量預(yù)測.

然而，神經(jīng)網(wǎng)絡(luò)本身存在一些不足，例如學(xué)習(xí)速度慢、容易陷入局部極值、學(xué)習(xí)及記憶不穩(wěn)定性.文化算法(CA)是一種包含信仰空間與種群空間的雙層進(jìn)化系統(tǒng)，種群空間的演化可以由信仰空間中保存的知識引導(dǎo)[12].螢火蟲算法(FA)是一種模擬自然界中成蟲發(fā)光的生物學(xué)特性優(yōu)化算法[13-16].螢火蟲算法與其他算法結(jié)合，能改善螢火蟲算法的性能.Rizk-Allah等[17]將蟻群算法和螢火蟲算法結(jié)合起來求解無約束優(yōu)化問題，該算法綜合了蟻群算法和螢火蟲算法的優(yōu)點，螞蟻的演化是通過整合蟻群算法和螢火蟲算法來完成的，螢火蟲算法作為局部搜索來改進(jìn)螞蟻發(fā)現(xiàn)的位置.Rahmani等[18]提出了一種混合進(jìn)化螢火蟲遺傳算法，受螢火蟲社會行為和生物發(fā)光通信現(xiàn)象的啟發(fā)，將離散螢火蟲算法與標(biāo)準(zhǔn)遺傳算法相結(jié)合，提高了新算法的適應(yīng)性，并應(yīng)用于求解容量化設(shè)施選址問題.Kora等[19]提出了一種使用混合螢火蟲和粒子群優(yōu)化技術(shù)，結(jié)合神經(jīng)網(wǎng)絡(luò)分類器來檢測束分支塊的技術(shù)，新算法的光強(qiáng)吸引力是模擬最優(yōu)解的主要控制力之一，光強(qiáng)吸引過程取決于搜索的隨機(jī)方向，局部搜索是通過粒子群優(yōu)化算法(PSO)操作修改后的光強(qiáng)吸引步驟執(zhí)行.本文設(shè)計一種新的文化螢火蟲算法，證明了算法的收斂性.利用新算法調(diào)整廣義回歸神經(jīng)網(wǎng)絡(luò)(GRNN)的權(quán)值和閾值，提高神經(jīng)網(wǎng)路的學(xué)習(xí)速度和性能，并將其應(yīng)用于船舶交通流量的預(yù)測.

1 基于自動識別系統(tǒng)的航道交通流量統(tǒng)計

1.1 船舶交通流量模型

船舶交通流量的計算公式為

δ=ρvw

(1)

式中：ρ為船舶交通流密度；v為船舶交通流速度；w為船舶交通流寬度.

采用自動識別系統(tǒng)(AIS)進(jìn)行船舶交通流量的測量統(tǒng)計，記錄某一段時間內(nèi)船舶通過航道某一斷面的交通流量.如圖1所示，P1P2為觀測斷面，Q1Q2為船舶相鄰兩次AIS采集的位置軌跡連線.

圖1 航行船舶航跡與航道觀察線Fig.1 Ship trajectory and channel observation line

兩條線段的位置關(guān)系有重合，相交及相離.按以下算法進(jìn)行判斷圖1中船舶航跡Q1Q2與觀察線P1P2是否相交，如圖2所示.

圖2 兩條線段相交關(guān)系判別Fig.2 Judge the intersection relationship of two line segments

(1) 快速排斥試驗：設(shè)以線段P1P2為對角線的矩形為R， 設(shè)以線段Q1Q2為對角線的矩形為T，如果R和T不相交，顯然兩線段不會相交.

(2) 跨立試驗：如果兩線段相交，則兩線段必然相互跨立對方.若P1P2跨立Q1Q2，則P1和P2兩點位于Q2和Q1的所在直線的兩側(cè)，用rxy表示x,y兩點間的矢量，則：

(rQ1P1×rQ1Q2)·(rQ1Q2×rQ1P2)>0

(2)

同理，判斷Q1Q2跨立P1P2的依據(jù)為

(rP1Q1×rP1P2)·(rP1P2×rP1Q2)>0

(3)

1.2 經(jīng)緯度與直角坐標(biāo)轉(zhuǎn)換

原始AIS數(shù)據(jù)只提供了船舶離散的地理坐標(biāo)點位置，屬于經(jīng)緯度類型的數(shù)據(jù)，而船舶的軌跡擬合基于船舶整體航行軌跡上進(jìn)行分析處理.在軌跡擬合的過程中，需要將離散的軌跡點在擬合程序中有序重組，轉(zhuǎn)化為船舶航行軌跡線，在此過程中需要使用的數(shù)據(jù)是平面直角坐標(biāo)系，所以在進(jìn)行擬合之前需要將從AIS系統(tǒng)中獲取的經(jīng)緯度數(shù)據(jù)轉(zhuǎn)換為墨卡托直角坐標(biāo)系數(shù)據(jù).

將經(jīng)緯度數(shù)據(jù)轉(zhuǎn)換為墨卡托平面坐標(biāo)系數(shù)據(jù)的方法如下[20].基準(zhǔn)緯度處橢球的卯酉圈曲率半徑

(4)

式中：a為地球橢球長半徑；e為橢球體第一偏心率；φ0為墨卡托投影的基準(zhǔn)緯度.基準(zhǔn)緯度的緯圈半徑為

r0=N0cosφ0

(5)

等量緯度為

(6)

式中：φ為緯度.墨卡托直角坐標(biāo)的橫坐標(biāo)為

xn=r0λ

(7)

λ為經(jīng)度.墨卡托直角坐標(biāo)的縱坐標(biāo)為

yn=r0q

(8)

2 廣義回歸神經(jīng)網(wǎng)絡(luò)基礎(chǔ)

(9)

式中：k為訓(xùn)練樣本數(shù)；Yi是y的第i個訓(xùn)練樣本值；Di為X與Xi間歐氏距離；σ為平滑參數(shù).模式層各單元傳遞函數(shù)為

(10)

求和層一個單元的傳遞函數(shù)為

(11)

其他單元傳遞函數(shù)為

(12)

式中：Yij式中為Yi的第j個元素.樣本預(yù)測值的第j個元素，即輸出層單元yj為

(13)

GRNN的結(jié)構(gòu)如圖3所示.

圖3 GRNN結(jié)構(gòu)Fig.3 Structure of GRNN

3 文化螢火蟲算法

結(jié)合文化算法與螢火蟲算法的優(yōu)點, 提出一種基于文化算法框架的螢火蟲優(yōu)化算法，即CFA.

3.1 CFA的模型

螢火蟲個體i相對個體j的熒光亮度為

Ii=I0ie-γ rij

(14)

式中：I0i為個體i的最大熒光亮度，與其目標(biāo)函數(shù)相關(guān)，目標(biāo)函數(shù)越優(yōu)則亮度越高；γ為光強(qiáng)吸收系數(shù)；rij為螢火蟲個體i和j之間的距離，

(15)

式中：d為問題的維數(shù)；xik為第i個螢火蟲的第k維分量.螢火蟲個體i和j之間的吸引力為

(16)

式中：β0為個體最大吸引力.若個體j的亮度高于個體i的亮度，即Ij>Ii，則個體j吸引個體i向自己方向移動，其位置x更新如下式：

(17)

式中：t為時間；α為步長；rand(0,1)為[0,1]之間均勻分布的隨機(jī)數(shù).

信仰空間由在進(jìn)化過程中獲取的經(jīng)驗和知識組成，根據(jù)現(xiàn)有的種群經(jīng)驗如最好個體以及新個體經(jīng)驗進(jìn)行更新.用當(dāng)前最優(yōu)個體xbest來更新信仰空間Belief的知識b，如下式：

(18)

信仰空間Belief通過下式規(guī)則來影響種群空間：

xik(t+1)=

(19)

式中：xij(t)為第t次迭代中第i個個體中第k維分量；ei為Belief中變量i可調(diào)整區(qū)間的長度；N(0,1)為服從標(biāo)準(zhǔn)正態(tài)分布的隨機(jī)數(shù).

3.2 CFA流程

CFA的偽代碼如下所示：

(1) CFA程序;

(2) 開始;

(3)t∶= 0;

(4) 初始化算法參數(shù)，例如β0，γ，α;

(5) 初始化種群空間，例如螢火蟲位置;

(6) 計算每個螢火蟲的適應(yīng)度，作為其最大亮度;

(7) 初始化信仰空間;

(8) 根據(jù)式(14)計算螢火蟲的相對亮度;

(9) 根據(jù)式(16)計算螢火蟲的吸引力;

(10) 根據(jù)式(4)更新螢火蟲的位置;

(11) 重新計算螢火蟲的亮度;

(12) 根據(jù)式(18)更新信仰空間;

(13) 根據(jù)式(19)影響種群空間;

(14)t∶=t+1;

(15) 當(dāng)t

(16) 輸出全局最優(yōu)值和最優(yōu)解;

(17) 結(jié)束CFA.

3.3 CFA的收斂性分析

先給出隨機(jī)優(yōu)化算法收斂于全局最優(yōu)的判據(jù)[22]，再證明CFA的收斂性.設(shè)A為優(yōu)化問題的可行解空間，f為目標(biāo)函數(shù)，ξ為算法D在迭代中搜索過的最優(yōu)解，第k次迭代結(jié)果生成的解為zk.

假設(shè)1若f(D(zk,ξ))≤f(zk)，且ξ∈A，則f(D(zk,ξ))≤f(ξ).

假設(shè)2對于A的任意Borei子集B，若其勒貝格測度ν(B)>0，μk(B)為由測度μk獲得B的概率，則有：

(20)

引理1假設(shè)目標(biāo)函數(shù)f為可測函數(shù)，可行解空間A為Rn上的可測子集，R為實數(shù)集，Rε為全局最優(yōu)解集，P(zk∈Rε)表示第k步算法生成的解zk屬于最優(yōu)解Rε的概率，若假設(shè)1和2成立，則有：

(21)

引理2有限齊次Markov鏈從任意非常返狀態(tài)出發(fā)依概率1要達(dá)到常返狀態(tài).

定理1FA所表示的Markov鏈?zhǔn)怯邢迺r齊的.

證明螢火蟲個體位置更新的轉(zhuǎn)移概率P為

Pi j{xt+1=j|xt=i,xt-1=it-1,…,x0=i0}=

P{xt+1=j|xt=i}

所以FA所表示的Markov鏈?zhǔn)怯邢迺r齊的.

定理2CFA中Belief的狀態(tài)的隨機(jī)變化屬于有限齊次Markov鏈.

證明CFA中Belief是用來記錄搜索過程中的進(jìn)化知識，CFA對目標(biāo)最優(yōu)值的搜索是在離散的、有限的空間中進(jìn)行的，因此用來記錄的Belief也是有限的.同時Belief每一次狀態(tài)轉(zhuǎn)移概率只和其當(dāng)前狀態(tài)有關(guān)，而與t無關(guān).因此，CFA的Belief的狀態(tài)的隨機(jī)變化屬于有限齊次Markov鏈.

定理3CFA依概率1全局收斂.

證明設(shè)信念空間Belief的全局知識的狀態(tài)空間為Ωg，局部狀態(tài)空間為Ωl，整個信念空間狀態(tài)空間為Π={(Ig,Is)}，其中Ig表示全局知識集，Is表示螢火蟲群體空間的局部集.

Π分為兩個子空間，即為常返狀態(tài)空間Π1和非常返狀態(tài)空間Π2，其中常返狀態(tài)空間可表示為

Π1={(Ig,Is)|Ig∈Ωg,Is∈Ωl}

則有Π2∩Π2=?，Π2∪Π2?Π.

若存在某知識狀態(tài)s∈Π1，則由全局知識的定義可知，將不再轉(zhuǎn)移到Π2，所以Π1為閉集.Π1中的狀態(tài)是相通的，所以Π1中狀態(tài)為常返狀態(tài).若s∈Π2，由于Π2中的狀態(tài)為非常返狀態(tài)，進(jìn)化過程中每一個體將在整個解空間范圍內(nèi)運動，s由Π2轉(zhuǎn)向Π1的轉(zhuǎn)移概率大于0，故進(jìn)化后的群體收斂于全局最優(yōu)解的概率大于0.由引理2知，信念空間知識的進(jìn)化由任意非常返狀態(tài)必然依概率1轉(zhuǎn)移到常返狀態(tài)，即依概率1收斂到全局最優(yōu)解.

3.4 基于CFA的神經(jīng)網(wǎng)絡(luò)預(yù)測船舶交通流量

CFA以神經(jīng)網(wǎng)絡(luò)的輸出均方差為適應(yīng)度函數(shù)，以神經(jīng)網(wǎng)絡(luò)的輸入層和隱含層中的權(quán)值、隱含層和輸出層中的權(quán)值、隱含層的閾值及輸出層的閾值為編碼，對神經(jīng)網(wǎng)絡(luò)進(jìn)行優(yōu)化，進(jìn)化目標(biāo)是得到最合適、最優(yōu)的神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu).神經(jīng)網(wǎng)絡(luò)用于預(yù)測前，需先進(jìn)行訓(xùn)練.預(yù)測流程如下：

步驟1初始化神經(jīng)網(wǎng)絡(luò)模型；

步驟2輸入訓(xùn)練樣本；

步驟3設(shè)置CFA參數(shù)和初始種群；

步驟4用CFA計算神經(jīng)網(wǎng)絡(luò)的最優(yōu)權(quán)值和閾值；

步驟5計算神經(jīng)網(wǎng)絡(luò)隱藏層；

步驟6計算神經(jīng)網(wǎng)絡(luò)輸出層；

步驟7計算神經(jīng)網(wǎng)絡(luò)的輸出誤差；

步驟8若誤差不滿足學(xué)習(xí)要求，轉(zhuǎn)步驟 4，直至訓(xùn)練結(jié)果達(dá)到要求為止；

步驟9用求得的最優(yōu)權(quán)值與閥值初始化神經(jīng)網(wǎng)絡(luò)；

步驟10輸入預(yù)測數(shù)據(jù)進(jìn)行預(yù)測，輸出預(yù)測結(jié)果.

4 實例與分析

4.1 實驗數(shù)據(jù)介紹

本文采用的AIS數(shù)據(jù)為舟山螺頭航道2016年1月的船舶交通流量數(shù)據(jù).下列是預(yù)測過程使用的數(shù)據(jù)，并將原始數(shù)據(jù)轉(zhuǎn)換為墨卡托直角坐標(biāo)系數(shù)據(jù).水上移動通信業(yè)務(wù)標(biāo)識碼(MMSI)為232971000的船在2016年1月28日的軌跡數(shù)據(jù)如表1所示.

我國學(xué)者提出考慮船舶尺寸大小和船舶總噸位的船舶換算系數(shù)，如表2所示，該系數(shù)能夠更準(zhǔn)確地反映水域船舶交通狀況和繁忙程度.

表1 船舶AIS數(shù)據(jù)與轉(zhuǎn)換結(jié)果Tab.1 Ship AIS data and conversion results

表2 船舶大小和船舶總噸位的船舶換算系數(shù)Tab.2 Conversion coefficient of ship size and gross tonnage

舟山螺頭航道2016年1月1日至6日的船舶交通流量按每小時統(tǒng)計，如表3所示.

表3 舟山螺頭航道船舶交通流量Tab.3 Ship traffic flow in Zhoushan Luotou Channel

圖4 CFA優(yōu)化的GRNN結(jié)構(gòu)收斂曲線Fig.4 Convergence curve of GRNN structure optimized by CFA

4.2 神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)優(yōu)化

用CFA優(yōu)化GRNN的結(jié)構(gòu)，進(jìn)化收斂曲線如圖4所示，圖中n為進(jìn)化代數(shù)，δ為適應(yīng)度.CFA種群規(guī)模設(shè)置為30，最大進(jìn)化迭代次數(shù)為30.

優(yōu)化后的神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)參數(shù)如表4所示.

表4 優(yōu)化后的神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)參數(shù)Tab.4 Optimized neural network structure parameters

4.3 船舶交通流量預(yù)測

本文采用前4 h的船舶交通流量數(shù)據(jù)，預(yù)測第5 h交通流量，如用2016年1月1日0∶00至4∶00的交通流量作為預(yù)測2016年1月1日4∶00至5∶00的船舶交通流量的訓(xùn)練樣本數(shù)據(jù)，以此類推，實現(xiàn)對舟山螺頭航道船舶交通流量的預(yù)測.表5為所預(yù)測的2016年1月7日0∶00至5∶00的交通流量.圖5為船舶交通流量預(yù)測結(jié)果曲線.

表5 優(yōu)化后的神經(jīng)網(wǎng)絡(luò)預(yù)測結(jié)果Tab.5 Prediction results of optimized neural network

圖5 船舶交通流量預(yù)測結(jié)果曲線Fig.5 Prediction result curve of ship traffic flow

表6為不同算法預(yù)測2016年1月1日4：00至5：00的船舶交通流量的預(yù)測結(jié)果比較，算法分別為支持向量機(jī)(SVR)，模糊控制(Fuzzy)，不經(jīng)過優(yōu)化的原始GRNN，只用FA優(yōu)化的GRNN(FA-GRNN)和用CFA優(yōu)化的GRNN(CFA-GRNN).

結(jié)果表明，CFA優(yōu)化神經(jīng)網(wǎng)絡(luò)后，相比支持向量機(jī)， 模糊控制， GRNN和FA-GRNN都能提高神經(jīng)網(wǎng)絡(luò)學(xué)習(xí)速度和預(yù)測準(zhǔn)確度，具有泛化性能好、不易陷入局部最優(yōu)等優(yōu)勢.因此CFA優(yōu)化神經(jīng)網(wǎng)絡(luò)能夠更好地預(yù)測海域的船舶交通流量.

4.4 潮汐變化因素分析

螺頭水道的船舶交通密度的變化受潮汐變化，晝夜變化的影響明顯，呈現(xiàn)出周期性的變化.舟山的潮汐情況如表7所示.

表7 舟山潮汐表Tab.7 Zhoushan tide

用高峰前10 h的船舶交通流量數(shù)預(yù)測高峰時段的船舶流量，如預(yù)測2016年1月10日 5∶00 至6∶00 的交通流量，2016年1月11日 6∶00 至 7∶00 的交通流量，2016年1月12日 7∶00 至 8∶00 的交通流量和2016年1月13日 8∶00 至 9∶00 的交通流量.預(yù)測結(jié)果如表8所示.

上表的結(jié)果進(jìn)一步檢驗了模型和算法的有效性和正確性.

4.5 船舶類型因素分析

螺頭水道通航的大型船舶和小型船舶航行規(guī)律明顯不同，對交通環(huán)境的影響不同，交通組織優(yōu)化的策略也不同.船長不小于200 m的船舶歸為大型船舶，船長小于120 m的船舶歸為小型船舶.流經(jīng)舟山螺頭水道的不同類型船舶舉例如表9所示，具體信息包括MMSI、船長、國際海事組織(IMO)、呼號和船籍國.

表8 優(yōu)化后的神經(jīng)網(wǎng)絡(luò)預(yù)測結(jié)果Tab.8 Prediction results of optimized neural network

表9 舟山螺頭航道不同類型船舶

Tab.9 Different types of ships in Zhoushan Luotou Channel

類型MMSI船長/mIMO呼號船籍國大船5667960003699631955S6LT9新加坡4772744003669467299VRJS5香港56574700033793425169VBN2新加坡小船41359113059703936YYYY中國413508170899359234BMFA中國6360099481139065429ELQE5中國

對大型船舶和小型船舶分組進(jìn)行預(yù)測.采用前10 h的船舶交通流量數(shù)據(jù)，預(yù)測接下來第11 h交通流量，如用2016年1月1日 0∶00 至 10∶00 的交通流量作為預(yù)測2016年1月1日 10∶00 至 11∶00 的船舶交通流量的訓(xùn)練樣本數(shù)據(jù)，以此類推，來實現(xiàn)對舟山螺頭航道船舶交通流量的預(yù)測.表10為所預(yù)測的2016年1月7日 0∶00 至 9∶00 的大型船舶和小型船舶的交通流量結(jié)果.

表10 不同類型船舶的預(yù)測結(jié)果Tab.10 Prediction results of different types of ships

從上表可以看出，總體而言，漲潮時較多船舶進(jìn)港，退流時較多船舶出港，通航密度大.大船的交通流量較穩(wěn)定，交通流量隨時間變化相對不明顯；小船的交通流相對較為復(fù)雜，穿越航道的時間段也不是很固定，交通流量隨時間變化較大.

5 結(jié)論

(1) 設(shè)計了新的文化螢火蟲算法，證明了其收斂性，能有效改進(jìn)GRNN的權(quán)值和閾值，泛化性能好，不易陷入局部最優(yōu)，提高預(yù)測性能.

(2) 應(yīng)用CFA-GRNN對舟山螺頭航道船舶交通流預(yù)測，實驗結(jié)果和誤差分析證明了CFA優(yōu)化廣義回歸神經(jīng)網(wǎng)絡(luò)能夠更好地預(yù)測海域的船舶交通流量.

(3) 用高峰前數(shù)小時的船舶交通流量數(shù)預(yù)測高峰時段的船舶流量，進(jìn)一步檢驗?zāi)Ｐ秃退惴ǖ挠行院驼_性，并對大型船舶和小型船舶分組進(jìn)行預(yù)測.