非線性模型精度評定的Bootstrap方法及其加權采樣改進方法

王樂洋，李志強

東華理工大學測繪工程學院，江西 南昌 330013

在大地測量數據處理領域，嚴格的線性模型并不多見[1-3]，而非線性特征一般更接近客觀事物的性質，已滲透于現代空間測量技術的各個領域，如地球重力場[4]、攝影測量與遙感、海洋測繪等領域，所涉及的函數模型普遍是非線性模型。傳統的獲取非線性模型的參數估值與精度信息的方法是最小二乘(least squares，LS)法，其函數模型為高斯-馬爾可夫(Gauss Markov，GM)模型，為線性模型[5-8]。隨著測量手段的多樣化和智能化，對平差結果的精度要求也越高，若采用傳統線性化的理論方法處理具有較高精度的觀測數據，可能會產生大于觀測誤差的模型誤差，進而導致平差結果的精度損失與特征轉變。因此，線性近似方法可能無法滿足部分高精度的要求，對非線性平差參數估計與精度評定方法的研究成為測繪領域所研究的熱點問題之一[9-11]，也是提升成果質量的重要需求。

已有能夠解決非線性模型中未知量最佳估值獲取問題的方法主要有線性化法、迭代解法及搜索算法3種[12-13]。線性化法是用線性模型的求解理論與方法近似求解，其弊端在于當模型的非線性程度較強時，將嚴重影響結果質量。迭代解法包含牛頓迭代法、高斯-牛頓法、信賴域法、擬牛頓法等，此類方法的局限性在于對迭代初值較為敏感，對較差的初值可能會出現迭代不收斂現象[9,13-14]。搜索算法包含粒子群算法、模擬退火法、遺傳算法等，該類方法無法獲取參數的解析表達式，通常以犧牲計算耗時來代替求導計算[15-17]。

已有能夠獲取參數估值二階精度信息的非線性精度評定方法主要有近似函數法和近似非線性函數的概率密度分布兩種方法。近似函數法利用泰勒公式對非線性模型展開并截取至二次項，根據誤差傳播定律獲取參數估值的方差-協方差陣。文獻[18]推導了非線性函數協方差傳播的一般公式和含有二次項的協方差傳播公式。文獻[19]將高階泰勒級數展開式方法用于GIS誤差傳播中，豐富了近似函數法在非線性精度評定的理論研究。文獻[20]將總體最小二乘(total least squares，TLS)視為一般的非線性模型，通過泰勒級數展開，導出了參數估值二階精度的協方差陣計算公式和改正數的偏差計算公式。近似函數法具有固定的方差陣解析表達式，但無法避免求導計算。當遇到未知參數向量與觀測值向量之間存在復雜非線性關系或觀測數據量較大的多維非線性函數時，求導計算復雜且獲取困難。

近些年，通過采樣來代替求導計算的近似非線性函數概率密度分布的精度評定方法流行起來，包括蒙特卡羅(Monte Carlo，MC)法[21-23]、Jackknife方法[24]等。蒙特卡羅法通過重復模擬隨機樣本來近似非線性函數的概率密度分布，進而輸出待統計量的精度信息。文獻[22]采用蒙特卡羅法計算近似單位權方差的偏差，并將偏差改正后的單位權方差估值用于求解參數估值的均方誤差。針對蒙特卡羅法模擬次數的選擇具有主觀性且無法對精度進行直接控制的問題，文獻[23]將自適應蒙特卡羅法用于非線性模型精度評定，采用對偶變量的思想，給出了獲取參數估值偏差的對偶自適應蒙特卡羅法的計算流程。當代計算機技術的發展使得蒙特卡羅方法的實現成為了可能，但該方法的模擬次數通常需要設置在105次以上，它往往以增加計算耗時來獲得精度信息，并且隨著采樣次數的變化其精度統計信息并不唯一[25-26]。Jackknife采樣方法按順序逐次減小特定數目的觀測值來生成一系列Jackknife樣本，從而獲取方差信息。文獻[24]將Jackknife方法用于TLS精度評定，給出了獲取參數估值精度信息的詳細采樣步驟，進一步擴展了非線性精度評定的近似函數概率分布法理論。但Jackknife方法需要保證足夠的觀測數據才能獲得合理的精度信息，并且隨著d值的增大，其計算量也會迅速增加[26]。

Bootstrap方法[27-29]稱自助法，是與Jackknife方法緊密相關的一種統計推斷方法，均用于估計或修正統計估計值的偏差或方差信息。與Jackknife方法一樣，Bootstrap方法通過檢查樣本數據內的變化，而不是通過參數假設來估計一個統計量的變異性，但Jackknife方法并沒有Bootstrap方法那么普及[30]。文獻[27]首次提出了通過重復采樣得到的自助世界的經驗分布來逼近總體真實分布的思想，并給出了其基本采樣策略和相關理論證明。由于Bootstrap方法通過重采樣代替求導計算，只需選擇采樣次數并結合有放回抽樣策略獲取自助樣本，最后解算每個樣本便可求得未知統計量的均值和協方差陣，因此歷年來有較多學者將該方法用于偏差分析和方差計算[31-34]。目前，Bootstrap方法的研究主要集中于該方法在統計學上的性質，將該方法用于解決大地測量領域非線性模型精度評定的問題有待研究；考慮到Bootstrap方法的原始采樣策略為等概率采樣，對于不等概率的采樣方式的適用性還有待驗證。

針對以上問題，本文從獲取更為合理的精度信息出發，針對現有非線性精度評定方法的不足，將Bootstrap方法用于獲取非線性模型參數估值的精度信息；顧及Bootstrap方法原始的等概率采樣方式，提出加權采樣策略的改進方法，以至獲取更加合理的方差估值；最后通過算例驗證本文非線性模型精度評定的Bootstrap方法的有效性及本文改進方法的正確性和優勢。

1 非線性精度評定及其二階近似方差-協方差

將非線性平差函數模型定義為

(1)

參照文獻[35]的推導思路，可得參數估值的二階近似方差-協方差陣為[35]

(2)

其中，J、H可以具體表示為

(3)

(4)

式中，

(5)

(6)

(7)

(8)

(9)

式中，ei=[0 0 … 1 … 0 0]T∈Rn×1，第i處為1，其余為零。

由式(2)可以看出，在非線性模型精度評定中，若采用近似函數法獲取非線性模型參數估值的二階近似協方差信息，其泰勒級數展開過程需要計算非線性模型的Jacobian矩陣和Hessian矩陣，且當遇到數據量較大的多維非線性函數或非線性程度強的模型時，偏導數的計算十分困難。因此，為避免復雜的求導計算并獲取更加合理的精度信息，本文將Bootstrap采樣方法及其加權采樣的改進方法引入非線性模型精度評定中，試圖有效解決近似函數法難以對非線性函數求導的問題。

2 非線性模型精度評定的Bootstrap方法

非線性平差模型的參數估計與精度評定問題可以理解為，利用從總體中抽取的固定原始觀測樣本對非線性模型中的未知統計量及其精度信息進行統計推斷。Bootstrap方法是將原始樣本作為總體，利用有放回隨機抽樣法從總體分布函數中得到一系列獨立樣本(稱為自助樣本)，并通過計算每個自助樣本來獲取未知統計量抽樣分布的經驗估計[27,30]，最后利用這個估計的抽樣分布來做總體推斷。該方法的優勢性主要表現在它不需要對未知模型及分布做任何假設，也無須推導估計量的精確表達式，只需通過檢驗樣本內統計量的變化便能夠估計未知參數的整個抽樣分布。

(10)

顧及非線性優化算法解算得到的參數估值不再是無偏估計量，而Bootstrap方法經有放回抽樣獲取的大量自助樣本能夠減小參數估值的偏差，在一定程度上能夠有效改善參數估值的質量[32]。因此，本文將式(10)的結果作為未知參數的自助估計。

(11)

Bootstrap方法的重采樣要點是采樣這個過程的隨機項。因此，根據重采樣對象的不同，Bootstrap方法又可分為：重采樣觀測值的Bootstrap方法和重采樣殘差的Bootstrap方法[30]。Bootstrap方法通過重采樣隨機項獲取自助樣本的方式代替復雜的或難以實現的求導計算，可以用于求取非線性平差參數估值的方差或協方差陣。因此，本文將這兩種采樣策略引入非線性模型的精度評定問題中，給出了獲取參數方差的具體步驟。

2.1 重采樣觀測值的Bootstrap方法

重采樣觀測值的Bootstrap方法，通過對觀測值的有放回隨機抽樣將原始觀測數據采樣成多個自助樣本，原始樣本集中的每個樣本點被采樣到的概率相同，每個自助樣本的樣本容量與原始樣本相同。對于不等精度的觀測數據，需要保證自助樣本中的觀測值與其權值相對應，因此在重采樣獲得觀測值自助樣本的過程中需要對觀測值的權Pi進行重采樣，兩者采樣方式一致，具體表現為以隨機數序列元素為自助樣本元素的下標，對觀測值及其權值進行賦值，進而獲取自助樣本及對應權陣。根據以上分析，總結得到重采樣觀測值的Bootstrap方法用于非線性精度評定的完整計算流程：

(1)假設觀測值N=(x1,x2,…,xn)為隨機項樣本X1，其權陣為P=diag(P1,P2,…,Pn)，將1/n的概率分別設置在(x1,x2,…,xn)各樣本值上。

(2)產生滿足均勻分布的n個隨機數(i1,i2,…,in)，其中1≤ij≤n,j=1,2,3,…,n。

(5)循環步驟(2)—步驟(4)M次，得到M組參數估計值。

(6)利用式(10)獲取參數均值，利用式(11)獲取參數估值的方差。

重采樣觀測值的Bootstrap方法的采樣過程將原始觀測值采樣成了多個自助樣本，盡管每個自助樣本的樣本容量與原始樣本相同，但并不是改變觀測值的先后排列順序，有放回隨機抽樣過程使得每個自助樣本中可能存在重復的原始數據點，而另外一些樣本點沒有出現。因此，每個自助樣本將隨機地異于原始樣本，導致每個自助樣本獲得的參數估值存在細微差異。

2.2 重采樣殘差的Bootstrap方法

文獻[27]在提出Bootstrap方法的同時，對抽取獨立同分布樣本點的采樣方式進行了延伸和拓展，結合多元回歸分析給出了一種對殘差進行抽樣的采樣策略。具體做法是，首先使用所有樣本點建立模型并獲取殘差，然后對殘差進行采樣并根據模型的方程重構觀測向量，最后計算參數并重復該步驟。因此，當模型的自變量為固定常數，因變量為自變量和一個隨機誤差項的函數時，基于模型的誤差結構，可將因變量的殘差設置為采樣過程的隨機項。

重采樣殘差的Bootstrap方法與重采樣觀測值的Bootstrap方法相比，主要的區別在于重采樣對象不同。采用重采樣殘差的Bootstrap方法解決精度評定問題的計算步驟為：

(7)重復步驟(4)—步驟(6)M次，共得到M組參數估值。

(8)通過式(10)獲取未知參數的自助估計值，利用式(11)計算參數估值的方差。

重采樣觀測值與重采樣殘差的樣本數據均來源于原始固定樣本，并未根據更多的觀測信息進行計算。因此，這兩種采樣方法均是利用相同的觀測數據最終得到未知參數估值及方差信息，僅是采樣過程中的隨機項不同，其采樣過程均不會產生額外的模型誤差，也不會改變模型態性；并且循環有放回隨機采樣過程獲得的大量自助樣本能夠減小迭代方法求解非線性參數估值的偏差，從而提高參數估值的精確度。

3 改進的Bootstrap方法

Bootstrap方法在重采樣隨機項元素構建自助樣本過程中，每個元素被采樣到的概率相同且均為1/n。可以看出，Bootstrap方法在采樣時將樣本數據視為等精度觀測，在構建自助樣本過程中忽略了數據的權值信息，使得精度不等的觀測數據出現在自助樣本中的概率卻相同，因此可能無法較好地逼近總體分布函數，在一定程度上有損理論嚴密性。

為使得Bootstrap方法能夠更加充分利用總體包含在樣本中的先驗信息和數據性質，使得自助樣本的經驗分布更加逼近總體分布函數，若利用觀測值的分布信息來構建經驗分布函數，通常需要借助分布假設，并且所得分布函數將嚴重依賴于所作的假設。當觀測樣本的分布假設與實際觀測樣本的分布存在較大偏差時，所得樣本經驗分布函數并不能很好地逼近總體分布函數，導致統計推斷結果與實際存在一定的偏差[36]。顧及觀測值的權的獲取方式較觀測值的分布信息獲取方式更為簡便，并且利用觀測權來構建經驗分布函數可以有效避免對觀測樣本分布的假設。因此，權作為比較觀測值之間精度高低的一種先驗信息，可用于調整隨機項的抽樣概率，以至獲取更加逼近總體分布函數的樣本經驗分布函數。因此，本文將原始觀測樣本的權值納入重采樣獲取自助樣本的過程。在該過程中，權值被用于構建一個自助樣本獲取準則，通過對隨機項元素進行采樣可能性的重新分配，使得每個元素被抽樣到的概率與其對應的觀測精度相匹配。

首先，對原始樣本所有元素的權值進行歸一化預處理，使得歸一化權值之和為1

(12)

式中，NORMi表示第i個觀測元素的歸一化權值；Pi是原始觀測樣本中第i個元素的權；n為原始樣本的樣本容量。

將NORMi作為隨機項中各元素xi的采樣概率，其概率值的大小取決于各元素所對應的權值與所有觀測權之和的比值大小。因此，可構建Bootstrap采樣過程中隨機變量的分布律函數

proxi=NORM{Y=xi}

(13)

式中，xi為隨機項X中的樣本值；Y為隨機項X中的隨機變量，其可能的取值為xi(i=1,2,…,n)；proxi表示事件{Y=xi}的概率。

(14)

由于Bootstrap法在采樣過程中隨機項的選擇不同，Bootstrap法可區分為重采樣觀測值的Bootstrap方法和重采樣殘差的Bootstrap方法，將以上加權采樣策略應用到該兩種方法中，分別給出將改進方法用于非線性模型精度評定中的詳細計算步驟。

3.1 改進的重采樣觀測值的Bootstrap方法

該方法首先對觀測數據的權值進行歸一化處理，計算每個樣本元素被采樣到自助樣本的概率值；構建采樣過程中的隨機變量的分布律函數；產生滿足該分布律函數的隨機數，通過賦值獲取自助樣本并計算參數；然后將生成自助樣本到獲取參數的整個過程循環M次；最后對未知統計量的均值和方差進行統計推斷。每個自助樣本的樣本容量與原始樣本相同，但原始樣本中的每個元素被采樣到的概率受其權值影響而不同，這取決于該樣本點的權值大小。根據以上分析，總結得到改進的重采樣觀測值的Bootstrap方法的計算流程：

(1)假設觀測值樣本N=(x1,x2,…,xn)為采樣過程中的隨機項X1，其權陣為P=diag(P1,P2,…,Pn)。

(2)利用式(12)對權陣進行歸一化處理，使得歸一化后的權值之和為1。

(4)通過式(13)、式(14)構建隨機變量的經驗分布函數。

(5)產生滿足該分布函數的n個隨機數(i1,i2,…,in)，其中1≤ij≤n,j=1,2,3,…,n。

(9)重復步驟(5)—步驟(8)M次，得到M組參數估計值。

(11)分別通過式(10)、式(11)獲取參數的均值和方差。

重采樣觀測值的Bootstrap方法及其改進方法均是將相同的原始觀測數據重采樣成了多個自助樣本，且獲取的自助樣本個數及每個自助樣本的樣本容量均相同，但改進方法由于采用了加權采樣的策略，將等概率重采樣轉變為不等概率采樣方式，使得權值大的原始樣本元素出現在自助樣本中的概率增大，從而使得自助樣本能夠更加貼合原始樣本。

3.2 改進的重采樣殘差的Bootstrap方法

與重采樣殘差的Bootstrap方法一致，改進的重采樣殘差的Bootstrap方法仍將因變量的殘差作為采樣過程中的隨機項，主要區別在于改進方法在采樣前對權值進行了歸一化預處理，通過獲取隨機項中隨機變量的分布函數，使得對殘差項的等概率采樣轉化為加權采樣。改進的重采樣殘差的Bootstrap精度評定方法的計算步驟總結為：

(1)從總體中獲取固定的觀測向量L及其權陣P，作為擬合該模型的原始樣本。

(5)對權陣P的對角元素Pi進行歸一化處理。

(7)通過式(13)、式(14)構建隨機變量的經驗分布律函數。

(8)產生滿足該分布律的隨機數(i1,i2,…,in)，其中1≤ij≤n,j=1,2,3,…,n。

(12)循環步驟(8)—步驟(11)M次，得到M組參數估值。

(14)分別通過式(10)、式(11)獲取未知參數的自助估計值及方差。

重采樣殘差方法和改進的重采樣殘差方法在自助樣本獲取過程中，均未對原始觀測樣本進行采樣，因此均利用了原始樣本中的所有觀測值，但對隨機項的有放回抽樣使得每個自助樣本包含的元素不同，進而使得每個自助樣本得到的估值存在差異。

根據以上4種精度評定流程可以發現，不管是采樣觀測值，還是采樣殘差，抑或是改進方法，重復M次的Bootstrap采樣過程僅由有放回抽樣的隨機性和滿足的抽樣分布所決定，即第i次的自助樣本獲取過程并不會影響第j次的重采樣，這使得自助樣本的獲取過程是獨立的。此外，由于自助樣本中的所有元素均來源于原始樣本，未利用其他的觀測值進行計算，并且由于受自助法自身有放回抽樣性質以及抽樣隨機性影響，使得自助樣本之間可能包含重復的隨機項元素，導致自助樣本存在一定的相關性。第i組參數估值僅僅由第i個自助樣本決定，因此每個自助樣本的解算過程是獨立的；但由于自助樣本i與自助樣本j可能包含相同元素，使得M組參數估值存在一定的相關性，使得各自助樣本的參數估值在自助均值左右浮動并且較為接近。但在采用矩法估計計算自助方差時，并不能簡單忽略M組參數估值的分布規律或者將其視為相同，這是因為自助重采樣的目的便是要得到一個較好的抽樣分布估計來對模型中的未知統計量進行統計推斷。

4 案例研究與分析

經試驗，Bootstrap采樣方法適用于線性平差模型的參數估計與精度評定問題，通過有放回隨機抽樣能夠在一定程度上提升參數質量，并且獲取合理的參數估值精度信息，表明在線性模型精度評定問題中，Bootstrap方法具有可靠性和優勢性。與最小二乘法相比，盡管實施Bootstrap采樣會增加精度評定過程計算耗時，但自助法具有操作簡單、易于編碼實現的優點，并且其結果存在一定程度上的質量提升，不失為一種較好的精度評定方法。在非線性模型精度評定問題中，本文第一節已證明已有的基于泰勒級數展開的近似函數法，在精度評定過程中需要計算非線性函數的一階或二階偏導數，當遇到高維且非線性程度較強的模型時，求導過程復雜難以獲取。而自助法用于解決非線性模型精度評定問題最大的優勢在于，獲取參數精度信息過程中能夠避免導數計算，該方法是通過對原始樣本的重采樣過程來代替復雜的Jacobian矩陣或Hessian矩陣計算過程，最后利用抽樣分布估計對未知參數進行統計推斷。因此，為充分體現Bootstrap方法在測量平差精度評定問題中的優勢性，本文考慮的是非線性模型情形。

為評估Bootstrap方法和本文改進方法在獲取非線性模型參數估值精度信息方面的性能，以及探討在大地測量領域中的應用，本文共設計了4個算例：兩個非線性強度不同的非線性回歸模型、邊角網平差模型、火山形變Mogi模型。通過算例1和算例2來驗證Bootstrap方法在非線性模型精度評定中的可行性及本文改進方法的正確性；通過算例3和算例4來展示Bootstrap方法在大地測量非線性精度評定中的應用價值，并驗證本文改進方法在參數方差獲取方面的優勢。對比試驗中，Monte Carlo方法的模擬次數取值為105。Bootstrap方法重采樣次數的選擇取決于統計量的具體研究問題并且依賴于需要做的檢驗[38]。文獻[38—41]指出，當重采樣次數為50～200時，Bootstrap方法能夠得到較為合理的未知參數方差或中誤差；當采樣次數超過200時，方差近似值的改進很小。經試驗，通過增加Bootstrap重采樣次數在一定程度上能夠提升精度評定結果質量，但隨著采樣次數的增加，中誤差結果趨于較為穩定的值。這與已有文獻的結論較為一致。因此，本文依據已有研究的經驗取值并顧及非線性平差精度評定問題，在保證較高的計算精度條件下，同時考慮較高的精度評定效率，將Bootstrap方法的采樣次數設置為103。案例研究中出現的字符及對應含義列于表1。

表1 字母縮寫及對應含義Tab.1 Abbreviations and their corresponding meanings

4.1 Bootstrap方法及本文改進方法的可行性及正確性驗證分析

針對非線性回歸模型的非線性強度不同，設計兩個仿真試驗，通過Bootstrap方法及加權采樣方法在參數估計中的應用，以及與AF法、JK法、MC法在參數標準差估計中的性能評估對比，以清楚地揭示本文精度評定方法的可行性和有效性。

4.1.1 算例1

算例1采用的非線性回歸模型函數形式如下

(15)

式中，xi為n維自變量向量的元素；yi為對應因變量值；εyi為因變量yi的隨機誤差；ξ1、ξ2為回歸參數。

4.1.2 算例2

為進一步探討AF法、JK法、MC法、Bootstrap方法及本文加權采樣方法在更復雜非線性情況下的方差估計性能和計算特性，采用比算例1非線性程度更強的回歸函數，其函數形式定義如下

(16)

式中，xi、yi分別為橫坐標x和縱坐標y的自變量；zi為n維自變量向量中的元素；εzi為zi的隨機誤差，ξ1、ξ2、ξ3、ξ4為回歸參數。

4.1.3 結果分析

在獲取模擬觀測數據后，分別采用LS法、GN法、BO法、MBO法、BR法及MBR法解算并獲取參數均值，其結果及差值二范數列于表2；再分別采用BO法、MBO法、BR法及MBR法對隨機項進行采樣，并通過高斯-牛頓解法解算自助樣本，獲取回歸系數的標準差及協方差信息，最后與FO法、SO法、JK法及MC法的結果進行對比，各方法的結果列于表3。

表2 非線性回歸模型參數估值及其與真值之間的二范數Tab.2 The estimated parameters and the norm of the nonlinear regression model

表3 非線性回歸模型參數估值的標準差和協方差計算結果Tab.3 The standard deviations and covariances of parameters estimations of the nonlinear regression model

由表2的參數估值及范數結果可以看出，在非線性強度不同的兩個仿真設置中，LS法的結果均為最差，表明LS法在將非線性模型線性近似的過程中，引入了較大的模型誤差，從而嚴重影響估計值的質量。而相比于LS法，GN法得到的參數估值較優，說明GN法能夠削弱非線性回歸模型線性近似所帶來的模型誤差的影響，其循環迭代過程可以不斷修正參數估值，在一定程度上能夠有效改善傳統LS法在對非線性模型線性化過程中所引起的精度損失，使其參數估值比傳統的線性近似方法更接近參數真值。BO法和BR法得到的參數估值與真值的范數均比LS法和GN法的結果更小，表明BO法和BR法通過實施Bootstrap并采用GN法解算自助樣本，在削弱線性化帶來的影響的基礎上，能夠在一定程度上減小參數估值的偏差，獲取更為精確的參數估值。與BO法相比，MBO法得到的參數估值與真值的范數更小，表明采用了加權采樣的重采樣觀測值方法在參數的精確度方面得到了較明顯的改善；而相比于BR法，MBR法也能夠提升參數估值的精確度。這表明MBO法及MBR法在獲取較高精度的參數估值方面比BO法和BR法更加有效。試驗結果表明，本文的加權采樣策略在非線性回歸模型的參數估計方面是可行且有效的，改進方法通過重采樣獲得的自助樣本能夠充分利用原始樣本的先驗信息，從而使得MBO法和MBR法獲得的估值較BO法和BR法更接近參數真值。

由表3可知，在兩個不同的仿真設置下，二階近似函數法獲取的參數估值標準差在數值上略大于一階近似函數法得到的標準差結果，且更接近于蒙特卡羅方法獲取的參數標準差。這是因為二階近似函數法顧及了非線性平差模型的二階項，說明由一階近似函數法得到的結果可能高估了非線性平差模型參數估值的精度。但由于二階近似函數法也僅僅包含至二階項，忽略了高階項，因此也說明了蒙特卡羅方法獲取的精度信息更加合理。并且試驗證明，當原始樣本量增大時，二階近似函數法所涉及的Hessian矩陣維數也將增大，也說明了無須求導計算的Bootstrap方法更加適用較大量的觀測數據。相比于蒙特卡羅方法的結果，Jackknife法獲得的標準差偏離較大，是因為Jackknife法易受原始樣本容量大小的影響且其重采樣過程獲取的Jackknife樣本較少，與文獻[24]中的結論一致。這表明Jackknife法往往會低估非線性模型參數估值的精度，因此該方法無法反映非線性模型參數估值的真實精度。

利用基于Bootstrap的方法對該非線性函數進行精度評定所獲得的標準差與近似函數法的結果在符號和量級上相同、在數值上相近，表明該4種方法均能夠對非線性平差模型進行有效的精度評定。從計算結果上看，BO法和BR法獲取的參數的標準差小于Jackknife法的結果且更接近于蒙特卡羅方法的計算精度，表明BO法和BR法較Jackknife法更加精確。相比于BO法，MBO法獲取的參數估值標準差在數值上更小，且與蒙特卡羅方法的結果更為接近；MBR法計算得到的參數標準差比BR法的結果更小并且同樣更靠近于蒙特卡羅方法獲取的標準差，說明MBO法和MBR法能夠獲取更為合理的參數標準差。結合表2中的參數估值結果，說明本文加權采樣的獲取自助樣本方式不僅適用于重采樣觀測值，也適用于重采樣殘差情形。試驗結果表明，將Bootstrap方法用于解決非線性參數估計與精度評定問題是可行且有效的，本文加權采樣方法獲取的更為精確的參數估值和精度信息，也證明了本文改進方法更具可靠性和優勢。

4.2 算例3：Bootstrap方法及本文改進方法在邊角網平差模型中的應用

將文獻[42]中的例3-2-1進行改化，獨立不等精度觀測如圖1所示，邊角網有12個角度和6條邊長，其中A、B、C為已知點，P、D為待定點。起算數據和觀測值分別列于表4和表5。

圖1 邊角網Fig.1 The sketch of triangulateration network

表4 起算數據Tab.4 The data of known control points m

在獲取角度和邊長觀測數據后，首先采用余切公式計算得到待定點坐標估值作為迭代初值，再分別采用BO法、MBO法、BR法及MBR法對隨機項采樣并通過高斯-牛頓迭代獲取未知點坐標，其參數結果及差值二范數列于表6。通過算例1和算例2的結果可知，Jackknife法的結果不足以反映參數估值的真實精度；而一階近似方法的結果可能會高估參數估值的精度；并且根據邊長觀測方程及測角觀測方程可以看出，邊角網平差模型中的未知點坐標與觀測值向量之間存在復雜的非線性關系，若采用泰勒級數展開并截取至二次項的方法獲取未知參數的方差信息，需要進行多次計算偏導，其求導過程十分復雜且較難獲得。針對以上問題，因此算例3不再使用一階近似函數法、二階近似函數法及Jackknife方法求解參數估值的精度信息，僅通過基于Bootstrap的4種方法與MC法進行對比，來驗證本文方法的正確性和優勢。待定點坐標標準差及協方差結果列于表7。

表5 角度和邊長觀測數據及其權值Tab.5 The data of angular and traverse measurement and corresponding weights

表6 邊角網平差模型坐標估值及其與真值之間的二范數Tab.6 The estimated coordinate parameters and the norm of the triangulateration network model

對比表6中BO法、MBO法、BR法及MBR法的結果可以看出，4種方法所得坐標估值均與坐標真值較為接近，再次驗證了Bootstrap方法及本文改進方法的正確性。同時，通過對比估值與真值的范數可知，MBO法所得估值較BO法的結果更接近于真值；與BR法的結果相比，MBR法所得結果的范數有較為明顯的改善。這表明MBO法及MBR法在獲取更高精度的邊角網坐標估值方面比BO法和BR法更加有效。本算例的參數估值結果證明了Bootstrap方法在邊角網平差模型中依然適用，也再次證明了本文改進方法的有效性及加權采樣的必要性。

由表7可知，在邊角網平差模型中，BO法、BR法、MBO法及MBR法計算得到的待定點坐標估值標準差在符號和量級上均與MC法的結果保持一致，同樣驗證了該4種方法在非線性平差精度評定問題中的正確性和有效性。由于邊角觀測方程的復雜非線性關系，使得二階近似函數法獲取參數中誤差較為困難，從側面體現了通過執行Bootstrap采樣來避免復雜求導運算仍能獲取合理精度信息的優勢，也說明了在面對高維復雜非線性模型的精度評定問題時，基于Bootstrap的方法具有更強的適用性。相比于BO法及BR法計算得到的標準差結果，MBO法和MBR法的結果在數值上更小，這表明通過充分利用觀測值先驗信息的加權采樣方法能夠在一定程度上提升參數精度的可靠性，也說明為更加合理反映參數估值的精度信息，采用加權采樣的抽樣方式更加有效。

表7 邊角網平差模型待定點坐標標準差和協方差計算結果Tab.7 The standard deviations and covariances of the estimated coordinate parameters in the triangu-lateration network model

4.3 算例4：Bootstrap方法以及本文改進方法在火山形變Mogi模型反演中的應用

為進一步驗證Bootstrap方法及本文改進方法在大地測量非線性精度評定中的廣泛適用性和應用價值，算例4采用文獻[43]首次提出的Mogi模型。該模型是一種研究由火山爆炸源所引起的地表形變位移及重力變化的地球物理模型，它能夠有效地用于模擬和解釋大量的火山區地表形變。應用火山區的地表形變觀測數據來反演巖漿壓力源參數，對火山的危險性評價有著重大的實際意義。地表位移又可分為垂直位移和水平位移，本文采用垂直位移單一反演的Mogi模型。

假設地面直角坐標系為(x,y,z)，巖漿壓力源的中心坐標表示為(x0,y0,D)，因此巖漿房膨脹所引起的地表垂直位移與等效壓力源參數之間的表達式可以表示為[44]

(17)

式中，Δh表示地表垂直位移；K表示為體積彈性模量；μ表示剪切模量；D表示為源的中心深度；r表示地表徑向距離且當r=0時垂直位移取到最大值；ΔV為巖漿房體增量。

當地殼為泊松介質時，取泊松比v=0.25且顧及巖漿房體積的膨脹為半徑為R的等效球體，因此體積彈性模量K與剪切模量μ存在以下關系式

(18)

將式(18)代入式(17)，地表垂直位移可以化簡為

(19)

假設地表存在n個觀測點，則第i個監測點到巖漿壓力源的地表徑向距離為

r2=(xi-x0)2+(yi-y0)2

(20)

式中，(x0,y0)為巖漿壓力源的中心在平面上的投影坐標。

將式(20)代入式(19)，地表垂直位移可表示為一個多維非線性函數

Δhi+ei=fi(ξ)=

(21)

式中，i=1,2,…,n；ei為第i個監測點的垂直位移Δhi的隨機誤差；f(·)表示未知參數向量ξ與觀測值之間的非線性映射關系；ξ=[ΔVDx0y0]T為待求的未知參數。

利用泰勒級數對式(21)在參數初值ξ0處展開并舍去二次及以上項，得到垂直位移反演壓力源參數的平差模型

(22)

通過仿真獲取某一火山區域的垂直位移觀測資料：觀測點的取值范圍為[x,y]∈[-5,5] km×[-5,5] km，相鄰垂直形變監測點的間隔為0.5 km，每個觀測點的權值在[0,25]之間隨機取值。壓力源參數真值設置為：體積增量ΔV=6000×103m3、源的中心深度D=4 km、膨脹源的中心坐標x0=0 km、y0=0 km；將地面監測點坐標及壓力源參數真值代入Mogi模型，正演得到地表垂直位移；最后對垂直位移模擬位移值施加均值為0、單位權方差為2 mm的隨機誤差，得到的垂直位移如圖2所示。

圖2 由火山Mogi模型正演得到的地表垂直形變Fig.2 The vertical displacement of forward the volcano Mogi model

由于垂直位移反演壓力源參數的平差模型是由泰勒級數對式(17)進行線性化處理后獲得，針對反演過程中可能出現的病態奇異情形[45]，算例4采用嶺估計法進行參數求解，分別采用BO法、MBO法、BR法和MBR法進行垂直位移單一反演，具體表現為在獲取自助樣本后，通過嶺估計方法來解決法矩陣求逆不穩定的問題并獲取巖漿壓力源參數[46-47]。4種方法反演得到的參數估值結果及差值二范數列于表8。

表8 火山Mogi模型壓力源參數反演結果及其與真值之間的二范數Tab.8 Inversion results of pressure source parameters and the norm of the volcano Mogi model

由表8中的參數估值可以看出，本文提出的4種方法獲取的估值與真值的差異主要表現在巖漿房體積增量和巖漿壓力源的中心深度。MBO法所得結果優于BO法，其參數估值與真值之間的范數減小了40.7%；與BR法相比，MBR法的范數減小了58.5%。結果表明，本文方法均能得到合理的參數估值，而采用了加權采樣的MBO法和MBR法在提升估值質量方面更具優勢。試驗結果說明，將Bootstrap方法用于法矩陣病態性嚴重的火山形變Mogi模型參數反演是可行且有效的，也說明了本文的加權采樣策略同樣適用于大地測量領域病態模型的參數估計問題。

由式(17)—式(22)可以看出，火山Mogi模型中的未知參數向量ξ與地表垂直位移觀測值向量Δh之間存在復雜的高維非線性關系，若采用對非線性模型泰勒級數展開的近似函數方法獲取未知參數的方差信息，需要進行多次計算偏導，使得精度評定過程十分復雜。因此，算例4僅采用本文提出的4種方法及MC法進行精度評定，各方法獲取的壓力源參數標準差結果列于表9。

由表9可知，BO法、BR法、MBO法及MBR法計算得到的參數估值標準差與MC法的結果在符號和量級上均保持一致，說明了本文Bootstrap方法精度評定的可行性和有效性。從估值的近似標準差結果來看，相比于BO法所獲得的參數標準差，MBO法計算得到的標準差比BO法的計算結果在數值上更小，且比BO法的結果更加接近于MC法所獲得的精度，說明MBO法能夠獲取比BO法更為合理的精度信息；與BR法的計算結果相比，MBR法的標準差有較為明顯的減小且更為靠近MC法的計算精度，結果同樣說明了MBR法計算得到的標準差比BR法更為精確。試驗結果表明，將Bootstrap方法用于獲取病態非線性平差模型的參數估值及精度信息依然是可靠且有效的，本文改進方法在集結Bootstrap方法的所有優點外，通過加權采樣的抽樣方式還能夠更加充分地利用權值信息，從而能夠獲得更加精確的估值與精度信息值。

表9 火山Mogi模型壓力源參數的標準差和協方差計算結果Tab.9 The standard deviations and covariances of the pressure source parameters of the volcano Mogi model

5 結 論

為獲得統計意義上更為合理的精度信息，若采用泰勒展開的近似函數法，則無法避免地需要計算非線性函數的Jacobian矩陣和Hessian矩陣；若采用蒙特卡羅方法，則通常需要105以上次的模擬計算。因此，提出無需求導計算且計算效率相對較高的非線性精度評定方法顯得尤為重要。本文在闡述非線性模型精度評定原理及Bootstrap方法重采樣原理的基礎上，給出了非線性模型精度評定重采樣觀測值Bootstrap方法及重采樣殘差的Bootstrap方法，并分別給出了具體的采樣步驟，試圖通過這兩種方法在避免求導計算的情形下依然能獲取合理的精度評定結果。針對Bootstrap方法的重采樣過程是對模型隨機項的等概率采樣，本文通過對原始樣本的權重信息進行歸一化處理并獲取隨機變量的經驗分布函數，將傳統自助法的等概率采樣轉化為加權采樣，嘗試能夠更加充分利用總體包含在原始樣本中的先驗信息和數據性質，以至于自助樣本的經驗分布更加逼近總體分布，進而獲取更為精確的參數方差估值。

對本文案例研究中的各方法進行對比分析可知，非線性模型精度評定的Bootstrap方法能夠獲取比近似函數法、Jackknife方法更為合理的參數標準差，這說明將Bootstrap方法用于解決非線性模型的精度評定問題是可行且有效的；同時，采用了本文給出的改進算法能夠得到比直接實施Bootstrap更加精確的精度信息，說明了本文加權采樣方法的可靠性與優勢性。因此，建議在獲取非線性模型參數估值的方差-協方差時，采用本文的加權采樣方法。

本文方法適用于觀測值相互獨立的情況，Bootstrap方法在相關觀測數據的非線性精度評定問題，以及該方法的進一步應用是接下來需要研究的內容。隨著Bootstrap方法引入非線性模型精度評定中，拓展該方法在大地測量領域中的應用，也將是非線性平差理論研究的一個重要補充。