基于大規模船舶軌跡數據的航道邊界提取方法

徐垚 李卓然 孟金龍 趙利坡 溫建新 王桂玲

摘 要：傳統的道路數據獲取方法成本高、更新慢等無法適用于海洋航道的獲取，從眾源軌跡數據中提取道路或航道信息具有成本低、更新快等特性，然而，由于船舶軌跡數據噪聲多、數據量大、不同區域分布不均使得航道邊界提取面臨較大挑戰。針對該問題，提出一種基于大規模船舶軌跡數據進行航道邊界提取的方法。首先對大規模的船舶軌跡數據進行并行化去噪、插值、軌跡分段;然后，基于并行化及基于Geohash編碼的空間聚類，將軌跡數據化簡為多個方形區域的點集數據;其次，對其進行窗口劃分，對傳統的NiBlack方法進行擴展，提出SpatialNiBlack算法，對方形區域進行航道識別;最后，提出一種新的提取算法del-alpha-shape，基于航道識別結果獲得航道邊界。理論分析與實驗結果表明，所提方法在最大密度值是200，最小密度值是10，窗口長和寬分別為5和5時，可同時達到86.7%的準確率和79.4%的召回率。實驗結果表明，該方法可以從大規模的軌跡數據中提取有價值的航道邊界，是一種有效的航道提取方法。

關鍵詞：軌跡數據;自動識別系統;時空大數據;Delaunay三角網;航道提取

中圖分類號： TP319

文獻標志碼：A

Abstract： The traditional road information extraction method is high-cost and slow-update. Compared with it， road or marine lane information extraction from crowdsourcing trajectory data is low-cost and easier to update. However， it is difficult to extract lane boundary due to vessel trajectory data with high noise， large data volume and uneven distribution across different regions. To solve this problem， an extraction method of marine lane boundary from exploiting trajectory big data was proposed. Firstly， the parallelized denoising， interpolation and trajectory segmentation for trajectory big data was conducted. Then， based on parallelization and Geohash-encoded spatial clustering， trajectory data was simplified into multiple square regions. The regions were divided and the NiBlack method was extended as SpatialNiBlack algorithm to recognize regions on lane. Finally， based on the filtering results， del-alpha-shape algorithm was proposed to construct a Delaunay triangulation network and obtain marine lane boundary. The theoretical analysis and experimental results show that the proposed method can achieve an accuracy of 86.7% and a recall rate of 79.4% when the maximum density value is 200， minimum density value is 10， length and width of window are 5 and 5 respectively. The experimental results show that the proposed method is effective to extract valuable marine lane boundaries from large-scale trajectory data.

Key words： trajectory data; Automatic Identification System （AIS）; spatio-temporal big data; Delaunay triangulation network; marine lane extraction

0 引言

由于傳統道路數據依靠人工測量、高分遙感影像等采集方式，具有價格昂貴、更新慢等特點，使得道路信息的快速獲取和更新成為一直以來亟待解決和廣受關注的問題。傳統道路信息提取的難點在于人工測量及高分遙感影像的獲取昂貴，以及基于圖像的準確道路檢測與識別技術。其中，人工測量在采集信息時耗時耗力，而且受到地理環境和自然條件等限制，會帶來安全問題和時間成本大問題。高分遙感影像通過衛星進行高清影像數據采集，需要對數以百萬量不同級別的移動物體長期且頻繁地進行數據采集和處理，由于成本高昂難以推廣使用。尤其是，基于圖像的識別檢測技術對已建設的陸地道路路網識別是一種常用的方式，但海洋航道不同于陸地道路，該方法無法基于海洋影像識別航道與非航道區域。近年來，隨著移動傳感器、物聯網、云計算、大數據等技術的發展，產生了海量的時空軌跡大數據，并能夠對其進行處理和分析。通過對眾多在道路上行駛的移動物體產生的時空軌跡大數據進行分析，提取其中的道路信息成為一種受關注的技術。與陸地上的道路相比，海上沒有人為建筑的航道，且受氣候、洋流、事故等各種外部因素影響動態變化，情形更加復雜。與陸地道路信息采集相比，利用傳統的人工測量等方法采集航道信息代價更為高昂，而利用來自船舶自動識別系統（Automatic Identification System， AIS）的眾源船舶軌跡數據提取航道的信息，如果可行其一定具有低廉、更新快等特性，不失為一種可以探索且較有前景的途徑;但是，眾源船舶軌跡數據具有海量、地理范圍大、數據質量差、數據密度在不同區域差異明顯等特征[1]，這給海上航道的獲取和更新提出了更大的挑戰，本文提出了基于眾源船舶軌跡數據進行航道信息提取的一整套方法，所瞄準的主要挑戰及本文貢獻包括：

1）針對眾源船舶軌跡數據規模大的特點，采用基于MapReduce的方法對眾源船舶軌跡數據進行并行化預處理，并基于并行化的Geohash編碼方法進行空間聚類，從而既保留了船舶軌跡數據中隱藏的航道信息，又將海量的眾源船舶軌跡數據進行了有效化簡[2-3]。

2）針對大范圍眾源船舶軌跡數據在不同區域分布不均的問題，提出SpatialNiBlack局部動態閾值過濾算法[4]，對軌跡數據的聚類結果進行局部閾值過濾，并提出航道邊界提取算法del-alpha-shape，能夠對大范圍的眾源軌跡數據進行航道平面邊界的準確提取[5]。

需要指出的是，航道本是具有寬度、水深等屬性的水域，本文主要關注航道的平面屬性，聚焦于如何基于海量的船舶軌跡數據提取航道的平面邊界信息。

1 相關工作

當前，利用船舶軌跡數據進行航道提取的同類研究還比較前沿，直接相關的工作還比較少。相關工作有兩類：一類工作通過對軌跡進行行為模式的統計，基于移動對象的行為模式與位置對比分析進一步建立航道[6-10];另外一類是利用K-Means、DBSCAN（Density-Based Spatial Clustering of Applications with Noise）等聚類算法挖掘航道的行為模式進一步提取航道[11-12]。雖然利用船舶軌跡數據進行航道提取的工作較少，但利用眾源軌跡進行道路邊界的提取是交通大數據領域的一個逐漸受到越來越多的研究者關注的研究問題，國內外學者主要圍繞道路邊界的提取和功能區域邊界的提取展開了一定的研究工作，已經有一些值得借鑒的相關工作[13-14]。這些工作主要通過約束構造平面點集形狀實現邊界提取。如文獻[15]運用約束德洛內（Delaunay）三角網從車輛軌跡線集中提取道路邊界[15]，主要通過三角形邊長度和Voronoi面積等幾何特征表達軌跡點分布的聚集性差異，然后將兩種不同維數的控制條件集成建立道路邊界識別模型，運用“種子點”擴展方法實現道路邊界的提取。也有一些學者提出雙閾值的Alpha Shapes算法[16]，利用簡單環的概念設計輪廓搜索算法，獲得既有較好完整性又有較高幾何精度的建筑物輪廓線。有的學者是基于凸殼的方法，采用凸殼的內縮特征獲取邊界的凹凸特征，這種方法一定程度上能處理空間密度分布不均的情況，但難以保證邊界的規則性，也無法識別空洞現象[17]。從內容上看，這些方法無法適應采樣稀疏、高噪聲的全球定位系統（Global Positioning System， GPS）船舶軌跡數據，針對大規模船舶軌跡數據難以取得良好的效果。

綜上分析，現有的平面點集形狀重構算法能夠通過多條件約束表征點集的形狀，但是對于大規模的船舶軌跡數據則存在以下兩點問題：1）由于采樣環境原因，原始的大規模船舶軌跡數據包含大量噪聲數據和丟失數據，增加了邊界提取的難度;2）海洋上的軌跡數據受船舶的行駛習慣和活動熱度影響，不同區域的船只由于活動程度不同，造成同一時間范圍下軌跡數據量的量級在不同區域存在差異，使提取航道邊界問題變得更加復雜。針對以上存在的問題，本文提出了一種新的算法模型，首先針對大規模的船舶軌跡數據，使用并行技術進行插值、去噪和Geohash聚類預處理;然后通過時空NiBlack局部動態閾值過濾和約束三角構網實現航道邊界的提取。

2 方法概述

本文提出的算法是一種運用空間NiBlack局部動態閾值過濾及約束三角網從海量船舶軌跡點集中提取航道邊界的過程。首先，對海量的船舶軌跡數據進行并行的去除異常軌跡點以及缺失數據補全預處理，提升數據質量;然后利用并行的Geohash編碼方法對其進行聚類，以在保障結果精度的前提下提高后續數據處理和分析的效率;其次，提出空間NiBlack算法，針對大范圍海量船舶軌跡數據在不同區域分布不均的問題，對聚類結果進行局部閾值過濾，對軌跡數據的預處理結果進行航道識別;最后，基于alpha shape邊界算法提出一種新的提取算法del-alpha-shape，利用上一步局部閾值過濾的結果基于Delaunay三角剖分算法構造三角網，并提取邊界輪廓，經過閾值過濾等處理后得到航道平面邊界。

3 并行化軌跡數據空間聚類

3.1 并行化預處理

由于測量誤差、傳感器誤差等因素，GPS軌跡數據存在數據缺失和噪聲的問題，這樣的數據在航道提取中很難發揮作用，所以在應用中需要對原始數據作預處理。本文軌跡預處理過程主要由過濾和插值兩部分組成，過濾是為了消除噪聲，插值是為了補全缺失數據。因為船舶在正常航行過程中不會頻繁轉向，產生的軌跡在一定范圍內近似于直線，所以為了提高程序效率，本文采用基于時間和空間閾值的過濾和插值算法，而不是其他復雜算法。過濾算法的原理是計算點到軌跡的距離，如果點到軌跡的距離大于預定閾值則刪除該點，從而消除孤立的噪聲點。插值算法的原理是計算兩個相鄰軌跡點的距離和整體軌跡相鄰兩點的平均距離，根據相鄰兩點的距離與平均距離的比值大小，在這兩個點之間插入不同數量的軌跡點，比值越大插入的點就越多。

本文采用某海洋船舶監測系統所采集的12個月的數據作為數據源，其中一個月約有60GB的數據，3億多條記錄。在使用CentOS 6.7單機實驗時，對于60GB的數據量預處理完成大概需要40h，因此需要采用基于集群的并行化處理技術進行效率的提升。這些數據以SequenceFile序列化文件類型保存在Hadoop分布式文件系統（Hadoop Distributed File System， HDFS）中。因為原始數據集中每條記錄均包含多個字段，本文對原始數據進行篩選，每條記錄處理為定義1的四個基本字段組成。

本文方法基于MapReduce對數據并行化預處理，算法流程如圖1所示。

1）數據分割，將HDFS中保存原始數據的SequenceFile文件劃分成m（mn，n表示數據節點數）個數據塊Split，每個塊交由一個數據節點Datanode處理。

2）Map階段，程序逐行讀取字段無缺失數據，提取每條數據的v、x、y、t四個屬性，以字段v為鍵，元組（x，y，t）為鍵值輸出，其形式為〈v，（x，y，t）〉。

3）在Reduce階段，每個Reduce處理具有相同鍵v的數據，操作如下：

第1步 將相同v的數據按時間t排序，設定時間閾值t_z和距離閾值d_z，軌跡閾值n。

第2步 將第i個軌跡點保存在數組list中，如果第i+1個軌跡點和第i個軌跡點之間的距離di小于d_z，且時間間隔Δti

第3步 計算數組list中軌跡點的個數N，如果N

第4步 計算相鄰兩點第j和第j+1之間的距離dj，如果dj小于d，或者兩點的經度之差大于300（分別在東經180°附近和西經180°附近的相鄰兩點之間不用插值，這個判斷的值可以是小于360°且盡可能大的值，設為300），則j=j+1，繼續執行此步;否則令s=dj/d+2，在第j和第j+1兩點之間插入s個點，并保存在list中， j=j+1，繼續執行此步，直到遍歷list，執行第5步。

第5步 保存list中的元素，并清空list，i=i+1，執行第2步，直到遍歷相同v的數據。將預處理后的數據以字段v為鍵、元組（x，y，t）為鍵值輸出，其形式為〈v，（x，y，z）〉。

4）將預處理后的數據按v保存在不同的文本文件中，最后得到預處理后的數據。在使用上述基于MapReduce的數據預處理之后，在表1所述的集群實驗條件下，整個過程從40h縮短為150min。

3.2 基于Geohash編碼的并行化空間聚類算法

軌跡聚類采用基于Geohash編碼的算法，當軌跡點較多時，對于軌跡點，并不需要遍歷計算數據集中所有軌跡點與當前點的距離，僅需關注小范圍內的相鄰軌跡點，因此考慮引入Geohash空間編碼技術并行化聚類GPS軌跡點。基于Geohash編碼空間聚類算法的原理是將地球看作一個平面，把地球平面劃分為大小相同的網格，所有的軌跡點都會對應一個柵格，從而將位于同一柵格的軌跡點聚類。在使用CentOS 6.7單機實驗時，對于60GB的數據量Geohash編碼完成大概需要18h，因此需要基于集群進行并行化處理。

本方法基于MapReduce對數據并行處理，算法流程如圖2所示。

程序的輸入數據是3.1節預處理后的數據，程序按照以下流程運行。

1）數據分割，將HDFS中以不同v命名的m個數據文件劃分成m′（m′≥n，n表示數據節點數）個數據塊Split，每個塊交由一個Datanode處理;

2）Map階段，抽取每條數據的經度x和緯度y并進行Geohash編碼，然后以編碼code為鍵，1為鍵值輸出，其形式為〈code，1〉;

3）Reduce階段，統計不同Geohash編碼code的密度dsy，計算其對應區域C的中心點經緯度C_x和C_y，以C_x、C_y和dsy為鍵，null為鍵值輸出，其形式為〈（C_x，C_y，dsy），null〉;

4）將所有數據按元組〈C_x，C_y，dsy〉保存。完成空間聚類。

在本文中，Geohash編碼位數為11位，大約等于11km見方的一個區域。在使用上述基于MapReduce的數據預處理之后，在表1所述的集群實驗條件下，整個過程從18h縮短為110min。

4 基于NiBlack局部動態閾值過濾及alpha shape的邊界提取算法

4.1 SpatialNiBlack局部動態閾值過濾算法

Niblack算法是一種比較常用、簡單、有效的局部動態閾值算法，通常應用于數字圖像分割、圖像二值化等領域，該算法根據局部平均值和局部標準差確定當前區域閾值。由于不同區域的軌跡點密度存在差異，固定閾值無法準確地過濾數據，因此本文基于NiBlack動態窗口過濾算法思想，提出適于軌跡數據分析的SpatialNiBlack算法，實現不同密度的軌跡點數據過濾，該方法有效克服了固定閾值的缺陷。

本文所用的數據是3.2節中聚類的結果，這些數據被保存在一個二維數組中，數組中的每一行用一個三元組〈C_x，C_y，dsy〉來表示一個柵格C，其中C_x和C_y分別表示柵格中心點的經度和緯度，dsy表示柵格C中軌跡點的密度。

本方法步驟如圖3所示，基于SpatialNiBlack局部動態閾值的過濾算法如下。

1）首先得到全球11位的柵格，設置大小為m×n的矩形窗口，因為數據從緯度變化方向來看，有密度差異，所以窗口設置成扁平的矩形，使窗口內大多的點的密度在同一個級別。窗口包含m×n個11位柵格，執行第2）步;

2）遍歷當前窗口中的柵格，統計窗口中所有柵格的dsy值，根據式（1）計算出閾值T（x，y），然后判斷窗口中心點所在柵格是否保留。執行下一步;

其中對于柵格中心點經緯度（x，y），T（x，y）為該位置的閾值，m（x，y）為該點所在窗口內軌跡點密度的均值，s（x，y）為該點所在窗口內軌跡點密度的標準差，k為修正系數（通常取-0.1）。

3）窗口中心柵格后移一個格，以該柵格為中心建立新的窗口，執行第2）步;直到窗口遍歷所有柵格。最后輸出過濾后的數據。

4.2 基于alpha shape的邊界提取

本文基于alpha shape邊界算法，提出一種新的提取算法del-alpha-shape，該算法能夠有效地對點集進行邊界提取。

alpha shape是一種經典的點集輪廓提取算法，即輸入點集，通過對任意兩點建立半徑為alpha的外接圓，如果沒有第三個點在該圓內，則這兩點連成的線段是邊界線段。通過對alpha值的調整可以得到不同精確程度的輪廓。如果兩個點連線距離大于2倍alpha閾值，那么無法構成半徑alpha的外接圓，對這兩點不進行判斷。

Delaunay三角剖分算法是將點集處理成由三角面構成的凸包平面圖，該平面圖內除了端點，不包含任何點集中的點。

定義9 三角形密度指數alpha。三角形的三個點可以構成唯一的一個外接圓，該外接圓半徑circum_r的倒數值，稱為三角形密度指數。

在Delaunay三角網圖（圖4）中，航道內的三角形的密度指數較大，航道外的三角形的密度指數較小。設定alpha閾值（用alpha_value表示），如果三角形alpha≥alpha_value，為航道三角形，保留該三角形的三條邊到Edges;否則，為空洞三角形，不保留邊。

該算法流程如圖5所示。

本文4.1節得到的柵格過濾結果是一個點集數據，以Points表示。該算法首先使用Delaunay方法對點集Points進行Delaunay三角化，得到三角面集合Triangles（第1）行）;然后遍歷Triangles，使用calCircum方法計算每個三角形的密度指數circum，如果circum不小于密度指數閾值，則將該三角形的邊加入邊集Edges（第2）～6）行）;下一步用polygonize方法對Edges進行多邊形化，構成多邊形集Polys（第7）行）;最后，遍歷Polys，判斷每個多邊形的面積是否大于面積閾值maxArea，如果大于，則將該多邊形的邊界坐標polygon.coords添加到Channels集合（第8）～11）行）。算法時間成本主要在第7）行多邊形生成過程，整體時間復雜度是O（N2）。

5 實驗分析

5.1 數據集介紹

船只航行軌跡數據是全球2016年內的所有船只軌跡。船只航行軌跡數據包括GPS時間戳、水上移動通信業務標識碼（MMSI）、GPS經度、GPS緯度等多項信息。軌跡點采樣間隔在10～120s，總數據量在1TB左右。

本實驗選取兩次不同的數據：一是2016年6月渤海區域的油輪航行軌跡數據;一是2016年6、7、8月“一帶一路”區域大宗貨物船只航行軌跡數據。本實驗的預處理步驟是在6臺搭建了Hadoop環境的機器下進行的，系統為CentOS release 7.0操作系統，JDK版本為1.7.0，Hadoop版本為2.6.0。集群具體配置如表1所示。

實驗合并過濾程序使用Java1.8編程語言實現，航道提取使用Python 3.5編程語言實現，在單臺物理機（64GB內存，系統CentOS release 7.0）運行。

5.2 實驗結果與分析

本文設計以下實驗驗證模型，第一組實驗使用數據是中國—東南亞—非洲連線區域2016年2—12月的大宗船只的航行數據。

圖6（a）是使用3.1和3.2節介紹預處理方法，預處理相關的參數取值為：距離閾值30km，時間閾值3600s，分段軌跡點數目閾值30，對中國—東南亞—非洲連線區域2016年2—12月的大宗船只原始軌跡數據預處理后的柵格化聚類結果上圖此句不通順，請作相應調整。。圖6（b）是使用4.1節介紹方法滑動窗口局部動態閾值過濾的結果。這里涉及到幾個參數設置，柵格位數11位、柵格最大閾值200、窗口寬50、窗口高3。圖6（c）是基于滑動窗口局部動態閾值過濾結果用4.2節介紹方法作航道提取的上圖圖6（c）是基于滑動窗口局部動態閾值過濾和4.2節介紹方法做航道提取的結果[18]此句不通順，請作相應調整。。

第二組實驗使用數據是渤海區域的2016年6—8月的油輪航行數據，主要為了進一步驗證較小地理空間范圍下模型提取航道和空洞的準確性。

圖7（a）是使用3.1和3.2節介紹的預處理方法，預處理相關的參數取值為：距離閾值30km，時間閾值3600s，分段軌跡點數目閾值30，對渤海區域2016年6—8月油輪原始軌跡數據預處理后的柵格化聚類結果如圖7。

圖7（b）是使用4.1節介紹方法滑動窗口局部動態閾值過濾的結果。這里參數設置為，柵格位數16位，柵格最大閾值500，窗口寬15，窗口高7。圖7（c）是基于滑動窗口局部動態閾值過濾結果用4.2節介紹方法作航道提取的上圖圖7（c）是基于滑動窗口局部動態閾值過濾和4.2節介紹方法做航道提取的結果。

本文以渤海某一區域的使用固定閾值進行航道提取的結果為參考，如圖8（a），然后使用本文提出的大規模窗口局部過濾法的航道結果與之對比分析。采用過濾結果疊置方法評價航道提取的有效性，并借鑒文獻的方法計算兩種常用的評價指標，即準確率（precision）和召回率（recall）來評價疊置結果。準確率反映實驗結果的準確度，召回率反映實驗結果的完整度。

本文準確率和召回率的計算如式（2）和式（3）：

其中：gridNumextstd表示新提取的航道柵格結果與標準航道柵格結果重合的數目，gridNumstd表示標準航道柵格結果的數目，gridNumext表示提取航道柵格結果的數目。

實驗數據是2016年6月至2016年—8月份油輪軌跡數據，在渤海某一區域內進行航道提取方法驗證，采取不同參數對比結果如圖8。

參數值含義

maxDsy100～∞這個數據項是何意？100后面加了“--”，何意？回復：表示取值范圍，可改為100-∞最大密度值

minDsy0～10最小密度值

w_size3×3、5×5、10×3滑動窗口大小

表3是選用不同參數進行多組實驗后，與標準航道結果疊置后求得的準確率和召回率結果。其中：precision表示準確率，recall表示召回率。從表中可以看出：

1）最大密度值maxDsy的限制能夠避免局部窗口內方差過大造成的閾值過大，圖8（b）導致窗口內所有值被過濾舍棄，造成召回率偏低，因此，設置最大密度值進行限制，并分多組實驗測試當該值為200時效果準確率和召回率都有較好的表現。

2）最低密度值minDsy可以過濾掉一些存在的噪聲數據，提升召回率，經過多組測試發現10是一個比較合適的值。

3）滑動窗口的大小選擇對結果有一定影響，當窗口較小時，準確率因閾值計算考慮范圍過小造成準確率下降，如圖8（k）;當窗口形狀不合適，可能造成閾值結果過大，過濾舍棄掉過多的值導致航道部分缺失，如圖8（l）。

實驗驗證了本文所提方法能夠在大范圍數據下保持較高的提取精確度，同時依然保持較高的完整度;同時，本文所提算法與其他相似的研究算法有所區別：

1）基于三角網的邊界提取。文獻[15]通過對軌跡數據構造三角網并約束條件提取道路，但其主要思路是通過自定義的長短邊概念找到邊界。本文主要通過三角網的外接圓半徑約束條件尋找邊界邊;同時，通過實現兩種算法進行實驗，發現本文算法在相同量的大規模數據下的效率更高，運行時間效率可提高約5倍此處表述不當，意義不對，“運行時間”應該是“減少幾分之幾”，請按照這個思路進行調整;或者改為“運行效率可提高約5倍”，這樣修改合適嗎？請明確。

2）基于GPS道路軌跡線提取。文獻[13]通過聚類方式改進算法，提取了不同級別道路的中心線結果。本文所提航道除了能夠表示路網結構外，還能夠看出航道的邊界大致范圍，由于航道在海面上沒有寬度限制，因此該范圍是一個不嚴格的大概范圍。

3）圖像輪廓線提取算法。文獻[16]借助圖像處理的研究中進行，其輸入大致是圖像數據，通過對像素值的識別計算，提取出包括空洞的輪廓線結果。本文算法是基于時空軌跡大數據進行的邊界提取，由于時空軌跡數據密度值相比0～255的像素值范圍具有更大的取值范圍，因此在算法思路上有所區別，本文在4.1節局部過濾算法內通過參數的限制有效解決了時空數據復雜的問題。

6 結語

本文提出一種基于大規模船舶軌跡數據進行航道邊界提取的方法，通過去噪、插值、分段預處理和Geohash編碼簡化處理數據;然后基于處理后的數據結果使用文中提出的SpatialNiBlack航道識別算法過濾航道;最后使用del-alpha-shape算法獲得有效的航道邊界。通過實驗可以發現：該方法能夠有效解決邊界提取的空洞問題;在較大的地理范圍下能夠進行航道邊界提取;通過改進的局部過濾算法，可以解決全局范圍下的數據因存在質量差異導致部分航道丟失的問題。進一步改進的方面：完善算法，使其能夠同時對近海區域的一些較為精細的航道進行更加有效的提取。

參考文獻 （References）

[1] WANG Y， ZHU Y， HE Z， et al. Challenges and opportunities in exploiting large-scale GPS probe data， HPL-2011-109 [R]. Palo Alto， CA： HP Laboratories， 2011.

[2] ZHANG W， GOERLANDT F， KUJALA P， et al. An advanced method for detecting possible near miss ship collisions from AIS data[J]. Ocean Engineering， 2015， 107（1）：60-69.

[3] WANG X， LIU X， LIU B， et al. Vessel route anomaly detection with Hadoop MapReduce[C]// Proceedings of the 2015 IEEE International Conference on Big Data. Piscataway， NJ： IEEE， 2015：25-30.

[4] NIBLACK W. An Introduction to Digital Image Processing[M]. Englewood Cliffs， NJ： Prentice-Hall， 1986：115-116.

[5] EDELSBRUNNER H， KIRKPATRICK D， SEIDEL R. On the shape of a set of points in the plane[J]. IEEE Transactions on Information Theory，1983，29（4） ：551-559.

[6] SUN F， DENG Y， DENG F， et al. Unsupervised maritime traffic pattern extraction from spatio-temporal data[C]// Proceedings of the 2015 International Conference on Natural Computation. Piscataway， NJ： IEEE， 2015：1218-1223.

[7] MAZZARELLA F， FERNANDEZ ARGUEDAS V， VESPE M. Knowledge-based vessel position prediction using historical AIS data[C]// Proceedings of the 2015 Sensor Data Fusion： Trends， Solutions， Applications. Piscataway， NJ： IEEE， 2015：1-6.

[8] FERNANDEZ ARGUEDAS V， PALLOTTA G， VESPE M. Maritime traffic networks： from historical positioning data to unsupervised maritime traffic monitoring[J]. IEEE Transactions on Intelligent Transportation Systems， 2017， 19（3）： 722-732.

[9] PALLOTTA G， VESPE M， BRYAN K. Vessel pattern knowledge discovery from AIS data： a framework for anomaly detection and route prediction[J]. Entropy， 2013， 15（6）：2218-2245.

[10] SPILIOPOULOS G， ZISSIS D， CHATZIKOKOLAKIS K. A big data driven approach to extracting global trade patterns[C]// Proceedings of the 2017 International Workshop on Mobility Analytics for Spatio-temporal and Social Data. Berlin： Springer， 2017：109-121.

[11] ETIENNE L， DEVOGELE T， BOUJU A. Spatio-temporal trajectory analysis of mobile objects following the same itinerary[C]// Proceedings of the 2010 International Symposium on Spatial Data Handling. Hong Kong： [s.n.]， 2010：86-91.

[12] 楊杰，李小平，陳湉.基于增量時空軌跡大數據的群體挖掘方法[J].計算機研究與發展，2014，51（s2）：76-85.（YANG J， LI X P， CHEN T. A group mining method for incremental spatio-temporal trajectory big data[J]. Journal of Computer Research and Development， 2014， 51（s2）：76-85.）

[13] 歐陽鴻，劉建勛，劉毅志，等.基于步行GPS軌跡的路網提取方法[J].計算機與現代化，2014（2）：124-128.（OUYANG H， LIU J X， LIU Y Z， et al. An extraction method of road network based on walking GPS trajectories[J]. Computer and Modernization， 2014，（2）：124-128.）

[14] 楊偉，艾廷華.基于車輛軌跡大數據的道路網更新方法研究[J].計算機研究與發展，2016，53（12）：2681-2693.（YANG W， AI T H. A method for road network updating based on vehicle trajectory big data[J]. Journal of Computer Research and Development， 2016，53（12）：2681-2693.）

[15] 楊偉，艾廷華.運用約束Delaunay三角網從眾源軌跡線提取道路邊界[J].測繪學報，2017，46（2）：237-245.（YANG W， AI T H. The extraction of road boundary from crowdsourcing trajectory using constrained Delaunay triangulation[J]. Acta Geodaetica et Cartographica Sinica， 2017， 46（2）： 237-245.）

[16] 李云帆，譚德寶，高廣，等.雙閾值Alpha Shapes算法提取點云建筑物輪廓研究[J].長江科學院院報，2016，33（11）：1-4.（LI Y F， TAN D B， GAO G， et al. Extraction of building contour from point clouds using dual threshold Alpha Shapes algorithm[J]. Journal of Yangtze River Scientific Research Institute， 2016，33（11）： 1-4.）

[17] 朱杰，孫毅中.多約束的平面點集形狀重構方法[J].測繪學報，2017，46（2）：253-264.（ZHU J， SUN Y Z. An efficient approach to shape reconstruction from planar point set based on multi-constraints[J]. Acta Geodaetica et Cartographica Sinica， 2017，46（2）：253-264.）

[18] WU L， XU Y， WANG Q， et al. Mapping global shipping density from AIS data[J]. Journal of Navigation， 2017， 70（1）：67-81.