考慮多維特征的船舶軌跡分層聚類算法

蘇俊杰， 蘭培真

(1. 集美大學海上交通安全研究所， 福建 廈門 361021；2. 交通安全應急信息技術國家工程實驗室， 福建 廈門 361021)

0 引 言

航運經濟的發展促使船舶數量日益增加。船舶自動識別系統(automatic identification system，AIS)作為船與船、船與岸之間船舶相關信息的交換設備在船舶監管中廣泛應用，產生了海量的船舶動靜態信息數據，而提取特征相似的船舶運動軌跡簇是分析船舶交通規律和船舶行為模式的基礎[1]，因此，急需開展船舶軌跡的聚類研究。早期軌跡聚類研究[2-3]以軌跡點為樣本，忽略了相鄰軌跡點的時空關聯，因此近期研究樣本大多為由連續點構成的軌跡段。陳錦陽等[4]通過設定Hausdorff距離閾值劃分軌跡類別，類似指標還有歐氏距離、最長公共子序列距離、動態時間規整(dynamic time warping，DTW)距離等，然而這些指標只考慮船位偏差，描述船舶相似運動的能力有限。劉磊等[5]考慮船位、航速和航向定義綜合距離，并以K最鄰近分類器劃分軌跡；王維剛等[6]以8個軌跡特征度量軌跡相似性，設計了加權樸素貝葉斯分類器，證明了考慮多維特征度量相似性能更好地捕捉不同軌跡的差異[7]：這些監督分類方法在航路固定時表現較好，但要求提供準確的先驗知識且難以發現隱蔽軌跡簇，無法滿足多數聚類情況的需要。因此，以具有噪聲的基于密度的空間聚類(density-based spatial clustering of applications with noise, DBSCAN)算法[8-11]為典型的無監督聚類方法的應用較多。然而，全局性的經驗參數難以在密度不均的船舶交通環境下對各局部進行合理的聚類，為此，趙梁濱等[12]通過層次聚類獲取局部密度環境的適應參數。此外，LI等[13]在K均值聚類算法的基礎上提出自適應K值的多步軌跡聚類算法，但該算法直接處理整條軌跡，時空復雜度較高且丟失了軌跡的局部相似性；馬文耀等[14]采用譜聚類(spectral clustering，SC)算法區分船舶在不同行為模式下的軌跡；牟軍敏等[15]提出船舶軌跡的自適應SC算法，但船位偏差仍是相似性度量的唯一依據；冷泳林等[16]以線型特點衡量軌跡的相似性，并以仿射傳播(affinity propagation，AP)算法進行聚類，但聚類后的船舶運動軌跡相似性較低。針對上述船舶軌跡聚類方法存在的不足，本文提出一種考慮多維特征的船舶軌跡分層聚類(multi-features hierarchical clustering, MHC)算法，主要由分層軌跡相似性度量和CFA-DBSCAN算法2部分構成，這里CFA(core firefly algorithm)指核心螢火蟲算法。

1 軌跡相似性度量

船舶運動慣性大，船舶行為模式具有區域性[10]，相似軌跡應為相似水域內軌跡線型特點相似且具有相似航行動態的軌跡，因此以地理網格和時段作為軌跡段的劃分依據。任意軌跡可表示為T=(p1,p2,…,pN)，其中：N為軌跡點數；pi為i時刻的船舶運動狀態，是由經度λi、緯度φi、航速vi、航向αi和時隙Δti構成的5維特征向量，即pi=(λi,φi,vi,αi,Δti)。軌跡T的線型特點及軌跡所處水域環境可由軌跡線型圖F1和水域信息圖F2描述，其中：F1為各軌跡點按順序連接而成的曲線圖；F2為軌跡所處地理網格范圍內的地圖圖像。以圖1為例說明水域信息圖的構建過程。為保證聚類結果的合理性，應盡量選擇噪聲點少的，即同類地理環境區分度低、異類地理環境區分度高的地圖圖像數據。因為遙感影像包含的噪聲點較多，容易降低不同水域環境的辨識敏感度，所以本文研究區域的水域信息圖采用如圖1a所示的公開地圖圖像數據。地理網格的設置應使得軌跡段的點數和水域信息圖采樣范圍合理，采用圖1b所示5×5地理網格將原始軌跡劃分到18個水域，即首先使由地理網格劃分的各軌跡段包含軌跡點數在[500，1 000]內，其次使地理網格內的地圖圖像能達到粗略區分港口、沿岸和開闊水域的目的，則軌跡所處水域信息可用該水域范圍內的地圖圖像表示。

a)研究區域的船舶軌跡

b)地理網格劃分的水域

將船舶的多維特征劃分為靜態圖像上層和動態數值下層，對不同層分別建立不同的軌跡相似性度量指標，經兩次遞進聚類逐步放大軌跡簇間的差異，最終實現隱蔽軌跡簇的辨識。MHC基本框架見圖2，其中，SIFT(scale-invariant feature transform)指尺度不變特征變換算法。

圖2 MHC基本框架

1.1 靜態圖像上層

上層聚類特征為軌跡線型圖和水域信息圖，軌跡線型圖和水域信息圖均為64×64的灰度圖像，各像素點取值區間為[0,255]。SIFT算法[17]魯棒性強，在光照、旋轉、視點和仿射變換條件下均能保持優良的穩定性，故采用SIFT算法建立上層相似性的度量指標。如圖3所示，SIFT算法本質上是將圖像相似性轉換為2張圖像所有匹配特征點的相似性，因此，需要檢測圖像包含的所有特征點，并為每個特征點生成描述符。

圖3 SIFT算法基本框架

特征點檢測階段需要確定特征點所在的尺度和位置。以圖4為例說明圖像特征點的檢測過程。給定圖像I(x,y)和尺度高斯變換函數G(x,y,σ)(σ為采樣尺度)，可構建圖像的尺度空間L(x,y,σ)。尺度空間是由不同分辨率圖像構成的高斯金字塔，用于描述圖像的尺度不變性。

(1)

L(x,y,σ)=G(x,y,σ)*I(x,y)

(2)

式中：σ取1.6；*為卷積運算。

高斯差分尺度空間D(x,y,σ)為相鄰尺度的高斯差分與原圖像的卷積，特征點由高斯差分尺度空間的局部極值點構成。為尋找局部極值點，一個像素點要與其所有的相鄰點比較，看該像素點是否比其圖像域和尺度域的相鄰點大或者小。如圖4c所示，檢測點需要與26個點(與檢測點同尺度圖像域的8個相鄰點+與檢測點上下相鄰尺度圖像域的9×2個點)進行比較，以確保在尺度空間和圖像空間都檢測到極值點。

D(x,y,σ)=

(G(x,y,kσ)-G(x,y,σ))*I(x,y)

(3)

a)高斯金字塔

b)高斯差分金字塔

c)特征點

特征點描述階段需要為每個特征點生成對應的特征描述符，即以特征點鄰域像素的梯度方向分布賦予特征點方向信息，L尺度下(x,y)處的梯度值m(x,y)和方向θ(x,y)分別為

m(x,y)=((L(x+1,y)-L(x-1,y))2+

(L(x,y+1)-L(x,y-1))2)1/2

(4)

(5)

以圖5為例說明特征描述符在實際應用中的構建。首先將特征點的鄰域像素劃分為2×2的區塊，每個區塊由4×4子域構成，計算每個子域的梯度并統計各區塊內8個方向(每個方向的范圍為0～45°)的梯度方向直方圖，其中，縱軸為每個梯度方向的梯度累計值，則每個區塊可描述為一個8維的向量，該特征點的特征描述符為一個2×2×8維的向量。

a)圖像梯度

b)統計直方圖

c)關鍵點描述

特征點匹配旨在尋找兩張圖像匹配的特征點對，匹配度量為特征描述符的余弦距離。匹配特征點對的查找方法：檢索某圖像特征點在另一圖像中的最鄰近和次鄰近特征點，若最近距離與次近距離的比值小于給定閾值(通常取0.6～0.8[18])，則該特征點與最鄰近特征點匹配。圖像距離為所有匹配特征點對余弦距離的平均值。上層船舶軌跡相似性可用軌跡線型圖和水域信息圖圖像距離的均值Dimg度量。

(6)

式中：d1為軌跡線型圖的圖像距離；d2為水域信息圖的圖像距離。

1.2 動態數值下層

下層聚類特征為船位、航速、航向和時隙，軌跡相似性即軌跡運動特征的相似性。船舶軌跡是不等長的時序數據，因此，采用DTW[10]算法建立下層軌跡相似性的度量指標。以圖6為例說明下層軌跡相似性的計算過程：給定任意兩條船舶軌跡TA和TB，構建初始距離矩陣，該矩陣的元素d(i,j)為對應兩軌跡點特征向量的余弦距離。余弦距離按式(7)計算，計算前需將各特征分量進行歸一化以消除不同量綱的影響。按式(8)計算初始距離矩陣的累積距離矩陣，即累積距離矩陣中Mc(i,j)的值等于其左、上和左上位置數值中的最小值與初始距離矩陣對應位置之和。下層軌跡的相似性可由遍歷軌跡TA和TB所有軌跡點的最小累積距離Mc衡量。

(7)

Mc(i,j)=d(i,j)+min(Mc(i-1,j-1),

Mc(i,j-1),Mc(i-1,j))

(8)

式中：xik和xjk分別為點pi和pj的第k個特征分量。通常d(i,j)越小，pi與pj越相似。假設pi=(1，1，0,0,1),pj=(1,1,1,0,1),則根據式(7)有

圖6 累積距離矩陣計算

2 CFA-DBSCAN算法

2.1 改進DBSCAN算法

DBSCAN算法是基于密度的無監督聚類算法[19]，辨識隱蔽簇的能力強且無須設定簇的數量，適用于船舶軌跡數據噪聲點多、船舶行為模式不確定性高的情況，相關定義如下：

定義1核心樣本：對某樣本li，若在給定距離半徑r的鄰域ε(li)中樣本的數量不低于給定閾值δ，則該樣本為核心樣本。

定義2密度直達：某核心樣本可密度直達其鄰域包含的所有樣本。

定義3密度可達：若li密度直達lk，lk密度直達lj，則li密度可達lj，即若某核心樣本li密度直達的樣本lk為核心樣本，則樣本li密度可達lk鄰域的所有樣本，密度可達具有傳遞性。

定義4密度相連：若核心樣本lk密度可達li和lj，則li與lj密度相連，即由密度可達的所有核心樣本的鄰域構成的區域內的樣本都密度相連。

DBSCAN算法旨在將密度相連的樣本聚成同簇，而在實際聚類過程中，樣本鄰域通常存在交集。為避免DBSCAN的重復查詢，僅聚類核心樣本，對同簇的核心樣本集對應的鄰域樣本集取并集來進行簇的擴充。

2.2 CFA

改進的DBSCAN算法同樣對參數的設置比較敏感，因此，基于螢火蟲算法(firefly algorithm，FA)[20]來優化參數。FA是新興的群智能優化算法，概念簡單、優化性能好，在參數調優問題上應用廣泛。FA基本步驟[21]如下：

步驟1初始化螢火蟲數量N，最大吸引度β0，光吸收系數γ，步長因子ρ，ρ∈[0,1]，最大迭代次數tmax和搜索精度ξ。

步驟2隨機初始化螢火蟲位置，計算目標函數值S作為螢火蟲各自的最大熒光亮度I0。

步驟3計算群體中螢火蟲的相對亮度I和吸引度β,相對亮度決定螢火蟲的移動方向。

(9)

式中：rij表示螢火蟲i與j的距離；β0為最大吸引度，即rij=0時的吸引度。

步驟4更新螢火蟲被吸引后的位置。

xi(t+1)=xi(t)+β(xj-xi)+ρ(μ-0.5)

(10)

式中：t為迭代次數；μ為[0,1]內均勻分布的隨機因子。

步驟5更新螢火蟲的亮度。

步驟6若達到最大迭代次數或滿足搜索精度，則轉下一步；否則重復步驟3～6。

步驟7輸出全局極值點和個體最優亮度值。

FA缺乏全局隨機性，易陷入局部最優[20]且螢火蟲群位置的頻繁更新使得算法時空復雜度較高，因此，提出CFA。區別于FA，CFA只選取亮度排在前20%的螢火蟲作為核心螢火蟲進行位置更新，并隨機分配其余螢火蟲的位置來確保全局隨機性；在沒有局部更亮的個體時，核心螢火蟲會向種群的全局最優點移動，最優位置即核心螢火蟲的收斂位置。為避免隨著迭代次數增加，核心螢火蟲圍繞最優位置振蕩，按式(11)分段設置步長因子。

(11)

式中：ρ0為初始步長因子，ρ0∈[0,1]；lij為螢火蟲i與j在不同維度上的距離。

2.3 CFA-DBSCAN算法

CFA-DBSCAN算法為改進DBSCAN算法與CFA參數尋優算法的結合，即采用CFA優化改進DBSCAN算法的參數r和δ以獲取當前最優的軌跡聚類效果。為尋求參數r和δ的最優值，定義參數優化的目標函數為軌跡聚類的總輪廓系數S[22]，即使得同簇內樣本相似度最大，不同簇間的相似度最小[23]。其中，S取值越接近1說明算法聚類性能越好。

式中：s(i)為軌跡的輪廓系數；a(i)為樣本i與同簇內其他樣本的平均距離；b(i)為樣本i與不同樣本的平均距離；NS為樣本總數。

3 檢驗分析

以廈門港及其附近水域的AIS數據來驗證MHC算法的有效性。AIS數據集包含船舶軌跡2 364條。為確保檢驗分析結果的準確性，參考ZHAO等[24]提供的方法對采集的AIS軌跡數據進行預處理，先用地理網格劃分船舶軌跡，再按小時時段劃分船舶軌跡，得到有效軌跡段共5 006條。該軌跡數據集的時間跨度為20190815T114425—20190820T114425，地理跨度為東經117°54′37″～118°36′39″，北緯24°06′38″～24°31′37″，航速范圍為0～14 kn，航向范圍為0°～359.9°。

為檢驗MHC算法的性能，以異常軌跡數、簇數、簇內DTW距離均值、簇間DTW距離均值和準確率為評價指標，對比分析MHC算法與CFA-DBSCAN、DBSCAN[25]、AP[16]、SC[14]、K均值聚類[13]等5種其他常用船舶軌跡聚類算法的評價指標值，見表1。5種其他常用算法均以船舶軌跡的DTW距離作為算法的樣本距離；CFA-DBSCAN算法的參數按螢火蟲種群規模N=100，最大迭代次數tmax=200，初始步長因子ρ0=0.8，光吸收系數γ=1設置；設定SC和K均值聚類算法的聚類簇數與MHC算法的聚類簇數相同。由表1可見，在無須事先設定簇數的聚類算法(MHC、DBSCAN和CFA-DBSCAN)中，MHC算法將軌跡聚為9簇，而DBSCAN算法和CFA-DBSCAN算法均只將軌跡聚為3簇，表明MHC算法能更好地辨識隱蔽軌跡簇。在將軌跡聚為3簇的算法中，CFA-DBSCAN算法的各指標均優于AP算法和DBSCAN算法，這表明根據經驗設置的DBSCAN算法參數不足以區分復雜的船舶軌跡，而通過CFA確定的最優參數能在一定程度上提高DBSCAN算法的性能。在將軌跡聚為9簇的算法中，MHC算法在各評價指標上都優于SC算法和K均值聚類算法，簇內DTW距離均值為5.199，簇間DTW距離均值為18.032，準確率為91.50%，表明MHC算法能更好地辨識軌跡間的異同，從而實現船舶軌跡的準確聚類。

表1 各軌跡聚類算法的評價指標值

為分析MHC算法各模塊在軌跡聚類過程中的影響，設置3組并行實驗：第一組作為對照實驗，僅考慮船位、航速和航向特征，并用DBSCAN算法對軌跡進行聚類；第二組實驗將DBSCAN算法替換為CFA-DBSCAN算法；第三組實驗在第二組實驗的基礎上引入水域環境、軌跡線型和時隙特征對船舶軌跡進行分層聚類。如圖7所示：DBSCAN算法將船舶軌跡劃分為3簇，基本對應著港區和外海的船舶，但局部軌跡分類明顯不符合實際，如外海船舶被歸到港區的船舶簇，且異常軌跡較多，聚類效果較差。這表明在較大水域環境范圍內船舶軌跡通常密度不均，DBSCAN算法全局性經驗參數的設置難以獲取局部密度環境下的合理聚類，且在區分船舶行為多變的軌跡方面存在局限性。CFA-DBSCAN算法同樣將船舶軌跡劃分為3簇，但異常軌跡數明顯減少，表明DBSCAN算法對參數的選取敏感，而通過CFA確定的參數能在一定程度上提高算法的質量。MHC算法通過上層聚類后得到3簇軌跡，基本對應著沿岸和外海曲度不同的船舶軌跡，異常軌跡相對較少，異常軌跡基本集中在水上交通環境復雜的港區；通過下層聚類后，船舶軌跡被劃分為9簇，基本對應船舶進出港區的軌跡匯聚和軌跡分流，與目前廈門的船舶交通狀況基本吻合，且異常軌跡明顯減少。這表明MHC算法通過上層和下層遞進聚類的方式能更好獲取局部密度環境的適應參數。

表2為算法的運行指標。由表2可見：CFA-DBSCAN算法的鄰域查詢次數比DBSCAN的少，表明鄰域查詢部分的改進減少了聚類過程中的重復操作；CFA-DBSCAN算法聚類時間比DBSCAN算法的多，這是由于DBSCAN算法參數是事先設定的，不需要通過多次迭代更新，節省了聚類時間，而CFA-DBSCAN算法需要通過多次迭代使得軌跡聚類結果的輪廓系數收斂于全局最優；MHC算法需要的聚類時間比前兩種算法的多，這主要是因為在傳統船舶軌跡相似性考慮標準的基礎上引入了水域環境、軌跡線型和時隙特征，增加了軌跡相似性度量復雜度，而圖像相似度的計算和上下層遞進聚類增加了軌跡相似性度量的時間和鄰域查詢次數，但也正是這些原因使得MHC算法更容易發現隱蔽軌跡簇，在聚類結果和準確度方面較優。

a)DBSCAN

b)CFA-DBSCAN

c)MHC上層聚類

d)MHC下層聚類

表2 軌跡聚類算法的運行指標

4 結 論

為準確聚類復雜的船舶軌跡和辨識隱蔽軌跡簇，本文提出考慮多維特征的船舶軌跡分層聚類(MHC)算法。在傳統船舶軌跡聚類特征的基礎上引入水域環境、軌跡線型和時隙特征分層建立軌跡相似性度量指標并以CFA-DBSCAN算法遞進聚類各層軌跡，以廈門港及其附近水域AIS數據驗證了MHC算法的有效性。檢驗結果表明，MHC算法相較于常用的聚類算法能更準確地區分船舶軌跡并發現隱蔽軌跡簇，聚類后同簇內軌跡的船舶運動特性更相似，而不同簇軌跡的船舶運動特性差異更明顯。聚類結果可為船舶交通流特性分析和船舶行為模式識別提供精確的樣本支持。MHC算法的聚類時間較長，因此，如何兼顧MHC算法聚類性能并優化聚類時間將是后續研究的重點。