孫同賀，羅志才，姚朝龍，宛家寬

?

病態(tài)整體最小二乘的迭代正則化算法

孫同賀1, 2，羅志才1, 3，姚朝龍1，宛家寬1

(1. 武漢大學(xué)測繪學(xué)院，武漢 430079；2. 內(nèi)蒙古科技大學(xué)礦業(yè)與煤炭學(xué)院，包頭014010；3. 武漢大學(xué)地球空間環(huán)境與大地測量教育部重點實驗室，武漢 430079)

當(dāng)誤差含變量(EIV)模型的設(shè)計矩陣病態(tài)時，采用普通整體最小二乘(TLS)算法得不到穩(wěn)定的數(shù)值解。為了減弱病態(tài)性，在整體最小二乘準(zhǔn)則的基礎(chǔ)上附加解的二次范數(shù)約束，組成拉格朗日目標(biāo)函數(shù)，推導(dǎo)EIV模型的正則化整體最小二乘解(RTLS)。然后將RTLS的求解轉(zhuǎn)換為矩陣特征向量問題，設(shè)計一個迭代方案逼近RTLS解。通過L曲線法求得正則化因子來確定正常數(shù)，從而避免人為選擇正常數(shù)的隨意性。數(shù)值實例表明，提出的迭代正則化算法是有效可行的。

EIV模型；病態(tài)問題；正則化整體最小二乘；L曲線法；正常數(shù)

在信號處理、計算機(jī)視覺、圖像處理、通信工程以及大地測量與攝影測量等相關(guān)領(lǐng)域中，待求參數(shù)一般通過建立觀測模型來求解。當(dāng)觀測系統(tǒng)的設(shè)計矩陣和觀測量都含有誤差時，便形成了誤差含變量(Error-in-variables, EIV)模型。整體最小二乘(Total least squares, TLS)方法是解決EIV模型的基本途徑。TLS的算法主要包括奇異值分解類算法和拉格朗日迭代算法[1?2]，并擴(kuò)展到加權(quán)TLS[3]以及附有等式和不等式約束的情況[4?7]。在TLS解的統(tǒng)計性質(zhì)方面，主要研究了TLS解與LS解的關(guān)系[8]，系數(shù)陣的隨機(jī)誤差對加權(quán)LS解[9]的影響。對TLS在測繪領(lǐng)域的典型應(yīng)用也進(jìn)行了相關(guān)研究[10]。

當(dāng)設(shè)計陣的列向量間呈近似線性相關(guān)時，觀測量的微小變化會引起最小二乘解的劇烈振蕩，測量上一般采用截斷奇異值或者Tikhonov正則化方法獲得穩(wěn)定的解[11?12]。根據(jù)奇異值分解，TLS解相較于LS解是一個降正則化過程，TLS解的條件數(shù)總是大于LS解的條件數(shù)[13]。GOLUB[14]提出了TLS的正則化解，給出了取不同的約束參數(shù)時，正則化整體最小二乘(Regularized total least squares, RTLS)解與TLS解的關(guān)系，并給出了一種直接算法。FIERRO等[15]給出了解RTLS的截斷奇異值分解算法。SCHAFFRIN等[16]基于Gauss-Helmert模型的混合近似解理論，提出了一種新的迭代計算方案，并運用到考古學(xué)中來確定賽馬場中賽道的原點和半徑。BECK等[17]將RTLS問題轉(zhuǎn)換成一個閉區(qū)間內(nèi)的單變量函數(shù)的最小化問題。袁振超等[18]推導(dǎo)了等權(quán)條件下病態(tài)總體最小二乘的正則化解法。王樂詳?shù)萚19]采用嶺估計方法解決加權(quán)總體最小二乘平差的病態(tài)性問題，給出了確定正則化參數(shù)的嶺跡法、廣義交叉核實法和L曲線法。又將參數(shù)作為虛擬觀測值，運用廣義最小二乘準(zhǔn)則，推導(dǎo)了病態(tài)總體最小二乘問題的虛擬觀測解法[20]。葛旭明等[21]提出了解RTLS的廣義正則化方法。然而，在應(yīng)用Tikhonov正則化方法解決病態(tài)問題時，在該方法的約束條件中需人為引入一個正常數(shù)，但是正常數(shù)的引入并沒有絕對的標(biāo)準(zhǔn)，也沒有一個客觀的參考值，因而其值的引入存在較強(qiáng)的主觀隨意性。針對此問題，引入了誤差約束的TLS正則化方法，在原有方程的基礎(chǔ)上加入2個不等式作為約束條件[22]。GUO等[23]直接給出了RTLS迭代算法的計算公式，而沒有給出詳細(xì)的推導(dǎo)過程。這些方法均是根據(jù)觀測量的精度及設(shè)計矩陣的形式人為選擇設(shè)計矩陣誤差上限值和正常數(shù)，因此具有隨機(jī)性。

本文作者基于文獻(xiàn)[23]中提出的RTLS的思想，詳細(xì)推導(dǎo)了病態(tài)EIV模型的正則化解及正常數(shù)與正則化因子的關(guān)系，然后將RTLS的求解轉(zhuǎn)換為求矩陣特征向量問題，設(shè)計了一個迭代方案逼近RTLS解。通過L曲線法求得正則化因子來確定正常數(shù)，從而避免了人為選擇正常數(shù)的隨意性。最后用兩個算例說明了本研究所提方法的可行性和有效性。

1 病態(tài)EIV模型的正則化方法

近年來，誤差含變量(EIV)模型在測量實際中廣泛出現(xiàn)。EIV模型采用奇異值分解算法求解時，發(fā)現(xiàn)它的整體最小二乘解相較于普通最小二乘解來說，法方程主對角線元素分別減去增廣矩陣的最小奇異值，從而是一個降正則化(De-regularization)過程。王樂洋[13]證明了EIV模型的條件數(shù)要大于相應(yīng)的線性G-M模型。因此，有必要研究EIV模型病態(tài)時的數(shù)據(jù)處理問題。

病態(tài)EIV模型可以表示為

在測量平差理論中，通過附加一些約束信息，可以削弱系統(tǒng)的病態(tài)性，保證估值的穩(wěn)定性。對于模型(1)所示的病態(tài)EIV模型，在觀測量和系數(shù)陣誤差殘差最小的基礎(chǔ)上，附加參數(shù)的二次范數(shù)約束，從而避免求解參數(shù)時因觀測量的微小變動而造成解的巨大波動。其平差準(zhǔn)則可表示成

式中：帽子符號表示估計值，波浪線符號表示預(yù)測值。根據(jù)式(4)和(5)可得

將式(9)、(10)代入式(8)，即

也可以表示成

整理上式得

將式(9)代入式(6)得

根據(jù)最優(yōu)化理論的Kuhn-Tucker條件和有效約束的概念，可以得到：若＞0，則，稱為有效約束；若，則＜為無效約束。由此，有

另外，由式(15)可以得到

令

則式(18)可以簡寫成：

由式(6)易知：

將式(25)代入式(23)，有

將式(17)代入式(27)，可得

聯(lián)立式(21)和(29)，可得如下的矩陣形式：

式中：

式(21)的求解也可轉(zhuǎn)換為求矩陣特征值特征向量對的問題，采用反冪法(Shifted inverse power)方法迭代求解，即求式(32)的特征向量：

定義第次迭代的殘差向量為

下面給出本研究所提方法的迭代步驟：

1) 給定約束陣，按照上述方法確定正則化因子，根據(jù)式(21)得到正則化整體最小二乘解；

2 算例分析

算例1：算例摘自文獻(xiàn)[24]，模型的設(shè)計矩陣和觀測值分別為

算例2：采用一組病態(tài)空間測邊網(wǎng)的算例，同時運用不同方法進(jìn)行解算，進(jìn)而探討本方法在測量數(shù)據(jù)處理中的可行性。

表1 模型參數(shù)X的計算結(jié)果

模擬的空間測邊網(wǎng)算例，1、2、…9為9個已知點，其坐標(biāo)具體數(shù)據(jù)可參考文獻(xiàn)[22]。9個已知點到兩個未知點10、11的觀測距離也給出，假設(shè)兩個未知點真值坐標(biāo)為(0，0，0)和(7，10，?5)，它們之間的觀測距離為10,11=13.1078?m，各距離為等精度觀測，測距中誤差為±0.001?m。要求根據(jù)19個邊長觀測值確定兩個未知點的坐標(biāo)。

經(jīng)計算，該測邊網(wǎng)觀測方程的系數(shù)矩陣嚴(yán)重病態(tài)，法方程的條件數(shù)為4.5886×103。在計算中，兩個未知點的坐標(biāo)近似值分別取為(0.01，?0.01，0.02)和(7.01，9.99，?5.01)。由文獻(xiàn)[22]可知，設(shè)計矩陣誤差上限可根據(jù)觀測量的精度確定，并在此基礎(chǔ)上略有加大。表2給出了各種算法的參數(shù)估計結(jié)果及與真值的差值范數(shù)，圖2給出了本方法得到的L曲線圖，橫縱坐標(biāo)軸同圖1。

通過公式推導(dǎo)及計算分析，可以得出如下結(jié)論：

2) 由表1和表2的統(tǒng)計結(jié)果可知，當(dāng)系數(shù)矩陣病態(tài)時，整體最小二乘解與真值的差值范數(shù)分別為6.7190和22.2787，較最小二乘解與真值的差值范數(shù)1.3088和10.4867要大，說明整體最小二乘解受病態(tài)的影響更為嚴(yán)重；

3) 采用文獻(xiàn)[22]方法，設(shè)計矩陣誤差上限取不同值，得到的解與真值的差值范數(shù)分別為(0.8253，0.8257，0.8262)和(0.0925，0.0932，0.0938)；采用本研究方法得到的解與真值的差值范數(shù)分別為0.8234和0.0932。當(dāng)選擇較為合適的正常數(shù)時，文獻(xiàn)[22]中的方法和本研究方法解算結(jié)果相當(dāng)，但如果正常數(shù)選擇不合理，結(jié)果將有較大的偏差；

表2 模型參數(shù)X的計算結(jié)果

圖1 L曲線(算例1)

圖2 L曲線(算例2)

4) 在正則化最小二乘算法和正則化整體最小二乘算法中，通常假設(shè)正常數(shù)足夠小，計算過程中忽略了的確切值。本研究在迭代過程中考慮了的確切值。從解與真值的差值范數(shù)看，正則化整體最小二乘算法較之本研究方法得到的結(jié)果分別有0.0056和0.0038的偏差，這主要是因為正常數(shù)帶來的影響；

5) 算例2中給定的初值并不理想，但利用本研究方法得到的解算結(jié)果要好于RLS和RTLS，由此說明本研究方法對初值精度的要求不高，推廣了該方法的使用。

3 結(jié)論

1) 當(dāng)EIV模型設(shè)計矩陣病態(tài)時，普通的整體最小二乘算法得不到穩(wěn)定的解。在觀測量和系數(shù)陣誤差殘差最小的基礎(chǔ)上，附加參數(shù)的二次范數(shù)約束，推導(dǎo)了病態(tài)整體最小二乘的迭代正則化算法，以及正常數(shù)與正則化因子的等式關(guān)系，采用shifted inverse power算法，將RTLS的求解轉(zhuǎn)換為求矩陣的特征向量問題。

2) 相比已有的正則化整體最小二乘算法，本研究提出的方法有兩個優(yōu)勢：首先，在迭代過程中考慮了的確切值；其次，通過L曲線法求得正則化因子來更合理地確定正常數(shù)，從而避免了不恰當(dāng)?shù)剡x取正常數(shù)，可能導(dǎo)致參數(shù)解存在較大偏差的情況。

REFERENCES：

[1] van HUFFEL S, VANDEWALLE J. The total least squares problem computational aspects and analysis[M]. Philadelphia: Society for Industrial and Applied Mathematics, 1991.

[2] GOLUB G H, van LOAN C F. An analysis of the total least squares problem[J]. SIAM Journal on numerical analysis, 1980, 17(6): 883?893.

[3] FANG Xing. Weighted total least squares solutions for applications in geodesy[D]. Hannover: Doctor’s Dissertation in Leibniz University, 2011.

[4] SCHAFFRIN B. A note on constrained total least-squares estimation[J]. Linear algebra and its applications, 2006, 417: 245?258.

[5] SCHAFFRIN B, FELUS Y A. An algorithmic approach to the total least-squares problem with linear and quadratic constraints[J]. Studia Geophysica et Geodaetica, 2009, 53: 1?16.

[6] ZHANG Song-lin, TONG Xiao-hua. A solution to EIV model with inequality constraints and its geodetic applications[J]. Journal of geodesy, 2013, 87: 23?28.

[7] FANG Xing. On non-combinatorial weighted total least squares with inequality constraints[J]. Journal of geodesy, 2014, 88: 805?816.

[8] van HUFFEL S, VANDEWALLE J. Algebraic connections between the least squares and total least squares problems[J]. Numerische Mathematik, 1989, 55: 431?449.

[9] XU Pei-liang, LIU Jing-nan. Effects of errors-in-variables on weighted least squares estimation[J]. Journal of Geodesy, 2014, 88: 705?716.

[10] 周擁軍, 鄧才華. 線性EIV模型的TLS估計及其典型應(yīng)用[J]. 中國有色金屬學(xué)報, 2012, 22(3): 948?953. ZHOU Yong-jun, DENG Cai-hua. Extended total least squares problems for linear errors-in-variables model and its typical applications[J]. The Chinese Journal of Nonferrous Metals, 2012, 22(3): 948?953.

[11] XU Pei-liang. Truncated SVD methods for discrete linear ill-posed problems[J]. Geophysical Journal International, 1998, 135: 505?514.

[12] TIKHONOV A N. Solution of incorrectly formulated problems and the regularization method[J]. Soviet Math Dokl, 1963, 4: 1035?1038.

[13] 王樂洋. 總體最小二乘解性質(zhì)研究[J]. 大地測量與地球動力學(xué), 2012, 32(5): 48?52. WANG Le-yang. Research on properties of total least squares estimation[J]. Journal of Geodesy and Geodynamics, 2012, 32(5): 48?52.

[14] GOLUB G H, HANSEN P C, O’LEARY D P. Tikhonov regularization and total least squares[J]. SIAM Journal on Matrix Analysis and Applications, 1999, 21(1): 185?194.

[15] FIERRO R D, GOLUB G H, HANSEN P C. Regularization by truncated total least squares[J]. SIAM Journal on Scientific Computing, 1997, 18(1): 1223?1241.

[16] SCHAFFRIN B, SNOW K. Total least-squares regularization of Tykhonov type and an ancient racetrack in Corinth[J]. Linear Algebra and Its Applications, 2010, 432: 2061?2076.

[17] BECK A, BEN-TAL A. On the solution of the Tikhonov regularization of the total least squares problem[J]. SIAM Journal on Optimization, 2006, 17(1): 98?118.

[18] 袁振超, 沈云中, 周澤波. 病態(tài)總體最小二乘模型的正則化算法[J]. 大地測量與地球動力學(xué), 2009, 2(2): 131?134. YUAN Zhen-chao, SHEN Yun-zhong, ZHOU Ze-bo. Regularized total least squares solution to ill-posed error-in- variable model[J]. Journal of Geodesy and Geodynamics, 2009, 29(2): 131?134.

[19] 王樂洋, 許才軍, 魯鐵定. 病態(tài)加權(quán)總體最小二乘平差的嶺估計方法[J]. 武漢大學(xué)學(xué)報信息科學(xué)版, 2010, 35(11): 1346?1350. WANG Le-yang, XU Cai-jun, LU Tie-ding. Ridge estimation method in ill-posed weighted total least squares adjustment[J]. Geomatics and Information Science of Wuhan University, 2010, 35(11): 1346?1350.

[20] 王樂洋, 于冬冬. 病態(tài)總體最小二乘問題的虛擬觀測解法[J]. 測繪學(xué)報, 2014, 43(6): 575?581. WANG Le-yang, YU Dong-dong. Virtual observation method to ill-posed total least squares problem[J]. Acta Geodaetica et Cartographica Sinica, 2014, 43(6): 575?581.

[21] 葛旭明, 武吉倉. 病態(tài)總體最小二乘問題的廣義正則化[J]. 測繪學(xué)報, 2012, 41(3): 372?377. GE Xu-ming, WU Ji-cang. Generalized regularization to ill-posed total least squares problem[J]. Acta Geodaetica et Cartographica Sinica, 2012, 41(3): 372?377.

[22] 葛旭明, 武吉倉. 誤差限的病態(tài)總體最小二乘解算[J]. 測繪學(xué)報, 2013, 42(2): 196?202. GE Xu-ming, WU Ji-cang. A regularization method to ill-posed total least squares with error limits[J]. Acta Geodaetica et Cartographica Sinica, 2013, 42(2): 192?202.

[23] GUO Hong-bin, RENAUT R A. A regularized total least squares algorithm[C]// Total Least Squares and Errors-in-Variables Modeling. Copenhagen: Kluwer Academic Publishers, 2002: 57?66.

[24] 魯鐵定, 寧津生.總體最小二乘平差理論及其應(yīng)用[M]. 北京: 中國科學(xué)技術(shù)出版社, 2011. LU Tie-ding, Ning Jin-sheng. Total least squares adjustment theory and its applications[M]. Beijing: China Science and Technology Press, 2011.

(編輯 王 超)

Iterative and regularized algorithm to ill-posed total least squares

SUN Tong-he1, 2, LUO Zhi-cai1, 3, YAO Chao-long1, YUAN Jia-kuan1

(1. School of Geodesy and Geomatics, Wuhan University, Wuhan 430079, China; 2. Mining and Coal Institute, Inner Mongolia University of Science and Technology, Baotou 014010, China; 3. Key Laboratory of Geospace Environment and Geodesy, Ministry of Education, Wuhan University, Wuhan 430079, China)

When the design matrix of errors-in-variables (EIV) model was ill-conditioned, the ordinary total least squares (TLS) solution was unstable. In order to weaken the ill-conditioning, an Euclid norm constraint of the solution was added to the TLS minimization rule. Then, the Lagrange objective function was formed and the regularized total least squares (RTLS) solution was deduced. Afterwards, the RTLS was transformed to a problem of looking for a matrix’s eigenvector. An iterative program was designed to approximate the solution. The L-curve method was used to choose the regularization factor to determine the positive constant, which can avoid the subjective decision. The simulations show the efficiency and feasibility of the algorithm.

EIV model; ill-posed problem; regularized total least squares; L-curve method; positive constant

Project(41474006) supported by the National Natural Science Foundation of China

2015-12-04; Accepted date:2016-04-10

SUN Tong-he; Tel: +86-472-5951692; E-mail: suntonghe202@163.com

1004-0609(2016)-10-2174-07

P207

A

國家自然科學(xué)基金資助項目(41474006)

2015-12-04；

2016-04-10

孫同賀，講師，博士；電話：0472-5951692；E-mail: suntonghe202@163.com