管制-飞行状态相依网络演化过程

2021-10-20李昂聂党民温祥西韩宝华曾裕景

航空学报 2021年9期

李昂，聂党民，温祥西,*，韩宝华，曾裕景

1. 空军工程大学空管领航学院，西安 710051

2. 中国人民解放军95429部队，昆明 650000

在管制系统的运行过程中，若能对管制系统未来一段时间内的运行态势的变化情况进行较为准确的预测，可以帮助管制人员提前规划，将矛盾端口前移，大大提高系统运行安全的稳定性。由于航空器的高速运动，运行态势瞬息万变，随着预测时间的增加，其发散性和波动性可能越大，如何分析与预测管制系统的演化过程成为近年来研究的热点。

国内外学者对空管领域的预测进行了大量研究，主要包括空中交通的预测和地面管制员有关方面的预测。张兆宁和张佳[1]提出了一种预测空中交通网络系统中大面积航班延误发生的方法；高旗等[2]提出了一种基于PSO-BP神经网络的多指标终端区拥堵等级预测模型；杨东玲[3]提出了一种基于自适应时空航迹聚类的航迹预测方法。不仅仅是空中交通的预测，国内外学者对地面管制员负荷预测、人因失误预测等方面也进行了一定的研究。刘继新等[4]筛选出了对管制员影响较大的7个显著相关指标，建立了基于Logistic的管制员应激程度预测模型；Vanderhaegen[5]提出了一种新颖的基于知识的人因失误预测学习方法；Corver等[6]研究了调查轨迹不确定性是否能调节交通冲突和工作负荷之间的关系。

广义的管制系统包括空域管理、空中交通服务和空中交通流量管理3个子系统，服务对象是飞行员。管制员依据这些子系统与飞行员进行直接的沟通，指挥航空器，提供服务。可以说管制员和空中的航空器是管制系统中最基础、最关键的组成部分。本文在分析管制系统时对其进行了简化，主要考虑管制员和空中航空器的关系，以它们的运行情况来代表整个管制系统的运行情况。

管制系统是一个由多个单元构成的复杂系统，复杂网络可以有效地对复杂系统进行分析，因而在其他的领域的预测研究中，复杂网络理论被广泛地运用。例如：在链路预测方面，刘大伟等[7]结合复杂网络理论提出了局部差异融合算法进行信息网络的链路预测；Liu等[8]提出了一个结合零模型的通用框架来量化加权复杂网络中拓扑结构、权重相关性和统计对链路预测的影响；Lü等[9]提出了一种无网络组织先验知识的通用结构一致性指标。在网络结构方面，Li等[10]提出了一个易感-感染-易感的复杂网络模型来预测网络流量；伍杰华等[11]提出了一个增强的朴素贝叶斯关系预测模型(ELNB)来建模局部子图结构。由于复杂网络理论在其他领域的成功应用，开始有学者将复杂网络理论运用到空管领域的相关研究中：王超和朱明[12]运用复杂网络方法对空中交通流时间序列的非线性特性进行了分析与预测；徐肖豪和李善梅[13]总结了近年来运用复杂网络在空中交通拥挤识别与预测领域进行的研究并指出了未来的研究方向。

而在现实生活中，每个网络都或多或少地与其他的网络间存在各种各样的关联，为了解决这类问题，Buldyrev等[14]首次提出了相依网络理论。同样的在管制系统中管制员需要对飞行员下达指令，飞行员在经过报告点时需要向管制员报告航空器的状态数据，两者间存在互相依赖的关系，所以本文准备在前期工作[15]的基础上用相依网络来表示管制系统。

在预测问题中，另一个重要的问题就是预测方法。对于同一个样本，使用不同的预测方法[16-19]往往会得到不同的预测结果。管制系统的演化有着较高的复杂性和不确定性，使用简单的机器学习算法很难对其演化过程进行较为准确的预测。为了使预测结果尽量接近实际情况，采用更加适用于复杂数据的深度学习算法[20-23]。本文首先使用管制-飞行状态相依网络代表管制系统，然后通过相依网络来分析管制系统的演化情况，最后使用深度学习算法对管制系统的演化过程进行预测。

1 基础理论

1.1 相依网络建模

本文对之前建立的管制-飞行状态相依网络模型进行完善。在文献[15]中，以航空器为节点，航空器之间的冲突关系为连边构建上层的飞行状态网络，以管制扇区的几何中心为节点，管制扇区之间的移交关系为连边构建下层的管制网络，在两层网络中间，则是以管制员对航空器的管制指挥情况建立相依边，以航空器节点的管制难度为相依边权重。在之前的模型中，飞行状态网络的层内边权只和航空器间的距离有关，而在实际中，两航空器间的相对速度也是评判其冲突情况的重要因素，所以本文在飞行状态网络的层内边权设置中加入了相对速度这一影响因素，使用迫近效应理论计算得到飞行状态网络的层内边权。

航空器对的迫近效应是指航空器对的汇聚/非汇聚态势对于管制难度的影响，而航空器位置、速度属性是分析这种影响的最基本要素[24]。航空器的位置和速度分别用P和V来表示，相对距离和相对速度分别用Dij和Vij来表示。

如图1所示，两机间的相对距离为

Dij=Pi-Pj=[xi-xj,yi-yj]

(1)

式中：xi和yi分别为航空器i的横坐标和纵坐标。

(2)

相对速度可以表示为

Vij=Vi-Vj

(3)

令Vij=[Δvx,Δvy]，由图2可得

图2 两机相对位置Fig.2 Relative position of two aircrafts

Δvx=Visinθi-Vjsinθj

(4)

Δvy=Vicosθi-Vjcosθj

(5)

则迫近率可用相对速度在两机间的连线上的分量表示为

(6)

由式(6)可看出，当(Vij·Dij)>0时,两机呈现发散态势，当(Vij·Dij)<0时,两机呈现汇聚态势。

飞行状态网络的层内边权设置同时考虑相对距离和相对速度，距离越近权重越大，迫近率越大权重越大，由于航空器间距离52 km时构成连边，而且航空器速度一般不会高于1 000 km/h，即迫近率最小为-2 000 km/h，为避免边权取值为负，加入了控制参数，层内边权ωij计算公式为

(7)

1.2 LSTM原理

长短期记忆(Long Short-Term Memory, LSTM)人工神经网络是一种具有记忆能力的循环神经网络(Recurrent Neural Network，RNN)，是深度学习常用的模型之一，它能够发现时间序列中的一些隐藏特征。LSTM神经网络在普通的RNN网络上进行改进，通过在隐藏层中引入记忆单元，并使用3个门控单元(遗忘门、输入门和输出门)来控制记忆单元的状态，从而解决了普通RNN因时间序列较长导致的学习效果较差的问题，使其具有长时间序列的记忆能力。LSTM神经元结构如图3所示。

图3 LSTM神经元结构Fig.3 Neuron structure of LSTM

LSTM的记忆单元与隐藏单元均记忆时间序列中的数据信息，而记忆单元中的数据受到3个门控单元的影响。

1) 遗忘门根据上一时刻的隐藏状态和当前输入对记忆单元中的部分无用信息进行删除，以减小记忆负荷，用公式表示为

ft=σ(WfH)+bf

(8)

H=[ht-1,xt]

(9)

式中：σ为sigmoid函数；Wf为遗忘门的权重；bf为偏置向量；ht-1为上一时刻的隐藏状态;xt为当前时刻的输入。

2) 输入门根据上一时刻的隐藏状态和当前时刻的输入控制当前的输入，将没有价值的信息过滤掉，选择性地向记忆单元中新增信息，用公式表示为

it=σ(WiH)+bi

(10)

(11)

当遗忘门和输入门计算结束后，由式(12)来更新记忆单元。

(12)

3) 输出门将上一时刻的隐藏状态、当前输入和更新后的记忆单元这三类信息汇总，得到当前的隐藏状态，用公式表示为

οt=σ(WοH)+bο

(13)

ht=οttanh(Ct)

(14)

式中：οt为输出值；Wο为输出门的权重；bο为偏置向量。

2 管制系统演化过程分析

2.1 指标选取

下面选取几个能够反应相依网络性能的拓扑指标，并进行简单介绍。

1) 平均节点度

航空器节点的度为与该航空器之间可能存在安全风险的航空器数量和与地面通信时可供其选择的管制扇区数量之和，管制节点的度为与该管制扇区存在移交关系的扇区数量和能够与其直接通信的航空器数量之和，而平均节点度则是所有节点的度的平均值，其数值可以反映管制系统中平均每架航空器周围与其存在冲突的航空器数量以及每位管制员需要监视的航空器数量。

(15)

2) 平均点强

加权的节点度之和即为节点点强si。节点的点强不仅能够反应与其相连节点数量的多少，还能反应其相邻节点对其造成的总影响。而平均点强则是所有节点点强的平均值，其数值可以反映管制系统中飞行员和管制员承受的平均压力，用公式表示为

(16)

3) 平均加权聚类系数

某一节点的所有邻居节点间实际相连的边数与理论上最多能够相连边数的比值叫做节点聚类系数。而加权聚类系数还考虑了节点间的权重，两节点间距离越近，权重越大，对加权聚类系数的贡献越大。航空器节点的加权聚类系数表示该航空器的周围航空器的聚集程度，管制节点的加权聚类系数表示该管制扇区以及相邻扇区内所有航空器之间的聚集程度，加权聚类系数越大，聚集程度越高，表达式为

式中：m、n分别为节点i的2个相邻节点。

则平均加权聚类系数则是所有节点加权聚类系数的平均值，其数值可以反映管制系统中航空器的聚集程度，表达式为

(17)

4) 网络密度

网络密度是网络中实际存在的边数与可容纳的边数上限的比值，在此处边数为层内连边与层间连边之和，其数值可以反映管制系统的饱和程度，用公式表示为

(18)

式中：ND为网络密度；L为网络中实际存在的连边数。

5) 网络效率

网络效率反映了网络的连通程度。任意2个节点间的效率表示为2个节点之间距离的倒数，而整个网络的效率为任意2个节点间效率的平均值，表示网络中任意一点联系到另一点需要的平均中转次数，其数值可以反映管制系统中管制员对空中航空器的整体管控力度，其表达式为

(19)

式中：NE为网络效率；dij为节点i和节点j间的最短路径。

这5项指标都能从不同的角度反映网络的整体性能，例如网络密度越大，说明网络中各节点之间的连边越多，节点间互相影响的程度越大，网络的复杂程度也就越高。

2.2 场景设置

由于目前实际运行中航空器飞行是以固定航路或航线为基础实施的，本文设置了固定航线飞行场景进行模拟。而为了提高空域资源的利用率，中国和一些欧美国家开始采用一些新技术，如基于航迹运行和自由飞行技术。在这种模式下，飞行员将不再需要按照划分好的高度层飞行，特别是在目前的军航训练中，在某一片空域内所有航空器都处于同一高度的情况是有可能存在的，因此本文还设置了自由飞行仿真场景。具体设定如下：

1) 固定航线飞行条件下的场景设置

以昆明长水机场附近的扇区空域为参照，构建的演化场景。该空域共分为6个扇区，包含26条航线，24个报告点，在实际中该空域横向最远距离约为1 500 km，纵向最远距离约为1 400 km。本文在仿真时对原始空域按1∶4的比例进行等比例缩小，初始航空器为28架，具体仿真场景如图4所示。图4中黑色数字表示各扇区编号，黑线实线代表扇区边界，黑色虚线代表进近管制区域边界，红线代表航线，黄线代表具有潜在冲突的航空器之间的连边，绿色方块代表航空器，黑色实心圆代表长水机场，蓝色三角形代表报告点。

图4 固定航线飞行场景设置Fig.4 Scene setting of flight on route

2) 自由飞行条件下的场景设置

通过MATLAB模拟一个300 km×300 km的空域并随机生成50架航空器。将演化场景设定如下：每隔30 s的时间间隔，1、3、7、9号4个管制扇区均有50%的概率进入一架航空器。初始的50架航空器和后来进入该空域的航空器，在刚出现在这片空域中时，航向任意，航速在600～800 km/h之间随机取值，演化场景如图5所示。

图5 自由飞行场景设置Fig.5 Scene setting of free flight

2.3 演化过程分析

由于自由飞行场景的演化过程更加复杂，本文以自由飞行的场景演化为例进行分析。

以该设置演化1 000 min(2 000次)，演化1 min 后、演化100 min后、演化500 min后和演化最终时刻的网络示意图如图6所示。

图6 自由飞行场景演化过程Fig.6 Evolution process of free flight scene

每次演化完成后，计算5项网络拓扑指标：平均节点度、平均点强、平均加权聚类系数、网络密度和网络效率，最终得到五项指标的时间序列，如表1所示。

表1 拓扑指标部分数据样本Table 1 Part of data samples of topological indicators

本文截取了5项指标400～600次的演化过程，如图7所示。由图7可知，5项指标的演化过程具有一定的相似性，大致上呈现出上升-下降-上升的变化趋势。为了更加具体地了解这5项指标之间的相关性，本文采用皮尔逊积矩相关系数(Pearson product-Moment Correlation Coefficient, PMCC)来判断各指标数据集之间是否为线性关系。该系数用ρ表示，取值范围在[-1,1]之间。若2个数据集之间的相关性越高，则|ρ|越接近1，反之，若2个数据集之间相关性越低，则|ρ|越接近0。

图7 5项指标数值变化Fig.7 Numerical changes of five indicators

(20)

式中：E(·)为数学期望；Cov(·)为协方差；Var(·)为方差。

由于各指标的数据存在数量级上的差异，对各指标数据进行标准化处理：

(21)

使用SPSS软件对上述指标数据进行相关性分析，计算各指标的PMCC值，分析各指标两两之间的相关性。相关性分析结果如表2所示。

由表2的相关性分析结果可知，这5项指标两两之间的相关性均较高，相对而言，平均节点度与其他4项指标之间的相关性较低。为了进一步分析这些指标中包含的信息量是否足够，判断其能否体现整个相依网络的综合性能，本文采用主成分分析法进行分析。主成分分析是一种在数据分析领域中广泛应用的降维技术，能够有效地将一组高维向量降维压缩成一组低维向量。本文利用SPSS软件对数据集进行主成分分析，通过比较第一主成分中的指标系数来判断各指标中所包含的信息量。将主成分数设置为1，各指标的分析结果如表3～表5所示。

表2 相关性分析结果Table 2 Results of correlation analysis

公因子方差表示主成分能够解释每个指标的程度，每个指标的初始值为1，提取值越大，主成分对指标的依赖性就越强。由表3可以看出，网络效率的提取值最高，平均加权聚类系数的提取值最低，但总的来看，5项指标的提取值都较高。

表3 公因子方差Table 3 Common factor variance

由主成分分析得到的各指标的特征值及贡献率如表4所示。本文取贡献率最高(64.236%)的第一主成分进行分析，得到第一主成分对应的成分矩阵，该矩阵可以反应出各指标在第一主成分中所占权重，如表5所示。

表4 主成分及其贡献率Table 4 Principal component and its contribution rate

综合表3和表5可以看出，这5项指标的提取值和成分值都处于一个较高的水平，包含的信息量充足，可以有效体现相依网络的整体性能，因此，本文用这5项指标来展现管制系统的演化过程是合理的。

表5 成分矩阵Table 5 Component matrix

3 管制系统演化过程预测

3.1 混沌性检验

在对时间序列进行预测之前，首先需要进行混沌性检验，判断其是否具有可预测性。对于每个时间序列，先由G-P算法计算其关联维d，再根据m≥2d+1得出嵌入维m，然后由快速傅立叶变换(FFT)计算平均周期T，最后根据T=(m-1)τ得到时间延迟τ，从而可由每个时间序列的嵌入维和时间延迟计算出最大李雅普诺夫指数，若最大李雅普诺夫指数大于0，则表明该时间序列具有混沌特征，即该时间序列具有可预测性[25]。对自由飞行和固定航线飞行场景下各指标时间序列的计算结果分别如表6和表7所示。

由表6和表7可知，在2种仿真场景下各指标时间序列的最大李雅普诺夫指数均大于0，说明在这2种仿真场景下各指标的时间序列均具有混沌性，是可预测的。

表6 嵌入维、时间延迟和最大李雅普诺夫指数(自由飞行)

表7 嵌入维、时间延迟和最大李雅普诺夫指数(固定航线飞行)

3.2 预测流程

管制系统演化过程预测流程如图8所示，主要步骤如下：

图8 预测流程图Fig.8 Prediction flow chart

步骤1在一段时间内，每隔0.5 min采集一次网络运行数据，分别进行相依网络建模。

步骤2计算不同时刻相依网络的各项指标，得到各项指标的时间序列。

步骤3将各指标的时间序列输入LSTM模型中进行训练，得到各自的预测模型。

步骤4将预测样本分别输入到各指标的预测模型中，输出各指标的预测结果。

4 仿真验证

4.1 自由飞行条件

将前1 900次的演化数据作为训练集对LSTM模型进行训练，后100次演化数据作为测试集输入训练好的预测模型中得到预测值和真实值，后100次演化数据的预测情况如图9所示。

图9 预测结果(自由飞行)Fig.9 Prediction results (free flight)

从整体上看，预测结果较为准确，尤其是图9(a)～图9(d)，预测结果与实际结果相差无几，依据预测结果，管制员可以对管制系统未来一段时间内的运行态势有一个大概的预估，根据实际情况可提前对部分航空器采取改航分流等调配措施，避免飞行事故征候的出现。

由图9(b)和图9(d)可知，在第47次演化时，相依网络的平均点强和网络密度这2个指标的真实值和预测值相差较大，且都是真实值高于预测值，而平均节点度、平均加权聚类系数和网络效率的值也都处于较高的水平。综合以上情况可知，当相依网络的复杂程度较高时，相依网络部分指标的预测值明显低于真实值，管制系统演化过程的预测精度会有所下降，且预测情况一般会好于实际情况。

为了进一步判断各项指标的预测精度，本文计算了各项指标的相对误差(误差值与真实值的比值)，如图10所示。

从图10中可以清楚地看到，平均节点度、平均加权聚类系数和网络效率这3项指标的相对误差一直在一个较小的范围内波动，不超过5%，而平均点强和网络密度这2项指标的相对误差变化幅度较大，但可以发现，只有演化次数在11次、21次、42次和47次时相对误差超过了20%，有96%的演化次数相对误差在20%以内。因此，从总体上看，这5项指标都具有较高的预测精度。

为突出LSTM算法相较其他算法的优势，另外采用贝叶斯和支持向量机算法对这5项指标的时间序列进行预测，计算每项指标的平均相对误差(相对误差绝对值的平均值)，将3种算法得到的结果进行对比，如表8所示。