整体最小二乘估计的深入研究.docx

整体最小二乘估计的深入研究 摘要: 整体最小二乘法是一种较为先进的最小二乘法结构，整体最小二乘法认为回归矩阵存在干扰，在计算最小二乘解时考虑了这个因素，而在一般最小二乘法时没有考虑该因素的影响。整体最小二乘法应用广泛，得到效果也比较好。本文主要讨论了整体最小二乘法的基本原理，给出了整体最小二乘的单位权中误差计算公式以及待估参数的近似精度评定公式。一、整体最小二乘的基本原理最小二乘法经历了百余年的发展考验，已经成为许多领域数据处理广泛应用的方法。测量数据的处理方法，通常是指按最小二乘法进行测量平差，它是测量数据处理中最基本、最广泛的应用方法，尤其是近几十年来得到了充分的发展和应用。最小二乘平差的基本思想是在最小二乘准则下进行测量数据的调整。 测量平差模型均可归结线性方程组 AX = L + 的求解问题。最小二乘准则要求残差的范数平方和极小，它主要是针对观测值中的偶然误差的。然而，实际问题中参数估计中的观测值和系数阵都可能存在误差，针对这种更复杂的情况，20 世纪 80提出了整体最小二乘法。先介绍整体最小二乘的基本思想：对于线性方程组 Ax = L，普通最小二乘的基本思想是在残差平方和极小的准则约束下求解最佳参数。这里有一个前提，系数矩阵 A 是没有误差的精确值，但是多数情况系数阵 A和观测向量 L 同时存在误差，若同时考虑二者的误差，此时，线性方程组可表示为 ( A +EA) x= L+EL其中A R ， LR m ， xR n ， rank (A ) = n , rank (A ) = n m； m 为观测值个数， n 为待估参数个数，EA为系数阵的噪声， EL 为观测噪声，误差矩阵EA EL 属于相互独立的白噪声误差。这一模型称为 EIV（ Errors-in-Variables）模型。解决这类问题的适宜方法是整体最小二乘法（ Total Least Squares, TLS）。对于线性方程组 Ax = L ，整体最小二乘问题就是在以下准则约束下minA;LRm(n-1)|A L-A L|F L R(A )寻求 A 、L，任何满足 Ax=L的 x均称为线性方程 A x = L的整体最小二乘解。E A-EL= A L-A L 为相应整体最小二乘改正数。式中， |MF|为 Frobenius范数，简称为 F 范数。整体最小二乘的求解是通过奇异值分解来实现的。将线性 Ax = L 改写为 A LxT-1T=0记增广矩阵 C = A L ，对增广矩阵 C 进行奇异值分解 C=UVT 其中 = diag(1,2,n,n+1) 12nn+10则整体最小二乘解可由增广矩阵右奇异向量的最后一列Vn+1 得到，即整体最小二乘解为 x=-1Vn+1,n+1V1,n+1Vn,n+1当 A 为列满秩时，整体最小二乘还有另一种解的形式 Xtls= ATA-n+12In-1ATL整体最小二乘的基本思想是同时考虑设计矩阵和观测向量的误差，而在许多情况下，设计矩阵的某一列或某几列是常数，如在直线拟合、曲面拟合、 GPS 非差定位等模型中都存在这种情况。因此，在这种情况下对 A 的不同列就应区别对待，与此相应的参数可分别采用最小二乘法和整体最小二乘法求解，简称为混合最小二乘将线性方程 Ax = L 表示为 A1 A2x1x2=L其中m 为观测值个数，n 为待估参数个数，n1 、 n2 分别为 A1 、A2 对应的参数个数， A1 的元素为常数。和整体最小二乘相比混合最小二乘问题就是 minA2;LRm(n2+1)|A2 L-A2 L|F准则下，寻求A2 L,任何满足 Ax =A1 x1 +A2 x2= L的x=x1T x2TT均称为混合最小二乘解。E A2 EL= A2 L-A2 L为相应的混合最小二乘改正量。混合最小二乘解的求解基本思路是首先采用 QR 分解法，或者约化的方法将系数矩阵分为常数部分和非常数部分， 后者采用整体最小二乘法求解， 后者采用普通最小二乘法求解。二 、 整体最小二乘求解附有限制条件的间接平差模型依据整体最小二乘原理的 附有限制条件的间接平差的误差方程为l+ V = ( B + VB) xCx W = 0 式中 ， l 为 n 1 的观测值向量; V 为 n 1 的观测误差向量; B 为 n m 的系数矩阵; VB 为 n m 的系数误差矩阵; x 为 m 1 的待估参数向量; C 为 c m 的限制系数矩阵; W 为 c 1 的向量， 满足 m c。 该问题依然可以用拉格朗日原理进行求解， 根据整体最小二乘原理， 建立拉格朗日目标函数如下 = VbT Vb + VTV + 2T( l Bx + V VBx) 2T ( W Cx)根据拉格朗日函数的必要条件， 分别对 V、 Vb、 、 x 求导， 经过转换可得 V+=0 VB xT = 0 l Bx + V VBx= 0 W + Cx = 0 BT + VBT CT= 0 再由上式可得l Bx= V + VBx = (1 + xTx)同时也可以得到误差改正数的计算公式为V = = ( l Bx) (1 + xTx)-1 VB= xT= ( l Bx) (1 + xTx)-1 xT根据上式，可得下列关系k = ( l Bx) T( l Bx) /(1 +xTx) = VbT Vb + VTV令 N = BTB， M = BTl，先将上式左右两边都乘以 BT 然后再加上 CT(1 + xTx) ， 联立式之前的公式有Nx M + CT(1 + xTx) = (BT CT) (1 +xTx) =VBT(1 + xTx) = xT( l Bx) =x( l Bx)T ( l Bx) /(1 + xTx) = xk最后可推出Nx M + CT1 +xTx=x k Cx W = 0根据整体最小二乘原理，所求得的最佳估计参数 x 使得 k 取得最小值。 对于参数的最优估计也可以采用迭代的方法进行。由于 (1 + xTx ) 形式复杂， 这里构建新的参数 t = (1 + xTx) ， k 一般为接近零的较小的值， 则先令 k = 0，得到初始化条件x1t(1)=NCTC0-1MW运用下面的式子来进行迭代计算k(i) = (l Bx (i) ) T(l Bx (i) ) /(1 + x (i)T x (i) )xi+1t(i+1)=NCTC0-1M+x (i)k(i)W式中，i1。求得待估参数 x 后，即可求得观测值的改正数 V 和系数矩阵的改正数 VB。三、精度评定单位权中误差计算公式为 0=+VbTVb+VTVf式中，f为自由度。无限制条件时， 单位权中方差的计算公式为0=+VbTVb+VTVn-m式中,n 为观测方程的个数; m 为待求参数个数。附有限制条件时， 单位权中方差的计算公式为0=+VbTVb+VTVn-m+c式中 ， n 为观测方程的个数; m 为待求参数个数; c 为限制条件方程个数。根据方差的定义,有 DX = E(X E(X) ) (X E(X) )T测量中的观测方程为 L = BX + d写成函数的形式为 F( L ， B， X ) = 0， 根据文献，误差传播在隐函数中,有 dX=XLdL+XBdB由误差传播率可得到 X 的中误差 X = KDKT通过对隐函数求导提取估计量对观测量的线性信息， 然后通过误差传播定律估计出待估参数的误差。整体最小二乘求解间接平差模型待估参数的近似方差计算公式如下DX 02(N klm)-1 N(N klm)-1 =(n m) -1k(N klm)-1 + k2(N klm) -1四、总结 整体最小二乘方法自 20 世纪 90 年代初正式提出以来，已在自动控制、信号处理、图像处理等许多领域取得了成功应用，作为一种新的数据处理方法是目前的一个研究热点之一。在测量数据处理中许多情况下系数矩阵和观测向量同时存在误差，如多元线性回归、GPS 高程拟合，图形图像纠正等学多情况都适于采用整体最小二乘法处理。运用整体最小二乘解算间接平差，增加了理论的严密性， 使得估计出的参数是最优的，且给出的整体最小二乘的迭代解法， 计算简便，易于编程实现。 实现了整体最小二乘求解附有限制条件的间接平差，有利于整体最小二乘法在测绘数据处理领域中的推广。 但从目前整体来看， 整体最小二乘方法在测绘领域的应用研究还比较片面，只是针对某些工程项目进行了应用研究的探索， 且整体最小二乘中定权策略、 待估参数的精确精度评定及可靠性理论等还有待作进一步研究。整体最小二乘的算法研究在理论上已经取得了较为丰富的成果，但由于整体最小二乘属于非线性估计，模型和算法的复杂性要远高于最小二乘估计，因此，在应用上收到一定的限制，如何进一步简化算法和提高算法的效率是今后整体最小二乘估计算法的重要目标。 参考文献：1丁克良,欧吉坤,陈义等.整体最小二乘法及其在测量数据处理中的应用C./中国测绘学会第九次全国会员代表大会论文集.2009:399-405. 2刘大杰，陶本藻.实用测量数据处理方法M.北京：测绘出版社，2000 3许超钤, 姚宜斌, 张豹,等. 基于整体最小二乘的参数估计新方法