張啟錢，徐禮鵬，張賽文

(南京航空航天大學 民航學院，江蘇 南京 211106)

0 引 言

在航班總量不斷增長的外部環(huán)境下，民航空域資源供給不足所引起的航班和空中交通管理問題日益凸顯。從空中交通管制扇區(qū)運行動態(tài)特征出發(fā)，對扇區(qū)展開分類研究，通過分析不同類別扇區(qū)的運行特性，可找出當前制約我國空管效能提升的瓶頸問題。

國內外空管業(yè)都在不斷探索空域和扇區(qū)的管理方法及相應的自動化技術，以期促進空中交通的安全、高效的運行。R.CHRISTIEN等[1]通過K-means方法根據交通流指標和空域結構指標進行扇區(qū)分類研究，改進了扇區(qū)容量評估模型，提高了管制員工作負荷評估和扇區(qū)容量預測的精度；B.HILBURN[2]在COCA (complexity and capacity)項目中同樣提出建立扇區(qū)容量模型對扇區(qū)進行分類的構想。NextGen的研究項目正在探索“通用空域”的概念[3]，為支撐“通用空域”的研究，有學者從扇區(qū)交通運行模式出發(fā)，識別不同的扇區(qū)，用于評估“通用空域”概念所需的運行相似性[4]。A.CHO等[5]基于扇區(qū)間共同結構特征，對360個高空扇區(qū)進行了整體式與分步式分類；S.MALAKIS等[6]通過決策樹和樣本場景的分類規(guī)則，對空中交通場景開展分類研究；高興[7]提出了面向動態(tài)容量的區(qū)域扇區(qū)分類，按靜態(tài)和動態(tài)因素進行扇區(qū)分類；王紅勇等[8]、溫瑞英等[9]運用K-means等方法對空中交通復雜程度進行評級，劃分出不同的扇區(qū)復雜性模式；董襄寧等[10-11]基于扇區(qū)結構和飛行流量等特征，發(fā)現區(qū)域和近扇區(qū)的復雜性指標中靜態(tài)指標的共性都強于動態(tài)指標；叢瑋等[12]從扇區(qū)航班分布和動態(tài)兩個維度著手，使用K-medoids算法進行了基于交通態(tài)勢的扇區(qū)聚類方法研究。

綜上，現有研究對扇區(qū)動態(tài)運行特征分析較少，且選用硬聚類算法進行類別劃分，存在一定程度的不合理性。基于此，筆者通過ADS-B(automatic dependent surveillance-broadcast)數據，提取具有時變性的動態(tài)因素，并歸納出適應扇區(qū)分類研究的普適性指標，采用改進的核模糊C均值聚類算法等機器學習算法，對選取的97個樣本扇區(qū)展開管制扇區(qū)分類研究。

1 指標構建

扇區(qū)是空中交通管理的基本單元。因扇區(qū)靜態(tài)物理結構因素與動態(tài)交通流因素不同，不同扇區(qū)運行特征存在差異，計算獲取能反映扇區(qū)運行特征和狀態(tài)的指標數據是扇區(qū)分類研究的基礎。

1.1 ADS-B航跡數據概述

航空器通過ADS-B系統(tǒng)向所有安裝ADS-B接收機的航空器或者地面臺分發(fā)實時空中交通數據，即ADS-B數據。數據中包含了時間戳記、航空器呼號標識、航空器類型、航空器三維位置信息(飛行高度、經度、維度)、航行方向、起飛機場和目的機場等信息。直接反映真實運行狀況的ADS-B數據能夠為基于動態(tài)運行特征的空中管制扇區(qū)分類研究提供科學的數據來源。

1.2 空中管制扇區(qū)分類指標

航空器數量是空中交通特征中最常被引用、研究和評估的一個指標。但航空器密度指標只能片面反映扇區(qū)運行情況和管制員工作負荷。因此，衍生了動態(tài)密度、內秉性等更多的復雜性相關研究，提出了更多的復雜性指標如：在不同飛行階段(平飛、上升、下降)的航空器比例、扇區(qū)內航向發(fā)生改變的航空器數量(比例)、管制飛行時長(里程)、交通混合程度、航空器沖突次數等。這些復雜性指標為空中交通管制環(huán)境提供了比單一的航空器密度指標更加全面的評估方法。

綜合扇區(qū)自身物理結構特點和動態(tài)交通流特征，通過對ADS-B數據解析計算，基于已有扇區(qū)分類相關研究成果[1,7,12-13]，考慮數據形式、可獲取性、用戶可接受性等，選取9個指標數據用于扇區(qū)分類研究，如表1。

表1 初始指標集Table 1 Initial index set

2 分析方法

在系統(tǒng)聚類法或K-means聚類法等傳統(tǒng)硬聚類算法中，每個樣本被嚴格劃分到某一類中，在處理一些復雜數據的類別劃分上存在不合理性。區(qū)別于硬聚類，模糊劃分允許每個樣本屬于具有不同隸屬度的多個聚類，隸屬度的取值范圍是[0,1]，一個樣本對于所有類別的值之和為1?；诖?，采用改進的核模糊C均值聚類算法對扇區(qū)分類進行研究。

2.1 核模糊C均值聚類(KFMC)算法

隨著模糊C均值聚類(fuzzy C-means，FCM)算法在各個研究領域的不斷應用，學者們發(fā)現FCM算法對噪聲和異常值敏感。此外，FCM算法采用歐式距離構造目標函數，聚類效果很大程度上依賴于聚類樣本的空間分布。針對這些問題，學者們提出了基于核函數改進的模糊C均值聚類KFCM(kernel fuzzy C-means)算法。KFCM算法的基本原理是將原空間的樣本，通過核函數映射到高位特征空間后，再利用FCM算法指導聚類[15]。

通過非線性映射Φ:x∈X?Rd→Φ(x)∈F?RH(d?H)，將樣本X從一個低維空間映射到高維空間F。對于所有的x,y∈X，X?Rd，若函數K滿足K(x,y)=〈Φ(x),Φ(y)〉=Φ(x)TΦ(y)，則稱函數K為核函數，其中Φ是從輸入空間X到特征空間F的映射， 〈Φ(x),Φ(y)〉 為Φ(x),Φ(y)的歐式內積。

KFCM算法的目標函數如下：

(1)

采用拉格朗日乘子法求解目標函數JKFCM的極小值，推導出KFCM算法的隸屬度和模糊聚類中心迭代式如下：

(2)

(3)

2.2 改進的核模糊C均值聚類(GA-KFCM)算法

KFCM算法通過將樣本映射到高維空間，突出樣本間的差異性，可獲得較好的聚類效果，使算法對噪聲和孤立點具有較好的魯棒性，克服了FCM算法不適合多種數據分布的缺陷。但KFCM算法仍然沒從根本上解決FCM算法對初始化聚類中心敏感的問題，聚類中心的隨機初始化會影響到算法聚類性能，不能保證聚類結果收斂到全局最優(yōu)解?；诖?，筆者基于遺傳算法改進KFCM算法得到GA-KFCM算法，結合遺傳算法的全局搜索尋優(yōu)能力在特征空間中搜索適當的聚類中心，彌補KFCM算法缺陷，獲得理想的聚類效果[16]。

2.2.1 GA-KFCM算法適應度函數

遺傳算法中的適應度函數可用于評價個體的性能優(yōu)劣，也被稱為評價函數。適應度函數通過對所求解問題的目標函數進行分析，根據優(yōu)化方向構造而成。聚類是將給定的樣本數據，按某種目標函數劃分成若干個類的過程[14]，對于KFCM算法而言，目標函數取得最小值時對應的聚類效果最好。根據目標函數值越小個體適應度越大原則，定義GA-KFCM算法適應度函數如下：

(4)

2.2.2 GA-KFCM算法基本算子

遺傳算法通過選擇，模擬自然界中的“優(yōu)勝劣汰”，適應度高的個體被遺傳到下一代的概率大。目前主要的選擇方法有：輪盤賭選擇法、排序選擇法和隨機聯(lián)賽選擇法等。采用非線性排序選擇法，按照適應度大小對個體進行排序，然后根據順序分配個體被選擇的概率。非線性選擇概率分布為：

f(z)=q(1-q)z-1,z=1,2,…,n

(5)

式中：q∈(0,1)為指定參數。

遺傳算法通過交叉操作可得到適應度較高的新個體，并通過變異操作保持種群多樣性。一般交叉概率建議取值范圍為0.40～0.99，變異概率建議取值范圍為0.001～0.100。但固定的概率取值會直接影響算法效果，因此采用自適應思想，動態(tài)確定交叉概率fcr(t)和變異概率fmr(t)為：

(6)

(7)

式中：pcr0為初始交叉概率；pcm0為初始變異概率；T為最大進化代數；t為當前進化代數。

2.2.3 GA-KFCM算法流程

GA-KFCM算法具體流程如下：

步驟1設定高斯核函數尺度參數σ，確定聚類數c(1

步驟2設置遺傳代數t=0，模糊聚類代數k=0。隨機產生n個聚類中心矩陣，形成初始種群。

步驟3計算種群中每個個體JKFCM，f(U,V)，f(z)，fcr(t)，fmr(t)。

步驟4對t代種群進行選擇、交叉、變異操作，形成第t+1代種群。

步驟5判斷是否滿足終止條件。不滿足則返回步驟3；如滿足，則解碼獲得最優(yōu)聚類中心。

步驟6根據最優(yōu)初始聚類中心對樣本數據進行聚類分析。

步驟7輸出結果。

2.3 GA-KFCM算法優(yōu)化驗證實驗

選取UCI機器學習測試數據庫中的IRIS和WINE數據集作為輸入樣本，針對FCM，KFCM和GA-KFCM 3種算法進行分析驗證實驗，通過對比不同算法結果，評價算法的聚類性能，檢驗算法的改進效果。設置GA-KFCM算法參數為：n=50，T=100，pcr0=0.6，pcm0=0.1，δ=0.001，σ=150，k=100，c=3，m=2，ε=0.000 01。將每種算法在IRIS和WINE數據集上運行30次，以平均結果進行對比。3種算法的正確率和算法迭代次數比對結果如圖1。

根據實驗結果得出如下結論：

1) GA-KFCM算法在兩套標準數據集上的聚類準確率較FCM算法提高了4.62%和26.96%。但由于GA-KFCM算法在聚類之前需要先通過遺傳算法進行最佳聚類中心尋優(yōu)，導致算法總體迭代次數較多。

圖1 FCM、KFCM、GA-KFCM算法性能對比結果Fig. 1 Comparison results of FCM、KFCM、GA-KFCM algorithmperformance

2)引入了核函數的GA-KFCM算法和KFCM算法在準確度上優(yōu)于FCM算法，在高維數據集上更加明顯。

3)KFCM算法在準確度上的表現與GA-KFCM算法相近，但實驗過程中發(fā)現KFCM算法輸出的聚類結果不穩(wěn)定，算法準確率存在一定浮動。而GA-KFCM算法根據遺傳算法輸出的搜索結果作為初始聚類中心，每次執(zhí)行后GA-KFCM算法的結果均是在最佳目標值附近，具有較高的魯棒性。

2.4 最佳聚類數確定

模糊聚類算法需要預先指定聚類數目c，不同c值對應的模糊劃分結果和質量不同。聚類有效性指標能夠正確評價聚類結果，幫助確定最佳聚類數。采用整體輪廓系數和核空間的VKXB指標來確定最佳聚類數。

輪廓分析可以用來研究聚類簇間的分離程度，輪廓系數定義為：

(8)

式中：ai是第i個點與同簇內的其他點之間的平均距離，用于量化簇內凝聚度；bi是第i個點與其他簇內各點之間的平均距離，用于量化簇間分離度；輪廓系數S(i)的值落在 [-1,1]的范圍之間。

S(i)值越大，說明第i個點分類越合理，負值表示樣本可能指定給了錯誤的簇。因此，可以通過計算樣本數據在不同的聚類數目下的聚類結果所對應的整體輪廓系數值，來反映當前聚類的效果，整體輪廓系數值定義為：

(9)

Xie-Beni有效性函數是目前應用最廣泛的模糊聚類有效性指標，定義為：

(10)

在GA-KFCM算法中，需要將VXB推廣到核空間，對應的核化Xie-Beni有效性函數VKXB為[17]：

(11)

當VKXB取得極小值時可得樣本數據最佳劃分。

3 扇區(qū)聚類實例分析

根據2019年民航空域發(fā)展報告，2019年日均流量前20位扇區(qū)有90%集中在華東和中南地區(qū)，約扇區(qū)94%呈增量趨勢。華東和中南地區(qū)是我國最繁忙的空域，具有很強的典型性，選取華東和中南地區(qū)107 個管制扇區(qū)內2019年4月11日—18日中00:00—24:00的ADS-B雷達航跡數據作為數據源。通過整理解析ADS-B數據，計算獲取得到能夠客觀刻畫管制扇區(qū)整體運行特征的動態(tài)交通流指標，最終獲得具有有效指標數據的97個管制扇區(qū)(89個區(qū)域扇區(qū)，8個進近扇區(qū))，采用GA-KFCM算法進行扇區(qū)分類研究。

3.1 聚類變量選取

基于ADS-B數據，提取管制扇區(qū)的9個特征因素作為扇區(qū)分類的初始變量。為提高數據質量，避免因選取的原始指標變量在數量級和單位不同，對結果可靠性造成負面影響，需采用 Z-Score規(guī)范化將原始數據作正態(tài)分布處理[18]。

通過主成分分析(PCA)，將多個變量轉換為少數幾個主成分(綜合變量)，用于反映原始變量的大部分信息[19]，消除指標信息重疊的影響，從而得到聚類算法輸入特征向量。

表2 主成分提取結果Table 2 Principal component extraction results

表2中的主成分能解釋的原樣本中約80% 的信息量，依據主成分選擇原則，選擇這3個主成分進行后續(xù)研究。

實驗得到的各主成分系數表達式數據見表3。根據表3中數據分析可得：主成分1在總占據計數、進入計數、離開計數和總飛行時間4個指標上有明顯且相近的正載荷，為0.948～0.981，定義為扇區(qū)運行交通流量成分；主成分2在巡航交通流比例上有較明顯的正載荷(0.916)，在爬升交通流比例及下降交通流比例上有明顯的負載荷，分別為-0.713和-0.621，定義為航空器飛行姿態(tài)成分；主成分3可定義為扇區(qū)的空域結構成分。

表3 主成分系數表達式矩陣Table 3 Principal component coefficient expression matrix

9個特征指標在3個主成分方向上的投影情況如圖2。

圖2 主成分空間投影圖Fig. 2 Principal component space projection

3.2 確定最佳聚類數

聚類數c取[3,10]區(qū)間內的整數，對樣本數據進行8組聚類分析實驗，并分析不同聚類數c下對應的聚類有效性評價指標數值，再結合相應的空管背景，綜合確定聚類最佳數目。不同聚類簇數取值下聚類結果有效性評價指標結果如圖3。由圖3可知：當c=5時，整體輪廓系數值最大為0.466 1，VKXB=0.215 4取得極小值為。根據最佳聚類簇數，確定基于動態(tài)運行特征的空中管制扇區(qū)分類研究的分類個數為5。

圖3 不同聚類數目對應的聚類有效性指標值Fig. 3 Cluster validity index values corresponding todifferent cluster numbers

3.3 扇區(qū)聚類結果

根據確定的最佳聚類數，結合聚類結果的可解釋性，確定c=5。在實際應用中，模糊指數m的最佳取值范圍為 [1.5, 2.5 ][20]，選擇常見取值m=2。算法經35次迭代至目標函數收斂，得到最終聚類結果，將97個樣本扇區(qū)劃分為運行特征鮮明的A、B、C、D、E類共5類。各類樣本扇區(qū)指標數據箱線分析如圖4。

圖4 扇區(qū)運行特征指標數據箱線分析Fig. 4 Box line analysis on sector operation characteristic index data

除箱線圖的標準元素外，樣本扇區(qū)總體的每個指標變量的中位數都以實心三角形符號表示，從而可以將簇中某個指標變量的分布與所有分析的樣本扇區(qū)的中位數進行比較；盒須外的數據被認為是異常值，用黑色圓圈表示。每個運行特征指標的標準化值已在y軸上顯示，消除了表示中的比例問題。

扇區(qū)總體運行特征體現出的復雜性從A到E逐漸降低，各類扇區(qū)詳細運行特征有：

A類扇區(qū)：交通流量大，總占據架次數遠大于所選取的樣本扇區(qū)中值，該類扇區(qū)的日均航空器進入架次為900～1 500架次，平均飛行時間中等，交叉匯聚點個數較多。扇區(qū)內航空器主要高度改變情況為垂直方向上有高度改變的航空器比例占進入架次的72%以上，屬于兼顧大流量和中高垂直方向復雜性的扇區(qū)，扇區(qū)整體運行狀況繁忙，呈復雜態(tài)勢，管制員負荷較高。

B類扇區(qū)：從交通流量維度看，該類扇區(qū)總占據架次數和進入架次數偏低。巡航交通流遠低于樣本扇區(qū)中值，意味著扇區(qū)內存在大量垂直方向上的高度調配，占到進入該扇區(qū)航空器約80%的比例，且集中于航空器的爬升。同時，扇區(qū)內航空器平均飛行時間較短，導致調配裕度偏小，空域環(huán)境復雜交叉匯聚點數量較多。根據分析得出，B類扇區(qū)的運行瓶頸并不是扇區(qū)內的交通量，而在于大量的垂直過渡帶來的垂直方向上調配難度和潛在沖突，呈現出的高復雜性態(tài)勢和管制壓力。

C類扇區(qū)：交通流量在所選樣本扇區(qū)里處于中等水平，進入該類扇區(qū)的航空器超過六成產生了垂直方向上的高度變化，且主要集中為下降交通流。平均飛行時間較大，且扇區(qū)內交叉匯聚點個數較少，管制調配裕度較大。整體運行狀態(tài)由于中等的流量和垂直方向復雜性，仍然呈現較復雜和繁忙態(tài)勢。

D類扇區(qū)：交通流量小，平均小時進入流量為15～25架次，瞬時管制壓力較??；在高度過渡方面，該類扇區(qū)所管制的航空器中，平飛航空器所占比例和垂直方向產生上高度改變的航空器所占比例持平。盡管扇區(qū)交通流量低，垂直方向復雜性一般，但扇區(qū)內航空器平均飛行時間較短，管制員存在一定調配壓力，扇區(qū)運行狀態(tài)呈中等復雜態(tài)勢。

E類扇區(qū)：在交通流量上分布較廣，扇區(qū)日均進入流量為450～1 300架次不等；從箱線圖巡航交通流的統(tǒng)計可以看出，航空器基本都是以平飛姿態(tài)通過該類扇區(qū)，產生高度過渡的航空器不足1/3；交叉點數量中等，平均飛行時間較長，調配難度小；運行壓力主要與交通流量相關，垂直方向上的復雜性較低，總體運行復雜態(tài)勢一般，管制壓力較小。

3.4 分類結果驗證

根據2019年空域發(fā)展報告統(tǒng)計數據，流量增幅最大的扇區(qū)分別是成都區(qū)域08扇區(qū)(26.51%)、合肥區(qū)域01扇區(qū)(11.59%)和合肥區(qū)域04扇區(qū)(11.22%)。

選取合肥的區(qū)域01、區(qū)域02、區(qū)域03和區(qū)域04號扇區(qū)中2019年4月11—18日的歷史運行數據進行仿真實驗。通過AirTop快時仿真軟件評估得到的管制員工作負荷和扇區(qū)容量結果，結合實際運行情況，對管制扇區(qū)分類研究結果進行驗證。

根據3.3節(jié)中扇區(qū)分類結果，合肥區(qū)域01、04號扇區(qū)分別屬于E類和A類扇區(qū)，合肥區(qū)域02、03號扇區(qū)同屬于D類扇區(qū)。圖5為合肥區(qū)域01號～04號扇區(qū)1 h管制員工作負荷與航空器進入架次關系擬合圖，圖中實線為擬合曲線，虛線為置信區(qū)間。

圖5 管制員工作負荷與航空器架次關系擬合圖Fig. 5 Fitting diagram of relationship between controller workloadand aircraft sorties

如圖5(a)可知：合肥區(qū)域01號扇區(qū)小時進入架次集中分布在25～40架次，根據管制員工作負荷閾值70%回歸得到扇區(qū)小時容量為53架次，且從圖上可以看出管制員工作負荷波動范圍較小。8 d歷史最高小時進入架次為45時，僅對應63%管制員工作負荷，表示扇區(qū)運行仍有潛力可挖。因此，該扇區(qū)符合E類扇區(qū)的復雜度一般，工作負荷主要來源于交通量的運行特征，扇區(qū)分類結果合理。

由圖5(b)、圖5(c)可知：根據管制員工作負荷閾值70%回歸得到的區(qū)域02、03號扇區(qū)小時容量分別為31、28架次。扇區(qū)歷史運行中存在著大量超過管制員工作負荷閾值70%的進入架次數據點，處于高負荷運行狀態(tài)；同一架次數下對應的管制員工作負荷區(qū)間波動明顯。由此可以得出，合肥區(qū)域02、03號扇區(qū)運行的復雜程度和管制壓力及調配難度與B類扇區(qū)一致，分類結果合理科學。

由圖5(d)可知：合肥區(qū)域04扇區(qū)小時進入流量集中分布在33～50架次之間，扇區(qū)交通流量較大。管制員工作負荷存在一定程度上的波動，說明在同一進入架次數值下對應的運行情況有較大不同。根據管制員工作負荷閾值70%回歸得到的扇區(qū)小時容量為43架次，對比歷史運行情況，扇區(qū)超容明顯。因此，將合肥區(qū)域04扇區(qū)分為大流量且運行較復雜的A類扇區(qū)是合理的。

4 結 論

基于扇區(qū)ADS-B數據，選取的我國較為繁忙具有典型性的中南地區(qū)和華東地區(qū)的區(qū)域和進近共97個管制扇區(qū)，采用基于遺傳算法改進的核模糊C均值聚類算法進行分類研究，研究結果表明：

1)不同交通流和物理結構反映出來的運行特征，對應著扇區(qū)內指揮調配的差異?；谶\行特征的扇區(qū)分類研究，有助于減少管制員在相似運行特征扇區(qū)之間轉換工作時所需的培訓。同時為管制扇區(qū)的分類專項管理提供理論基礎。

2)結合聚類有效性指標和聚類的可解釋性，將97個管制扇區(qū)劃分為5種不同類型，扇區(qū)運行特征在不同扇區(qū)間存在明顯差異，表明了模糊聚類算法相比于傳統(tǒng)硬聚類算法在處理扇區(qū)分類這種現實問題方面的有效性。同時采用遺傳算法尋找最優(yōu)聚類中心的方法，對傳統(tǒng)的核模糊聚類算法進行了優(yōu)化，提高了算法劃分的準確程度。

3)選取合肥4個區(qū)域扇區(qū)進行扇區(qū)管制員工作負荷與進入架次的擬合分析，各扇區(qū)的運行特征與擬合結果的匹配程度驗證了扇區(qū)分類結果的合理性。各類扇區(qū)不同的交通流反映出來的運行動態(tài)特征，是將航空器空間行為和扇區(qū)繁雜程度及管制員工作負荷進行關聯(lián)比較的有效途徑。結合管制運行實際情況，對運行態(tài)勢復雜的高負荷運行扇區(qū)可以進一步優(yōu)化其空域結構、科學調整扇區(qū)數量和使用方式，以期增加調配裕度，降低管制難度，使各扇區(qū)間工作負荷更趨均衡。同時對尚未飽和的扇區(qū)，深挖內部潛力，充分利用空中交通管制資源。以求實現容量、流量、人員之間的最佳匹配，全面提升扇區(qū)空域運行品質。