李 楠，強(qiáng)懿耕，樊 瑞

(中國(guó)民航大學(xué) 空中交通管理學(xué)院，天津 300300)

0 引 言

近年來(lái)隨著我國(guó)民航業(yè)的高速發(fā)展，復(fù)雜機(jī)場(chǎng)已陸續(xù)實(shí)施基于性能導(dǎo)航(PBN)的導(dǎo)航技術(shù)來(lái)提高運(yùn)行和指揮效率，因此需要設(shè)計(jì)新的PBN程序來(lái)簡(jiǎn)化傳統(tǒng)的標(biāo)準(zhǔn)飛行程序。筆者針對(duì)飛行密度大、飛行情況比較復(fù)雜的終端區(qū)空域，通過(guò)聚類算法分析航空器歷史軌跡，挖掘航空器潛在的飛行模式，其中所蘊(yùn)含的管制員指揮模式可作為評(píng)估、優(yōu)化進(jìn)離場(chǎng)程序的重要參考。

軌跡聚類分析是目標(biāo)軌跡模式挖掘過(guò)程中常用的方法，用以探究隱藏在海量軌跡信息中目標(biāo)的運(yùn)動(dòng)特征與規(guī)律，在探索人類的移動(dòng)模式、颶風(fēng)的運(yùn)動(dòng)模式和城市交通結(jié)構(gòu)智能劃分等[1-3]諸多領(lǐng)域已經(jīng)得到了廣泛的應(yīng)用。目前，在航空領(lǐng)域中，國(guó)內(nèi)外學(xué)者對(duì)航空器飛行軌跡的聚類分析多集中于兩個(gè)方面。一是借助經(jīng)典的聚類方法(K-means，F(xiàn)CM等)來(lái)分析整條軌跡或軌跡段，挖掘軌跡的運(yùn)行模式，如C. WANG等[4]針對(duì)具體的一類航跡簇，采用基于K-Means算法對(duì)航跡相似性矩陣進(jìn)行聚類分析得到該類簇的中心航跡，計(jì)算中心軌跡與標(biāo)準(zhǔn)進(jìn)場(chǎng)飛行程序(STAR)的偏差來(lái)為重新設(shè)計(jì)飛行程序提供參考；R. FRANK等[5]采用軌跡點(diǎn)比對(duì)的方法計(jì)算出軌跡間的相似性距離，再針對(duì)相似性距離矩陣進(jìn)行聚類分析；Q. Y. YU等[6]提出了一種基于多特征軌跡相似性度量的新軌跡聚類算法，可以最大化同一聚類軌跡的相似性；A. SALAPOUR等[7]采用時(shí)間扭轉(zhuǎn)匹配的方法來(lái)度量軌跡間的相似性，提出了一種基于方向的軌跡聚類相似性度量；李楠等[8]建立了基于速度修正系數(shù)的軌跡相似性模型，再應(yīng)用譜聚類實(shí)現(xiàn)了航空器飛行軌跡的合理分類；郭威等[9]沿用傳統(tǒng)的最小描述長(zhǎng)度方法劃分軌跡，再結(jié)合密度聚類和掃描線算法得出表征目標(biāo)物體運(yùn)動(dòng)的特征軌跡；徐濤等[10]引入懲罰弧形面積改進(jìn)了在計(jì)算軌跡間相似度時(shí)軌跡航向的影響，再利用聚類方法準(zhǔn)確地挖掘出了飛行程序軌跡。二是通過(guò)借助幾何模型或其他機(jī)器學(xué)習(xí)算法對(duì)軌跡重采樣，使重采樣后每條軌跡的軌跡點(diǎn)數(shù)和維度均保持一致，直接進(jìn)行聚類分析，避免了使用相似性度量方法計(jì)算軌跡間的距離，如趙元棣等[11]在保持較低失真率的前提下使用主成分分析法降低數(shù)據(jù)維度,再運(yùn)用均值漂移聚類算法實(shí)現(xiàn)飛行軌跡聚類分析并提取出異常軌跡；王濤波等[12]在采用模糊減法聚類提取軌跡特征點(diǎn)后，運(yùn)用改進(jìn)后的模糊聚類方法直接對(duì)航跡特征點(diǎn)進(jìn)行聚類并構(gòu)建相應(yīng)的中心航跡，更好的反映了真實(shí)飛行中航空器的飛行模式。以上研究存在有待改進(jìn)的地方，僅僅考慮軌跡間的空間距離對(duì)航空器進(jìn)行聚類劃分是不夠的，時(shí)間是影響聚類結(jié)果的一個(gè)重要因素；由于軌跡點(diǎn)自身包含位置、速度和航向等多維特征，建立相似性模型時(shí)可以在考慮位置特征的前提下引入速度和航向等動(dòng)態(tài)特征；當(dāng)軌跡規(guī)模較大且噪聲點(diǎn)較多時(shí)，常用的聚類算法容易受初始點(diǎn)對(duì)噪聲敏感、多個(gè)參數(shù)輸入且需要先驗(yàn)知識(shí)等因素影響，陷入局部最優(yōu)搜索，導(dǎo)致聚類效果變差。

1 航空器飛行軌跡數(shù)據(jù)處理

1.1 航空器飛行基準(zhǔn)軌跡

在監(jiān)視數(shù)據(jù)庫(kù)中，航空器飛行軌跡是由離散的數(shù)據(jù)點(diǎn)依時(shí)間順序連接而成的點(diǎn)序列。采用某機(jī)場(chǎng)終端區(qū)ADS-B監(jiān)視數(shù)據(jù)作為原始數(shù)據(jù)?？紤]到原始軌跡中包含大量不完整軌跡、飛越航班和臨近機(jī)場(chǎng)起降航班的軌跡，需要進(jìn)行軌跡預(yù)處理。

根據(jù)航空器飛行特點(diǎn)篩選出進(jìn)場(chǎng)航班，刪除軌跡點(diǎn)嚴(yán)重缺失的軌跡，采用時(shí)間間隔為3 s的線性插值更新剩余軌跡。鑒于ADS-B接收機(jī)受地形和傳輸損耗的影響，導(dǎo)致測(cè)量值與真實(shí)值存在一定的誤差，所以結(jié)合CV卡爾曼濾波算法平滑處理軌跡，最大限度地接近飛機(jī)的真實(shí)位置，建立航空器飛行基準(zhǔn)軌跡。

1.2 基于時(shí)間比的自頂向下飛行軌跡壓縮算法

經(jīng)線性插值和平滑處理后的軌跡包含大量的軌跡點(diǎn)，采用Top-Down Time Ratio算法壓縮軌跡，以減小計(jì)算開銷，從根本上提高軌跡大數(shù)據(jù)分析能力。Top-Down Time Ratio算法是Douglas-Pecuker(DP)算法的擴(kuò)展。與傳統(tǒng)DP算法不同的是，DP算法在考慮空間特性的前提下還考慮了時(shí)間特性，判斷是否為軌跡關(guān)鍵點(diǎn)的依據(jù)也由DP算法中的垂直距離變成了同步歐氏距離 (SED)[13]。同步歐氏距離指原始軌跡中真實(shí)軌跡點(diǎn)與其在軌跡首尾連接線上按時(shí)間比例所求得的歐氏距離。

如圖1，A到K為原始軌跡中相鄰的軌跡點(diǎn)，具體的軌跡壓縮步驟為：

圖1 Top-Down Time Ratio軌跡壓縮算法Fig. 1 Top-Down Time Ratio trajectory compression algorithm

步驟1：將軌跡首尾兩點(diǎn)A,K連接成一條直線，該直線稱為軌跡的弦。

步驟2：應(yīng)用同步歐氏距離(如圖1中DO)，計(jì)算離弦距離最遠(yuǎn)的點(diǎn)及其距離，比較該軌跡點(diǎn)距離與誤差閾值距離D的大小。

步驟3：若該軌跡點(diǎn)距離小于誤差閾值距離D，則舍棄該點(diǎn)，以弦作為該軌跡的近似軌跡；若大于誤差閾值距離，則以此關(guān)鍵點(diǎn)將曲線劃分成兩段，并對(duì)分割后的兩段軌跡分別重復(fù)步驟1～3操作。

如此循環(huán)，直到所有距離均小于誤差閾值距離，則保存這些關(guān)鍵點(diǎn)構(gòu)成軌跡，即為壓縮后的軌跡。軌跡壓縮后的壓縮率(compression ratio，RC)計(jì)算公式如式(1)：

(1)

式中：M0和M1分別為原始軌跡壓縮前后的軌跡點(diǎn)數(shù)目。同時(shí)，為了表征原始軌跡與壓縮后軌跡之間的差距，引入平均距離誤差S來(lái)衡量，計(jì)算公式如式(2)：

(2)

式中：dist(pi)為原始軌跡中的第i個(gè)點(diǎn)距其壓縮后對(duì)應(yīng)的關(guān)鍵點(diǎn)的歐氏距離；n為原始軌跡的軌跡點(diǎn)數(shù)目。

2 航空器飛行軌跡聚類方法的改進(jìn)

2.1 建立多維特征的軌跡相似性模型

設(shè)航空器飛行軌跡可表示為TR，則某條飛行軌跡i可表示TRi={TRi1,TRi2,…,TRij,…,TRin}，其中n為該條軌跡的軌跡點(diǎn)數(shù)，且對(duì)于不同的飛行軌跡，壓縮后每條軌跡的軌跡點(diǎn)數(shù)不一定相同?？紤]到兩軌跡點(diǎn)之間的位置、速度和航向?qū)傩?，定義軌跡TRx內(nèi)任一點(diǎn)TRxa到另一條軌跡TRy內(nèi)任一點(diǎn)TRyb的多因素距離如式(3)：

multidist(TRxa,TRyb)=wd·dist(TRxa,TRyb)+wv·dist(Vxa,Vyb)+wθ·dist(θxa,θyb)

(3)

式中：Vxa,Vyb,θxa和θyb分別為軌跡點(diǎn)TRxa和TRyb的速度和航向；dist(Vxa,Vyb)和dist(θxa,θyb)分別為兩點(diǎn)之間速度屬性和航向?qū)傩缘臍W式距離；wd,wv和wθ分別為位置、速度和航向?qū)傩缘臋?quán)重，權(quán)重因子的選擇主要依賴于各機(jī)場(chǎng)的不同標(biāo)準(zhǔn)進(jìn)離場(chǎng)程序，滿足wd≥0,wv≥0,wθ≥0和wd+wv+wθ=1。為了消除各特征屬性之間，不同量綱會(huì)對(duì)多因素距離結(jié)果產(chǎn)生影響，采用歸一化處理各特征屬性，再進(jìn)行距離的度量。

同樣，定義基于Hausdorff距離的軌跡相似性度量方法,計(jì)算兩條軌跡之間的距離,計(jì)算公式如式(4)：

(4)

同理，求出r(TRy,TRx)，再利用公式R(TRx,TRy)=max{r(TRx,TRy),r(TRy,TRx)}，得到兩條軌跡間的相似性距離,最后求出所有軌跡間的相似性矩陣。

2.2 模糊C-均值聚類的實(shí)現(xiàn)

模糊C-均值聚類算法(FCM)[14]是用隸屬度確定每個(gè)數(shù)據(jù)點(diǎn)屬于某個(gè)聚類程度的一種聚類算法，其聚類準(zhǔn)則是求U，V,使得J(U,V)取得最小值。應(yīng)用FCM將給定的樣本空間X={x1,x2,…,xn}劃分為c類,分類結(jié)果用隸屬度矩陣U=(μij)c×n表示；μij表示第i個(gè)樣本屬于第j類的隸屬度，μij滿足式(5)：

(5)

定義目標(biāo)函數(shù)：

(6)

式中：m為冪指數(shù)，且m>1；vj(j=1,2,…,c)是聚類中心。

2.3 改進(jìn)的模糊C-均值聚類

由于FCM采用的梯度法搜索方向總是沿著能量減小的方向，使得算法本身對(duì)初始值敏感。針對(duì)FCM存在的問(wèn)題，采用將粒子群(PSO)算法和禁忌搜索(TS)算法相結(jié)合的組合優(yōu)化算法，通過(guò)PSO算法的全局尋優(yōu)能力和TS算法跳出局部極值的能力得到最優(yōu)初始聚類中心，從而更加快速、有效地收斂到全局最優(yōu)解。具體步驟如下：

步驟1：輸入待聚類數(shù)據(jù)，給定FCM算法和組合優(yōu)化算法的初始參數(shù)，如冪指數(shù)m，最大迭代次數(shù)Tmax，學(xué)習(xí)因子c1、c2，速度上下界vmax、vmin，慣性權(quán)重w，鄰域解、候選解個(gè)數(shù)及禁忌表長(zhǎng)。

步驟2： 初始化粒子群規(guī)模n，隨機(jī)選擇s個(gè)聚類中心組成n×s個(gè)粒子。設(shè)初始位置為Xi，并隨機(jī)產(chǎn)生此n個(gè)粒子的隨機(jī)速度，t=0。

步驟3： 按照式(6)計(jì)算每個(gè)粒子個(gè)體的適應(yīng)度值fit(i)，并找出個(gè)體極值pbest(i)和全局極值gbest(i)。再比較當(dāng)前粒子群的局部和全局最優(yōu)位置和每個(gè)粒子的適應(yīng)度值，最后決定是否替換。

步驟4：根據(jù)式(7)和式(8)更新每個(gè)粒子的速度和位置，如果滿足誤差閾值要求或算法的迭代次數(shù)則繼續(xù)下一步，否則返回步驟3。

V′i=w·Vi+c1·r1·[pbest(i)-Xi]+c2·r2·(gbest-Xi)

(7)

X′i=Xi+Vi

(8)

步驟5：以PSO算法的輸出解作為TS算法的初始解，置空禁忌表。

步驟6：如果滿足收斂準(zhǔn)則直接輸出最優(yōu)初始中心；如果不滿足，根據(jù)領(lǐng)域函數(shù)的定義得到當(dāng)前解的所有(或若干)鄰域解，并從中確定若干候選解。

步驟7：根據(jù)藐視規(guī)則來(lái)判斷每個(gè)候選解是否均滿足，若滿足，則將此候選解認(rèn)定為當(dāng)前解，同時(shí)將最早放入禁忌表的對(duì)象用當(dāng)前解的對(duì)象替換；若不滿足，則判斷候選解的禁忌屬性，此時(shí)當(dāng)前解為非禁忌對(duì)象對(duì)應(yīng)的最佳解，同時(shí)將最早放入禁忌表的對(duì)象用與之對(duì)應(yīng)的禁忌對(duì)象替換。

步驟8：判斷是否滿足收斂準(zhǔn)則，如果滿足則直接輸出FCM的最優(yōu)初始中心；如果不滿足，返回步驟6。

3 算例仿真

3.1 飛行軌跡的壓縮處理

以華北地區(qū)某機(jī)場(chǎng)終端區(qū)第16R號(hào)跑道的679條進(jìn)場(chǎng)歷史軌跡為研究對(duì)象，共有694 239個(gè)原始軌跡點(diǎn)。應(yīng)用Top-Down Time Ratio算法壓縮飛行軌跡，誤差閾值距離D是一個(gè)很關(guān)鍵的參數(shù)指標(biāo)，直接決定了是否保留某些關(guān)鍵特征點(diǎn)。試驗(yàn)取一條軌跡點(diǎn)數(shù)目為941個(gè)的飛行軌跡為研究對(duì)象，選擇5個(gè)不同的誤差閾值距離來(lái)衡量其對(duì)軌跡特征點(diǎn)提取的影響，從而確定可以更好反映軌跡分布的軌跡特征點(diǎn)數(shù)目，結(jié)果如表1和圖2。

表1 軌跡壓縮算法的誤差閾值對(duì)比Table 1 Error threshold comparison of trajectory compression algorithm

圖2 不同誤差閾值距離對(duì)應(yīng)的軌跡壓縮Fig. 2 Trajectory compression corresponding to different error threshold distances

由表1可以看出，隨著誤差閾值距離D的不斷增加，軌跡壓縮效果越來(lái)越好，但平均誤差距離也進(jìn)一步增大。同時(shí)由圖2明顯地看到，軌跡的關(guān)鍵點(diǎn)多集中在曲線轉(zhuǎn)彎處，近似直線段用較少的點(diǎn)來(lái)表征，且隨著誤差閾值距離的減小，更多的特征點(diǎn)被提取出來(lái)用來(lái)描述近似直線段特征。同時(shí)為了檢驗(yàn)筆者算法在軌跡壓縮時(shí)的精確性和穩(wěn)定性，同現(xiàn)有的DP算法進(jìn)行對(duì)比，以平均距離誤差S作為誤差指標(biāo)，結(jié)果如圖3。

由圖3可知，在相同的壓縮率下，Top-Down Time Ratio算法的平均距離誤差更小些，能較好地還原原始軌跡。此外，隨著壓縮率的不斷增加，兩者的誤差均呈上升趨勢(shì)。因此，在保證飛行特征和有效降低數(shù)據(jù)規(guī)模的前提下，此處暫取誤差閾值距離為30 m來(lái)壓縮軌跡。

圖3 不同壓縮率下的平均距離誤差Fig. 3 Average distance error under different compression ratios

3.2 飛行軌跡聚類分析

在不同的飛行階段，要均衡體現(xiàn)出位置、速度和航向?qū)壽E的影響，需要選取不同的權(quán)重值。經(jīng)試驗(yàn)，取w1=0.7,w2=0.2,w3=0.1。利用式(3)計(jì)算出所有軌跡間的相似性距離，得到軌跡間的相似性矩陣。

設(shè)置FCM算法參數(shù)。取聚類數(shù)c=60，即飛行軌跡的聚類中心數(shù)；冪指數(shù)m=2；目標(biāo)函數(shù)的終止容限為1×10-6；最大迭代次數(shù)Tmax=100。

設(shè)置PSO算法參數(shù)。取粒子群規(guī)模N=50；學(xué)習(xí)因子c1=2,c2=2;最大速度Vmax=2,最小速度Vmin=-2;慣性權(quán)重w=0.6；最大迭代次數(shù)Tmax=100。

設(shè)置TS算法參數(shù)。取禁忌表長(zhǎng)為25；鄰域解和候選解均取20個(gè)；終止判據(jù)是最大迭代次數(shù)為100。

應(yīng)用FCM和TSPSO_FCM算法對(duì)軌跡間的相似性距離進(jìn)行聚類，已知該機(jī)場(chǎng)至少有3條傳統(tǒng)的標(biāo)準(zhǔn)進(jìn)場(chǎng)程序，此處暫取聚類數(shù)為3類，聚類結(jié)果如圖4。

圖4 TSPSO_FCM算法的飛行軌跡聚類Fig. 4 Flight trajectory clustering of TSPSO_FCM

由于聚類算法是以最終保證組內(nèi)相似度最高，組間相似度最低為聚類原則，結(jié)合FCM聚類原理及式(7)知可用目標(biāo)函數(shù)值來(lái)衡量聚類效果，F(xiàn)CM與TSPSO_FCM算法各自的聚類性能結(jié)果如圖5。

圖5 基于目標(biāo)函數(shù)值的聚類性能比較Fig. 5 Comparison of clustering performance based on objective function values

由圖4可以看出，TSPSO_FCM 算法逐漸減小并且趨于穩(wěn)定。單獨(dú) FCM 算法的目標(biāo)函數(shù)值在逐步減小的過(guò)程中沒(méi)有出現(xiàn)反復(fù)情形，這也說(shuō)明基于梯度下降的 FCM 算法在樣本空間進(jìn)行聚類時(shí)就有可能陷入局部極小且目標(biāo)函數(shù)值較大。

3.3 參數(shù)變化對(duì)聚類效果的影響

針對(duì)在飛行軌跡壓縮過(guò)程中不同的誤差閾值、聚類數(shù)會(huì)直接影響最終的聚類結(jié)果。因此，在其他各類參數(shù)不變的情況下，取不同的誤差閾值距離和聚類數(shù)進(jìn)行試驗(yàn)，結(jié)果如表2。

表2 誤差閾值距離和聚類數(shù)的影響Table 2 Influence of error threshold distance and number of clusters

由表2可以看出，當(dāng)誤差閾值距離逐步增大時(shí)且對(duì)比未進(jìn)行軌跡壓縮時(shí)，使用 TSPSO_FCM算法得到的目標(biāo)函數(shù)值相比FCM算法更優(yōu)，聚類效果更好。將飛行軌跡簇k取4類和5類進(jìn)行聚類時(shí)，TSPSO_FCM算法得到的的目標(biāo)函數(shù)值均小于k取3類時(shí)，這是因?yàn)殡S著聚類個(gè)數(shù)的增多，各類間的緊湊性變大，導(dǎo)致目標(biāo)函數(shù)值變小。由圖6可知，當(dāng)軌跡簇分為4類時(shí)，從東南方向進(jìn)場(chǎng)的航空器飛行軌跡分為兩簇，符合該終端區(qū)的標(biāo)準(zhǔn)進(jìn)場(chǎng)程序；當(dāng)軌跡簇分為5類時(shí)，很明顯地看到識(shí)別出了一種不同于機(jī)場(chǎng)傳統(tǒng)進(jìn)場(chǎng)程序的新的飛行模式。經(jīng)驗(yàn)證可知，該飛行模式下的航空器飛行為繁忙時(shí)段需調(diào)配航空器安全間隔，由管制員直接引導(dǎo)航空器進(jìn)場(chǎng)所形成。

圖6 不同聚類數(shù)對(duì)應(yīng)的飛行軌跡簇Fig. 6 Flight trajectory cluster corresponding to different cluster numbers

因此，該方法通過(guò)分析參數(shù)的變化對(duì)聚類結(jié)果的影響，有效地挖掘出了該機(jī)場(chǎng)終端區(qū)的航空器進(jìn)場(chǎng)飛行模式，為該機(jī)場(chǎng)由傳統(tǒng)程序向基于性能導(dǎo)航(PBN)程序轉(zhuǎn)變提供重要參考。

4 結(jié) 論

1)在保持較低失真率的情況下，采用基于時(shí)間比的自頂向下算法壓縮飛行軌跡降低了計(jì)算開銷，提高了軌跡數(shù)據(jù)分析能力。

2)建立基于多維特征屬性的相似性模型，并結(jié)合組合優(yōu)化算法改進(jìn)FCM算法，避免了單獨(dú)使用FCM聚類易于陷入局部最優(yōu)的缺陷。試驗(yàn)結(jié)果表明，改進(jìn)的FCM算法具有很強(qiáng)的全局搜索能力，降低了FCM算法對(duì)初始值的敏感度，可以得到較優(yōu)的滿意解。

3)飛行軌跡聚類時(shí)參數(shù)的變化會(huì)直接影響聚類算法對(duì)航空器飛行模式的挖掘能力，對(duì)于重新設(shè)計(jì)和優(yōu)化傳統(tǒng)飛行程序具有重要參考價(jià)值。