-
尽管测量手段越来越丰富,测量多余观测也空前丰富,但是异常误差的影响和探测依然是测量数据处理研究的热点之一。实践中,常将粗差如Baarda[1]纳入函数模型,提出的粗差数据探测法[1],该方法将粗差的观测数据看成是具有不同数学期望,但方差相同的母体;另一种常用方法是将粗差纳入随机模型,引进统计学中的抗差估计[2-3],该方法将含粗差的观测数据看成是期望不变,但方差异常大的另一个母体。基于统计学原理对粗差进行探测和定位,并将其剔除,其研究成果主要包括数据探测法[1]、多维粗差同时定位与定值法(simaltaneous locating and evaluating multiple gross errors, LEGE)[4]、拟准检定法(quasi-accurate detection of gross errors, QUAD)[5-6]、部分最小二乘法(partly least squares, PLS)或预测残差法[7-8]等。针对上述几种粗差探测和定位方法的异同,文献[9]对均值漂移模型粗差探测法和LEGE法的原理、方法和粗差探测过程进行了比较,证明两者在原理上基本等价,探测结果基本相近。文献[10]对独立等精度情形下部分最小二乘法和QUAD法进行了比较,得出其结果具有相同效果。文献[11]研究表明,独立等精度下,QUAD、LEGE法和部分最小二乘法在粗差估值方面具有等价性,而在相关观测情形下存在一定差异。
可靠性是研究平差系统发现粗差的能力和不可发现的粗差对平差结果的影响程度。Baarda最早提出测量系统的可靠性理论,并推导了有关计算公式[1]。Teunissen[12]、Förstner[13]、Schaffrin[14]、李德仁等[15]、王金岭等[16-17]、欧吉坤等[18]对可靠性理论进行了深入的研究,并做了重要发展。基于均值漂移模式的可靠性理论包括最小可探测偏差(minimal detectable bias, MDB)、内部可靠性和外部可靠性等3部分[18]。
虽然已有的文献对上述方法进行了比较,但没有对不同方法的可靠性进行比较。本文推导了几种方法的MDB的计算公式,通过算例分析了几种方法的MDB指标的异同,进一步厘清几种方法的本质和联系。
HTML
-
设线性化误差方程式为:
式中,V为n×1维残差向量;A为秩rank(A)=u的n×u阶设计矩阵;x为u×1维的未知参数向量;l为n×1维观测值向量;P为n×n阶观测值权阵,且P=Qu-1=σ02(D(l))-1。
按最小二乘法求得参数估值:
改正数平方和为:
改正数为:
改正数的协因数阵为:
-
模型误差的存在将会反映到平差结果的残差上,因此可以通过检验验前和验后单位权方差是否一致来判断模型误差的正确性[19]。
用来检验的原假设和备选假设为:
式中,$\hat \sigma _0^2$为后验单位权方差因子;σ02为先验单位权方差因子;f为自由度(多余观测数),其满足关系式:
检验量为:
如果拒绝原假设,则表示观测值含有粗差(此处仅考虑粗差的情形),表明检验量不服从中心卡方分布,其非中心化参数为[17]:
式中,z为粗差向量的真值;H为相应的粗差向量的系数矩阵。
-
为了检验备选假设,现将粗差归入函数模型,设有n2观测量含有粗差,其取值范围为[20]:
按照均值漂移模型原理,将观测值l分为l1和l2两组,第一组观测值l1不含粗差,其维数为n1,第二组观测值l2可能含有粗差∇S,其维数为n2。∇S视为观测值l2的均值漂移向量,即[15]
包含未知参数∇S的扩展模型为:
设H=[0 E]T,E为n2×n2单位阵,式(12)可表示为:
根据最小二乘原理得到估值[15]:
式中,R=I-A(ATPA)-1ATP [21]。
顾及粗差的未知参数x的估值为[15]:
式(14)可以表示为[22]:
式中,${\mathit{\boldsymbol{\hat l}}_1} = {\mathit{\boldsymbol{A}}_1}{\mathit{\boldsymbol{\hat x}}_\nabla }, {\mathit{\boldsymbol{\hat l}}_2} = {\mathit{\boldsymbol{A}}_2}{\mathit{\boldsymbol{\hat x}}_\nabla }$。其改正数平方和为:
若已知单位权方差σ02,则可构造统计量[15]:
非中心化参数为:
其原假设和备选假设为:
如果原假设成立,则统计量T1为中心化卡方分布,式(21)中的非中心化参数λ=0。若给定显著性水平α0,便可对统计量进行检验。若T1>χ1-α0, n22,则拒绝原假设H0。
-
内部可靠性是在一定的假设检验条件下,系统发现模型误差(包括系统误差和粗差)的能力,能以一定的检验功效β0,通过显著水平为α0的统计量检验可发现的粗差的下界值。根据α0、β0可以确定统计量的非中心参数λ0。
-
若只含有一个粗差的单个备选假设情形,其MDB为:
式中,${\mathit{\boldsymbol{h}}_i} = {\left[{\begin{array}{*{20}{l}} 0& \cdots &1&0& \cdots &0 \end{array}} \right]^{\text{T}}}$为单位化向量,第i个元素为1。
单个粗差的可控性数值为:
式中,qii为观测值协因数阵Qll相应的对角线元素。
-
如果存在多个粗差的情形,给定非中心化参数λ0的情况下,不能从式(21)得到MDB。若MDB向量表示为模∇0(S)和单位方向矢量∇u的形式,则向量的下界域为:
对于多个粗差可发现的下界值在不同方向上是不同的,因此可以通过计算矩阵PSS的最小和最大特征值及其特征向量,从而得知在什么方向上多个粗差最难和最容易发现[15]。
根据文献[17]给出了第i个方向的最大MDB值,根据广义Rayleigh-Ritz理论,其值满足关系式:
式中,cn2i为1×n2的向量,其元素对应于第i个粗差为1,其他为0,所以cTn2icn2i为n2×n2方阵,对应于第i个粗差的对角线元素为1,其余均为0。式(26)的特征值和特征向量可以通过下式解算得到:
式中,λ为特征值;u为特征值对应的特征向量。
因此向量∇0S在第i个观测值方向的误差下界域为:
相应特征值对应的特征向量为:
其可控数值为[17]:
其中,
1.1. 高斯马尔科夫线性模型最小二乘估计
1.2. 平差模型的检验
1.3. 粗差估值与假设检验
1.4. 内部可靠性理论
1.4.1. 一个粗差单个备选假设下的内部可靠性
1.4.2. 多个粗差单个备选假设下的内部可靠性
-
文献[7]将部分最小二乘法应用于粗差的定位和定值,仅采用未受污染的观测值根据最小二乘准则计算待估参数,从而得到受污染的观测值l2对应的预测残差为[7]:
V2对应的协因数矩阵为:
式中,M′=A2(A1TQ11-1A1)-1A1TQ11-1。
原假设和备选假设为:
若已知单位权方差σ02,则可以构建下列统计量:
-
若只含有一个粗差的单个备选假设情形,统计量为:
在原假设条件下, wi为标准正态分布,在备选假设下,wi有如下平移参数:
若只含有一个粗差的单个备选假设情形,其MDB为:
式中,δ0为正态分布的非中心参数,在给定显著水平α和检验功效1-β的情况下可以计算得到。
-
在多个粗差情形下,其粗差向量的下界域满足如下关系:
式中,V20为粗差向量的下界域;λ0为自由度n2卡方分布非中心化参数。
给出了第i个方向的最大MDB值,根据广义Rayleigh-Ritz理论,其值满足关系式:
式中,cn2i的含义同式(26)。式(39)的特征值和特征向量可以通过下式解算得到:
式中,λ为特征值;u为特征值对应的特征向量。
因此向量V20在第i个观测值方向的误差下界域为:
其特征值对应的特征向量为:
2.1. 一个粗差单个备选假设下的内部可靠性
2.2. 多个粗差单个备选假设下的内部可靠性
-
为了求解秩亏方程,需要选择r(r≥u)个拟准观测值,将第一组选择为拟准观测值,则有:
式中,${\mathit{\boldsymbol{\bar P}}_{11}} = \mathit{\boldsymbol{Q}}_{11}^{-1}$为l1对应的权矩阵。式(43)等价于:
其解为:
式中,$\mathit{\boldsymbol{G}} = \left( {\mathit{\boldsymbol{A}}_1^{\text{T}}{{\mathit{\boldsymbol{\bar P}}}_{11}}\;\;\;\;\;\;{\bf{0}}} \right)$; ${\mathit{\boldsymbol{ \boldsymbol{\hat \varDelta} }}_2}$为受粗差污染观测值的真误差估值[6]。
对式(46)进一步分析可得到:
从式(47)可知,用QUAD法解算的粗差估值和用部分最小二乘解算的结果一致。
应用协因数传播律可知,${\mathit{\boldsymbol{ \boldsymbol{\hat \varDelta} }}_2}$对应的协因数矩阵与部分最小二乘法结果相同,即为:
式中,M′=A2(A1TQ11-1A1)-1A1TQ11-1。
因此其后续的假设检验及可靠性理论与部分最小二乘法一致。
-
由于l2为受粗差污染的观测值,所以R2Δ2对于观测值残差V的影响远远大于R1Δ1[4]。于是得到近似关系[4]:
按照最小二乘准则,其估值为[4]:
式(50)为受粗差污染观测值的真误差估值公式,表现为粗差的影响[4]。
根据协因数传播律,得到${\mathit{\boldsymbol{ \boldsymbol{\hat \varDelta} }}_2}$的协因数阵为:
式中,M″=(R2TR2)-1R2TR1。于是粗差统计量为:
非中心化参数为:
原假设和备选假设为
若原假设成立,则非中心化参数为0;若拒绝原假设,则其非中心化参数不为0。
-
若只含有一个粗差的单个备选假设情形,统计量为:
若原假设成立, 则正态分布的非中心化参数为0;若接受备选假设,则正态分布的非中心化参数为:
若只含有一个粗差的单个备选假设情形,其MDB为:
式中,δ0为正态分布的非中心化参数,在给定显著水平α和检验功效1-β的情况下可以计算得到。
-
在多个粗差情形下,其粗差向量的下界域满足如下关系:
式中,Δ20为粗差向量的下界域;λ0为自由度n2卡方分布非中心化参数。
给出了第i个方向的最大MDB值,根据广义Rayleigh-Ritz理论,其值满足如下关系式:
式(58)的特征值和特征向量可以通过下式解算得到:
式中,λ为特征值;u为特征值对应的特征向量。
因此向量Δ20在第i个观测值方向的误差下界域为:
其特征值对应的特征向量为:
4.1. 一个粗差单个备选假设下的内部可靠性
4.2. 多个粗差单个备选假设下的内部可靠性
-
利用文献[17]中水准网,如图 1所示,设观测值为6个观测高差,已知点高程为1 000 m,观测值的协因数矩阵为:
方案Ⅰ:假定为单个粗差的情形,设显著水平α0=0.1%,检验功效β0=80%,其卡方分布的非中心化参数λ0=17.07,分别设:①观测量的协因数矩阵为QlI1=I;②观测量协因数矩阵QlI2=diag(Ql);③观测量协因数阵QlI3=Ql。不同方法粗差可发现下界域MDB计算结果比较如表 1所示。
方案 方法 1 2 3 4 5 6 σli 2.35 1.97 0.89 2.32 0.45 1.18 方案Ⅰ①MDBi/m 数据探测法 5.630 4 6.662 0 6.662 0 5.266 8 5.266 8 6.081 5 PLS法 5.630 4 6.662 0 6.662 0 5.266 8 5.266 8 6.081 5 QUAD法 5.630 4 6.662 0 6.662 0 5.266 8 5.266 8 6.081 5 LEGE法 5.630 4 6.662 0 6.662 0 5.266 8 5.266 8 6.081 5 方案Ⅰ②MDBi/m 数据探测法 10.686 1 10.065 4 10.065 4 9.768 7 6.501 2 8.393 5 PLS法 10.686 1 10.065 4 10.065 4 9.768 7 6.501 2 8.393 5 QUAD法 10.686 1 10.065 4 10.065 4 9.768 7 6.501 2 8.393 5 LEGE法 10.723 1 10.400 4 10.400 4 9.771 6 6.793 3 8.695 5 方案Ⅰ③MDBi/m 数据探测法 2.979 5 10.346 1 10.346 1 2.595 6 1.321 7 2.591 7 PLS法 9.397 1 11.512 9 10.704 8 9.530 5 2.020 0 5.032 8 QUAD法 9.397 1 11.512 9 10.704 8 9.530 5 2.020 0 5.032 8 LEGE法 9.161 6 13.575 2 13.575 2 6.110 2 3.465 3 8.249 8 Table 1. Comparison of MDBs of Different Methods
方案Ⅱ:假定为两个粗差的情形,设显著水平α0=0.1%,检验功效β0=80%,其卡方分布的非中心化参数λ0=19.67,分别设:①观测量的协因数矩阵为QlII=I;②观测量协因数阵矩QlII=diag(Ql);③观测量协因数阵矩QlII=Ql。不同方法MDB计算结果比较如表 2-4所示。
MDBi/m MDBi/m i j 数据探测法 PLS法或QUAD法 LEGE法 i j 数据探测法 PLS法或QUAD法 LEGE法 1 2 7.012 5 7.012 5 7.012 5 4 1 5.867 1 5.867 1 5.867 1 1 3 7.012 5 7.012 5 7.012 5 4 2 5.725 7 5.725 7 5.725 7 1 4 6.272 2 6.272 2 6.272 2 4 3 5.725 7 5.725 7 5.725 7 1 5 6.272 2 6.272 2 6.272 2 4 5 7.242 5 7.242 5 7.242 5 1 6 7.681 8 7.681 8 7.681 8 4 6 6.272 2 6.272 2 6.272 2 2 1 8.297 3 8.297 3 8.297 3 5 1 5.867 1 5.867 1 5.867 1 2 3 ∞ ∞ ∞ 5 2 5.725 7 5.725 7 5.725 7 2 4 7.242 5 7.242 5 7.242 5 5 3 5.725 7 5.725 7 5.725 7 2 5 7.242 5 7.242 5 7.242 5 5 4 7.242 5 7.242 5 7.242 5 2 6 7.681 8 7.681 8 7.681 8 5 6 6.272 2 6.272 2 6.272 2 3 1 8.297 3 8.297 3 8.297 3 6 1 8.297 3 8.297 3 8.297 3 3 2 ∞ ∞ ∞ 6 2 7.012 5 7.012 5 7.012 5 3 4 7.242 5 7.242 5 7.242 5 6 3 7.012 5 7.012 5 7.012 5 3 5 7.242 5 7.242 5 7.242 5 6 4 7.242 5 7.242 5 7.242 5 3 6 7.681 8 7.681 8 7.681 8 6 5 7.242 5 7.242 5 7.242 5 Table 2. The MDB for Two Outliers when Observations are Independent and Equal Precision
MDBi/m MDBi/m i j 数据探测法 PLS法或QUAD法 LEGE法 i j 数据探测法 PLS法或QUAD法 LEGE法 1 2 11.811 7 11.811 7 11.823 4 4 1 10.489 4 10.489 4 10.491 2 1 3 11.811 7 11.811 7 11.823 4 4 2 10.490 1 10.490 1 10.495 4 1 4 11.474 5 11.474 5 11.512 1 4 3 10.490 1 10.490 1 10.495 4 1 5 12.761 3 12.761 3 12.799 2 4 5 13.550 1 13.550 1 13.738 1 1 6 14.164 5 14.164 5 14.165 7 4 6 10.495 3 10.495 3 10.495 9 2 1 11.125 7 11.125 7 11.148 1 5 1 7.763 7 7.763 7 8.121 2 2 3 ∞ ∞ ∞ 5 2 7.969 9 7.969 9 8.008 1 2 4 10.808 7 10.808 7 11.167 1 5 3 7.969 9 7.969 9 8.008 1 2 5 12.339 3 12.339 3 12.421 8 5 4 9.017 9 9.017 9 9.325 5 2 6 14.164 5 14.164 5 14.165 7 5 6 10.495 3 10.495 3 10.495 9 3 1 11.125 7 11.125 7 11.148 1 6 1 11.125 7 11.125 7 11.148 1 3 2 ∞ ∞ ∞ 6 2 11.811 7 11.811 7 11.823 4 3 4 10.808 7 10.808 7 11.167 1 6 3 11.811 7 11.811 7 11.823 4 3 5 12.339 3 12.339 3 12.421 8 6 4 9.017 9 9.017 9 9.325 5 3 6 14.164 5 14.164 5 14.165 7 6 5 13.550 1 13.550 1 13.738 1 Table 3. The MDB for Two Outliers when Observations are Uncorrelated
i j MDBi/m i j MDBi/m 数据探测法 PLS法或QUAD法 LEGE法 数据探测法 PLS法或QUAD法 LEGE法 1 2 3.514 6 8.136 9 8.625 5 4 1 9.834 4 9.839 6 9.838 6 1 3 3.514 6 8.136 9 8.625 5 4 2 2.991 9 9.917 2 8.274 2 1 4 11.289 2 13.300 5 13.746 2 4 3 2.991 9 9.917 2 8.274 2 1 5 18.465 2 18.488 7 18.584 8 4 5 14.420 5 15.474 7 16.774 1 1 6 14.030 6 15.833 9 14.428 6 4 6 7.354 0 9.251 9 7.413 1 2 1 12.204 1 12.516 1 13.011 3 5 1 8.190 9 8.320 7 9.276 0 2 3 ∞ ∞ ∞ 5 2 1.635 2 1.637 8 2.908 9 2 4 11.925 7 12.595 3 17.975 9 5 3 1.635 2 1.637 8 2.908 9 2 5 12.800 0 13.271 6 14.790 7 5 4 7.343 0 7.416 2 10.593 0 2 6 14.030 6 14.302 7 14.428 6 5 6 7.354 0 7.427 1 7.413 1 3 1 12.204 1 12.473 6 13.011 3 6 1 12.204 1 12.725 6 13.011 3 3 2 ∞ ∞ ∞ 6 2 3.514 6 5.009 9 8.625 5 3 4 11.925 7 12.096 6 17.975 9 6 3 3.514 6 6.340 8 8.625 5 3 5 12.800 0 13.006 2 14.790 7 6 4 7.343 0 7.345 9 10.593 0 3 6 14.030 6 14.749 6 14.428 6 6 5 14.420 5 14.422 0 16.774 1 Table 4. The MDB for Two Outliers when Observations are Correlated
-
从表 1可以看出,对于单个粗差单个备选假设条件下,若观测值独立等精度时,数据探测法、QUAD法、PLS法和LEGE法的MDB相同;若观测值独立不等精度时,数据探测法、QUAD法和PLS法的MDB相同,LEGE法的MDB略大;若观测值相关时,QUAD法和PLS法的MDB相同,与数据探测法、LEGE法的MDB都不同,而数据探测法的MDB最小。
从表 2、表 3、表 4可以看出,对于两个粗差单个备选条件下,若观测值独立等精度时,数据探测法、QUAD法、PLS法和LEGE法的MDB相同;若观测值独立不等精度时,数据探测法、QUAD法和PLS法的MDB相同,LEGE法的MDB略大;若观测值相关时,QUAD法和PLS法的MDB相同,且与数据探测法、LEGE法的MDB都不同,数据探测法的MDB最小。
5.1. 最小可探测偏差MDB比较
5.2. 结果分析
-
本文对基于均值漂移模式的几种粗差定位和定值方法进行了深入分析,给出了几种方法的MDB计算公式;通过数据分析,当观测值独立等精度时,几种方法的MDB值相同,说明几种方法可能探测的最小边界误差效果相同,当观测值相关或观测精度不等时,几种方法的MDB不同,其中数据探测法的MDB值最小,说明其可能探测出的最小边界误差最小。
本文仅限于对单个备选假设下的内部可靠性进行比较,对于从两个多维备选假设下的内部可靠性理论,其可发现且可与其他模型误差相区分的模型误差下界值,以及不可发现和区分的模型误差对平差结果的影响,需进行更深入的比较分析和研究。