-
测量数据处理中,解决粗差问题是不可避免的环节。处理粗差目前有粗差探测与抗差估计两种方法[1-12]。粗差探测把粗差归结为函数模型的不完善,而抗差估计则通过方差膨胀的方式削弱甚至消除粗差的不良影响。选择合适的统计检验量以准确地探测和识别异常观测是使用这两种方法的前提条件。
由于尺度因子通常是未知的,可基于最小二乘平差残差的加权平方和进行估计,然后构造出服从
分布的统计量或者服从学生氏分布的 统计量[8-10]。但该方法几乎不具有抵御粗差影响的能力,因此需要做进一步的工作。尺度因子的绝对中位差(median absolute deviation,MAD)估计由于其简单易操作、抗差性强等特点,是目前粗差探测与抗差估计领域应用最为广泛的一种估计尺度因子的方法[1,13-14]。 本文首先对样本尺度因子的MAD估计进行了讨论,详细探讨了Fisher一致性调节因子的确定。在此基础上,论述了平差模型中尺度因子的MAD估计及其在测量数据处理中的应用。最后结合全球导航卫星系统(global navigation satellite systems,GNSS)网平差的数值计算结果讨论了测量平差成果的精度评定问题。
-
设一元随机变量
服从期望为 、方差为 的分布,对其进行 次独立观测,得到概率分布相同的样本序列为 。如果总体的期望 已知,基于样本序列可对方差 做如下估计: σ ¯ 2 = 1 n ∑ i = 1 n ]]>Y i - μ 2 由于算术平均值(mean)运算不具有抗差性,上述方差估计方法对粗差非常敏感,中位数(median)运算的崩溃污染率却达到50%,因此以中位数代替式(1)中的平均值运算,得到:
σ ˜ = k ⋅ m e d i ]]>Y i - μ 式中,med表示中位数运算;参数
为保证 满足Fisher一致性的调节因子。 与 的概率分布有关,不同的概率分布对应的调节因子一般不同[1,13-14]。 记
,则根据中位数的定义可得: P Y - μ < m = P Y - μ ≥ m = 0.5 ]]>如果随机变量
的概率密度函数关于 对称,则有: P Y - μ ≥ m = P ]]>Y - μ ≤ - m 从而得到:
P Y - μ σ ≥ m σ = 0.25 ]]>这表明
是随机变量 的累积分布函数的上0.25分位点。进而可证明,当 服从正态分布时, 服从标准正态分布,其上0.25分位点为: Φ - 1 1 - 0.25 = 0.674 5 ]]>此时,
σ ˜ = m e d i Y i - μ 0.674 5 = 1.482 6 ⋅ m e d i Y i - μ ]]> 特别地,当
时,式(7)变为: σ ˜ = 1.482 6 ⋅ m e d i ]]>Y i 如果总体的期望
未知,可先求样本均值,再得到样本方差的估值如下: σ ¯ 2 = 1 n ∑ i = 1 n ]]>Y i - Y ¯ 2 式中,
。根据文献[4],式(9)得到的估计是有偏的,为保证其无偏,还需要乘以因子 。 由式(9)可以得到尺度因子
的估计如下: σ ̂ = 1 n ∑ i = 1 n ]]>Y i - Y ¯ 由于估计量
的崩溃点是0,同时其影响函数是无界的,因此不能错误地认为式(10)是尺度因子的抗差估计[14]。 为提高式(10)的抗差性,将其中的两次平均值运算都以中位数运算代替,得到:
M A D = m e d i ]]>Y i - m e d j { Y j } 式中,MAD表示尺度因子
的(关于中位数的)绝对偏差的中位数估计。 若以MAD作为尺度因子
的估计,同样必须使其满足Fisher一致性,可直接以式(2)中的调节因子来代替。当 服从正态分布时,得到尺度因子 的如下估计[1,13-14]: σ ˜ = 1.482 6 ⋅ M A D ]]>尺度因子的所有估计形式里,基于MAD的估计具有最高的崩溃污染率,然而MAD估计的抗差性是以牺牲统计量的渐进效率为代价的,其原因在于大多数观测量对MAD估计没有做出直接贡献。在统计学领域,文献[14]提出了几种替代MAD估计的方法,但MAD估计由于其简单易操作、抗差性强等优点,仍然是目前应用最为广泛的一种方法。
-
参考文献[4]中给出的线性平差模型:
E ( L ) = A X ]]>式中,
为服从正态分布的 观测向量; 为 列满秩的设计阵; 为 的未知参数向量; 为单位权方差因子;权阵 满足对称正定性。 基于极大似然估计原理,得到[4]:
]]>X ̂ = ( A T P A ) - 1 A T P L σ ¯ 2 = 1 n V T P V 其中,
V = L - A ]]>X ̂ 称为残差向量。
由此易证明
是 的有偏估计,为使其满足无偏性质,需做如下修正[4]: σ ̂ 2 = 1 n - u V T P V ]]>然而,此时这两种估计都不能有效抵御粗差的影响。如果假定残差向量的期望为零,则对于残差序列
可采用式(8)对尺度因子进行估计,即绝对残差的中位数(median absolute residual,MAR)估计,具体计算如下: σ ˜ = 1.482 6 ⋅ m e d i ]]>v i 采用式(8)估计尺度因子需符合以下条件:(1)观测序列应该是独立的;(2)序列应服从相同的概率分布;(3)观测序列的期望值为零。当观测量受到粗差污染时,残差序列或将不满足以上3个条件中的任何一条,因此使用式(8)存在一定的风险。为此,本文首先对残差进行标准化,得到均服从正态分布
的序列 。这里符号 表示第i个标准化残差,计算如下: v ˜ i = ]]>v i V a r v i 进而可得到尺度因子的MAD估计如下:
σ ˜ v = 1.482 6 ⋅ m e d i ]]>σ v ˜ i - m e d j { σ v ˜ j } 以
记第 个元素为1其余均为0的 维单位向量,则 为第 个一致最大功效(uniformly most powerful,UMP)统计检验量[16]。UMP统计量在独立观测情况下退化为标准化残差,但在相关观测情况下,二者是完全不同的统计量[15-17]。构造如下具有相同概率分布 的序列 。类似于标准化残差的情形,可得到尺度因子的MAD估计如下: σ ˜ w = 1.482 6 ⋅ m e d i ]]>σ w i - m e d j { σ w j } 式(17)、式(19)和式(20)分别给出了平差模型中尺度因子3种不同形式的抗差估计。
-
粗差探测与抗差估计都必须面对统计量的选择问题。虽然可供选择的统计量多种多样,但由于尺度因子在许多情况下是未知的,因此不可避免地首先要对尺度因子进行估计。一种简单直接的方法是基于式(16)给出的方差的无偏估计构造出统计量,比如
分布和 分布[8-10]。但由于式(16)几乎不具有抗差性, 分布和 分布使用起来均潜在风险,尤其在只有一个模型冗余度的极端条件下,根本无法定义 统计量,而此时无论观测量受到粗差污染与否,所有 统计检验量的绝对值始终恒等于1,因此同样无法用于粗差探测与识别[17]。 如果以服从正态分布的标准化残差或UMP统计量为出发点,并将涉及的标准差
分别用其抗差MAD估计(式(19)或式(20))代替,可得到两类相关观测条件下近似服从正态分布的抗差统计检验量。由于UMP统计量在同类统计检验量中具有最大功效,因此使用抗差UMP统计量可以更为准确地进行粗差的探测与识别[15-17]。 平差的目的主要是对未知参数进行估计,以及对得到的未知参数估值进行精度评定。使用基于MAD的尺度因子可以使统计检验量无量纲化,进而可将离群的统计量对应的观测值标记为可疑观测值,对可疑观测值进行剔除或降权处理后,最终将得到可靠的测量成果。
§2给出的式(17)、式(19)和式(20)都是对样本的尺度因子进行估计。在对平差成果进行精度评定时,需要估计总体的尺度因子。为保证估值的渐进无偏性,综合式(9)和式(16),本文提出尺度因子的近似估计如下:
σ ̂ w = n n - 1 1.482 6 m e d i ]]>σ w i - m e d j { σ w j } 如果使用的是粗差剔除方法或者是包含拒绝域的抗差估计方法,式中的
还应扣除那些被剔除或者权降为零的观测量的数目。显然,当观测量个数很大时,式(21)与式(20)仅有细微差异。 -
本文实验中的GNSS网包含一个控制点和5个未知点,共观测了45条独立基线的数据[18]。借助于已有的设计阵和待定点坐标参考值,以及观测量的方差-协方差阵,可以模拟出服从正态分布的观测量[19]。在此基础上,随机指定1~6个观测量分别独立加入大小、符号不同的粗差,具体可参看文献[16],虽然基于UMP统计量和基于标准化残差的抗差估计方案都能够有效抵御粗差的不良影响,但基于UMP统计量的抗差估计在抗差性和计算效率两个方面都优于基于标准化残差的抗差估计。
为进一步研究MAD估计在精度评定中的应用,本文进行了如下模拟计算:(1)方案1:对没有植入粗差的数据采用最小二乘平差进行参数估计,最终依据式(16)对尺度因子进行估计;(2)方案2:同方案1,但采用式(21)估计尺度因子;(3)方案3:对受粗差污染后的数据采用基于UMP统计量的抗差估计进行参数解算,根据式(16)估计尺度因子并用于精度评定;(4)方案4:同方案3,但采用式(21)估计尺度因子。
为保证评价的客观性,上述模拟进行了10 000次,图 1统计了4种方案10 000次GNSS网平差模拟数据得到的尺度因子估值序列,统计结果见表 1。
图 1 4种方案10 000次GNSS网平差模拟数据得到的尺度因子估值序列
Figure 1. Estimated Scale Factors of Four Schemes over 10 000 Repetitions in GNSS Network Adjustment
表 1 4种方案得到的尺度因子估值序列的统计结果/mm
Table 1. Statistics of Estimated Scale Factor of Four Schemes/mm
方案 最大值 平均值 中位数 方案1 1.28 1.00 1.00 方案2 1.52 1.00 1.00 方案3 2.09 0.84 0.85 方案4 1.56 0.89 0.90 由图 1、表 1可知,没有粗差污染的情况下,使用残差加权平方和方法是最好的选择,此时MAD估值序列呈现出一定的波动性;对于受污染的数据,采用抗差估计尽管可以在一定程度上抵御了粗差的影响,但基于残差加权平方和方法得到的尺度因子明显偏小,而且序列中存在“突刺”现象;相对而言,基于MAD的尺度因子估值表现更为稳健。
-
粗差探测与抗差估计是目前处理粗差的两种常用方法。选择合适的统计量是探测和识别异常观测的关键。本文讨论了样本尺度因子的MAD估计,以及Fisher一致性调节因子的确定方法,并对平差模型中尺度因子的MAD估计及其在粗差处理中的应用进行了探讨。GNSS网平差的数值例子表明,由于具有良好的抗差性,尺度因子的MAD估计不但可用于构造粗差探测与识别的统计检验量,还可用于平差成果的精度评定。
-
摘要: 观测数据受到粗差污染时,平差结果往往失真,不可采用。选择合适的统计检验量是对测量数据进行粗差处理的关键一环,而构造统计检验量必须对尺度因子做出估计。首先对样本尺度因子的绝对中位差(median absolute deviation,MAD)估计进行了讨论,并详细探讨了其中涉及到的Fisher一致性调节因子的确定;然后在此基础上,分别提出了平差模型中基于标准化残差和一致最大功效(uniformly most powerful,UMP)统计检验量序列的尺度因子的两种抗差估计方法,尺度因子的两种抗差估计都可用于构造相应的统计检验量以识别和定位可疑观测量;最后对全球导航卫星系统(global navigation satellite systems,GNSS)网平差进行具体数值计算,结果表明,尺度因子的MAD估计具有良好的抗差性,不但可用于粗差处理,还可用于平差成果的精度评定。Abstract:
Objectives The least-squares method is very sensitive to outliers, and the adjustment outputs will usually be unacceptable when some of the observations are contaminated. Selection of appropriate statistical tests plays a pivotal role both in robust estimation and conventional outlier detection procedures. Methods The MAD (median absolute deviation) estimate of scale factor in the univariate case is discussed firstly. Determination of the Fisher-consistency factor is described for Gaussian normal distribution. Robust estimates of scale factor in linear adjustment model are addressed based on standardized least-squares residuals and the uniformly most powerful test statistics, respectively. Both of them can be used for constructing statistical tests, to identify the potential outlying observations, and therefore their deterioration effect will be mitigated. For illustrative purpose, Monte Carlo simulations in GPS network adjustment scenario are performed. Results Numerical results show that the MAD-based estimate of scale factor is robust and works well in accuracy assessment for adjustment outputs. Conclusions Explicit formula for estimating the scale factor, the MAD is a very robust scale estimator and has low computation complexity. It is therefore appropriate to use the MAD for adjustment computations and accuracy assessment when outliers are present. -
Key words:
- scale factor /
- outlier detection /
- robust estimation /
- MAD estimate
-
表 1 4种方案得到的尺度因子估值序列的统计结果/mm
Table 1. Statistics of Estimated Scale Factor of Four Schemes/mm
方案 最大值 平均值 中位数 方案1 1.28 1.00 1.00 方案2 1.52 1.00 1.00 方案3 2.09 0.84 0.85 方案4 1.56 0.89 0.90 -
[1] Huber P J. Robust Statistics[M]. New York: Wiley, 1981 [2] 周江文, 黄幼才, 杨元喜, 等. 抗差最小二乘法[M]. 武汉: 华中理工大学出版社, 1997 Zhou Jiangwen, Huang Youcai, Yang Yuanxi, et al. Robust Least Squares Method[M]. Wuhan: Huazhong University of Science and Technology Press, 1997 [3] Wolf P R, Ghilani C D. Adjustment Computations: Statistics and Least Squares in Surveying and GIS [M]. New York: Wiley, 1997 [4] Koch K R. Parameter Estimation and Hypothesis Testing in Linear Models[M]. Berlin: SpringerVerlag, 1999 [5] Leick A, Rapoport L, Tatarnikov D. GPS Satellite Surveying[M]. New York: Wiley, 2015 [6] Barnett V, Lewis T. Outliers in Statistical Data [M]. New York: Wiley, 1994 [7] Baarda W. A Testing Procedure for Use in Geodetic Networks[J]. Netherlands Geodetic Commission, Publications on Geodesy, 1968, 2(5): 1-97 [8] Pope A J. The Statistics of Residuals and the Detection of Outliers[R]. Rockville: NOAA Technical Report, 1976 [9] Chatterjee S, Hadi A S. Sensitivity Analysis in Linear Regression[M]. New York: Wiley, 1988 [10] Kargoll B. On the Theory and Application of Model Misspecification Tests in Geodesy[D]. Bonn: University of Bonn, 2007 [11] 鲁铁定, 杨元喜, 周世健. 均值漂移模式几种粗差探测法的MDB比较[J]. 武汉大学学报·信息科学版, 2019, 44(2): 185-192 doi: 10.13203/j.whugis20140330 Lu Tieding, Yang Yuanxi, Zhou Shijian. Comparative Analysis of MDB for Different Outliers Detection Methods[J]. Geomatics and Information Science of Wuhan University, 2019, 44(2): 185- 192 doi: 10.13203/j.whugis20140330 [12] 王海涛, 欧吉坤, 袁运斌, 等. 估计观测值粗差三种方法的等价性讨论[J]. 武汉大学学报·信息科学版, 2013, 38(2): 162-166 http://ch.whu.edu.cn/article/id/6096 Wang Haitao, Ou Jikun, Yuan Yunbin, etal. On Equivalence of Three Estimators for Outliers in Linear Model[J]. Geomatics and Information Science of Wuhan University, 2013, 38(2): 162- 166 http://ch.whu.edu.cn/article/id/6096 [13] Hampel F R. The Influence Curve and its Role in Robust Estimation[J]. Journal of the American Statistical Association, 1974, 69(346): 383-393 doi: 10.1080/01621459.1974.10482962 [14] Rousseeuw P J, Croux C. Alternatives to the Median Absolute Deviation[J]. Journal of the American Statistical Association, 1993, 88(424): 1 273-1 283 doi: 10.1080/01621459.1993.10476408 [15] 郭建锋. 模型误差理论若干问题研究及其在GPS数据处理中的应用[D]. 武汉: 中科院测量与地球物理研究所, 2007 Guo J. Theory of Model Errors and its Applications in GPS Data Processing[D]. Wuhan: Institute of Geodesy and Geophysics of Chinese Academy of Sciences, 2007 [16] Guo J, Ou J, Wang H. Robust Estimation for Correlated Observations: Two Local Sensitivity-based Downweighting Strategies[J]. Journal of Geodesy, 2010, 84(4): 243-250 doi: 10.1007/s00190-009-0361-y [17] Guo J. A Note on The Conventional Outlier Detection Test Procedures[J]. Boletim Ciencias Geodesicas, 2015, 21(2): 433-440 doi: 10.1590/S1982-21702015000200024 [18] Snow K B, Schaffrin B. Three- Dimensional Outlier Detection for GPS Networks and Their Densification via the BLIMPBE Approach[J]. GPS Solutions, 2003, 7(2): 130-139 doi: 10.1007/s10291-003-0058-2 [19] Yang Y, Song L, Xu T. Robust Estimator for Correlated Observations Based on Bifactor Equivalent Weights[J]. Journal of Geodesy, 2002, 76(6): 353-358 doi: 10.1007%2Fs00190-002-0256-7.pdf -