加權遞推最小二乘在絕對重力儀數據處理中的應用

羅棋，龐聰，李查瑋

（中國地震局地震研究所（地震大地測量重點實驗室），湖北武漢430071）

加權遞推最小二乘在絕對重力儀數據處理中的應用

羅棋，龐聰，李查瑋

（中國地震局地震研究所（地震大地測量重點實驗室），湖北武漢430071）

針對絕對重力儀高精度的需求，在抗差參數估計的基礎上提出了加權遞推算法及其改進算法，加權遞推和改進算法改良抗差權重因子并且避免求復雜多維的逆矩陣。通過仿真表明，該算法對含有白噪聲或者異常數據都起到很好的仰制作用并提高了重力測量精度。

絕對重力儀；抗差估計；加權遞推；最小二乘

絕對重力測量在計量學，地球物理和大地測量有著重要的作用，其精度的提高對地震活動、地質勘測都有不可估量的作用[1]。一般測量絕對重力加速度g的方法是基于落體每下落半個波長距離產生一個干涉條紋的規律，獲得物體下落的位移信息，對比銣原子鐘給出的時間，得到時間距離對，擬合得到重力加速度值[2]。

目前普遍采用最小二乘和非線性擬合時間距離對[3]，與非線性相比，最小二乘原理簡單，易于理解和掌握，且在編程中易于實現，但是最小二乘算法適用于近似服從正態分布的測量數據，參數擬合結果對非正態分布的觀測粗差極為敏感[4]，沒有考慮觀測誤差的不同對算法解的影響[5]。用雙采樣過零檢測，時刻數據將越來越密集，時間間隔越來越短，測量的誤差影響將會放大，因此需要對最小二乘進行改進。文獻[4]提出了抗差估計，但是涉及到求解復雜多維的逆矩陣。基于上述問題，為了進一步提供算法的可行性和精度要求，提出加權遞推形式，避免了復雜逆矩陣的求解。文中基于文獻[4]中的迭代算法，采用文獻[5]的加權陣，對擬合的時間距離對進行加權處理，擬合誤差方差越大的賦予的權重越小，誤差方差小的賦予的權重大，這也充分體現了加權遞推最小二乘算法的優越性[5]。

1 加權遞推最小二乘算法

1.1 加權遞推最小二乘

絕對重力儀研究的對象是做自由落體運動的落體，針對落體運動，絕對重力儀數據符合二次多項式曲線（包含重力梯度和光速有限性請參考文獻[6]中的算法，在擬合公式中增加正弦項和余弦項以仰制激光器調制頻率對輸出的影響參考文獻[4]，將其整理后依舊可以帶入最小二乘擬合參數，并不影響算法的通用性），通常假設時間距離對數據為(ti,si)(i=1,2,…m)，最小二乘的基本思想即選擇一個n次多項式s(t)用一種使均方誤差極小的方式來逼近(ti,si)，其中n遠小于m。

寫成向量模型，模型偏差為：

其中ε為偏差，。

有m個等式，將其合成矩陣向量形式：

選取適當的對角加權矩陣W，得到加權非遞推最小二乘解：

初始化：P(0)=aIm，a為很大的實數，Im為m×m單位矩陣，（0）=0 。

圖1給出了對應的加權遞推算法計算框圖，絕對重力儀數據處理系統每讀取一個時間距離對(ti,si)，根據式（2）算出偏差ε(ti)，取權值帶入加權遞推最小二乘中，通過迭代，擬合出高精度的參數估計。對式（4）取加權矩陣W(m)=Im，則上述遞推公式則為遞推最小二乘。

圖1 加權遞推算法計算框圖

1.2 加權遞推算法改進

在上述遞推公式中，每次讀取到新的時間距離對數據，立刻參與計算得到參數，更新的參數取代舊參數參與下一步的擬合，即可以立刻判斷參數的優異程度，但是數據精度與初始化參數有關，并且需要大量擬合數據來消去初始化參數帶來的影響。實際運用中，先存取所有的時間距離對，再取出合適數據來參與計算，很明顯這種分段式處理數據，其精度將得到提高。文中將結合這種方式，對加權遞推方式做出改進。

將公式（5）做變形處理

上式中s0,v0,g0分別為初始位移、初始速度、當地重力加速度，在實際中只需要求出重力加速度參數g值，因此對式（8）做變形寫成矩陣形式：

式中δi=w(i)為權值 ，當δi=1 時，上式（9）（10）為改進最小二乘遞推。明顯看出改進之后的算法避免了求逆矩陣，由于時間距離對的選取可以是非連續的，在數據選取和處理上是非常容易實現的。

2 仿真實驗

為驗證算法在絕對重力儀數據處理中的可行性，文中基于MATLAB軟件分別針對最小二乘，加權最小二乘遞推，改進最小二乘遞推和改進加權最小二乘遞推算法這4種情況進行了模擬仿真。對文獻[7]中的掃頻正弦信號模擬干涉條紋進行雙采樣提取時間距離對，往信號中添加18%的幅值噪聲。

圖2給出了這4種算法的g值仿真結果，其中直線代表真實值g=9.8 m/s2，曲線代表估計值，用上述4種算法及其抗差最小二乘得到的數據處理結果見表1。

表1 4種算法及抗差加權仿真結果

如圖2和表1可以看出，加權后的算法精度和準度比未加權的高出兩個量級，同時改進后的算法比未改進的算法更快速的收斂于真實值，即改進算法只需要相對少量數據便可得到高精度的估計值。在干擾嚴重時，無論是抗差還是加權遞推都有效地提高了精度，當只加入幅值白噪聲時，抗差效果比最小二乘提高不大，但是加權遞推依舊提高了精度，這充分的說明了本文提出的算法能夠滿足準確性。

圖2 4種算法參數估計收斂圖

3 結束語

文中借鑒文獻[4]中抗差的思想和迭代流程，就算法的改進提出了加權遞推最小二乘算法，采用平差準則作為加權陣，經過改進算得到的估計g值不僅具有抗差參數估計的優越性，還具有運算簡單方便，易于實現、移植的特征。結合計算機運算特點，改進后的算法通過仿真實驗驗證，具有快速收斂特性，不依賴初始值設定的特點，經過多次實驗，改進算法非常穩定，具有工程應用的可行性。該算法還可以在參數細節上進行修改，以達到實際絕對重力儀數據處理中最佳的估計值。

[1]Svitlov S.Frequency-domain analysis of absolute gravimeters[J].Metrologia，2012，49（6）:706-726.

[2]吳瓊.高精度絕對重力儀關鍵技術研究[D].北京：中國地震局地球物理研究所，2011.

[3]Svitlov S，Maslyk P，Rothleithner C，et al.Comparison of Three Digital Fringe Signal Processing Methods in a Ballistic Free-Fall Absolute Gravimeter[J].Metrologia，2010，47（6）:677-689.

[4]胡明，張為民.抗差估計在自由落體式絕對重力儀中的應用[J].大地測量與地球動力學，2016，9（9）:822-825.

[5]劉謝進.遞推加權最小二乘算法的研究[J].系統仿真學報，2007，21（14）:4248-4250.

[6]吳書清，吉望西，劉達倫，等.絕對重力儀數值計算及不確定度評定方法[J].計量學報，2009，30（3）：212—215.

[7]楊厚麗，鄒彤，郭唐永，等.復外差法處理絕對重力儀數字干涉信號[J].大地測量與地球動力學，2015，4（2）:346-348.

[8]周姍姍，柴金廣，李丹.實時遞推的最小二乘預測跟蹤算法[J].電子設計工程，2011，19（10）：53-55.

[9]強明輝，張京娥.基于MATLAB的遞推最小二乘法識別與仿真[J].自動化與儀器儀表，2008（6）:4-5，39.

[10]賈沛璋.誤差分析與數據處理[M].北京：國防工業出版社，1992.

[11]秦延，陳宗海，李衍杰.遞推最小二乘算法的補充證明[J].系統仿真學報，2004，16（10）:2159-2160.

[12]Orlob M，Braun A.Impact Estimation and Filtering ofDisturbancesin FG5 AbsoluteGravimeter Observations[J].International Journal of Geosciences，2013，4（2）:302-308.

[13]鄧自立，王欣，高媛.建模與估計[M].北京：科學出版社，2007.

[14]邵建新.最小二乘法線性擬合中參數的確定問題[J].大學物理，2003，22（1）.

[15]卓金武.MATLAB在數學建模中的應用[M].北京：北京航空航天大學出版社，2014.

Application of weighted recursive least square to data processing of absolute gravimeter

LUO Qi，PANG Cong，LI Cha?wei
（Institute of Seismology，CEA，Wuhan430071，China）

In view of the high precision demand of absolute gravimeter，the weighted recursive algorithm and its improved algorithm are proposed based on the robust parameter estimation，weighted recursive algorithm and improved robust weighting factor and inverse matrix for complex multidimensional.Simulation results show that the proposed algorithm has a good effect on the white noise or abnormal data and improves the accuracy of gravity measurement.

absolute gravimeter；robust estimation；weighted recursive；least square

P315

A

1674-6236（2017）22-0065-04

2016-09-21稿件編號：201609188

羅棋（1990—），男，湖北天門人，碩士。研究方向：絕對重力儀落體控制。