雷前坤,張俊
(貴州大學(xué)礦業(yè)學(xué)院,貴州 貴陽(yáng) 550025)
在水準(zhǔn)測(cè)量中,受觀(guān)測(cè)條件的影響,觀(guān)測(cè)值中常常含有系統(tǒng)誤差。當(dāng)系統(tǒng)誤差小于或至多等于偶然誤差的量級(jí)時(shí),通常不考慮系統(tǒng)誤差對(duì)觀(guān)測(cè)值的影響[1]。目前主要采用經(jīng)典最小二乘對(duì)水準(zhǔn)網(wǎng)進(jìn)行平差處理,利用該方法進(jìn)行水準(zhǔn)網(wǎng)平差時(shí)只能顧及偶然誤差對(duì)觀(guān)測(cè)值平差結(jié)果的影響,而無(wú)法顧及系統(tǒng)誤差的影響。當(dāng)系統(tǒng)誤差較大時(shí),繼續(xù)采用經(jīng)典最小二乘進(jìn)行平差,將得不到可靠的參數(shù)估值,甚至?xí)a(chǎn)生錯(cuò)誤結(jié)論[2]。半?yún)?shù)模型是統(tǒng)計(jì)學(xué)界在20世紀(jì)80年代提出的,并被引用到測(cè)量領(lǐng)域進(jìn)行測(cè)量數(shù)據(jù)處理,在模型精化、系統(tǒng)誤差分離和減弱方面得到了廣泛應(yīng)用[3]。文獻(xiàn)[2~7]對(duì)半?yún)?shù)模型從多個(gè)方面進(jìn)行了討論,證明了半?yún)?shù)模型在某些測(cè)量數(shù)據(jù)處理方面相對(duì)于經(jīng)典最小二乘具有更高的理論精度。本文嘗試?yán)冒雲(yún)?shù)模型來(lái)對(duì)水準(zhǔn)網(wǎng)進(jìn)行平差,并將平差結(jié)果與經(jīng)典最小二乘的平差結(jié)果進(jìn)行對(duì)比分析,得出有益結(jié)論。
半?yún)?shù)平差模型可以表示為[2]:
L=BX+S+△
(1)
半?yún)?shù)平差模型中既含有參數(shù)分量,又含有非參數(shù)分量,參數(shù)分量用以表示參數(shù)與觀(guān)測(cè)值之間確定的函數(shù)關(guān)系,而非參數(shù)分量則表示系統(tǒng)誤差或模型誤差與觀(guān)測(cè)值之間不確定的函數(shù)關(guān)系。若系數(shù)矩陣B=0,則半?yún)?shù)平差模型變?yōu)榉菂?shù)模型;若模型誤差S=0,則半?yún)?shù)平差模型變?yōu)榫€(xiàn)性參數(shù)模型。因此,半?yún)?shù)平差模型是參數(shù)模型與非參數(shù)模型的混合模型[4,5]。
根據(jù)式(1)可以寫(xiě)出半?yún)?shù)平差模型的誤差方程:
(2)
根據(jù)最小二乘原理有
(3)
式中,P為正定方陣,是觀(guān)測(cè)值L的權(quán)陣;R是一個(gè)適當(dāng)給定的正定矩陣,稱(chēng)為正規(guī)化矩陣,R的選擇與具體問(wèn)題有關(guān);α是一個(gè)給定的純量因子,稱(chēng)為平滑因子(或稱(chēng)平滑參數(shù)),極小化過(guò)程中在S和V之間起到平滑的作用。
根據(jù)式(2)、式(3)通過(guò)構(gòu)造拉格朗日函數(shù)解算可得半?yún)?shù)平差模型的法方程:
(4)
將法方程式(4)寫(xiě)成下列形式
(5)
(6)
令N=BTPB,由于N可逆,于是有
(7)
將式(7)代入式(6)經(jīng)過(guò)整理得:
(8)
令M=P+αR-PBN-1BTP,則可得到
(9)
(10)
式中,H(α)=s+(I-s)B[BTP(I-s)B]-1BTP(I-s),其中s=(P+αR)-1P。
從上述半?yún)?shù)平差模型的解算過(guò)程中不難看出,解算半?yún)?shù)平差模型的關(guān)鍵和核心是如何合理的構(gòu)造正規(guī)化矩陣R和選擇平滑因子α。正規(guī)化矩陣R的構(gòu)造方法歸納起來(lái)主要有構(gòu)造矩陣法、樣條函數(shù)法和時(shí)間序列法等,平滑因子α的選擇方法主要有廣義交叉核實(shí)分值函數(shù)法、信噪比值法、L-曲線(xiàn)法等[7]。本文中正規(guī)化矩陣R采用時(shí)間序列法來(lái)構(gòu)造,平滑因子α采用L-曲線(xiàn)法中的最小距離法來(lái)求解。
(11)
則可?。?/p>
R=GTG
(12)
NN(α)=VT(α)PV(α)
從上列兩加權(quán)范數(shù)可以看出,當(dāng)平滑因子α取不同的值時(shí),可以得到不同的SN(α)和NN(α),以SN(α)為橫坐標(biāo),NN(α)為縱坐標(biāo)畫(huà)圖,可以得到許多不同的點(diǎn)對(duì)(SN(α),NN(α))。將這些點(diǎn)進(jìn)行曲線(xiàn)擬合,可以得到一條形狀像“L”的光滑曲線(xiàn),稱(chēng)之為“L-曲線(xiàn)”。用該曲線(xiàn)來(lái)求解平滑因子α的方法叫作“L-曲線(xiàn)法”。用“L-曲線(xiàn)法”來(lái)求解平滑因子有兩種方法,最短距離法和最大曲率法,即曲線(xiàn)上曲率最大的點(diǎn)所對(duì)應(yīng)的α或曲線(xiàn)上到坐標(biāo)原點(diǎn)最近的點(diǎn)所對(duì)應(yīng)的α即為所要求的平滑因子的值,本文中采用最小距離法來(lái)求解平滑因子α。用最短距離法求解平滑因子α的步驟如下:
(1)通過(guò)黃金搜索法來(lái)確定平滑因子α的取值范圍;
(2)在平滑因子α的取值范圍內(nèi),每取一個(gè)α值,進(jìn)行一次平差計(jì)算,求出一組(SN(α),NN(α)),通過(guò)取不同的α值,分別進(jìn)行一次平差計(jì)算,得到許多組(SN(α),NN(α)),將這些點(diǎn)進(jìn)行曲線(xiàn)擬合,則曲線(xiàn)上距離原點(diǎn)最近的點(diǎn)對(duì)應(yīng)的α即確定為最終的平滑因子。
圖1 水準(zhǔn)網(wǎng)
水準(zhǔn)路線(xiàn)觀(guān)測(cè)值 表1
采用經(jīng)典最小二乘和半?yún)?shù)模型對(duì)高差觀(guān)測(cè)值進(jìn)行平差。
(1)列誤差方程
(2)計(jì)算各高差觀(guān)測(cè)值的權(quán)
取10 km的觀(guān)測(cè)高差為單位權(quán)觀(guān)測(cè),即按:
定權(quán),得觀(guān)測(cè)值的權(quán)陣為:
經(jīng)典最小二乘平差的計(jì)算結(jié)果:
半?yún)?shù)平差模型的計(jì)算結(jié)果:
精度統(tǒng)計(jì) 表2
文中針對(duì)經(jīng)典最小二乘平差無(wú)法顧及水準(zhǔn)網(wǎng)中系統(tǒng)誤差這一問(wèn)題,采用半?yún)?shù)模型來(lái)對(duì)水準(zhǔn)網(wǎng)進(jìn)行平差并獲得了優(yōu)于經(jīng)典最小二乘的平差結(jié)果。半?yún)?shù)模型中含有描述系統(tǒng)誤差的非參數(shù)分量,在平差過(guò)程中顧及到了系統(tǒng)誤差對(duì)觀(guān)測(cè)值平差結(jié)果的影響,因此所得平差結(jié)果的精度要優(yōu)于經(jīng)典最小二乘平差。