張 明，陶延林，山中雪

(青海省第一地質(zhì)勘查院，青海 海東 810699)

為了獲得精確、實(shí)時(shí)的礦區(qū)空間地理信息，測(cè)量人員需要采用多種方法對(duì)測(cè)量數(shù)據(jù)函數(shù)模型進(jìn)行研究，隨著測(cè)量技術(shù)及數(shù)據(jù)處理理論的不斷發(fā)展，測(cè)量數(shù)據(jù)的處理精度也在不斷提高，目前，采用總體最小二乘算法來(lái)解決變量中的誤差以及模型參數(shù)估值問(wèn)題受到科技工作者的廣泛關(guān)注[1]。

國(guó)內(nèi)外學(xué)者對(duì)最小二乘算法在礦山測(cè)量數(shù)據(jù)中的應(yīng)用進(jìn)行了大量的研究，王樂(lè)洋采用總體最小二乘算法對(duì)大地測(cè)量數(shù)據(jù)進(jìn)行了反演理論及應(yīng)用研究[2]；汪奇生等研究了總體最小二乘算法的線(xiàn)性回歸迭代算法，針對(duì)線(xiàn)性回歸模型的特點(diǎn)，以一元線(xiàn)性回歸為例，運(yùn)用拉格朗日原理推導(dǎo)了線(xiàn)性回歸迭代算法，該算法推導(dǎo)過(guò)程簡(jiǎn)單且易于實(shí)現(xiàn)，通過(guò)算例驗(yàn)證了該算法的正確性[3]。

1 總體最小二乘算法函數(shù)模型求解

本文主要運(yùn)用基于拉格朗日函數(shù)與奇異值分解函數(shù)兩種算法對(duì)總體最小二乘算法函數(shù)模型參數(shù)進(jìn)行求解，在無(wú)加權(quán)條件下，采用奇異值分解函數(shù)算法對(duì)等權(quán)條件下的總體最小二乘算法對(duì)測(cè)量數(shù)據(jù)進(jìn)行求解；在加權(quán)條件下，則采用拉格朗日函數(shù)算法更適合對(duì)測(cè)量數(shù)據(jù)進(jìn)行求解，并對(duì)測(cè)量數(shù)據(jù)進(jìn)行迭代處理[4]。

1.1 基于奇異值分解算法

采用矩陣奇異值算法對(duì)解決最小二乘問(wèn)題以及最優(yōu)化問(wèn)題起到重要作用，奇異值分解定義如下[5]：

(1)

BCrm×n

(2)

式中:σi為矩陣B的奇異值，r為矩陣B的秩，BTB的特征值為λ1≥λ2≥…≥λr=…=λn。

根據(jù)推導(dǎo)可得出基于奇異值分解的參數(shù)向量總體最小二乘解為：

(3)

(4)

1.2 基于拉格朗日函數(shù)的迭代算法

由于等權(quán)條件下利用奇異值分解函數(shù)的總體最小二乘算法能夠有效解決礦山測(cè)量數(shù)據(jù)中模型參數(shù)的求解問(wèn)題，但是在大多數(shù)情況下，該方法不能滿(mǎn)足測(cè)量精度的要求，因此，為了解決不等精度條件下模型參數(shù)的求解問(wèn)題，研究了基于拉格朗日函數(shù)的迭代算法，該算法模型[6]：

f(v,b,λ,x)=vTQ-1v+bTQB-1b+2λT(l+v-Bx-EBx)

(5)

式中：λ為拉格朗日乘數(shù)；Q為權(quán)矩陣P的協(xié)因數(shù)陣。

經(jīng)過(guò)推導(dǎo)可得出：

(6)

對(duì)式(6)整理可得到各個(gè)變量的函數(shù)值：

(7)

(8)

(9)

式中：P與PB的關(guān)系為：vTPv+bTPBb=min，r為多余觀測(cè)，通過(guò)采用改正后的系數(shù)矩陣迭代和觀測(cè)向量，從實(shí)現(xiàn)對(duì)模型參數(shù)的求解。

2 模型參數(shù)求解計(jì)算

以某礦區(qū)為例，對(duì)其高程點(diǎn)進(jìn)行高程異常擬合試驗(yàn)，在試驗(yàn)區(qū)內(nèi)布置控制點(diǎn)27個(gè)，采用三等高程控制測(cè)量方法得到控制點(diǎn)高程(h)，通過(guò)使用全球定位系統(tǒng)C級(jí)控制網(wǎng)，采用三維無(wú)約束平差方法得到大地高(H)與平面坐標(biāo)(x,y)，其中大地高(H)是地面點(diǎn)沿參考橢球面法線(xiàn)到參考橢球面的距離。試驗(yàn)采用WGS-84坐標(biāo)系統(tǒng)的參考橢球面作為投影面。高程異常δ的取值：

δ≈H-h

(10)

控制點(diǎn)點(diǎn)位分布如圖1所示，采用二次多項(xiàng)式模型對(duì)描述的高程異常與模型變量進(jìn)行統(tǒng)計(jì)性規(guī)律[7]。

圖1 控制點(diǎn)的點(diǎn)位分布

把觀測(cè)值看作為等權(quán)觀測(cè)，采用總體最小二乘算法(基于奇異值分解)、最小二乘算法，對(duì)模型參數(shù)進(jìn)行求解，選取控制點(diǎn)(2,4,5,7,9,11,12,15,16,18,19,20,22,23,25)作為已知點(diǎn)，把其余控制點(diǎn)作為檢核點(diǎn)，對(duì)擬合模型參數(shù)(α0,α1,α2,α3,α4,α5)進(jìn)行求解，總體最小二乘算法與最小二乘算法求解的參數(shù)見(jiàn)表1。

表1 總體最小二乘算法與最小二乘算法求解的參數(shù)

通過(guò)使用總體最小二乘算法與最小二乘算法，對(duì)已知點(diǎn)的擬合殘差與檢核點(diǎn)的擬合殘差進(jìn)行求解，結(jié)果見(jiàn)表2。

表2 模型擬合殘差 mm

已知點(diǎn)的擬合殘差與檢核點(diǎn)的擬合殘差如圖2所示。

圖2 已知點(diǎn)與檢核點(diǎn)的擬合殘差

根據(jù)上述分析的擬合殘差，對(duì)外符合精度與內(nèi)符合精度分別進(jìn)行了計(jì)算模擬，把單位權(quán)重誤差δ0作為衡量模型內(nèi)符合精度的依據(jù)，把檢核點(diǎn)擬合的殘差值中的誤差作為衡量模型外符合精度的依據(jù)[8]：

(11)

表3 模型擬合的精度

由表3可知：當(dāng)采用總體最小二乘法進(jìn)行模擬參數(shù)推估時(shí)，其擬合精度并不是顯著高于最小二乘法，但是算法的幾何意義與理論意義表明，總體最小二乘算法處理數(shù)據(jù)的精度要遠(yuǎn)遠(yuǎn)高于最小二乘法[9]。推斷原因可能為：礦區(qū)測(cè)量數(shù)據(jù)的精度受數(shù)據(jù)自身特征以及模型的適應(yīng)性影響，因此對(duì)礦區(qū)高程異常擬合進(jìn)行加權(quán)下的研究，對(duì)模型的觀測(cè)量向量進(jìn)行定權(quán)：

(12)

式中：di為控制點(diǎn)間的平面距離；Pi為元素的權(quán)值，[]表示取不大于元素的整數(shù)，根據(jù)元素的權(quán)值對(duì)角權(quán)矩陣P進(jìn)行確認(rèn)，利用加權(quán)最小二乘算法對(duì)模型參數(shù)進(jìn)行推估；采用總體最小二乘法，對(duì)模型參數(shù)進(jìn)行4次迭代，得到的參數(shù)估值見(jiàn)表4。

通過(guò)加權(quán)條件下的模型對(duì)模型參數(shù)進(jìn)行求解，擬合點(diǎn)的殘差見(jiàn)表5。

通過(guò)使用檢驗(yàn)點(diǎn)擬合的殘差來(lái)實(shí)現(xiàn)對(duì)總體最小二乘算法和最小二乘短發(fā)的方差，模型擬合精度見(jiàn)表6。

從表6分析可知：在加權(quán)條件下，采用總體最小二乘算法計(jì)算出的模型參數(shù)的精度要高于最小二乘算法計(jì)算出的模型參數(shù)的精度，而且在加權(quán)條件下，計(jì)算出的模型參數(shù)的精度較高，但是并不是很顯著。

當(dāng)對(duì)礦區(qū)高程異常進(jìn)行擬合時(shí)，根據(jù)現(xiàn)在礦區(qū)測(cè)量要求(對(duì)觀測(cè)值的精度)，采用加權(quán)條件下的總體最小二乘算法，對(duì)擬合模型參數(shù)求解出的數(shù)據(jù)在精度上是遠(yuǎn)遠(yuǎn)滿(mǎn)足要求的。

表4 加權(quán)條件下總體最小二乘法和最小二乘算法模型參數(shù)

表5 加權(quán)條件下模型擬合殘差

圖3 加權(quán)條件下已知點(diǎn)與檢核點(diǎn)的擬合殘差

表6 模型擬合的精度

3 結(jié) 語(yǔ)

本文總結(jié)了基于拉格朗日函數(shù)算法與奇異值分解函數(shù)的總體最小二乘算法的函數(shù)模型，分析了函數(shù)本身的特點(diǎn)，然后在某礦區(qū)進(jìn)行應(yīng)用與對(duì)比，應(yīng)用得出，當(dāng)對(duì)礦區(qū)高程異常進(jìn)行擬合時(shí)，根據(jù)現(xiàn)在礦區(qū)測(cè)量要求，采用加權(quán)條件下的總體最小二乘算法對(duì)擬合模型參數(shù)的求解出的數(shù)據(jù)在精度上是遠(yuǎn)遠(yuǎn)滿(mǎn)足要求的，研究為解決礦山測(cè)量數(shù)據(jù)存在的問(wèn)題提供一定理論依據(jù)。