最近几年，作为变量误差模型(errors-in-variables,EIV)的参数估计方法，整体最小二乘(total least-squares,TLS)或加权整体最小二乘(weighted TLS,WTLS)方法在大地测量领域获得了广泛关注^{[1, 2]}，但大多数学者的工作主要集中于WTLS估计的算法研究^[3-8]。WTLS估计是在一般的变量误差模型中引入观测误差和系数矩阵误差的方差-协方差矩阵的参数估计方法^[3]。事实上，变量误差模型中并不是所有的系数矩阵元素都是随机的。例如，在坐标变换中，一些元素不含误差，可以认为是固定的，而另外一些元素则是由观测值所组成的，可以认为是随机的。顾及这一事实，通常考虑的变量误差模型被推广到更一般的部分变量误差模型^[9]。部分变量误差模型的WTLS估计方法，其优势在于不用拉格朗日乘子且法方程未知参数的个数大大减少，因而获得了更多关注^{[10, 11]}。但是，WTLS估计与经典的LS估计一样，对粗差很敏感，不具有抵御粗差影响的能力^[12-17]。

如何有效地处理观测数据中的粗差一直是测量领域关注和研究的热点^{[18, 19]}，抗差估计是一种处理粗差的有效方法，并且针对Gauss-Markov模型的抗差估计成果较多^[20-26]。虽然针对变量误差模型，文献[14-17]提出了为数不多的几个抗差估计，但这些估计均采用残差和验后单位权方差确定降权因子，易受粗差的影响，进而影响其抗差功效。由于部分变量误差模型的系数矩阵也可能含有粗差，所以经典的抗差估计方法并不能直接应用于部分变量误差模型，更为有效的方法值得研究。为此，本文将在提出部分变量误差模型WTLS估计的两步迭代解法的基础上，运用抗差LS估计的等价权方法，提出一种整体抗差最小二乘(TRLS)估计方法，并构建一致最大功效统计量确定降权因子。针对WTLS估计两步迭代解法的特点，设计两个不同的降权方案：第1方案是在估计系数矩阵元素时，不对观测值降权，仅对系数矩阵降权;第2方案是在估计系数矩阵元素时，既对系数矩阵降权，同时也对观测值降权。

1 WTLS估计的两步迭代解法

考虑如下部分变量误差模型^[9]

式中，L是n×1观测向量;X是t×1未知参数向量;I_n是n×n单位阵;Δ是n×1随机误差向量;h是由系数矩阵A中零元素和常数元素所组成的nt×1向量;B是给定的nt×s矩阵,s表示A中不同随机元素的个数;a是由A中不同随机元素所组成的s×1向量;表示a的真值;e=a-。随机模型为

式中，σ²为未知的单位权方差。

下面提出部分变量误差模型(1)的WTLS估计的两步迭代解法。首先，将任意给定的X的初值X⁽⁰⁾代入式(1)，得到

令，则式(3)可表示为

再令

并应用LS原理^[27],可得向量的估计为

对应的残差向量为

然后，将回代入式(1)的第1个方程，得到

应用数学函数invvec，将记为A^*。其中，invvec是一个将nt×1向量变换为n×t矩阵的函数，例如

式中，β_i(i=1,2,…,n)均为t×1列向量。于是式(8)可写为

再次应用LS原理^[19]，可得未知参数向量X的估计为

对应的残差向量为

最后，可得单位权方差σ²的估计为

值得指出的是，上述WTLS估计的两步迭代解法本质上与现有算法^[9]等价。鉴于上述公式是在LS框架下推导得到的，所以经典LS估计和抗差LS估计的相关理论和方法^[21-24]都可以直接应用于部分变量误差模型。

2 整体抗差最小二乘(TRLS)方法

如上所述，当部分变量误差模型(1)中的观测值和系数矩阵元素同时含有粗差时，对式(6)和式(11)运用抗差LS估计的等价权方法^[21-24]，可得向量和未知参数向量X的整体抗差LS(TRLS)解分别为

和

式中，和为相应的等价权阵。

根据基于敏感度分析的抗差估计理论^[26]，本文将按下列方式确定式(14)、(15)中的等价权阵

式中，ω_i为降权因子，可通过Huber权函数^[20]进行确定

式中，为残差或者标准化残差;c为临界值，其取值区间为[1.0,2.0]^[20]。降权因子ω_i也可通过IGGⅢ权函数^[22]进行确定

式中，c₀和c₁的取值区间分别为[1.0,2.0]和[3,5,6.0]。然而，在现有变量误差模型的抗差估计方法^[14-16]中，其降权因子是由残差v_i和单位权方差的验后估计²共同决定的

式中，m为给定的常数。

众所周知，残差和单位权方差的验后估计易受粗差的影响^[26];因此，按式(19)确定降权因子，会造成抗差功效的损失。基于敏感度分析的抗差估计理论表明，基于一致最大功效统计量的抗差估计优于基于残差或标准化残差的抗差估计^[26]。令人遗憾的是，通过现有的WTLS估计方法均无法构造此一致最大功效统计量，为此本文利用第1节提出的WTLS估计的两步迭代解法，构造如下统计量，以将基于一致最大功效统计量的抗差估计拓展到部分变量误差模型

式中，

相应的，式(14)、(15)中的等价权阵和分别构造如下

对应的降权因子分别为

式中，，其中未知的单位权方差σ²可利用具有高崩溃污染率的最小二乘中位数(least median square,LMS)方法^{[18, 28]}估计

注意到式(14)中既包含系数矩阵的权，也包含观测值的权，为此提出如下两个降权方案，以便作比较分析。

方案1，在式(14)中仅对P_e进行降权。

方案2，在式(14)中对进行降权，即同时对P_e和P_Δ进行降权。

综上所述，计算TRLS估计的迭代算法设计如下：

步骤1，给定初值

步骤2，构造等价权矩阵和

步骤3，计算向量和未知参数向量估值

步骤4,如果或者迭代次数k超过给定的常数，则迭代结束，否则，转向步骤2。

3 算例及其分析

算例1，2D仿射变换的数学模型为

式中，(x_s,y_s)和(x_t,y_t)分别为同一点在旧坐标系和新坐标系中的坐标;a_i、b_i和c_i(i=1,2)为坐标变换参数。本次试验采用文献[15]提供的部分数据。

假定变换参数的参考值为

而无误差的观测数据如表 1所示。

根据Monte Carlo策略^[15]，将协方差矩阵为

的随机误差分别加在13个点的无误差的原始观测值上，得到模拟观测值，并且在第4个点的旧坐标分量x_s上加入大小为2 m的粗差。临界值ε和m分别选为0.000 01和1.5，最大迭代次数为200，权函数取Huber权函数。为消除随机因素的影响，随机模拟运行了200次。单个粗差情形下，由各方案计算得到的变换参数估值与参考值之间的欧氏距离如图 1所示，图 2则给出了单位权方差的估值。

与此同时，对上述200次模拟结果进行了统计分析，其结果如表 2、表 3和表 4所示。从表中不难看出，当系数矩阵元素受到粗差污染时，WTLS估计扭曲得十分厉害，其抗差性很差; 而基于方案1的TRLS估计则较为有效地削弱了粗差对未知参数估值的影响，获得了可靠的参数估值;但方案2并没有达到预期的抗差效果，其稳定性还不如WTLS估计，这说明在估计系数矩阵时，无须进行相应的降权。如果将TRLS估计抵御粗差的成功率定义为

则两种方案的成功率分别为100%和29.5%，这表明方案1的TRLS估计更具稳定性和可靠性。

下面考虑多个粗差的情形。为此，继续在第4个点的旧坐标分量y_s上加入大小为2 m的粗差，同时在第2个点的新坐标两分量x_t和y_t上也都加入大小为2 m的粗差。模拟计算200次，其变换参数估值与参考值之间的欧氏距离以及单位权方差估计分别如图 3和图 4所示，而变换参数估值与参考值之间的平均距离如表 5所示。不难看出，方案2的TRLS估计基本上不具有抵御粗差的能力，而方案1的TRLS估计则克服了粗差的影响，并且成功率达到100%。

综上可知，不论是单个粗差还是多个粗差，方案2的TRLS估计均不具有较强的抗差能力，而方案1的TRLS估计则能够有效削弱单个或多个粗差对未知参数估值的不良影响，这说明在处理部分变量误差模型中的粗差时，应采用方案1的TRLS估计方法。为了显示本文算法的优越性，在相同模拟条件下，仍然应用前面的相关数据，将本文算法与基于残差和验后单位权方差确定降权因子(式(19))的抗差估计进行比较，其计算结果如表 6和表 7所示。相关结果表明本文算法优于现有算法。

算例2，考虑如表 8所示的实测数据^[29]，该数据是关于20个孩子的血液样本中血清中卡那霉素的水平。其中一项数据是通过足跟方式获取，而另外一项则是用导管的方式获取。文献[29]指出表 8中的第2个观测值含有粗差，据此本文将去掉该观测值所得到的WTLS估值作为未知参数的参考值。在计算过程中，迭代终止的条件为相邻迭代估值之间的距离小于0.000 01，或者是迭代次数超过200次，权函数取Huber权函数。分别采用WTLS估计、方案1和方案2的TRLS估计进行解算，并与文献[15]中的结果进行对比，如表 9所示。表 10则给出了由各种方法得到的单位权方差估值。图 5表示的是不同参数估值对应的拟合结果。

从表 9可以看出，方案1和方案2的TRLS方法均具有一定的抵御粗差的功效，但就未知参数估值与参考值之间的欧氏距离来看，方案1明显优于方案2。与文献[15]中的结果进行比较发现，本文的TRLS估值更稳健。从表 10可以看出，本文算法的计算精度高于文献[15]，这也有力地验证了本文方法的有效性和优良性。

4 结 论

(1) 部分变量误差模型的WTLS估计同LS估计一样，对粗差很敏感，即抗差性能很差。

(2) 提出了部分变量误差模型WTLS估计的一种两步迭代解法，该解法的计算公式均是在LS框架下推导得到的，具有易于理解和可操作性强等特点，为构建部分变量误差模型的整体抗差估计奠定了基础。

(3) 利用部分变量误差模型的两步迭代解法，构建了与抗差LS估计相类似的TRLS解式，保持了传统的数据处理方式，并构造了一致最大功效统计量以确定降权因子，其中未知的单位权方差可采用具有高崩溃污染率的LMS方法进行估计。

(4) 结合部分变量误差模型两步迭代解法的特点，设计了两个降权方案：第1个方案是在估计系数矩阵元素时，不对观测值进行降权，仅对系数矩阵降权，而第2个方案则是在估计系数矩阵元素时对观测值和系数矩阵同时降权。

(5) 2D相似变换模拟算例的计算结果表明，单个粗差情形下两种不同降权方案的TRLS估计均具有一定的抗差能力，但方案1的抗差功效更优。多个粗差情形下，方案2的TRLS估计基本上不具有抗差性。因此，处理部分变量误差模型中的粗差时，应采用方案1的TRLS估计。线性拟合实测数据的计算结果表明本文方法优于现有变量误差模型的抗差估计方法。

(6) 需要注意的是，统计量(式(20)、式(21))并非严格意义上的一致最大功效统计量，因此，部分变量误差模型的整体抗差估计及其算法值得进一步研究。