病態(tài)不確定性平差模型的嶺估計(jì)算法

魯鐵定，吳光明，周世健

1. 東華理工大學(xué)測(cè)繪工程學(xué)院,江西 南昌 330013； 2. 流域生態(tài)與地理環(huán)境監(jiān)測(cè)國(guó)家測(cè)繪地理信息局重點(diǎn)實(shí)驗(yàn)室,江西 南昌 330013； 3. 江西省數(shù)字國(guó)土重點(diǎn)實(shí)驗(yàn)室,江西 南昌 330013； 4. 南昌航空大學(xué)，江西 南昌 330063

測(cè)繪數(shù)據(jù)常存在復(fù)雜的不確定性，不確定性是不精確性、模糊性、不明確性等概念的總稱。不確定性與誤差意義接近，但它涵蓋的內(nèi)容更廣，如屬性不確定性、模糊不確定性等[1]。不確定度是不確定性的度量，是不確定性的一種指標(biāo)[2]。不確定度與測(cè)量界的精度度量方式幾乎完全一致，不確定度可以用方差、均方差、誤差區(qū)間、誤差橢圓、誤差橢球表示[1]。測(cè)量數(shù)據(jù)的不確定性不再是一個(gè)具體的數(shù)值，有時(shí)僅知道它們各自在一定的區(qū)間內(nèi)變動(dòng)，有時(shí)僅是一個(gè)模糊數(shù)，沿用隨機(jī)誤差的分布限定會(huì)給測(cè)量平差數(shù)據(jù)處理帶來(lái)困難[3]。在測(cè)繪數(shù)據(jù)處理領(lǐng)域，應(yīng)用不確定度理論，研究不確定度評(píng)定方法，研究有效降低不確定性的影響等成為研究熱點(diǎn)[4-7]。文獻(xiàn)[8]在不確定度理論下構(gòu)建海底數(shù)字高程模型；文獻(xiàn)[9]從幾何的角度分析了最大不確定度對(duì)數(shù)據(jù)的影響，并提出不確定性平差模型解算方法；文獻(xiàn)[10—11]應(yīng)用廣義拉格朗日法構(gòu)造目標(biāo)函數(shù)，更為具體地分析了不確定度的各種情況，并推導(dǎo)了相應(yīng)的求解公式。在傳統(tǒng)最小二乘方法中，由于系數(shù)矩陣的病態(tài)性進(jìn)而會(huì)導(dǎo)致因觀測(cè)量微小波動(dòng)而造成解算結(jié)果產(chǎn)生巨大波動(dòng)，參數(shù)估值已嚴(yán)重失真，不是最優(yōu)解[12]，很多學(xué)者對(duì)此展開(kāi)了研究，諸如提出了Tikhonov正則化法[13-14]、嶺估計(jì)法[15-16]等解決病態(tài)問(wèn)題的方法。對(duì)于病態(tài)情形下的EIV模型應(yīng)用總體最小二乘進(jìn)行解算，解算受系數(shù)矩陣誤差和觀測(cè)值誤差的影響將更加嚴(yán)重[17]，學(xué)者也提出了病態(tài)總體最小二乘的正則化法[18-20]、病態(tài)加權(quán)總體最小二乘的嶺估計(jì)解法[21]，得到參數(shù)估值更加穩(wěn)定。

當(dāng)不確定性平差模型出現(xiàn)病態(tài)，如何解決病態(tài)問(wèn)題的相關(guān)研究較少?？紤]到不確定性平差模型中系數(shù)矩陣可能出現(xiàn)接近于0的奇異值，如何處理此模型的病態(tài)問(wèn)題是本文研究重點(diǎn)?；趯?duì)病態(tài)G-M模型、EIV模型下的嶺估計(jì)法和不確定性平差模型的平差準(zhǔn)則分析，以及文獻(xiàn)[10]對(duì)不確定性平差模型迭代解算算法，本文建立了相應(yīng)的病態(tài)不確定性平差模型及平差準(zhǔn)則，該準(zhǔn)則加入了穩(wěn)定泛函，對(duì)法矩陣的奇異值進(jìn)行修正，將模型由嚴(yán)重病態(tài)變成病態(tài)性較弱或無(wú)病態(tài)，使得法矩陣求逆變得穩(wěn)定，并推導(dǎo)了病態(tài)不確定性平差嶺估計(jì)法的迭代解算公式，通過(guò)算例對(duì)算法進(jìn)行驗(yàn)證及討論。

1 不確定性平差模型及平差準(zhǔn)則

不確定平差模型為

(1)

文獻(xiàn)[10]總結(jié)出不確定性平差模型的不同情形并對(duì)其加以討論，從而建立的不確定性最小二乘(uncertainty least squares，ULS)平差準(zhǔn)則

(2)

構(gòu)造廣義拉格朗日目標(biāo)函數(shù)

(3)

式中,λ、μ、u是拉格朗日乘子，ΔA、ΔL的不確定度用代數(shù)式表示，整理得到法方程

(4)

(5)

當(dāng)不確定性平差模型病態(tài)時(shí)，法方程系數(shù)矩陣ATA求逆會(huì)變得極不穩(wěn)定，在以均方誤差作為估值參考依據(jù)，當(dāng)法矩陣存在特征值接近于零時(shí)，方差將會(huì)非常大，導(dǎo)致求解出的參數(shù)估值不可靠。

2 病態(tài)不確定性平差模型的嶺估計(jì)

根據(jù)不確定性平差模型，當(dāng)系數(shù)矩陣A接近0的奇異值，模型的計(jì)算公式(5)中，法矩陣求逆將極不穩(wěn)定，導(dǎo)致求解參數(shù)估值可靠性降低。為了降低病態(tài)對(duì)不確定性平差模型平差結(jié)果的影響，在ULS平差準(zhǔn)則式(2)中加入穩(wěn)定泛函，即

(6)

式中，α是嶺參數(shù)。構(gòu)建廣義拉格朗日目標(biāo)函數(shù)

(7)

求一階偏導(dǎo)得

(8a)

(8b)

(8c)

(8d)

(8e)

(8f)

(9)

根據(jù)式(9)可以求得

(10)

法方程式

(11)

從式(8d)可知

(12)

將式(12)代入式(11)得

(13)

再把式(10)代入式(13)得

(14)

(15)

可以看出式(15)為病態(tài)總體最小二乘嶺估計(jì)解算方法，與文獻(xiàn)[17]相一致，因此，病態(tài)總體最小二乘嶺估計(jì)是本文病態(tài)不確定性平差模型嶺估計(jì)的特例。

3 病態(tài)不確定性平差模型嶺估計(jì)解算步驟

解算步驟如下：

(1) 在計(jì)算時(shí)先給出系數(shù)矩陣A、觀測(cè)向量L、不確定度φ、β。

(3) 在一范圍[a，b]內(nèi)按一定步長(zhǎng)Δd選擇α。

(4) 確定u、μ，確定公式參考文獻(xiàn)[10]的方法，計(jì)算公式為

(16a)

(16b)

通過(guò)式(16)求出u、μ，要保證迭代值u(i)、μ(i)為非負(fù)數(shù)，u(i)、μ(i)求解參考文獻(xiàn)[10]中的解算步驟部分：

(a) 當(dāng)u、μ均大于0時(shí)

u(i)=u、μ(i)=μ；

(b) 當(dāng)u、μ均小于等于0時(shí)

u(i)、μ(i)均為0；

(c) 當(dāng)u>0、μ≤0時(shí)

(d) 當(dāng)u≤0、μ>0時(shí)

(5) 根據(jù)式(14)進(jìn)行迭代計(jì)算

則

(17)

4 算例及分析

4.1 算例1

采用文獻(xiàn)[17]中的模擬病態(tài)問(wèn)題算例，法矩陣條件數(shù)為2.083 8×104，嚴(yán)重病態(tài)。未知參數(shù)有5個(gè)，X=[x1x2x3x4x5]T，真值為X=[11111]T。為驗(yàn)證病態(tài)不確定最小二乘嶺估計(jì)法，分別用LS、TLS、ULS、嶺估計(jì)LS(R-LS)、嶺估計(jì)TLS(R-TLS)和嶺估計(jì)ULS(R-ULS)進(jìn)行解算。由于該數(shù)據(jù)的不確定度未知，計(jì)算時(shí)φ、β從不同區(qū)間((0，5]、(0，2.5]、(0，1]、(0，0.5]、(0，0.1])隨機(jī)取值，分別重復(fù)計(jì)算1000次，比較不同不確定度對(duì)結(jié)果的影響，結(jié)果如表1所示。

表1 不同方法參數(shù)估計(jì)結(jié)果

圖1 各方法L-曲線圖Fig.1 The L-curve of each methods

圖2 不同區(qū)間結(jié)果Fig.2 Different interval results

圖3 不同區(qū)間結(jié)果Fig.3 Different interval results

4.2 算例2

采用文獻(xiàn)[23]空間測(cè)邊網(wǎng)算例。P1、P2、…、P10為10個(gè)已知點(diǎn)，其坐標(biāo)具體數(shù)據(jù)略去。10個(gè)已知點(diǎn)到3個(gè)未知點(diǎn)P11、P12、P13(假設(shè)模擬坐標(biāo)真值分別為(0，0，0)、(68，-26，9)和(14，41，-11))的距離，以及3個(gè)未知點(diǎn)間的距離假定已通過(guò)測(cè)量得到。設(shè)各距離為等精度觀測(cè)，中誤差為±0.01 m。根據(jù)33個(gè)距離觀測(cè)值確定3個(gè)未知點(diǎn)坐標(biāo)。計(jì)算中3個(gè)未知點(diǎn)坐標(biāo)近似值分別取(0.03 m，-0.025 m，0.01 m)、(68.03 m，-25.97 m，8.98 m)和(14.04 m，40.97 m，-11.04 m)。該測(cè)邊網(wǎng)所建立觀測(cè)方程的系數(shù)陣A嚴(yán)重病態(tài)，法矩陣條件數(shù)為89 543。為驗(yàn)證病態(tài)不確定最小二乘嶺估計(jì)法，分別用LS、TLS、ULS、R-LS、R-TLS和R-ULS進(jìn)行解算，由于該數(shù)據(jù)的不確定度未知，計(jì)算時(shí)φ、β也從不同區(qū)間((0，5]、(0，2.5]、(0，1]、(0，0.5]、(0，0.1])隨機(jī)取值重復(fù)計(jì)算1000次，比較不同不確定度對(duì)結(jié)果的影響，結(jié)果見(jiàn)表2。

圖4 φ、β(0，0.1]結(jié)果Fig.4 φ、β(0,0.1) results

表2 不同方法參數(shù)估計(jì)結(jié)果

應(yīng)用L-曲線法確定嶺參數(shù)，嶺參數(shù)變化如圖5所示，確定嶺參數(shù)后，φ、β隨機(jī)取值計(jì)算中發(fā)現(xiàn)R-ULS最低的差值范數(shù)為0.112 2，但出現(xiàn)的次數(shù)太少，而出現(xiàn)最多的差值范數(shù)是0.116 1，所以把差值范數(shù)為0.116 1時(shí)的參數(shù)估值作為最優(yōu)估值。為驗(yàn)證不同不確定度時(shí)R-ULS的可行性，本文在給定5個(gè)取值區(qū)間再次分別重復(fù)計(jì)算1000次，并將所得差值范數(shù)減去0.112 2并作出相應(yīng)的散點(diǎn)圖(橫坐標(biāo)為計(jì)算次數(shù)，縱坐標(biāo)為差值范數(shù)減0.112 2)，縱坐標(biāo)為0時(shí)差值范數(shù)則為0.112 2，結(jié)果如圖6、圖7和圖8所示。

圖5 各方法L-曲線圖Fig.5 The L-curve of each method

圖6 不同區(qū)間結(jié)果Fig.6 Different interval results

圖7 不同區(qū)間結(jié)果Fig.7 Different interval results

圖8 φ、β(0，0.1]結(jié)果Fig.8 Results of φ、β(0,0.1)

4.3 算例分析

病態(tài)情況是法矩陣出現(xiàn)幾個(gè)接近于零的特征值，從而法矩陣求逆將變得極不穩(wěn)定，導(dǎo)致參數(shù)估值與真值偏差較大。在兩算例中，嶺參數(shù)用L-曲線法確定。L-曲線圖如圖1、圖5所示，嶺參數(shù)則為圖中的拐點(diǎn)處曲率最大的點(diǎn)對(duì)應(yīng)的值。兩算例LS估計(jì)差值范數(shù)均比TLS、ULS估計(jì)低，說(shuō)明TLS和ULS均比LS受病態(tài)情況影響更為嚴(yán)重，對(duì)病態(tài)性更敏感。

嶺估計(jì)削弱或消除法矩陣的病態(tài)性，使得法矩陣求逆變得穩(wěn)定，有效地抑制住病態(tài)帶來(lái)的影響。將嶺估計(jì)運(yùn)用于病態(tài)情況下的LS估計(jì)和TLS估計(jì)，算例1得到參數(shù)估值的差值范數(shù)分別是0.854 7、0.846 8，算例2得到結(jié)果分別是0.711 6、0.204 8，均優(yōu)于相應(yīng)的LS、TLS、ULS方法得到結(jié)果。病態(tài)不確定性平差模型應(yīng)用嶺估計(jì)法可以有效地提高了參數(shù)估計(jì)解算結(jié)果的穩(wěn)定性。兩算例中R-ULS最優(yōu)的差值范數(shù)分別為0.838 9、0.116 1，較R-LS、R-TLS解算結(jié)果更優(yōu)，說(shuō)明嶺估計(jì)可用于病態(tài)不確定性平差模型解算。

由于不確定度φ、β是未知的，在試驗(yàn)過(guò)程中不確定度是上述5個(gè)取值區(qū)間的隨機(jī)數(shù)，分別重復(fù)計(jì)算1000次，發(fā)現(xiàn)R-ULS存在最優(yōu)解。根據(jù)兩算例的散點(diǎn)圖，在(0，5]區(qū)間中，算例1的0.838 9出現(xiàn)次數(shù)最多且最低(縱坐標(biāo)為參數(shù)差值范數(shù)減0.838 9)；算例2的0.116 1出現(xiàn)次數(shù)最多(縱坐標(biāo)為參數(shù)差值范數(shù)減0.112 2)，存在最低的0.112 2，但出現(xiàn)次數(shù)太少。在參數(shù)真值已知情況下，可根據(jù)差值范數(shù)大小確定最優(yōu)估值。隨著區(qū)間進(jìn)一步縮小，最優(yōu)值出現(xiàn)概率在降低，不同的不確定度得出不同的結(jié)果，可從(0，1]、(0，0.5]、(0，0.1]3個(gè)取值區(qū)間相應(yīng)圖形看出最優(yōu)值出現(xiàn)概率在降低。此外，從兩算例的散點(diǎn)圖發(fā)現(xiàn)不同不確定度參數(shù)估值的差值范數(shù)是在一定區(qū)間內(nèi)變化，算例1是在區(qū)間范圍[0.838 9，0.938 9]，算例2是區(qū)間范圍[0.112 2，0.182 2]。說(shuō)明病態(tài)不確定性平差模型嶺估計(jì)解受不確定性有界的約束，在不確定度較小時(shí)，不確定性約束更為明顯。同時(shí)進(jìn)一步說(shuō)明隨著不確定度增大，病態(tài)不確定性平差模型嶺估計(jì)解對(duì)不確定度的敏感程度在降低。

5 結(jié) 論

本文分析了當(dāng)不確定性平差模型出現(xiàn)病態(tài)時(shí)，ULS受病態(tài)情形的影響，提出基于嶺估計(jì)，同時(shí)顧及系數(shù)矩陣和觀測(cè)向量出現(xiàn)不確定性誤差，推導(dǎo)了病態(tài)不確定性平差模型的嶺估計(jì)平差準(zhǔn)則，推導(dǎo)了迭代算法，以提高解的穩(wěn)定性。通過(guò)算例結(jié)果表明，R-ULS能夠有效地抑制病態(tài)的影響和降低差值范數(shù)，說(shuō)明了提出的病態(tài)不確定性平差模型嶺估計(jì)法具有一定的有效性。同時(shí)從算例得出，病態(tài)不確定性模型的嶺估計(jì)解受不確定度影響，影響程度隨著不確定度增大而降低。