超定测距定位方程参数估计的松弛重心迭代法

袁兴明，孙振

(1.山东工业职业学院 建筑与信息工程学院，山东 淄博 256414；2.山东理工大学 测绘工程系，山东 淄博 255049)

0 引言

大地测量与导航定位中，测距定位方程为超定非线性方程，可采用非线性最小二乘理论进行求解.该方程传统上依据泰勒级数展开取至一阶项，采用线性平差估计求解，这会导致信息量的缺失和模型特征的改变[1-2].当非线性模型强度和残差较大时，传统的线性化平差估计并不能取得有效的结果[3-4].若此时伴随不适定问题，会进一步降低参数估计解的精度和可靠性.固有曲率和参数效应是刻画非线性模型线性化的数量指标，能够评估推断效果的优劣程度、适应条件和容许误差，但无法评估非线性平差的收敛平稳情况[5-6].研究表明，非线性扰动主要来源于线性近似时系数矩阵的扰动、附加的截断误差及正交过程[7-8].

迭代估计是求解非线性无约束最优化方法常用的数值方法，即通过构造一定的迭代序列使其满足一定的收敛条件.虽然牛顿法能够处理非线性最小二乘问题，但由于其储存成本等因素，更倾向于采用高斯-牛顿法[9-10].研究表明，高斯-牛顿法的适应条件为残差小和非线性强度弱的非线性方程，可是当非线性方程含有不适定问题时，传统的线性化平差估计和经典数值方法如高斯-牛顿法会由于观测数据误差扰动而产生强烈的不稳定特征，导致参数估计解不稳定，甚至求解失效[11-14].由此可知，根据实际问题的线性化程度、秩亏或者病态程度，选择合适的非线性平差模型就显得尤为重要.

对于测距定位问题，非线性最小二乘解是观测向量末端以观测权为质量质点系的重心.针对测距定位方程最小二乘解性质，来构造非线性数值算法的研究较为困难.在这方面，Xue等[15]作出了一定的突破，他基于非线性最小二乘解性质发展了简单，且无需矩阵求逆的稳定数值方法，即重心迭代法，然而，该方法由于线性搜索因子1/tr(P)过于保守，收敛速度可能非常缓慢，具有较低的收敛效率，实用价值不高.为增加重心迭代法的实用价值，本文基于残差最小步长准则提出了一种改进的重心迭代法，即松弛重心迭代法来提高重心迭代法的收敛效率.最后采用全球卫星导航系统（GNSS）数据和水下定位数据，验证了该方法的主要结果.

1 超定测距定位方程的松弛重心迭代算法

设测距定位观测方程的向量表达式为

式中：L为m×n的观测向量；d(x)=∥xi−x∥2=为已知坐标至未知坐标的欧式距离；x=[x1,x2,···,xm]∈Rm为未知点构成的向量集；xi=[xi1,xi2,···,xim]∈Rm为已知坐标点构成的向量集；ε为观测误差.

非线性最小二乘平差实质是求泛函极小值的最优化问题，常用来处理式(1)，则非线性最小二乘参数估计解为

式中，

为残差向量.运用无约束非线性最小二乘求解式(1)，非线性最小二乘目标函数即式(3)满足如下正交条件方程：

由式(4)和式(5)联立可得

根据式（6），文献[15]导出了测距定位方程参数估计的重心迭代法，其迭代公式为

式（8）为重心迭代公式，式（7）实际是一种最速下降法，即选取重心迭代法不依赖于初始值的精度，具有计算简单，无需矩阵取逆、计算海森矩阵的优点.重心迭代法实质上是一种最速下降法，收敛速度主要依赖于线性搜索因子1/tr(P)的影响.1/tr(P)是重心法求解的保守策略，易受到观测数据个数n的影响，即随着观测数目增多，其线性搜索策略可能越保守，比如当权阵为单位阵，观测数目n为20时，重心迭代法的步长为0.05，这时重心迭代法具有较低的收敛效率.为增加重心法的实用性，必须提高其收敛效率.为此，本文将对重心迭代法进行如下改进

式(9)为松弛重心迭代公式，其中，ωk为松弛参数，主要依据残差性质来自适应调整迭代步长，提高数值解的整体效果.

文献[16]提出了残差最小准则确定步长的思想，验证了方法的计算效率.为此，本文采用了残差最小准则来确定松弛因子.残差最小准则方法充分采用目标函数的导数信息和函数值信息，相比于传统线搜索方法如黄金分割法和抛物线法具有一定优势.下面给出了该方法确定松弛参数的严密推导公式.设第k+1次的拟合残差为松弛因子的函数，即

式中，d(xk+1)为xk+1处的函数计算值，可将其在xk依据泰勒级数展开得

考虑到xk+1=xk+ωkd(xk)，联合式(1)和式(2)，可将松弛因子函数化简为

将式(9)代入式(12)联合式(4)可得

将式(13)关于ω求其偏导并进行化简可得

若要k+1次残差最小，需满足 ∂R(ω)/∂ω=0，可得松弛因子的确定公式.考虑到迭代过程，可将步长记为如下形式

2 实验分析

为了测试本文提出的新方法在观测病态和良态的表现效果，分别选取GNSS观测数据和南海水下定位数据进行测试.本文设置迭代收敛条件为该条件所消耗成本远高于重心法和松弛重心法的迭代成本，虽然本文用该条件来验证该方法.

2.1 BDS伪距单点定位

本文采用2017年4月28日采集的GNSS观测数据（该数据文件包括GPS、BDS及GLONASS的三系统观测值），采样间隔为1 s，选取其中1个历元北斗卫星导航系统（BDS）观测数据进行伪距单点定位实验.表1给出了原始伪距观测值和误差改正后的伪距观测值（误差改正包括：电离层改正、对流层改正和卫星钟差改正等）.迭代初值按照文献[15]给出的方法进行计算.

表1 单历元的BDS观测数据

分别采用高斯-牛顿迭代法、重心迭代法和松弛重心迭代法来处理水下实测数据，解算结果如表2所示.由表可知，高斯-牛顿法、重心法和松弛重心法的数值收敛解相同但耗时最短.该算例表明，高斯-牛顿法相对于重心法和松弛重心法具有更好的局部收敛性质，这是因为BDS卫星距离位置待测点相差3万多千米，非线性强度弱.松弛重心迭代法依据残差准则自适应确定松弛因子，提高了重心迭代法的收敛效率，节省了计算成本.三种数值相比较而言，高斯-牛顿法计算成本较低，建议当距离观测方程良好且距离观测量非常大时，采用高斯-牛顿法.

表2 不同算法的解算结果

图1给出了重心法和松弛重心法的点位迭代序列图.图中k表示迭代次数.由图可知，松弛重心迭代法依据残差准则自适应确定松弛因子，不会对最终数值收敛解和迭代序列稳定性产生影响，而明显提高了重心迭代法的收敛，节省了运算成本.以上表明，松弛重心迭代法在继承了重心迭代法迭代序列稳定，无需计算海森矩阵的优点之外，相对于重心迭代法具有更好的局部收敛性质，实用价值更高.

图1 重心迭代算法和松弛重心迭代算法的点位迭代序列图

2.2 水下应答器绝对位置计算

本文采用水下定位实测数据进行验证，数据主要是测量船围绕应答器航行获取.在1圈数据中选取15个相邻数据进行验证，数据坐标如图2所示.假设观测距离为等精度观测，要求根据15个已知点坐标和观测距离求待测点的坐标.

以1圈数据的单点定位结果作为真值.数据采集过程中，测量船在水平面航行，Z方向上具有较强的病态性；此外，测量船航行过程中，采集数据相隔时间较短，X方向和Y方向也会产生微弱的病态性.迭代初值按照文献[15]给出的方法进行计算.经计算初始设计矩阵条件数为cond(N(x0))=3.97×1013.分别采用高斯-牛顿法、重心迭代法和松弛重心迭代法来处理水下实测数据，解算结果如表3所示.

图2 水下实测数据平面坐标变化图

表3 不同算法的解算结果

由表3可知，高斯-牛顿法易受到线性初值和方程组模型态性影响，迭代序列受到观测数据误差扰动而无法收敛，建议处理病态问题时一般不采用高斯-牛顿法.重心迭代法和松弛重心迭代法无需矩阵取逆，能够避免观测误差扰动对取逆算子扰动的影响.松弛重心迭代法根据残差准则来自适应确定松弛因子，主要是用来提高重心迭代法的收敛效率，并不会对数值收敛结果产生影响，其结果与重心迭代法的结果相同.由表进一步能够看出，重心迭代法具有较低的收敛效率，其消耗时间为132.4420 s，迭代次数为28280次；松弛重心迭代法明显提高了重心迭代法的收敛效率，时间大约可节省130 s，迭代次数大约是重心迭代法的1/8.

图3给出了松弛重心迭代法和重心迭代法的点位迭代序列图.由图可知，重心迭代法和松弛重心迭代法在X方向和Y方向点位迭代序列变化较为平缓；在Z方向上，重心迭代法和松弛重心迭代法点位迭代序列整体变化较为平缓，可是分别在迭代1425次和3410次都产生了较大的扰动，这是由于短程水下测距定位方程病态性主要集中于Z方向，由于观测数据误差而导致迭代序列产生不稳定特征.以上结果表明，松弛重心迭代法在继承重心迭代法迭代序列稳定、无需矩阵求逆和计算海森矩阵的优点之外，具有更好的局部收敛性质，实用价值更高.

图3 松弛重心迭代算法和重心迭代算法的点位迭代序列图

3 结束语

无论处理病态问题和良态问题，重心迭代法由于线性搜索因子 1/tr(P)过于保守，具有较低的收敛效率.本文提出的方法在继承了重心迭代法不依赖于初始值的精度，具有计算简单，无需矩阵求逆、计算海森矩阵的优点外，依据残差最小准则确定松弛参数来自适应更新迭代步长，明显提高重心迭代法的收敛速度.算例结果表明：松弛重心迭代法无论处理良态问题和不适定问题，相比于重心迭代法具有更好的局部收敛性质，应用价值更高.

高斯-牛顿法由于忽视了距离观测方程的二阶项信息，适合处理残差小且观测态性良态测距方程.但是当处理不适定测距方程时，高斯-牛顿法由于设计矩阵秩亏或者病态而解算失败.本文提出的新方法在处理不适定问题时，相比于高斯-牛顿法有更好的收敛效果.为此，我们应该充分考虑松弛重心迭代法和高斯-牛顿法的优点来处理尽可能多的问题.