-
卫星钟差数据的质量直接影响星载原子钟性能分析结果的可靠性和钟差预报的精度[1]。星载原子钟钟差数据的预处理是原子钟稳定性分析、钟差预报的前提和基础,因此欲获得有意义的分析结果,在进行原子钟性能评估之前,需进行数据预处理[2]。在实际应用中,原始钟差数据常会出现数据间断、跳变和粗差等异常情况[3],处理不同类型的异常值需要针对性的方法。已有不少学者对于异常值处理开展了研究[1-12]。对于数据跳变和数据间断问题,郭海荣[2]提出一种利用移动窗口寻找钟差数据的跳变点的算法;郭吉省[4]提出一种基于希尔伯特-黄变换的钟差跳变检测算法;Hackman等[5]给出了非等间隔数据的处理方法。对于粗差问题,冯遂亮[3]利用小波分析多分辨率的特点对该方法在原子钟数据粗差处理方面的应用进行了初步的探讨;方书山[6]研究了多元统计方法中的聚类分析在卫星钟差粗差探测中的应用;魏道坤[7]针对一次差分数据呈现线性趋势的情况提出了一种中值线粗差探测法;张倩倩等[8]提出一种基于Bayesian原理的星地时间同步钟差粗差探测方法。而最常用的粗差探测传统方法是中位数绝对偏差(median absolute deviation,MAD)[2-3,9-11],该方法原理简单,计算复杂度低,实际应用效果好。但传统的MAD方法在粗差探测存在不足:模型表达式存在歧义;对趋势特征明显的数据探测效果有限。
本文通过对传统MAD方法的深入分析,提出一种改进的MAD钟差粗差探测方法,相较于传统MAD方法,其粗差探测的准确率和查全率有了显著的提升。
-
已知钟差数据
,对相邻历元间的钟差数据进行一次差分,得到钟差数据对应的一次差分数据 ,其中 , 。传统MAD通常应用于钟差一次差分数据,其基本原理为:将钟差一次差分数据 与一次差分数据的中位数 及MAD数倍之和进行比较: 式中,
,其中 表示中位数; 。若 满足式(1)时(常数n取值根据需要确定),就认为 为粗差点。 在粗差探测出来之后,一般通过将粗差值设为0或其他的值,完成钟差异常值的预处理,本文主要探讨钟差的粗差探测方法。
-
本文使用的是国际全球卫星导航系统服务组织提供的全球定位系统(global positioning system,GPS)精密钟差数据,图 1为PRN-01卫星(2017-01-01—2017-01-31,5 min采样间隔)的钟差数据和一次差分数据。钟差数据的有效位较多,对异常值不敏感(图 1(a)),而一次差分数据异常值特征明显(图 1(b)),故一般采用钟差一次差分数据进行异常值处理。传统MAD粗差探测方法原理简单,探测效率高,但仍存在需改进的方面,本文对此进行了分析,提出了改进方法。
式(1)中,
是钟差一次差分数据,左边项 表示钟差一次差分数据的绝对值,右边项中MAD是正值,系数 一般取3~5[10], 可视为恒为正值的变量,当右边项 时,则 均判断为粗差,由此可知,传统MAD模型数学表达存在歧义。MAD方法本质上是判断每个 偏离中数的距离是否超限,当左项为 时,其数学含义为 偏离中数 的距离;相应的右项为 时,其数学含义为判断 是否为粗差的门限距离值,因此传统的MAD方法可修改为: 式中,
值为常数,一般情况下可以准确表达一次差分数据的分布特征,但对于趋势特征变化明显的数据,理想状态下, 值应该是随趋势变化而不断变化的,采用固定 值,会使得偏离 值的 难以通过距离约束条件探测粗差。 图 2为在同等参数条件下(
),分别采用式(1)和式(2)对PRN-01卫星一次差分数据进行粗差探测的结果(粗差均以中位数替换,下同)。可见,采用式(2)的探测结果(见图 2(b))优于式(1)的探测结果(见图 2(a)),式(1)中多个粗差值未被探测出来。对于大部分一次差分数据,其趋势变化特征不明显。本文提出基于岭回归[13-16]的基本原理,生成由 值组成的一次差分趋势线,得到动态 ,可有效克服部分具有显著趋势变化特征的钟差一次差分数据对粗差探测的干扰。 岭回归是一种有偏估计回归方法,实质上是一种改良的最小二乘估计法。包含粗差的钟差一次差分数据是一种“病态”数据,采用经典的最小二乘估计法会产生过拟合问题,而岭回归通过损失部分精度,可以得到更可靠的拟合结果。已知岭回归的数学模型为:
式中,
为钟差一次差分数据的时间序列; 为钟差一次差分数据序列; 是单位矩阵; 为正则化系数(通常经多次实验确定); 为回归值。式(3)联合式(2)可得动态MAD方法,数学模型为: 在同等参数条件下(
),分别采用式(4)和式(2)对GPS系统PRN-32卫星的一次差分数据进行粗差探测的结果如图 3所示。图 3(a)为式(4)探测的结果,图 3(b)为式(2)探测的结果,图 3中红色实线为 值,显而易见,图 3(a)的红色实线相较于图 3(b)更能客观地反映数据的趋势变化特征。在同等参数条件下,采用式(4)的探测结果优于式(2),式(2)中多个粗差值未被探测出来,需进一步减小 的取值,但过小的取值意味着判断粗差的距离约束进一步减小,容易导致部分正常值被判断为粗差。 如图 4所示,本文方法的改进主要体现在两方面:(1)针对传统MAD方法在数学表达方面存在的歧义,调整了MAD方法的结构,消除了模型表达上的歧义;(2)传统MAD方法在面对趋势性特征明显的数据时,粗差探测的效果不好,本文采用有偏估计的方法(岭回归),重新定义了
的取值,消除了趋势性的影响。 -
本文采用国际全球卫星导航系统服务组织的GPS系统(2017-01-01—2017-01-31)的5 min采样间隔的精密钟差数据进行实验,验证本文方法的有效性。GPS系统在轨运行32颗卫星,由于PRN-03、PRN-06、PRN-08、PRN-21和PRN-23星数据存在间断,PRN-04星无数据,因此上述卫星不参与实验。另外,本文方法主要与传统MAD方法进行对比,利用准确率(P)和查全率(R)两个指标对不同方法的粗差探测能力进行量化评估。
假设钟差数据中存在粗差的个体总数为C,探测出粗差的个体总数为B,正确探测的粗差个体总数为A,则准确率和召回率计算公式为:
式中,P为准确率;R为召回率。
-
分别采用改进MAD方法(式(4)和式(2))和传统MAD方法(式1)进行了对比实验,本文参数均通过多次实验获得,为保证实验条件的一致性,对于式(1)和式(2)一般取同一个
值,对于式(2)和式(4)一般通过优化 的取值,使得约束条件保持同一性。本文实验结果如图 5和图 6所示,图 5和图 6中灰色实线为原始数据,黑色实线为处理后数据,红色实心圆为粗差值。图 5(c)在 , 参数条件下,准确探测出所有的粗差值,准确率和召回率均为100%(见表 1)。同样,图 5(b)在 参数条件下,准确探测出所有的粗差值,准确率和召回率均为100%,这是由于PRN-02星的钟差一次差分数据的趋势变化特征不明显,没影响粗差探测,但对于图 5(a),在 参数条件下,式(1)右边项 为负值,导致所有的钟差一次差分数据均被探测为粗差,表 1中准确率和召回率的结果没有实际意义,若继续增加 的取值,使得 为正,即增加了参数寻优的计算成本,也不适合实际应用。 表 1 实验结果统计
Table 1. Statistics of Experiment Results
卫星 式(1) 式(2) 式(4) n A B C P/% R/% n A B C P/% R/% n A B C P/% R/% PRN-01 4 67 67 77 100 87.01 4 77 77 77 100 100 8 77 77 77 100 100 PRN-02 5 68 8 927 68 0.76 100 5 68 68 68 100 100 5 68 68 68 100 100 PRN-05 8 55 55 59 100 93.22 8 58 58 59 100 98.31 8 59 59 59 100 100 PRN-07 8 59 8 927 59 0.66 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-09 10 52 55 75 94.55 69.33 10 75 77 75 97.40 100 10 75 75 75 100 100 PRN-10 2 63 69 76 91.30 82.89 2 73 82 76 89.02 96.05 10 76 76 76 100 100 PRN-11 10 54 54 58 100 93.10 10 54 54 58 100 93.10 10 58 58 58 100 100 PRN-12 8 59 59 59 100 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-13 8 42 4 080 63 1.03 66.67 8 63 63 63 100 100 8 63 63 63 100 100 PRN-14 8 26 900 59 2.89 44.67 8 59 59 59 100 100 8 59 59 59 100 100 PRN-15 8 59 74 59 79.73 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-16 8 58 58 59 100 98.31 8 59 59 59 100 100 8 59 59 59 100 100 PRN-17 8 57 57 60 100 95.00 8 60 60 60 100 100 8 60 60 60 100 100 PRN-18 8 55 55 59 100 93.22 8 59 59 59 100 100 8 59 59 59 100 100 PRN-19 8 58 58 60 100 96.67 8 59 59 60 100 98.33 8 60 60 60 100 100 PRN-20 8 50 50 57 100 87.72 8 51 51 57 100 89.47 8 57 57 57 100 100 PRN-22 8 59 8 927 59 0.66 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-24 8 54 56 54 96.43 100 8 54 54 54 100 100 8 54 54 54 100 100 PRN-25 10 86 8 927 86 0.96 100 10 86 88 86 97.72 100 10 86 86 86 100 100 PRN-26 3 72 8 927 72 0.81 100 3 69 69 72 100 95.83 15 72 72 72 100 100 PRN-27 8 56 56 75 100 74.67 8 74 75 75 98.67 100 8 75 75 75 100 100 PRN-28 8 55 59 59 93.22 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-29 8 31 5 187 59 5.98 52.54 8 59 59 59 100 100 8 59 59 59 100 100 PRN-30 8 72 72 75 100 96.00 8 75 75 75 100 100 8 75 75 75 100 100 PRN-31 8 59 7 853 59 0.75 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-32 3 75 8 927 75 0.84 100 3 72 72 75 100 96.00 8 75 75 75 100 100 图 6(c)在
, 参数条件下,准确探测出所有的粗差值,准确率和召回率均为100%,图 6(b)在 的参数条件下,准确率为89.02%,召回率为96.05%,其中,探测出粗差的个体总数B=82,实际的真实值为76,表明其中至少有6个错误探测值,准确探测的个数为73,其中有3个粗差未被探测出来(图 6(b)绿色实线区域),这进一步说明错误探测的个数为9(与图6(b)蓝色实线区域红色实心点个数一致)。分析错误粗差探测点的分布特征可知,由于PRN-10星钟差一次差分数据趋势性特征明显,错误的探测点集中分布在数据首尾部分,由此表明,对于趋势特征明显的钟差一次差分数据,采用动态MAD的方法能取得更好的探测效果。图 6(a)在 的参数条件下,准确率为91.30%,召回率为82.89%,其中,探测出粗差的个体总数B=69,实际的真实值为76,正确探测出的个数A=63,其中有6个错误探测值(图 6(a)蓝色实线区域),这说明有13个粗差未被探测出来(图 6(a)绿色实线区域),表明传统MAD方法由于模型结构的不合理,导致探测能力不及本文方法。 特别指出,在图 6的实验中,公共参数
设置存在差异(图 6(c)中, ;图6 (a)、6(b)中, ),这是由于图 6(c)采用的是式(4)模型,动态MAD方法准确地表达了趋势性的变化特征,使得距离约束条件 的取值不会偏大,而采用固定 值会使 的取值偏大,相应的 的取值偏小,通过对参数的不断优化,动态MAD方法可以最大限度地准确识别出所有的粗差。 实验结果表明,本文方法通过优化MAD模型结构和改进动态MAD,使粗差探测的准确率和召回率大幅提高。对GPS系统26颗可用卫星的钟差数据进行了实验,说明本文方法的普适性。式(1)、式(2)和式(4)的公共参数
的取值存在差异,为保证实验环境的同一性,本文利用式(4)模型,首先,通过多次重复实验取得最优参数(是指在该参数条件下,粗差识别的准确度和召回率达到最优);然后,以最优参数为参考,对式(1)、式(2)模型采用同一参数,当一次差分数据趋势变化特征明显时,式(2)采用固定 值,导致距离约束条件偏大,故参数设置小于式(4)的参数,与式(1)保持一致,但距离约束条件 相等,可以认为实验环境具备同一性。表 1列出了GPS系统26颗数据可用卫星的钟差一次差分数据粗差探测结果。 分析表 1可得,实验结果分为3类:
1)表 1中加粗字体为趋势特征明显的一次差分数据的实验结果,包括PRN-01、PRN-10、PRN-26和PRN-32星。其中,PRN-01(图 1)存在一个明显的频跳,但幅度有限,式(4)中的
变化范围不大,与式(2)的固定值差异较小,所以式(2)和式(4)的探测结果保持一致,均优于式(1)的结果;PRN-10(图 6)趋势特征显著,式(4)的结果优于式(2)和式(1),同理PRN-26和PRN-32星的实验结果与PRN-10星一致,但式(1)的结果由于 值为负导致所有数据被判定为粗差,实验结果无意义。 2)表 1中PRN-13、PRN-14、PRN-29和PRN-31星采用式(1)的探测结果是部分正常值被探测为粗差,这是由于
为负值,但 值为正,远小于式(2)中的右边项 ,在同一参数条件下,式(1)的距离约束过小,部分正常值被探测为粗差,对于PRN-02、PRN-07、PRN-22、PRN-25、PRN-26和PRN-32星,被探测出粗差的个体总数B为8 927(样本总数),所有的数据均被探测为粗差,这是由于式(1)中的 为负,且 也为负,导致结果无意义。 3)表 1中剩余卫星钟差一次差分趋势稳定,
值为正,所以在同一参数条件下,3种模型的实验结果较好,但式(4)的实验结果整体上仍优于其他两种模型。 由上述分析可知,本文模型结构合理,不仅可有效避免
值为负导致大量的错误探测结果,而且显著的趋势变化特征对粗差探测结果影响不大,具有较强的普适性。 -
本文通过分析传统MAD粗差探测方法的工作原理,针对其在模型结构方面存在的不足,提出了一种改进的MAD钟差粗差探测方法,具有以下特点:(1)本文方法优化了传统MAD的模型结构,修正了数学模型表达中的歧义,避免了大量错误探测结果的出现;(2)本文方法基于岭回归的基本原理,优化
值的选取,形成动态MAD,有效克服了趋势变化特征显著的钟差一次差分数据对粗差探测的干扰,提高了粗差探测的准确率和召回率,有效控制了算法的复杂度,具备实际应用价值。
-
摘要: 针对传统中位数绝对偏差(median absolute deviation,MAD)方法在探测钟差粗差方面的不足,提出一种改进的MAD钟差粗差探测方法。以一次差分数据为研究对象,通过深入分析传统MAD方法的工作原理,发现传统MAD方法的数学模型表达存在歧义。在优化探测模型结构的基础上,基于岭回归的基本原理,生成动态MAD,有效克服了部分具有显著趋势变化特征的钟差一次差分数据对粗差探测的干扰。使用国际全球卫星导航系统服务组织提供的GPS精密钟差数据进行了实验,并与传统MAD方法进行比较,结果表明,所提方法相较于传统的MAD方法在粗差探测准确率和查全率方面具有明显优势,实际应用价值较高。Abstract:
Objectives Satellite clock difference (SCB) is affected by many factors, such as data discontinuity, jump and gross error, among which the gross error is the most representative. Traditional SCB gross error detection method is the median absolute deviation (MAD) method Methods With the single difference of SCB as the research object, we analyzed the working principle of traditional MAD method. In view of the deficiency of traditional MAD in SCB gross error detection, an improved gross error detection method was proposed. Firstly, the ambiguity in mathematical expression was eliminated by optimizing the model structure, and then dynamic MAD was generated based on the basic principle of ridge regression, which effectively overcame the interference of some single difference data with significant trend changes on gross error detection. Results We used International Global Navigation Satellite System Service precision SCB as reference for experiments. Compared with traditional MAD, the detection accuracy and recall rate of the proposed method were better than traditional MAD, especially for the data with obvious change of trend characteristics. Conclusions The proposed method has obvious advantages in gross error detection accuracy and recall rate and has high application value. -
表 1 实验结果统计
Table 1. Statistics of Experiment Results
卫星 式(1) 式(2) 式(4) n A B C P/% R/% n A B C P/% R/% n A B C P/% R/% PRN-01 4 67 67 77 100 87.01 4 77 77 77 100 100 8 77 77 77 100 100 PRN-02 5 68 8 927 68 0.76 100 5 68 68 68 100 100 5 68 68 68 100 100 PRN-05 8 55 55 59 100 93.22 8 58 58 59 100 98.31 8 59 59 59 100 100 PRN-07 8 59 8 927 59 0.66 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-09 10 52 55 75 94.55 69.33 10 75 77 75 97.40 100 10 75 75 75 100 100 PRN-10 2 63 69 76 91.30 82.89 2 73 82 76 89.02 96.05 10 76 76 76 100 100 PRN-11 10 54 54 58 100 93.10 10 54 54 58 100 93.10 10 58 58 58 100 100 PRN-12 8 59 59 59 100 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-13 8 42 4 080 63 1.03 66.67 8 63 63 63 100 100 8 63 63 63 100 100 PRN-14 8 26 900 59 2.89 44.67 8 59 59 59 100 100 8 59 59 59 100 100 PRN-15 8 59 74 59 79.73 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-16 8 58 58 59 100 98.31 8 59 59 59 100 100 8 59 59 59 100 100 PRN-17 8 57 57 60 100 95.00 8 60 60 60 100 100 8 60 60 60 100 100 PRN-18 8 55 55 59 100 93.22 8 59 59 59 100 100 8 59 59 59 100 100 PRN-19 8 58 58 60 100 96.67 8 59 59 60 100 98.33 8 60 60 60 100 100 PRN-20 8 50 50 57 100 87.72 8 51 51 57 100 89.47 8 57 57 57 100 100 PRN-22 8 59 8 927 59 0.66 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-24 8 54 56 54 96.43 100 8 54 54 54 100 100 8 54 54 54 100 100 PRN-25 10 86 8 927 86 0.96 100 10 86 88 86 97.72 100 10 86 86 86 100 100 PRN-26 3 72 8 927 72 0.81 100 3 69 69 72 100 95.83 15 72 72 72 100 100 PRN-27 8 56 56 75 100 74.67 8 74 75 75 98.67 100 8 75 75 75 100 100 PRN-28 8 55 59 59 93.22 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-29 8 31 5 187 59 5.98 52.54 8 59 59 59 100 100 8 59 59 59 100 100 PRN-30 8 72 72 75 100 96.00 8 75 75 75 100 100 8 75 75 75 100 100 PRN-31 8 59 7 853 59 0.75 100 8 59 59 59 100 100 8 59 59 59 100 100 PRN-32 3 75 8 927 75 0.84 100 3 72 72 75 100 96.00 8 75 75 75 100 100 -
[1] 王宁. BDS星载原子钟钟差数据预处理与钟性能分析研究[D]. 郑州: 信息工程大学, 2017 Wang Ning. Study on the Methods of Clock Bias Data Preprocessing and Performance Analysis for On-Board Atomic Clocks of BDS[D]. Zhengzhou: Information Engineering University, 2017 [2] 郭海荣. 导航卫星原子钟时频特性分析理论与方法研究[D]. 郑州: 信息工程大学, 2006 Guo Hairong. Study on the Analysis Theories and Algorithms of the Time and Frequency Characteriztion for Atomic Clocks of Navigation Satellite[D]. Zhengzhou: Information Engineering University, 2006 [3] 冯遂亮. 原子钟数据预处理与钟性能分析方法研究[D]. 郑州: 信息工程大学, 2009 Feng Suiliang. Study on the Methods of Data Preprcessing and Performance Analysis for Atomic Clocks [D]. Zhengzhou: Information Engineering Universty, 2009 [4] 郭吉省. UTC(NIM)原子时标驾驭研究[D]. 北京: 北京工业大学, 2013 Guo Jixing. Time Scale Steering in UTC(NIM) [D]. Beijing: Beijing University of Technology, 2013 [5] Hackman C, Parker T E. Noise Analysis of Unevely Spaced Time Series Data[J]. Metrologia, 1996, 33(5): 457-466 doi: 10.1088/0026-1394/33/5/4 [6] 方书山. GPS卫星钟差完备性监测研究与算法实现[D]. 阜新: 辽宁工程技术大学, 2011 Fang Shushan. Research of GPS Satellite Clock Eror Integrity Monitoring and Algorithm Implementtion[D]. Fuxin: Liaoning Technical University, 2011 [7] 魏道坤. 卫星钟差预报模型研究[D]. 西安: 长安大学, 2014 Wei Daokun. Study on the Satellite Clock Bias Forcast Model[D]. Xi'an: Chang'an University, 2014 [8] 张倩倩, 韩松辉, 杜兰, 等. 星地时间同步钟差异常处理的Bayesian方法[J]. 武汉大学学报·信息科学版, 2016, 41(6): 772-777 doi: 10.13203/j.whugis20140666 Zhang Qianqian, Han Songhui, Du Lan, et al. Bayesian Methods for Outliers Detection and Estimtion in Clock Offset Measurements of Satellite - Ground Time Transfer[J]. Geomatics and Informtion Science of Wuhan University, 2016, 41(6): 772-777 doi: 10.13203/j.whugis20140666 [9] 黄观文. GNSS星载原子钟质量评价及精密钟差算法研究[D]. 西安: 长安大学, 2012 Huang Guanwen. Research on Algorithms of Prcise Clock Offset and Quality Evaluation of GNSS Satellite Clock[D]. Xi'an: Chang'an University, 2012 [10] 王宇谱, 吕志平, 陈正生, 等. 一种新的钟差预处理方法及在WNN钟差中长期预报中的应用[J]. 武汉大学学报·信息科学版, 2016, 41(3): 373-379 doi: 10.13203/j.whugis20140216 Wang Yupu, Lü Zhiping, Chen Zhengsheng, et al. A New Data Preprocessing Method for Satellite Clock Bias and Its Application in WNN to Predict Medium-Term and Long-Term Clock Bias[J]. Gematics and Information Science of Wuhan Universty, 2016, 41(3): 373-379 doi: 10.13203/j.whugis20140216 [11] Wang Y P, Lu Z P, Qu Y Y, et al. Improving Prdiction Performance of GPS Satellite Clock Bias Based on Wavelet Neural Network[J]. GPS Soltions, 2017, 21(2): 523-534 doi: 10.1007/s10291-016-0543-z [12] 马朝忠, 朱建青, 韩松辉. 基于ARIMA模型的卫星钟差异常值探测的模型选择方法[J]. 武汉大学学报·信息科学版, 2020, 45(2): 167-172 doi: 10.13203/j.whugis20180230 Ma Chaozhong, Zhu Jianqing, Han Songhui. Model Selection Method Based on ARIMA Model in Ouliers Detection of Satellite Clock Offset[J]. Gematics and Information Science of Wuhan Universty, 2020, 45(2): 167-172 doi: 10.13203/j.whugis20180230 [13] Hoerl A E, Kennard R W. Ridge Regression: Biased Estimation for Nonorthogonal Problems[J]. Technmetrics, 1970, 12(1): 55-67 doi: 10.1080/00401706.1970.10488634 [14] Dwivedi T D, Srivastava V K, Hall R L. Finite Sample Properties of Ridge Estimators[J]. Technmetrics, 1980, 22(2): 205-212 doi: 10.1080/00401706.1980.10486136 [15] Golub G H, Heath M, Wahba G. Generalized CrosValidation as a Method for Choosing a Good Ridge Parameter[J]. Technometrics, 1979, 21(2): 215- 223 doi: 10.1080/00401706.1979.10489751 [16] Hastie T. Ridge Regularization: An Essential Cocept in Data Science[J]. Technometrics, 2020, 62 (4): 426-433 doi: 10.1080/00401706.2020.1791959 -