李善梅，王雨鑫，雷青蕾，宋思霓，王超

(1.中國民航大學(xué)，空中交通管理學(xué)院，天津 300300；2.浙江省通用航空產(chǎn)業(yè)發(fā)展有限公司，低空飛行服務(wù)中心，杭州 311600)

0 引言

目前，空中交通擁堵相態(tài)的識別主要依靠管制員的經(jīng)驗認知，復(fù)雜交通情況下?lián)矶伦R別的準確性難以保障。如何從科學(xué)角度，利用實際交通數(shù)據(jù)和管制員認知，準確識別空中交通擁堵，并在此基礎(chǔ)上，挖掘復(fù)雜空中交通現(xiàn)象背后隱含的交通擁堵時空演化規(guī)律，成為亟需解決的問題。終端區(qū)是空中交通運行的主要瓶頸區(qū)域，終端區(qū)內(nèi)的交通流具有高密度性和強匯聚性等特點，易引發(fā)交通擁堵，嚴重制約著整個空中交通運行的順暢。因此，準確識別終端區(qū)交通擁堵狀態(tài)，是制定有效擁堵管理措施的前提。

以往有關(guān)終端區(qū)交通流特性的研究多集中于理論建模與仿真分析。陳勇等[1]提出終端區(qū)容量建模的基本假設(shè)，并據(jù)此推導(dǎo)出終端區(qū)的基本容量模型。黎新華等[2]建立終端區(qū)交通跟馳模型，并提出基于跟馳穩(wěn)定的終端區(qū)容量評估方法。張洪海[3]提出基于航跡運行(Trajectory Based Operation,TBO)模式下截點直飛方式與融合點方式進場交通流優(yōu)化模型。MENON等[4]建模描述了空中交通環(huán)境中自定義的流量計數(shù)空域系統(tǒng)要素，建立空中交通流歐拉模型。張洪海等[5]運用Netlogo 仿真系統(tǒng)推演空中交通流基本參數(shù)變化趨勢。IDRISSI 等[6]利用時間分離、進近階段速度調(diào)整和點合并系統(tǒng)(PMS)保證交通流暢通。由于空中交通系統(tǒng)是一個復(fù)雜的開放系統(tǒng)，交通運行受天氣、軍航活動及通信導(dǎo)航監(jiān)視(Communication,Navigation,Surveillance,CNS)設(shè)備等不同因素的影響，建模仿真方法難以將所有因素均加以考慮，并且，在建模仿真過程中需設(shè)置理想的假設(shè)條件，導(dǎo)致模型仿真結(jié)果與實際運行存在一定偏差，影響研究結(jié)果的準確性，難以體現(xiàn)真實交通流的運行狀況。

隨著空中交通數(shù)據(jù)采集越來越便捷，越來越多的航跡數(shù)據(jù)可以獲取得到，而這些數(shù)據(jù)正是實際交通系統(tǒng)在多種影響因素綜合作用下的真實反映。因此，基于數(shù)據(jù)驅(qū)動的空中交通流特性研究，越來越受到人們的關(guān)注，并取得了一定的研究成果。TAN 等[7]從復(fù)雜網(wǎng)絡(luò)角度提出一種深度學(xué)習(xí)方法，用于識別擁塞情況。許炎等[8]結(jié)合交通流參數(shù)關(guān)系基本圖與終端區(qū)空中交通運行方式和管制規(guī)則劃分交通流相態(tài)，并建立交通流模型。ZHANG 等[9]利用測地線距離度量飛行軌跡相似度，采用改進的譜聚類算法對軌跡樣本進行分類，并利用基于最小生成樹的改進算法提取典型交通流。HU等[10]計算離散時間損耗隊列的流量保護方程的凈量變化，進一步確定流量與密度之間的關(guān)系。JIANG 等[11]從復(fù)雜網(wǎng)絡(luò)的角度描述空中交通狀況，使用獨立分量分析(ICA)在線監(jiān)測空中交通擁堵。

現(xiàn)有數(shù)據(jù)驅(qū)動的空中交通擁堵識別方法主要為聚類方法，此類方法完全從航跡數(shù)據(jù)出發(fā)劃分交通狀態(tài)，未考慮管制員的經(jīng)驗知識，因此，識別結(jié)果存在不被管制員接受的情況。近年來，核密度估計方法越來越廣泛地應(yīng)用于各領(lǐng)域的類別劃分問題中。通過采用核密度估計將航跡數(shù)據(jù)和管制員經(jīng)驗相結(jié)合，不僅考慮數(shù)據(jù)分布情況，同時，將管制員對于交通態(tài)勢總體認知經(jīng)驗考慮進來，使識別結(jié)果更容易被管制員接受。

鑒于此，本文將交通數(shù)據(jù)與管制經(jīng)驗相結(jié)合，提出一套基于交通流基本圖理論與核密度估計算法的空中交通流擁堵識別方法，在此基礎(chǔ)上，挖掘終端區(qū)交通擁堵的時空分布規(guī)律。首先，基于終端區(qū)實際航跡數(shù)據(jù)識別空中交通流；接著，基于基本圖理論挖掘表征交通流狀態(tài)的特征參數(shù)，分析各參數(shù)之間相互關(guān)系；然后，將核密度估計方法應(yīng)用于交通流基本圖當中，采用數(shù)據(jù)與經(jīng)驗雙重驅(qū)動，獲取交通流相態(tài)劃分結(jié)果；最后，挖掘不同相態(tài)下的航跡分布特征及交通擁堵的時空分布規(guī)律。與現(xiàn)有空中管理人員憑借經(jīng)驗對空中交通相態(tài)進行分類相比，本文工作有助于提升空管人員對終端區(qū)交通相態(tài)掌握的科學(xué)性和準確性，以便采取更有效的空中交通擁擠管理措施，降低管制員工作負荷，提升空中交通運行效率。

1 空中交通流識別

終端區(qū)為航空器起降的繁忙區(qū)域，其內(nèi)部交通具有交通結(jié)構(gòu)復(fù)雜和航空器密度高等特點。受不確定因素的影響，航空器在飛行過程中，管制員會對其速度、高度及航向進行調(diào)整，導(dǎo)致其飛行航跡并非嚴格按照終端區(qū)內(nèi)部事先劃定的進離場飛行程序飛行，航跡空間分布較為散亂，導(dǎo)致無法直接提取終端區(qū)內(nèi)各條交通流的航跡數(shù)據(jù)。因此，需要識別終端區(qū)交通流，以便后續(xù)繼續(xù)挖掘空中交通流的演化特性。

航跡聚類是目前交通流識別的主要手段，其目的為從混亂的軌跡數(shù)據(jù)中將相似航跡歸為一類，其中，軌跡的相似性度量為航跡聚類的關(guān)鍵基礎(chǔ)。常見的度量方法有：歐式距離和hausdorff距離等。不同數(shù)據(jù)結(jié)構(gòu)需選擇相適應(yīng)的相似度度量滿足聚類要求。由于終端區(qū)交通密度大，航跡具有強機動性與差異化等特征，基于歐式距離的度量方法可以較好地適應(yīng)軌跡數(shù)據(jù)結(jié)構(gòu)。同時，DBSCAN(Density-Based Spatial Clustering of Applications with Noise)聚類算法可以對任意形狀的數(shù)據(jù)集進行聚類，而K-Means 等聚類算法一般只適用于凸數(shù)據(jù)集；且DBSCAN聚類結(jié)果沒有偏倚，而K-Means等聚類結(jié)果對初始值較為敏感。故本文選用基于歐氏距離的DBSCAN聚類方法對重采樣后的航跡樣本進行交通流聚類劃分。依據(jù)航跡聚類的劃分結(jié)果，可將樣本中所有航跡劃分為t個類別，每一類別視為1條交通流。

采用2019年10月13日～26日北京終端區(qū)進場航跡數(shù)據(jù)，對航跡數(shù)據(jù)進行插值及重采樣等預(yù)處理工作，通過線性插值將原航跡數(shù)據(jù)填充為1 s 時間間隔的航跡點序列，通過重采樣工作將各條航跡轉(zhuǎn)變?yōu)榫哂邢嗤桔E點數(shù)據(jù)的航跡，便于后續(xù)航跡聚類時計算各航跡之間的距離。在此基礎(chǔ)上，使用基于歐式距離的DBSCAN 聚類方法，對終端區(qū)內(nèi)部約9483 條重采樣后的航跡P={p1,p2,…,p9483} 進行聚類操作。DBSCAN算法步驟如下。

Step 1 初始化參數(shù)，設(shè)定鄰域半徑(ε)和最小鄰居數(shù)(κ)。

Step 2 選擇一個未被訪問過的航跡pi。

Step 3 查找pi的ε-鄰域內(nèi)的所有航跡，即與pi歐式距離小于等于ε的航跡，若ε-鄰域內(nèi)的航跡數(shù)量大于等于κ，則將pi標記為中心航跡，并將其ε-鄰域內(nèi)的航跡添加到一個新的類中；否則，將pi記為噪聲數(shù)據(jù)。

Step 4 對于新簇中的每條航跡，如果其鄰域內(nèi)還有未被分配到任何類的航跡，則將其添加到當前類中。

Step 5 重復(fù)Step 2～Step 4，直到所有航跡均被訪問過。

北京終端區(qū)進場程序如圖1 所示，有5 條進場程序，故設(shè)定航跡聚類數(shù)目為5，聚類結(jié)果如圖2所示。得到5 條交通流，其走勢和程序圖基本一致，不同顏色代表不同交通流。

圖1 北京終端區(qū)飛行進場程序Fig.1 Beijing terminal area flight approach procedure

圖2 航跡聚類結(jié)果Fig.2 Trajectory clustering effect

2 基于基本圖的空中交通流三參數(shù)挖掘

2.1 空中交通流特征參數(shù)提取

交通流特征是構(gòu)成交通流的各條航跡基本屬性的集體表現(xiàn)。當進入終端區(qū)的航空器數(shù)量較多，交通處于擁堵情況下，終端區(qū)管制員經(jīng)常會減慢航空器的下降速度，實行雷達引導(dǎo)，拉開航空器之間的間隔(導(dǎo)致航空器在終端區(qū)內(nèi)飛行距離變長，從而相對速度降低)?；谏鲜龉苤菩袨樘攸c，考慮采用流入交通量、相對速度和高度變化率表征空中交通流狀態(tài)。

首先，挖掘單條航跡的基本屬性，具體如下。

(1)航跡時長tji

tji為第j條交通流中航空器i，自進入終端區(qū)的起始航跡點至飛離終端區(qū)的末端航跡點之間的飛行時間。

(2)航跡距離dji

dji為第j條交通流中航空器i，自進入終端區(qū)的起始航跡點至飛離終端區(qū)的末端航跡點之間的飛行距離。

(3)相對速度vji

vji為第j條交通流標稱航跡距離(航空器所屬交通流中心航跡的長度)與該交通流中航空器i航跡時長的比值。

(4)高度變化率hji

hji為第j條交通流中航空器i，單位時間內(nèi)上升/下降的距離。

然后，基于上述單個航跡的基本屬性，定義交通流的特征向量ej(t)=[Nj(t),Hj(t),Vj(t)]，具體含義如下。

(1)Nj為統(tǒng)計間隔t內(nèi)第j條交通流中流入終端區(qū)的交通流量，即

式中：Ii為0-1 變量，航空器i在時段t之前不在終端區(qū)內(nèi)，在時段t內(nèi)進入終端區(qū)時，Ii=1；否則，Ii=0。

(2)Hj為統(tǒng)計間隔t內(nèi)第j條交通流中航空器的平均高度變化率，即

(3)Vj為統(tǒng)計間隔t內(nèi)第j條交通流中航空器的平均相對速度。

2.2 空中交通流基本圖分析

基于第1節(jié)交通流聚類結(jié)果，以終端區(qū)西南方向DUGEB 入口交通流為例，進行空中交通流的基本圖分析?；谠摲较蚪煌?788 條航跡數(shù)據(jù)，統(tǒng)計每30 min該交通流的流量、平均相對速度和平均高度變化率，繪制三參數(shù)時序圖，如圖3 所示?？梢园l(fā)現(xiàn)，平均高度變化率與平均相對速度具有較強的波動性，而流量特征序列呈現(xiàn)出較強的周期性。下面進一步分析不同特征間的相關(guān)關(guān)系。

圖3 三參數(shù)時序圖Fig.3 Timing diagram for three parameters

繪制交通流三參數(shù)的基本圖，如圖4 所示，圖4(a)～圖4(c)為流量、相對速度與高度變化率的關(guān)系圖。為便于對空中交通流基本圖的解釋，圖4(d)給出了經(jīng)典的交通流基本圖。

圖4 空中交通流基本圖Fig.4 Fundamental diagram of air traffic flow

(1)從圖4(a)可以看出：流量和相對速度近似服從拋物線關(guān)系，隨著流量的增加，平均相對速度整體呈下降趨勢，交通趨于擁擠，數(shù)據(jù)點密度增加；當流量較小時，相對速度較大，數(shù)據(jù)點分布較分散。

(2) 從圖4(b)可以看出：流量和高度變化率也近似服從拋物線關(guān)系，隨著流量的增加，高度變化率呈上升趨勢，交通趨于擁擠，數(shù)據(jù)點密度增加；當流量較小時，高度變化率絕對值較小，且數(shù)據(jù)點分散?？梢?，對于空中交通流基本圖，較難發(fā)現(xiàn)如同交通流基本圖呈現(xiàn)出的拋物線右半部分非常擁堵情況。這是由于終端區(qū)內(nèi)部航空器的運行完全由管制員進行全局性質(zhì)的管理，在扇區(qū)內(nèi)，航空器過多或飛行態(tài)勢過于復(fù)雜時，管制員往往提前采取限流措施，避免嚴重擁堵情況的發(fā)生。

(3)從圖4(c)可以看出：平均高度變化率與平均相對速度呈現(xiàn)較為嚴格的線性關(guān)系，隨著交通擁堵程度的增加，平均相對速度和高度變化率均呈現(xiàn)下降的趨勢，且中間部分數(shù)據(jù)點分布集中，兩端較為分散。

(4) 從圖4(d)可以看出：地面交通領(lǐng)域常用流量、密度和速度作為交通流的三參數(shù)，并根據(jù)其相互關(guān)系表征交通擁堵狀態(tài)。由于地面交通中車輛主要在二維空間沿道路行駛，而終端區(qū)內(nèi)的航空器主要在三維空間運行，具有雷達引導(dǎo)、變高度和變速度等行為。因此，空中交通流的運行特點與地面交通有所不同，不能直接套用交通流基本圖理論分析空中交通。通過擬合空中交通流高度變化率、流量和相對速度三參數(shù)兩兩特征間的關(guān)系發(fā)現(xiàn)，其呈現(xiàn)出與交通流基本圖相一致的分布特征，因此，本文采用這3個參數(shù)描述空中交通流。

3 基于核密度估計的空中交通擁堵識別

空中交通流相態(tài)劃分具有一定的模糊性，空管實際運行中，基于經(jīng)驗的相態(tài)劃分方法具有較大的主觀性。現(xiàn)在流行的聚類方法主要從數(shù)據(jù)驅(qū)動角度進行相態(tài)劃分，識別結(jié)果較難被管制員認可。因此，有必要將數(shù)據(jù)與管制員的經(jīng)驗相結(jié)合，即考慮管制員對交通態(tài)勢的總體認知，提高相態(tài)識別結(jié)果的可用性。空中交通流相對速度與高度變化率的線性相關(guān)性為采用核密度估計方法進行交通相態(tài)劃分提供了條件。因此，本文擬采用高斯核密度估計方法劃分交通流相態(tài)。

3.1 高斯核密度估計原理

高斯核密度估計方法屬于數(shù)據(jù)驅(qū)動方法，避免了主觀經(jīng)驗影響，無需概率密度分布的假設(shè)。該方法基本思想為采用平滑的峰值函數(shù)(“核”)擬合觀察到的數(shù)據(jù)點，從而模擬真實的概率分布曲線。該方法認為，當已知某事件的概率分布時，若某數(shù)值在觀察中出現(xiàn)，則可以認為這個數(shù)的概率密度較大，這個數(shù)值周圍數(shù)值的概率密度也較大，與該數(shù)值相差較多的數(shù)值的概率密度則較小。

假設(shè)x1,x2,…,xn為n個離散隨機樣本，xi為第i個觀測數(shù)據(jù)。x為一個具有d維特征的樣本點，x∈?d，其概率密度函數(shù)f(x)未知，則此概率密度函數(shù)的高斯核密度估計為

式中：h為帶寬(窗寬或平滑系數(shù))；n為樣本容量；K(·)為核函數(shù)，并且滿足條件為

式中：u為K(u)中的自變量，即K的參數(shù)；μ2為核函數(shù)K(u)的二階矩。

求解MISE 最小化有兩種方法。第1 種，采用交叉驗證方法，使用部分數(shù)據(jù)作為訓(xùn)練集估計fh(x)，剩余數(shù)據(jù)作為測試集，用于以某種方式估計MISE。第2種，當N→∞的極限下，尋找最佳帶寬的解析公式，MISE為

式中：o為所需的函數(shù)項的階數(shù)小于或等于其他已知函數(shù)項的階數(shù)。

略去式(8)中的高階項，則得到漸進積分均方誤差 AMISE(Asymptotic Mean Integrated Square Error)為

對式(10)求偏導(dǎo)數(shù)可得到使AMISE 取最小值時的帶寬，此帶寬h 為最優(yōu)帶寬。

3.2 基于高斯核密度估計的空中交通流擁堵識別流程

基于高斯核密度估計原理，建立基于高斯核密度估計的交通流相態(tài)劃分方法。首先，基于核密度估計獲得平均高度變化率H和平均相對速度V聯(lián)合分布的核密度函數(shù)；然后，尋找水平曲線輪廓，將H-V散點劃分為兩種狀態(tài)，曲線內(nèi)部的點為典型狀態(tài)，曲線外部的點為非典型狀態(tài)，即分別對應(yīng)交通流的典型行為和非典型行為。具體步驟如下。

Step 1 核密度函數(shù)計算

利用高斯核密度估計算法，計算平均高度變化率與平均相對速度的聯(lián)合概率密度函數(shù)。

Step 2 曲面切割閾值的確定

Step 3 水平曲線分割

基于求解得到的z*值，在z*值處切割曲面(H,V)，從而定義一個水平曲線c(H,V)，該曲線內(nèi)點的總數(shù)占總體的概率1-α部分，視為典型狀態(tài)下的點，曲線外的點視為非典型狀態(tài)下的點。

Step 4 相態(tài)劃分閾值確定

曲線將數(shù)據(jù)點劃分成不同部分，各部分對應(yīng)不同相態(tài)。由于平均變化率與平均相對速度呈線性關(guān)系，因此，可基于曲線劃分結(jié)果確定平均相對速度閾值，用于后續(xù)交通流相態(tài)的識別。提取Step 3得到的水平曲線兩端端點，分別提取各端點周圍最近的n個點(本文以n=10 為例)，計算平均相對速度的平均值，作為劃分該端點兩端相態(tài)的閾值。獲得平均相對速度的兩個閾值，用來將交通流劃分為3個相態(tài)。

3.3 北京終端區(qū)交通擁堵識別算例分析

基于上述方法，對北京終端區(qū)DUGEB 方向交通流進行高斯核密度估計。采用該交通流的平均相對速度和平均高度變化率數(shù)據(jù)，繪制其3D 渲染圖，如圖5 所示?？梢姡麄€分布區(qū)域主要集中在兩種模式上，即高密度模式和低密度模式，大部分區(qū)域?qū)儆诟呙芏饶Ｊ?，表明該交通流處在高密度模式下的時間較多。

圖5 變化率-相對速度聯(lián)合分布核密度估計Fig.5 Decline rate-relative velocity joint distribution kernel density estimation diagram

依據(jù)交通流擁堵識別方法，分別令α等于0.1，0.2，0.3，0.4，0.5，識別交通流狀態(tài)，并統(tǒng)計不同狀態(tài)下的平均相對速度，如表1 所示?？梢钥闯?，當α越大，自由態(tài)平均相對速度越小，擁堵態(tài)平均相對速度越大，平穩(wěn)態(tài)平均相對速度范圍則越小。

表1 不同α 下各狀態(tài)的平均相對速度閾值Table 1 Average relative speed threshold of each state under different α

通過咨詢北京終端區(qū)資深管制員，采用專家打分法讓專家基于經(jīng)驗給出不同相態(tài)之間的平均相對速度閾值，統(tǒng)計各專家對速度閾值打分的平均值，得到綜合評分結(jié)果。將專家評分結(jié)果與表1中各α值下各相態(tài)的相對速度統(tǒng)計結(jié)果進行比較發(fā)現(xiàn)，當α=0.1 時，聚類出來的各相態(tài)平均相對速度分布情況與專家打分最為接近。因此，選用α=0.1的聚類結(jié)果進行后續(xù)研究?？梢?，相較于典型統(tǒng)計學(xué)確定α閾值方法，結(jié)合專家經(jīng)驗的方法使計算結(jié)果更好地符合管制員對交通態(tài)勢的實際認知，聚類結(jié)果更易被管制員所接受。

當α=0.1 時，相態(tài)劃分結(jié)果如圖6所示。水平曲線將數(shù)據(jù)點大體分為兩類：一類是，圈內(nèi)的點，相對速度和高度變化率的值處于中等；另一類是，圈外的點，主要分布于圈的左上側(cè)和右下側(cè)，左上側(cè)的相對速度和高度變化率較大，右下側(cè)的相對速度和高度變化率較小。可見，核密度估計法成功地將交通流劃分為3個相態(tài)：相對速度和高度變化率較高的自由態(tài)、相對速度和高度變化率中等的平穩(wěn)態(tài)與相對速度和高度變化率均較低的擁堵態(tài)。其中，位于右上方圈外的點，可視為異常點，相態(tài)分析時不做考慮。

圖6 基于核密度估計的劃分結(jié)果Fig.6 Partitioning results based on kernel density estimation

進一步，確定各相鄰相態(tài)之間平均相對速度的閾值。自由態(tài)與平穩(wěn)態(tài)之間的閾值為8.9 km·min-1，平穩(wěn)態(tài)與擁堵態(tài)之間的閾值為6.5 km·min-1。因此，可得到交通流相態(tài)判別標準為

4 終端區(qū)不同交通相態(tài)的航跡分布特性

為挖掘不同相態(tài)的交通特征，同時，驗證基于核密度估計的相態(tài)劃分方法的有效性，本文從飛行航跡和管制行為的角度，進一步挖掘不同相態(tài)下交通流的微觀特征。

提取不同交通流相態(tài)下的航跡數(shù)據(jù)，繪制水平剖面和垂直剖面航跡圖，如圖7所示?？梢钥闯?，3種相態(tài)下的航空器航跡具有如下特征。

圖7 3種相態(tài)下航空器水平剖面和垂直剖面航跡Fig.7 Flight path of horizontal section and vertical section of aircraft in three states

(1)自由態(tài)

航空器水平和垂直剖面的航跡分布均較為整齊，航空器數(shù)量較少，管制員大多指揮航空器直飛到跑道入口。航跡相似性較高，復(fù)雜性和沖突態(tài)勢較低，可實現(xiàn)連續(xù)下降。

(2)平穩(wěn)態(tài)

航空器水平和垂直剖面的航跡分布整齊性較自由態(tài)有所下降，航空器數(shù)量適中，管制員大多指揮航空器有一個較小范圍的雷達引導(dǎo)飛行，直飛情況較少，航跡相似性較自由態(tài)有所下降，航跡復(fù)雜性有所增強。

(3)擁堵態(tài)

航空器水平和垂直剖面的航跡分布整齊性在3種相態(tài)中最差，航空器數(shù)量最多，管制員往往指揮航空器有一個較大范圍的雷達引導(dǎo)飛行，基本沒有直飛情況。航跡差異性很大，部分航空器存在繞飛的情況，高度剖面上存在多階段分布下降情況，航空器之間的復(fù)雜性最高。

為進一步從定量的角度挖掘3 種交通流相態(tài)的航跡特點，分別統(tǒng)計不同相態(tài)下航跡的飛行時長、飛行距離和轉(zhuǎn)向次數(shù)，并進行對比分析。

(1)飛行時間對比分析

不同相態(tài)下，航跡的平均飛行時長、最小飛行時長、最大飛行時長及其方差如表2 所示，并繪制不同相態(tài)下飛行時長的箱型圖，如圖8(a)所示?？梢钥闯?，自由態(tài)下航空器的平均飛行時間最短，擁堵態(tài)的平均飛行時間最長。隨著交通流擁堵的嚴重性不斷增強，航空器的平均飛行時間隨之增長，自由態(tài)的平均飛行時間為1240.9 s，且各航空器飛行時間分布較為集中；平穩(wěn)態(tài)的平均飛行時間為1469.8 s，分布范圍較暢通態(tài)有所擴大；擁堵態(tài)的平均飛行時間為1827.5 s，大約是自由態(tài)的1.5 倍，且各航空器飛行時間分布較為發(fā)散，方差最大。

表2 3種相態(tài)飛行時間對比Table 2 Comparison of flight time in three states

圖8 航空器飛行特征箱體圖Fig.8 Aircraft flight characteristics box diagram

(2)飛行距離對比分析

不同相態(tài)下，航跡的平均飛行距離、最小飛行距離、最大飛行距離及其方差如表3 所示，并繪制不同相態(tài)下飛行距離的箱型圖，如圖8(b)所示。

表3 3種相態(tài)飛行距離對比Table 3 Comparison of flight distance in three states

可以看出，自由態(tài)下，航空器的飛行距離最短，擁堵態(tài)的飛行距離最長。與圖7 中呈現(xiàn)出的結(jié)論一致，自由態(tài)下，航空器直飛現(xiàn)象較為明顯；平穩(wěn)態(tài)下，航空器沿進場飛行程序飛行較多；擁堵態(tài)下，航空器繞飛現(xiàn)象較為嚴重。隨著交通流擁堵的嚴重性不斷增強，航空器的飛行距離隨之增長，自由態(tài)的平均飛行距離為167595.9 m，且各航空器飛行時間分布較為集中；平穩(wěn)態(tài)的平均飛行距離為183705.5 m，分布范圍較暢通態(tài)有所擴大；擁堵態(tài)的平均飛行距離為211845.9 m，大約是自由態(tài)的1.3 倍，且各航空器飛行距離分布較為發(fā)散，方差最大。

(3)轉(zhuǎn)向次數(shù)對比分析

統(tǒng)計不同交通流相態(tài)下航空器的平均轉(zhuǎn)向次數(shù)依次為：自由態(tài)，1.2593次，大部分航空器只在由進場方向轉(zhuǎn)向跑道方向時實施1次機動轉(zhuǎn)彎；平穩(wěn)態(tài)，3.5231次，多數(shù)航空器實施標準進場程序，在對準跑道方向前進行3次機動轉(zhuǎn)彎；擁堵態(tài)，5.4545次，大多航空器的飛行軌跡較為復(fù)雜，存在“繞飛”情況，在對準跑道方向前多次實施機動轉(zhuǎn)彎飛行。進一步，繪制平均轉(zhuǎn)向次數(shù)箱型圖如圖8(c)所示，可以發(fā)現(xiàn)，自由態(tài)下轉(zhuǎn)向次數(shù)最低，與圖7中呈現(xiàn)的該相態(tài)下管制員大多指揮航空器進行直飛相一致；擁堵態(tài)下，管制員給航空器發(fā)布的雷達引導(dǎo)指令較多，且引導(dǎo)范圍較大，因此，該相態(tài)下的平均轉(zhuǎn)向次數(shù)最多。從箱型圖可以看出，平穩(wěn)態(tài)下轉(zhuǎn)向次數(shù)的分布最為集中，主要是由于該相態(tài)下，管制員大多引導(dǎo)航空器沿進場程序飛行，轉(zhuǎn)向次數(shù)較為固定。

5 終端區(qū)交通擁堵的時空分布特性

基于交通流相態(tài)的劃分結(jié)果，從時間和空間兩個維度研究終端區(qū)交通擁堵的分布規(guī)律。

5.1 網(wǎng)格劃分及參數(shù)統(tǒng)計

首先，對終端區(qū)進行網(wǎng)格劃分，為便于分析，基于航跡分布情況，將終端區(qū)人為劃分成4×5 的矩形，如圖9(a)所示。對劃分所得的網(wǎng)格進行編號，從左到右，從上到下，網(wǎng)格編號從1 開始遞增，其中，沒有航跡信息或者航跡信息較少的網(wǎng)格，屬于無效網(wǎng)絡(luò)，無需對其編號，中間部分(例如，網(wǎng)格3和網(wǎng)格8中間的網(wǎng)格)的航空器趨于落地，基本沿進近程序運行，管制員不再對其有較多調(diào)整，大部分航班已移交給塔臺管制員管理，故此部分不是終端區(qū)管制員的重點關(guān)注區(qū)域，如圖9(b)所示，編號網(wǎng)格區(qū)域涵蓋了終端區(qū)管制員關(guān)注的主要區(qū)域。

圖9 終端區(qū)劃分結(jié)果Fig.9 Terminal area meshing results

然后，計算各網(wǎng)格每30 min為一時間片的平均相對速度V，并根據(jù)式(13)得到各網(wǎng)格不同時段的交通相態(tài)。

5.2 交通擁堵的時間分布特征

基于各網(wǎng)格擁堵相態(tài)的識別結(jié)果，統(tǒng)計終端區(qū)各時段不同相態(tài)網(wǎng)格的比例，并繪制分布圖，如圖10 所示，發(fā)現(xiàn)終端區(qū)交通相態(tài)具有時間不均衡特性?？梢钥闯?，0:00-6:00終端區(qū)內(nèi)各網(wǎng)格均處于自由或平穩(wěn)狀態(tài)，并以自由態(tài)為主；6:00-9:00 平穩(wěn)態(tài)比例逐漸上升，交通態(tài)勢由自由態(tài)開始向擁堵態(tài)過渡；9:00開始，擁堵態(tài)比例開始上升，平穩(wěn)態(tài)和自由態(tài)比例均有所下降，出現(xiàn)較多擁堵區(qū)域，在14:00時，擁堵最為嚴重。終端區(qū)整體交通狀況在20:00開始逐漸好轉(zhuǎn)。

圖10 終端區(qū)交通各狀態(tài)時間分布Fig.10 Terminal area traffic states time distribution diagram

進一步挖掘各個網(wǎng)格交通狀態(tài)隨時間的變化情況，以網(wǎng)格3、網(wǎng)格5、網(wǎng)格7 及網(wǎng)格8 為例，給出其交通狀態(tài)分布情況，如圖11所示?？梢钥闯?，網(wǎng)格3和網(wǎng)格5的擁堵情況較網(wǎng)格7和網(wǎng)格8嚴重，基本上，5:00-23:00 均有擁堵現(xiàn)象發(fā)生。若按擁堵比例超30%時段視為擁堵時段，則網(wǎng)格3比較嚴重的時段主要集中在10:30-12:00、14:30-15:30 和16:30-17:00，共3 h。網(wǎng)格5 比較嚴重的時段集中在11:30-12:00、14:00-14:30、18:30-19:00和19:30-20:30，共2.5 h。因此，當處于15:00-17:00時段時，管制員應(yīng)重點關(guān)注網(wǎng)格3的交通態(tài)勢；當處于18:00-21:00時段時，管制員應(yīng)重點關(guān)注網(wǎng)格5 的交通態(tài)勢；當處于10:00-12:00 和14:00-15:00 時段時，管制員應(yīng)同時關(guān)注網(wǎng)格3和網(wǎng)格5的交通態(tài)勢。

圖11 典型網(wǎng)格交通狀態(tài)時間分布Fig.11 Time distribution of typical grid traffic conditions

由于經(jīng)過網(wǎng)格7 和網(wǎng)格8 的航班量相對較少，其擁堵態(tài)勢沒有網(wǎng)格3和網(wǎng)格5嚴重。網(wǎng)格7和網(wǎng)格8 在0:00-8:00 大部分時間均處于自由態(tài)和平穩(wěn)態(tài)，從8:00 開始有少部分時間出現(xiàn)擁堵現(xiàn)象，8:00-14:00自由態(tài)較少，平穩(wěn)態(tài)較多，擁堵態(tài)較少，14:00-19:00擁堵情況比較頻繁，擁堵情況在19:00以后開始緩解。

5.3 交通擁堵的空間分布特征

基于各網(wǎng)格擁堵相態(tài)的識別結(jié)果，統(tǒng)計終端區(qū)交通相態(tài)在空間上的分布情況發(fā)現(xiàn)，其具有空間不均衡特性。部分時段各個網(wǎng)格交通狀態(tài)分布情況如圖12 所示?？梢钥闯?，交通擁堵相態(tài)主要集中在網(wǎng)格2、網(wǎng)格3、網(wǎng)格5、網(wǎng)格7、網(wǎng)格8及網(wǎng)格9，這些區(qū)域也是終端區(qū)進場交通流匯聚之處，交通態(tài)勢較為復(fù)雜。另外，隨著時間的變化，擁堵網(wǎng)格區(qū)域也發(fā)生著變化。例如，13:30-14:00 時段，網(wǎng)格3 和網(wǎng)格5 為擁堵態(tài)，網(wǎng)格7、網(wǎng)格8 及網(wǎng)格10 為平穩(wěn)態(tài)，其余網(wǎng)格為自由態(tài)；14:00-14:30時段，網(wǎng)格3和網(wǎng)格5 已由擁堵態(tài)轉(zhuǎn)變?yōu)槠椒€(wěn)態(tài)，網(wǎng)格8 和網(wǎng)格10由平穩(wěn)態(tài)轉(zhuǎn)為自由態(tài)，網(wǎng)格13 由自由態(tài)轉(zhuǎn)為平穩(wěn)態(tài)；14:30-15:00時段，網(wǎng)格3和網(wǎng)格7由平穩(wěn)態(tài)轉(zhuǎn)為擁堵態(tài)，網(wǎng)格8 和網(wǎng)格9 由自由態(tài)轉(zhuǎn)為擁堵態(tài)?？梢姡苤茊T應(yīng)重點關(guān)注網(wǎng)格2、網(wǎng)格3、網(wǎng)格5、網(wǎng)格7、網(wǎng)格8及網(wǎng)格9區(qū)域，及時發(fā)現(xiàn)擁堵的產(chǎn)生，并對其進行疏導(dǎo)，以期提高終端區(qū)整體的運行效率。

圖12 終端區(qū)交通狀態(tài)空間分布Fig.12 Spatial distribution of traffic state in terminal area

6 結(jié)論

本文提出一種新的空中交通流擁堵識別方法及其時空演化分析思路。將核密度估計方法和交通流基本圖理論與管制員的經(jīng)驗知識相結(jié)合。傳統(tǒng)方法僅依靠聚類算法劃分交通狀態(tài)，忽略了管制員對交通態(tài)勢的整體認知，導(dǎo)致識別結(jié)果容易受到質(zhì)疑。本文通過將數(shù)據(jù)分布情況與管制員的經(jīng)驗相結(jié)合，使識別結(jié)果更符合實際，易被管制員所接受。此外，本文從航跡分布特征和基于網(wǎng)絡(luò)劃分的擁堵時空分布特征等角度進行分析，揭示了終端區(qū)交通擁堵的時空分布規(guī)律，提高管制員對交通態(tài)勢把握的準確性。

研究結(jié)論如下。

(1)考慮終端區(qū)管制員雷達引導(dǎo)、調(diào)速度及調(diào)高度等管制行為特點，基于基本圖理論挖掘能夠表征空中交通流相態(tài)的交通流的特征參數(shù)：交通量、平均相對速度和平均變化率，且平均相對速度與平均變化率呈現(xiàn)較強的線性關(guān)系，可用于識別交通擁堵態(tài)勢。

(2)將航跡數(shù)據(jù)與管制知識相結(jié)合，基于高斯核密度估計方法識別終端區(qū)交通流擁堵狀態(tài)，得到劃分交通相態(tài)的相對速度閾值。擁堵相態(tài)下的航跡存在雷達引導(dǎo)等復(fù)雜管制行為，飛行軌跡、飛行時長、飛行距離和轉(zhuǎn)向次數(shù)均較高，說明本文擁堵識別方法的有效性，識別結(jié)果更符合管制員認知。

(3)北京終端區(qū)交通擁堵具有時間和空間分布的不均衡特性。時間上，該終端區(qū)擁堵時段主要分布于9:00、14:00和19:00左右；空間上，擁堵常發(fā)區(qū)域主要位于中部區(qū)域，并挖掘出不同區(qū)域易發(fā)生交通擁堵的時段。這些擁堵時空分布規(guī)律可提高管制員對擁堵態(tài)勢把握的準確性與科學(xué)性。