基于改进KNN检测的ADS-B轨迹插补研究.pdf

基于改进 KNN 检测的 ADS-B 轨迹插补研究靳慧斌a，冯朝辉a，张召悦b，王志森b（中国民航大学 a.安全科学与工程学院；b.空中交通管理学院，天津300300）摘要院 针对广播式自动相关监视渊ADS-B袁automatic dependent surveillance broadcast冤地面接收站数据存在缺失点尧离群点等异常点的问题袁在建立多约束条件判别模型的基础上筛选 ADS-B 轨迹缺失点袁并对航迹进行插值袁采用标准化欧氏度量的改进邻近算法对整条轨迹进行异常点数据检测袁最后应用多重插补法对异常轨迹点进行修复遥 以实际运行数据为例袁验证对比了改进的 K 最近邻渊KNN袁K-nearest neighbor冤算法用于检测 ADS-B 异常轨迹点的准确率袁该方法可有效解决轨迹异常点检测疏漏问题袁为航空器 ADS-B 原始轨迹数据挖掘提供支撑遥关键词院 ADS-B 数据筛选曰异常点检测曰标准化欧氏度量曰改进邻近算法曰多重插补方法中图分类号院 V355文献标志码院 A文章编号院 1674原5590渊圆园23冤04原园园23原园6Research on ADS-B trajectory interpolation based on improved KNN detectionJIN Huibina,FENG Chaohuia,ZHANG Zhaoyueb,WANG Zhisenb(a.College of Safety Science and Engineering;b.College of Air Traffic Management,CAUC,Tianjin 300300,China)Abstract:Aiming at the problem of anomalous points such as missing points and outliers in the data of automatic dependentsurveillance broadcast(ADS-B)ground receiving station,the missing points of ADS-B trajectory are screened onthe basis of establishing a multi-constraint discriminant model,and the trajectory is interpolated.The improvedneighbor algorithm with standardized Euclidean metric is used to detect anomalous point data of the wholetrajectory.Finally,the multiple interpolation method is applied to repair the anomalous trajectory points.Takingactual operational data as an example,the accuracy of the improved K-nearest neighbor(KNN)for detectingADS-B anomalous trajectory points is verified and compared.This method can effectively solve the problem ofdetection omission of trajectory anomalous point and provide support for mining original trajectory data of aircraftADS-B.Key words:ADS-B data screening;anomalous point;standardized Euclidean metric;improved neighbor algorithm;multi鄄ple interpolation method收稿日期院 2021-09-06；修回日期院 2021-12-29基金项目院 中国民航大学研究生科研创新基金项目（10502756）作者简介院 靳慧斌（1976），男，天津人，副研究员，博士，研究方向为安全风险分析.广播式自动相关监视（ADS-B，automatic depende-nt surveillance broadcast）是航空器周期性广播自身位置及相关信息，供空管对其进行监视的一种监视技术1。ADS-B 数据传输过程中可能会因为定位设备故障、恶劣天气、飞行姿态变化等因素导致部分报文数据缺失或数据跳点，使地面站接收的 ADS-B 轨迹数据出现缺失点、离群点等异常点2。异常点的产生不利于系统和相关用户对轨迹位置数据进行应用分析，并有可能产生轨迹点间隙过大、错误位置告警、图形轨迹误判等影响航空安全水平的事件。为解决数据缺失问题，国内外学者从不同方向展开了研究。鞠时光等3从连续性及保形性两方面提出了用 姿 参数三次样条函数对等高线轨迹数据进行修补论证，有效地解决了绘制地形图过程中出现的等高线相交或缺失问题，但该方法只适用于线性函数。针对非线性数据插补修复，Zamani4采用二次样条函数处理函数边界非线性问题；Farhan 等5建立了一种随机多重插补（RMI，random multiple interpolation）方法来克服非线性轨迹不确定性问题。但两种方法对于处理大量数据具有局限性。为解决大数据插补局限性问第 41 卷第 4 期圆园23 年 8 月中 国 民 航 大 学 学 报JOURNAL OF CIVIL AVIATION UNIVERSITY OF CHINAVol.41 No.4August20232023 年 8 月中 国 民 航 大 学 学 报题，Peeters 等6评估了处理丢失数据的不同技术，得到了对于大数据插补问题使用预测均值匹配的多重插补（MI，multiple interpolation）是最佳方法的结论。为提升 ADS-B 数据完整性，邹文华等7建立了一种引入时间序列的基于贝叶斯常均值模型的数据增广算法，解决 ADS-B 轨迹点回跳现象。此外为确保插补数据精度，Sun 等8构建了一种使用最小二乘样条逼近（LSB，least square b-spline）的轨迹数据重建算法，能够避免因数据差值而导致插补误差过大问题。上述研究多侧重从单一角度对数据缺失问题进行插补分析，或以局部数据是否间断为依据进行插补验证，未考虑数据可能存在异常点的问题。针对这一情况，本文采用多约束条件筛选出航空器轨迹缺失点并将其清洗插入，然后利用改进后的 K 最近邻（KNN，K-nearestneighbor）方法检测 ADS-B 异常点，解决了ADS-B 样本数据量过大带来的聚类点权重分配不平衡而无法有效检测异常点的问题，并对轨迹进行多重回归插补修复，进一步提高轨迹数据完整性，从而保障航空安全监管水平。1基于 KNN 的轨迹数据异常点检测1.1缺失点筛选为确保经度、纬度和高度使用同一量纲（m），需要利用米勒投影完成经纬度和平面坐标的相互转换9。转换过程如图1 所示。图 1 中：ARP（aerodrome reference point）为机场基准点；L 为地球周长；W 为地球周长的平面展开，并将周长视为 X 轴；H 为 Y 轴，约等于地球周长一半（L/2）；mill 为米勒投影常数。ADS-B 地面接收机在同一时间接收不同航空器机载 ADS-B OUT 发出的报文信息是有限的10，中国民用航空局规定接收机标准是 600 架飞机的解码11，一台 ADS-B OUT 发射全类型报文为 4 次/s，则 ADS-B地面接收机每秒最少接收 2 400 条全类型报文。假设一台 ADS-B OUT 单位时间 T 内发射报文次数为 C，则每条报文时间间隔To=TC（1）从航路进入终端区的进近航路相对地面高度一般为 6 000 m 左右，无线电传播速度为 299 792 458 m/s，因此在终端区范围内，因距离而导致的接收延迟问题可以忽略。一条全类型报文解码平均时间长度为 Ta。地面接收机接收报文并完成解码，传入终端显示界面的时间间隔为 Tu，则同一航班的两条报文时间间隔应为Top=To+Ta+Tu（2）一条轨迹点的报文信息包括经度、纬度、高度、地速、航向角、时间，所以设轨迹点 P=Plo，Pla，Pal，Psp，Pha，Pt，轨迹点速度均为瞬时速度，因此将轨迹点之间的时间段位移作为平均速度位移，通过距离和时间差判断是否有缺失点，其判断公式如下Pt，i-Pt，i-1 Top异常Psp，i-1驻t+0.5 伊 驻v驻t屹Px缺失扇墒设设设设设缮设设设设设（3）式中：时间差 驻t=Pt，i-Pt，i-1；速度差 驻v=Psp，i-Psp，i-1；Px为两点之间距离；i 为轨迹点序号。其缺失点判断如图2 所示，其中 Pvt为同一架飞机 t 时间段内飞行的距离。当两点之间出现缺失数据，应插入空白轨迹点，其插入轨迹点数目满足1 臆 NwPt，i-Pt，i-1Top（4）1.2基于标准欧氏度量的离群点筛选KNN 算法具有最近邻特性，对于非线性数据筛选可避免偏离数据与正常数据出现较高拟合12-13，因此适用于终端区 ADS-B 数据检测。1.2.1数据处理在 1.1 节中已经完成缺失点的筛选和空白轨迹点将弧度转为实际距离将纬度从度数转为弧度配置参数 ARP、L、W、H、mill输入经纬度将经度从度数转为弧度输出直角坐标系 X、Y图 1经纬度米勒投影转换Fig.1Miller projection conversion of longitude and latitude图 2缺失点判断示意图Fig.2Diagram of missing point judgmentPxPiPvtPi-1PiPi-1Pt，i-Pt，i-1 Top缺失区域24-第 41 卷 第 4 期插入。但在计算轨迹点相似度时，需要赋予空白轨迹点标签，避免出现因无标签而导致分类错误的问题，通过下三角距离矩阵对空白轨迹点进行赋值，其矩阵如下W=0（t2-t1）3-Top33姨00（t3-t2）3-Top33姨00（tn-tn-1）3-Top33姨0杉删山山山山山山山山山山山山山山山山煽闪衫衫衫衫衫衫衫衫衫衫衫衫衫衫衫衫（5）式中 t 为轨迹点时间，当矩阵中元素 Wi，j 0 时，对应 Pt，i=0。因 ADS-B 轨迹在界面只显示二维坐标轨迹，所以只需要插补经度和纬度。1.2.2距离计算KNN 距离度量有两种，一种是欧氏度量（Eu原clidean metric），另一种是曼哈顿度量。欧氏度量用于衡量各点之间的绝对距离。如三维坐标点 M1（X1，Y1，Z1）和 M2（X2，Y2，Z2），则 M1和 M2两点之间的距离为D=（X2-X1）2+（Y2-Y1）2+（Z2-Z1）2姨（6）在 n 维空间中，两个 n 维向量 A=（X1，X2，Xn）、B=（Y1，Y2，Yn）间的欧式度量为D=ni=1移（Yi-Xi）2姨（7）曼哈顿度量是向量之间距离的范数。如向量 A=（X1，X2，Xn）和 B=（Y1，Y2，Yn）之间的曼哈顿距离为ni=1移Xi-Yi（8）选择欧氏度量或曼哈顿度量主要取决数据维度特征13。对于高维向量，曼哈顿度量比欧氏度量计算距离的效果更好；对于低维向量，欧氏度量要优于曼哈顿度量。因轨迹点 P经过筛选，只考虑经度、纬度和时间，所以本文选择欧式度量。欧式度量适用于样本向量的各个分量度量标准统一的情形。ADS-B 轨迹样本各分量权重相同，当时间与经纬度波动范围量纲值差距较大时，可能会引起各个分量分类不平衡问题14。因此可以通过标准化欧氏度量（standardized Euclidean metric）解决 Pt和 Plo、Pla分量差值问题。标准化欧氏度量是针对欧氏度量缺点的一种改进。标准化欧氏度量的思路为：既然数据各维分量的分布不一样，那先将各个分量都“标准化”到均值、方差相等15。假设样本集 X 的均值为 m，标准差为 s，X 的“标准化变量”表示为X忆=X-ms（9）标准化欧氏度量公式为d=ni=1移Xi-Yisi蓸蔀2姨（10）由式（10）可以得到不同轨迹点之间标准化欧氏度量为d=乙嗓ni=1移zla，nzt，n（Plo，i-Plo，i-1）2+zlo，nzt，n（Pla，i-Pla，i-1）2+zlo，nzla，n（Pt，i-Pt，i-1）2/zlo，nzla，nzt，n乙瑟12（11）式中 z 为分量度量，因 KNN 遵循最近邻原则，为避免同一航班轨迹点之间不同分量差值较大，导致线性方差误差过大，因此需要限制轨迹点方差样本数目i-n1臆i臆i+n2，n1为第 i 点前一个不相邻异常点位置，n2为第 i 点后一个不相邻异常点位置，将式（11）的变量进行如下改进，即zla，n=z1a，（i-n1，i+n2）zt，n=zt，（i-n1，i+n2）zlo，n=zlo，（i-n1，i+n2）扇墒设设设设设设缮设设设设设设（12）通过改进后的 KNN，检测出缺失点和其他异常点的总集合O1，O2，On，然后通过程序遍历总集合，将异常点清洗为空白轨迹点，便于后续插补。2基于多重插补的 ADS-B 轨迹补全多重回归插补理论定义16为：设 茁 是从所要修复样本数据 X 中估计得到的参数向量，茁 的先验分布为f（茁），使用贝叶斯准则可得到其后验概率分布为f（茁|X）邑Q（X|茁）f（茁）（13）式中 Q 是分析样本数据 X 得到似然函数。设 Xp为 X 中缺失信息的变量。pobs和 pmis分别为 p中的已存元素和缺失元素，并假设 p 中的信息缺失机制为随机丢失（MAR，missing at random）。缺失数据样本的参数向量期望值 茁|Xp，pobs为E（茁|Xp，pobs）=pmis乙E（茁|Xp，pobs，pmis）dP（pmis|Xp，pobs）=pmis乙E（茁|Xd）dP（pmis|Xp，pobs）（14）ADS-B 轨迹的多重插补方法等同于通过蒙特卡靳慧斌，冯朝辉，张召悦,等：基于改进 KNN 检测的 ADS-B 轨迹插补研究25-2023 年 8 月中 国 民 航 大 学 学 报罗积分法计算 E（茁|Xp，pobs）。每随机抽取一次缺失值，pgmis（坌g=1，2，.，G）从其分布 pgmis|Xp，pobs中产生一个向量 pg以pobspgmisT中的输入值。Xg以Xp，pg表示具有输入值的完整数据集。Bernstein-von Mises 定理17将贝叶斯估计与经典的最大似然估计联系起来，总结这种关系如下E（茁|Xg）p寅茁赞g（15）（茁-E（茁|Xg）g寅Normal（0，-Hg-1N）（16）式中：茁赞g是具有输入值的第 g 个完整数据集的最大似然估计；Hg是相关联的费舍尔信息矩阵；N 是轨迹样本数目大小。将在每次推算中计算得到的最大似然估计代入式（15）中，并转化为式（16），由此可得E（茁|Xp，pobs）抑茁=1GGg=1移茁赞g（17）式中茁 为茁赞g的平均值。利用条件方差公式，计算随机参数向量 茁|Xp，pobs方差如下Var（茁|Xp，pobs）=E（Var（茁|Xg）+Var（茁赞g）（18）根据式（18），在二维多重插补随机修复下（G=1，2），轨迹样本数据的 茁|Xp，pobs计算公式如下Var（茁|Xp，pobs）抑1GGg=1移Var（茁|Xg）+1G-1Gg=1移（茁赞g-茁）（茁赞g-茁）T（19）式（19）可以用于计算轨迹样本的最终方差。通过从轨迹数据集中抽取多个随机样本，可以推断出总体参数的均值和方差的无条件估计18。均值和方差-协方差矩阵计算如下兹軈=Ex（茁|X）（20）专赞=Varx（E（茁|X）+Ex（Var（E（茁|X）（21）式中：Ex（茁|X）为样本数据 X 中参数向量的期望；Varx（E（茁|X）为 Ex（茁|X）期望的方差。轨迹缺失数据多重插补的步骤如下：步骤 1从要修复的轨迹数据随机抽取部分样本数据 xg，g=1，2，G；步骤 2引入误差项 e 的方差，计算为滓赞2g=移eTexg步骤 3从 N（0，滓赞2g）计算得到随机向量着赞g；步骤 4从具备均值兹軈和方差-协方差矩阵专赞且符合多元正态分布的轨迹样本中计算得到兹赞g；步骤 5计算插入缺失值向量 pgmis，如下ZZ=Gg=1移（茁赞g-茁軍）（茁赞g-茁軍）T着赞g随机误差项pgmis=Z兹赞g+着赞gg=1，2，G扇墒设设设设设设设设缮设设设设设设设设（22）式中 Z为数据插补间可变性系数，在多次插补步骤结束时，可得到 G 个完整的数据集。因航空器性能的影响，轨迹数据插补不应以整个轨迹数据方差为系数，为得到更加准确的插补值，采用轨迹分隔和轨迹排序判断对多重插补方法进行改进，其过程如图 3 所示。因软件界面只显示航空器 ADS-B 二维坐标轨迹，所以只需要插补经度和纬度即可。3算例分析3.1数据处理本文选择天气状况良好的时间段，排除天气干扰因素，筛选出某机场终端区 4 天内轨迹坐标数据，如图 4 所示。通过匹配机型，从中随机挑选 2 个不同机型的航班，其轨迹信息如表 1 所示。样本轨迹数据通过米勒投影转换和缺失点筛选图 3轨迹修复流程图Fig.3Flow chart of trajectory repair子轨迹插补修复以不相邻空白点为特征点对轨迹进行分割提取经纬度数据输入异常点轨迹数据集子轨迹空白点 3pg，minobs pgmis pg，maxobs输出插补值否否是是26-第 41 卷 第 4 期图 5算法结果对比Fig.5Comparison of algorithm results1.000.980.960.940.920.900.9660.9330.9530.9300.9740.9510.9430.9810.9620.9830.9510.967欧氏度量曼哈顿度量标准化欧式度量（a）二维数据（b）三维数据图 4终端区不同维度轨迹显示Fig.4Trajectory display of different dimensions in terminal area表 1样本轨迹信息Tab.1Information of sample trajectory航班号 机型起始位置点/（毅）终止位置点/（毅）轨迹点总数780310 A321 E116.758 8，N39.014 13 E117.402 1，N39.077 02704781181 B739 E116.858 5，N38.937 87 E117.399 8，N39.080 88660航班号缺失点数目占比/%最小间隔/s平均间隔/s780310537.00.4073.059781181456.40.4062.401表 2缺失点信息Tab.2Information of missing points得到结果如表 2 所示。每条全类型报文平均解码时间 Ta=1.12 伊 10-4s，每条全类型报文发射间隔 To=0.25 s，数据经过处理传入终端界面时间 Tu 0.1，从表 2 可以得到最小间隔Tmin Ta+To+Tu均成立，说明缺失点判别模型算法正确。通过距离矩阵对缺失点赋值，便于改进后的 KNN算法检测分类。3.2算法对比在轨迹点分类的具体应用中，选择适用的距离计算方法是决定 KNN 分类结果质量的关键性因素。常见的空间距离算法有欧氏度量和曼哈顿度量，如 1.2节所述。本节将欧氏度量、曼哈顿度量和标准化欧氏度量进行实验对比，通过对比 3 种距离算法检测ADS-B 异常轨迹点的结果，验证 3 种不同算法检测ADS-B 异常点的准确率。所有样本训练集与测试集比例设置为 6 颐 4，每个样本设置不同 k 值（样本包括780310、781181），交叉验证采用不同距离算法的KNN检测ADS-B 异常轨迹点，并计算每种算法检测出的异常点数量与每个样本异常点总数之间占比，其结果如图 5 所示。从图 5 可得标准化欧氏度量的精确度均高于欧氏度量和曼哈顿度量。通过对比不同 k 值各算法准确率，当 k=2 时为最佳选择。通过同一终端区的不同航班轨迹、不同 k 值和不同距离算法的交叉验证，验证了改进后的 KNN 算法检测 ADS-B 异常点是完全可行的。3.3轨迹修复在进行插补修复前，需要将检测出的异常值的经度和纬度归 0，并赋值序号标签，方便插补算法识别需插补位置，以免忽略。完成赋值后，筛选出样本中纬度和经度将其作为修补序列 CX，Y输入多重插补算法中，其对比结果如图 6 所示。由图 6 可知，ADS-B 轨迹包含线性轨迹与非线性轨迹，经过缺失点筛选和改进后的 KNN 检测出异常点，将其清洗后通过多重插补能够有效完善轨迹图像，进一步提高数据的完整性。4结语航空器 ADS-B 技术的普及应用，使得 ADS-B数据成为当前航空器交通流分析、轨迹监管和空域运行靳慧斌，冯朝辉，张召悦,等：基于改进 KNN 检测的 ADS-B 轨迹插补研究图 6实验结果对比Fig.6Comparison of experimental results117.45117.40117.35117.30117.25（a）780310 插补前纬度/（毅）38.8539.0539.00117.45117.40117.35117.30117.25（b）780310 插补后纬度/（毅）38.8539.0539.00117.45117.40117.35117.30117.25（c）781181 插补前纬度/（毅）38.9439.02117.45117.40117.35117.30117.25（d）781181 插补后纬度/（毅）38.9838.9439.0238.9827-2023 年 8 月中 国 民 航 大 学 学 报态势分析的重要数据来源之一。然而 ADS-B 数据的缺失及异常点的出现给相关系统及用户的应用分析造成一定困难。本文方法不仅可以有效判别是否存在缺失点，而且可以有效检测出除缺失点外的其他异常点，能够有效插补航空器 ADS-B 轨迹，提升轨迹图像连续性，降低因异常点或图形问题而造成的不安全事件发生概率。针对 ADS-B 数据缺失水平问题，应进一步研究 ADS-B 数据质量分析指标体系，并根据用户需求和实际应用场景不断改善和发展 ADS-B 技术，提高数据链传输质量和地面站接收能力，从根本上解决数据缺失和异常点出现问题。参考文献:1张召悦援 空管监视技术M援 北京:国防工业出版社,2017援2曹娜援 基于海量实测的 ADS-B 数据质量分析D援 天津:中国民航大学,2017援3鞠时光,郭伟刚援 实用三次样条插值函数J援 小型微型计算机系统,1992,13(9):53-60援4ZAMANI M援 A simple 2D interpolation model for analysis of nonlineardataJ援 Natural Science,2010,2(6):641-645援5FARHAN J,FWA T F援 Airport pavement missing data management andimputation with stochastic multiple imputation modelJ援 TransportationResearch Record:Journal of the Transportation Research Board,2013,2336(1):43-54.6PEETERS M,ZONDERVAN-ZWIJNENBURG M,VINK G,et al援 Howto handle missing data:a comparison of different approachesJ援 Euro原pean Journal of Developmental Psychology,2015,12(4):377-394援7邹文华,左谛,金开研援 一种结合多重插补的 ADS-B 数据筛选算法J援 电讯技术,2021,61(4):461-467援8SUN J Y,WANG Y Q,SHEN Y Y,et al援 High-precision trajectory datareconstruction for TT&C systems using LS B-spline approximationJ援IEEE Signal Processing Letters,2020,27:895-899援9张瑜援 某型地波雷达终端显示软件的设计与实现D援 西安:西安电子科技大学,2017援10 RTCA.Minimum operational performance standards(MOPS)for 1090MHz extended squitter automatic dependent surveillance-broadcast(ADS-B)and traffic information services-broadcast(TIS-B):DO-260BS.Washington DC:RTCA,2009.11 中国民用航空局.一般运行和飞行规则:CCAR-91-R1S.北京:中国民用航空总局飞行标准司,2007.12 黄建理,杜金燃,谢家全,等援 一种基于改进 KNN 的大数据离群点检测算法J援 计算机与现代化,2017(5):67-70,75援13 邓楠洁援 基于数据依赖的高维大数据相似性度量方法研究D援 北京:北京邮电大学,2018援14 丁义,杨建援 欧式距离与标准化欧式距离在 k 近邻算法中的比较J援 软件,2020,41(10):135-136援15 CHIKODILI N B,ABDULMALIK M D,ABISOYE O A,et al援 Outlier de原tection in multivariate time series data using a fusion of K-medoid,standardized euclidean distance and Z-scoreC/International Confer原ence on Information and Communication Technology and Applications援Cham:Springer,2021:259-271援16 LITTLE R J A,RUBIN D B援 Statistical analysis with missing dataM援2nd ed援 Hoboken,NJ:John Wiley&Sons,Inc.,2002.17 TRAIN K援 Discrete choice methods with simulationM援 2nd ed援 Cam原bridge:Cambridge University Press,2009.18 BODNER T E援 What improves with increased missing data imputations?J援 Structural Equation Modeling:A Multidisciplinary Journal,2008,15(4):651-675援（责任编辑：明月）中国民航大学学报 投稿须知本刊投稿采用网上投稿，不接受电子邮件等其他方式投稿，投稿网址 https:/ 中国民航大学学报 并惠赐佳作。28-