陳 帥,詹本勇
(1.天津市普迅電力信息技術(shù)有限公司,天津 300000;2.華能云南滇東能源有限責(zé)任公司白龍山煤礦,云南 曲靖 655508)
最小二乘法(又稱最小平方法)是一種數(shù)學(xué)優(yōu)化方法,采用最小二乘原理進(jìn)行的最小二乘估計(jì)(LS估計(jì))具有良好的統(tǒng)計(jì)性質(zhì),但在當(dāng)今工程測(cè)量及數(shù)據(jù)處理過(guò)程中,當(dāng)自變量較多時(shí),往往會(huì)使法方程的系數(shù)陣病態(tài),使得采用LS估計(jì)求得的參數(shù)估值與真值相差比較大。再者,由于觀測(cè)值受到異常污染,使得觀測(cè)值中含有粗差時(shí),因LS估計(jì)不具有抗粗差干擾的能力,也會(huì)導(dǎo)致最終求取的參數(shù)估值與真值偏差較大,甚至完全失真。當(dāng)存在上述問(wèn)題時(shí),需要采用抗差有偏估計(jì)的方法進(jìn)行處理[1]。
現(xiàn)有的抗差有偏估計(jì)方法主要有抗差嶺估計(jì)、抗差泛嶺估計(jì)等,上述方法在對(duì)受到系數(shù)陣病態(tài)及觀測(cè)粗差影響的參數(shù)估計(jì)取得了一定的成果,但大多數(shù)在迭代過(guò)程中采用的嶺參數(shù)與最小二乘嶺估計(jì)中的相同,而實(shí)際上在迭代過(guò)程中等價(jià)權(quán)是在不斷變化的,所以在每次迭代計(jì)算中應(yīng)選擇最合適的嶺參數(shù)。嶺參數(shù)確定是影響估計(jì)結(jié)果非常重要的因素,因此,如何準(zhǔn)確的選取嶺參數(shù)是抗差有偏估計(jì)方法的關(guān)鍵。
L曲線法[2]是確定嶺參數(shù)的一種有效方法,且在最小二乘嶺估計(jì)中已成功得到應(yīng)用。本文采用中位數(shù)和標(biāo)準(zhǔn)化殘差構(gòu)建抗差嶺估計(jì)模型[3],使用L曲線法進(jìn)行嶺參數(shù)求解,通過(guò)迭代過(guò)程,最終選取精度最高、最合適的嶺參數(shù)。
設(shè)參數(shù)平差的函數(shù)模型[4]為
(1)
誤差方程為
(2)
由最小二乘原理可知,式(1)的最小二乘解為
(3)
=min,
求得
(4)
(5)
采用雙因子等價(jià)權(quán)模型[5],取
(6)
式中:wjj和wii表示自適應(yīng)降權(quán)因子。由IGGⅢ計(jì)算得
wjj=
(7)
由式(2)和式(4)推導(dǎo)Qvv的計(jì)算公式為
Qvv=NggPNgg-2Ngg+P-1,
(8)
其中,Ngg=A(ATPA+αI)-1AT.
σ0為單位權(quán)中誤差,由中位數(shù)計(jì)算,計(jì)算公式為[6]
(9)
導(dǎo)出抗差嶺估計(jì)的迭代解式:
(10)
以對(duì)數(shù)形式推導(dǎo),設(shè)
(11)
兩邊同時(shí)取對(duì)數(shù),有
(12)
L曲線上各點(diǎn)的曲率可由下式計(jì)算[7]:
(13)
設(shè)計(jì)病態(tài)矩陣為
對(duì)數(shù)據(jù)做如下方案處理:
方案一:無(wú)粗差數(shù)據(jù);方案二:將L1加入-5.5的模擬粗差;方案三:將L1、L14分別加入-5.5和-4.0的模擬粗差;方案四:將L1、L5和L7分別加入-5.5、-2.0和+4.5的模擬粗差;方案五:將L4、L12和L13分別加入+5.0、-3.0及+2.5的模擬粗差。
對(duì)上述無(wú)粗差數(shù)據(jù)及加入不同數(shù)量粗差的數(shù)據(jù)分別進(jìn)行最小二乘估計(jì)、最小二乘嶺估計(jì)以及抗差嶺估計(jì),通過(guò)三種不同方法求得的參數(shù)估計(jì)值的結(jié)果如表1所示。
表1 五種方案的參數(shù)估計(jì)值對(duì)比表
表2示出了加入粗差模擬的數(shù)據(jù)經(jīng)最小二乘估計(jì)、最小二乘嶺估計(jì)以及抗差嶺估計(jì)后的殘差及權(quán)因子為零的觀測(cè)值序號(hào)。
表2加入粗差后的殘差值及相應(yīng)觀測(cè)值序號(hào)統(tǒng)計(jì)表
方案零等價(jià)權(quán)觀測(cè)值序號(hào)LS估計(jì)殘差LS嶺估計(jì)殘差抗差嶺估計(jì)殘差模擬粗差數(shù)值 二1+4.076 2+4.139 7+5.722 7-5.5 三114+3.300 2+0.460 7+3.445 0+1.945 6+5.797 4+4.304 5-5.5 -4.0 1+3.027 8+3.192 4+5.761 3-5.5 四5+0.568 4+0.684 7+2.209 4-2.0 7-1.338 5-2.478 2-4.602 8+4.5 4-2.784 7-2.742 5-5.005 9+5.0 五12+1.131 6+2.153 9+3.147 0-3.013-1.547 9-1.700 8-2.303 3+2.5
通過(guò)以上計(jì)算結(jié)果的對(duì)比分析,可以看出:
若在抗差迭代過(guò)程中始終使用LS嶺估計(jì)所采用的嶺參數(shù)則難以達(dá)到這樣的計(jì)算效果,這是由于在迭代過(guò)程中使用L曲線法選取適合每次迭代計(jì)算的嶺參數(shù)在理論上更加合理。
4)由于抗差嶺估計(jì)模型在迭代過(guò)程中采用的是L曲線法來(lái)確定嶺參數(shù)的值,與LS嶺估計(jì)模型中的嶺參數(shù)比較是變化的。圖1是采用L曲線法求解抗差嶺估計(jì)模型最終確定的嶺參數(shù)和LS嶺估計(jì)模型中的嶺參數(shù)的示意圖。
圖1 LS嶺估計(jì)與抗差嶺估計(jì)確定的嶺參數(shù)示意圖
由于受到系數(shù)陣病態(tài)以及觀測(cè)粗差的影響,采用經(jīng)典最小二乘估計(jì)求取參數(shù)估值是極其不穩(wěn)定的,導(dǎo)致最終的計(jì)算結(jié)果與真值偏差較大,這時(shí)需要采用抗差有偏估計(jì)的方法進(jìn)行處理。文章介紹了抗差嶺估計(jì)模型,給出了關(guān)鍵公式的推導(dǎo),并說(shuō)明利用L曲線法求解抗差嶺估計(jì)模型中嶺參數(shù)的原理。最后結(jié)合實(shí)際算例作出五種方案的計(jì)算比較,通過(guò)對(duì)無(wú)粗差數(shù)據(jù)以及加入不同數(shù)量粗差的數(shù)據(jù)分別進(jìn)行LS估計(jì)、LS嶺估計(jì)、抗差嶺估計(jì),并對(duì)結(jié)果就行了詳細(xì)的分析,結(jié)果表明:LS估計(jì)不能夠克服系數(shù)陣病態(tài)的影響,抗粗差能力也非常弱;LS嶺估計(jì)雖然能夠有效的抵制系數(shù)陣病態(tài)的影響,但抗粗差能力也比較弱;與前兩種方法相比較,基于L曲線法的抗差嶺估計(jì)模型能夠有效的改善和抵制系數(shù)陣病態(tài)和觀測(cè)粗差帶來(lái)的影響。
[1]GUI Q, J. ZHANG J. Robust biased estimation and its applications in geodetic adjustments[J]. Journal of Geodesy, 1998(72):430-435.
[2]HANSEN P C. Analysis of discrete ill-posed problems by means of the L-curve [J]. SIAM Review, 1992, 34(4):561-580.
[3]王 彬,高井祥,王 堅(jiān),等.一種高精度GPS基線網(wǎng)抗差估計(jì)方法[J].全球定位系統(tǒng),2011,36(3):39-42.
[4]隋立芬.抗差嶺估計(jì)原理及其應(yīng)用[J].測(cè)繪通報(bào),1994(1):7-8.
[5]YANG Y, SONG L,XU L. Robust estimation for correlated observations based on bifactor equivalent weights [J]. Journal of Geodesy, 2002(76):353-358.
[6]宋力杰.測(cè)量平差程序設(shè)計(jì)[M].北京:國(guó)防工業(yè)出版社,2009.
[7]王振杰,歐吉坤.用L-曲線法確定嶺估計(jì)中的嶺參數(shù)[J].武漢大學(xué)學(xué)報(bào)·信息科學(xué)版,2004,29(3):235-238.