仝海波,張國柱

國防科學技術大學電子科學與工程學院,湖南長沙 410073

改進M估計的抗多個粗差定位解算方法

仝海波,張國柱

國防科學技術大學電子科學與工程學院,湖南長沙 410073

隨著導航衛(wèi)星數量的增多,觀測數據中出現多個粗差的概率顯著增大,基于單個粗差假設的RAIM算法不能保證多個粗差的有效抑制。抗差估計在定位可靠性要求高的場合受到了廣泛關注。針對傳統M估計受初值誤差影響的問題,提出一種基于改進M估計的抗差定位解算方法。該算法采用S估計方法計算初值,根據可用衛(wèi)星數實時調整S估計中的參數使得初值能夠最大限度抑制粗差。GPS觀測數據處理結果表明改進的M估計能夠有效抑制多個粗差。

全球衛(wèi)星導航系統;定位;抗差估計;粗差探測;M估計;S估計

1 引 言

隨著北斗、伽利略以及GPS現代化的推進,未來的用戶機將擁有更多、更精確的偽距觀測值。由于可見星的增多,偽距中粗差發(fā)生的概率增大,但常用的最小二乘定位結果容易受到粗差的影響。一種簡單且成熟的方法是采用RAIM (receiver autonomous integrity monitoring)算法進行粗差的探測和排除。由于RAIM算法設計之初是為了處理GPS應用中可能出現的單個較大的粗差,所以在處理多個粗差以及單個較小的粗差時暴露出諸多不足[1]。同時,以M估計為代表的抗差(穩(wěn)健)估計[1-4]在多衛(wèi)星導航系統組合定位中的應用得到了越來越多的關注。

目前,基于M估計的衛(wèi)星導航定位算法研究主要集中在兩個方面:一方面是針對特定的觀測噪聲模型,構建合適的目標函數對單歷元的觀測數據進行處理。當觀測噪聲互相關時,文獻[5—6]提出了IGGIII算法以及改進后的雙因子方法抑制粗差,其采用三段式ρ函數與傳統M估計的Huber函數不同。類似的,通過對Danish函數的改進,文獻[7]給出了更一般的表達式來應對相關的觀測噪聲。當觀測噪聲為理想高斯模型時,經典的Huber函數可有效地抵制單個粗差,基于牛頓迭代的快速計算方法[8]使其工程應用成為可能。當觀測噪聲為混合模型時,文獻[9]提出了一種自適應的抗差估計方法,通過誤差統計量的假設檢驗,自適應調整ρ函數。另一方面充分利用多個歷元的觀測數據,將M估計應用于時域濾波過程中以提高抗差性,文獻[10]解決了抗差估計在卡爾曼濾波應用中的秩虧問題。自適應的抗差卡爾曼濾波[11]不僅能有效抑制粗差,而且適應觀測噪聲模型的變化,進一步優(yōu)化了定位結果。

以上兩個方面的研究和探討均以迭代計算的初值沒有嚴重偏離真值為前提條件。當粗差污染較嚴重且接收機運動狀態(tài)沒有先驗信息時,傳統的最小二乘初值計算結果作為初值時會含較大的偏差,從而使M估計抗差性受影響[12-13]。通過預測殘差的方法[14],雖然可以改善卡爾曼濾波的抗差性,但對用戶的動態(tài)范圍有一定的限制且不適用于單歷元的抗粗差定位。通過改進S估計的目標函數,解決了多個粗差引起的初值偏離問題,使得基于M估計的定位算法能夠有效抑制多個粗差。

2 傳統M估計抗差定位的局限性

為減少觀測數據的粗差對定位結果的影響,可采用抗差M估計。該估計方法與最小二乘估計不同之處主要是目標函數的選取。下面簡單介紹基于M估計的定位原理,然后分析與最小二乘估計的不同之處以及該估計在抗差中的局限性。

式中,Δy是n×1矢量,其元素為各衛(wèi)星觀測偽距與線性化點預測的偽距之差,n是可見衛(wèi)星數; Δx是n×1矢量,其元素是相對于線性化點的偏差;H是Δx和Δy的關聯矩陣,也稱觀測矩陣;ε是n×1的測量誤差矢量,其元素一般假設為相互獨立的高斯隨機過程。

為求解線性方程(1),常見的最小二乘估計的目標函數如下

當測量噪聲為高斯分布時,最小二乘估計為最優(yōu)估計,但實際觀測量中常常含有粗差,由式(2)和式(4)可知,粗差會直接影響定位結果,從而使得觀測矩陣H偏離真值,導致粗差向偽距殘差投影時產生偏差。因此,采用最小二乘準則解算后,偽距殘差的大小不能正確反映粗差的大小,整個定位結果會受粗差影響而遠遠偏離真值。

式中,函數ρ(·)為一階可導的偶函數,由M估計根據具體情況可設計不同的加權函數;^σ是Δy方差的估計值,通常由下式計算[9]式中,med(·)表示取向量元素的中值。

M估計主要通過選擇合適的函數ρ,自動為含粗差觀測值分配小的加權值,從而達到抑制粗差的目的。根據觀測噪聲模型的不同,國內外學者對多種有效的ρ函數進行了深入的研究[15-16]。從抑制幅值較大的粗差方面考慮,本文采用文獻[16]中的雙邊加權函數

式中,K是根據M估計方差性能確定的常數。此時,M估計結果的迭代式如下

式中,和的上標均為迭代次數,在實際計算過程中,為避免加權引起秩虧,上式中0常以較小的實數代替。由于式(3)和式(8)相似的特點,M估計又可看做是選擇權迭代最小二乘估計。文中稱該估計為傳統M估計,其計算流程總結如下:

(1)采用最小二乘計算結果作為初值。

(2)根據式(4)計算偽距殘差矢量并估計方差,再由式(9)得到等價權矩陣。

(3)由M估計的迭代式(8)計算第k＋1次估計結果,其中k＝0,1,2,…。

由上述計算流程可知,傳統基于M估計的抗差前提條件是初值不能偏離真值太遠。由于偽距觀測值與用戶位置、鐘差之間是非線性關系,而H矩陣是觀測方程在初值點展開的雅可比矩陣。當粗差的幅值較大時,最小二乘估計結果會出現遠遠偏離真值的情況,此時矩陣H不再是偽距和用戶位置線性關系的有效近似。受到衛(wèi)星幾何結構的影響,每個偽距測量值精度對定位結果的影響各不相同。當粗差存在于影響最大的偽距上時,定位結果會產生明顯偏移,這種情況的定量分析可參考接收機自主完好性算法。由此可見,當初值遠遠偏離真值時,偽距殘差不能有效反映偽距測量精度,那么基于偽距殘差加權的抗差M估計性能將會受到嚴重影響。

3 改進后的M估計定位解算方法

由前一節(jié)的分析可知,用戶位置和鐘差的初值誤差會直接影響到抗差定位的性能。當可用的觀測數據量趨于無窮大時,文獻[17]給出了一種抗差估計方法即S估計,該估計方法能有效抑制多個粗差。S估計的目標函數如下

在樣本數趨于無窮多時,b通常取0．5,相應的雙邊加權函數ρ(t)中的K＝1．547,而此時S估計的效率約為0．287?？梢奡估計抗差性能的提高是以方差增大為代價的。本文將利用S估計的這個特點,把其估計結果作為M估計迭代開始的初值,從而進一步提高定位結果的精度和可靠性。

GPS用戶通?？梢娦l(wèi)星數為7～12顆。當偽距觀測值數量有限時,S估計的抗差性能會有所下降[18],崩潰點(breakdown point,BP)描述了該估計所能容忍的含粗差樣本比例上限[19]。樣本數為n時,S估計的BP可通過式(12)計算

式中,符號[·]表示取整數;p為待估參數的數量。在本文中,n表示可用的衛(wèi)星數,待估參數為接收機的X、Y、Z坐標和鐘差,即p＝4。當含粗差觀測值的數量增多到滿足式(12)時,S估計的抗差性將無法保證。那么S估計用于定位初值解算時最多能抑制的粗差數m滿足

由式(13)可知,當n＝10時,m＝2;當n＝12時,m＝3。即10顆衛(wèi)星可用時最多能有效抑制2個粗差,12顆衛(wèi)星可用時最多能有效抑制3個粗差。隨著可用衛(wèi)星數的增多,抗差性能也不斷增強。

直接將S估計用于求解方程(1)得到定位結果不是高斯噪聲下的最優(yōu)估計,且定位精度下降明顯。比較式(11)和式(12),發(fā)現有限個觀測值無法達到接近0．5抗差極限,而此時雙邊加權函數的K仍取1．547不僅不能得到預期的抗差性能,而且使得估計結果方差過大。為了讓估計值的方差減小,同時不影響S估計的抗差性能,令b＝,展開后得到

式中,n為可見衛(wèi)星數;K為待確定的參數。通過數值計算方法,得到K與n的對應關系如表1所示。

表1 衛(wèi)星數n與K的關系表Tab．1 Look-up table ofnandK

在利用S估計計算初值前,將加權函數中的系數K根據可用觀測偽距的數量進行合理選擇,避免由于固定K引起的發(fā)散,從而在不降低抗差性能的前提下提高估計效率。

綜合以上分析,基于改進后M估計的定位解算流程可總結:首先,根據可用衛(wèi)星數n查表選擇K;其次,結合所選K,聯立式(7)、式(10)和式(11)求解抗差初值;最后,以抗差初值代替迭代最小二乘法得到的初值進行M估計的定位解算,基于M估計的迭代過程見上一節(jié)流程步驟(2)—步驟(4)。通過M估計的迭代計算將進一步減小估計結果的方差,提高定位精度。

4 實測GPS數據處理結果分析

GPS實測數據采用美國聯邦航空局提供的參考站觀測數據,文件為“Acv_EPak_1330_1616_ 06”,其觀測值更新頻率為1次/s,持續(xù)時間為1 h,可見衛(wèi)星數為10～12顆,圖2(c)中給出了可見衛(wèi)星數隨時間的變化情況。采用斯坦福大學開發(fā)的Matlab工具包SGMP[20]進行LS定位解算。目前,該工具包僅提供最小二乘的定位解算,將基于Huber函數的M估計[8]和本文改進后的抗差估計分別替換RAIM得到傳統的M估計和本文改進后的M估計定位結果。為比較3種定位算法的抗差性,設計了4個場景:一是不含粗差,即粗差為0;二是指定PRN3衛(wèi)星含粗差且粗差隨時間增大,最大幅值為300 m;三是指定PRN3和PRN6兩顆衛(wèi)星含粗差,兩個粗差隨機生成且互不相關,服從0到300的均勻分布;四是指定PRN3、PRN6和PRN9共3顆衛(wèi)星同時含有粗差且3個粗差之間相互獨立,服從場景三中的均勻分布。

將RAIM、傳統M估計和改進后的M估計3種抗差定位算法分別對粗差污染程度不同的4種場景數據進行處理。3種定位算法在不同場景下的三維位置(地心地固系)的定位精度比較詳見表2。從表中可以看出,在無粗差情況下,3種定位解算結果一致;改進前后的M估計結果與RAIM相比,Z軸方向精度均略有下降。隨著粗差個數的增加,RAIM算法定位誤差明顯增大。傳統M估計可有效抵制一個粗差,由于較小的粗差的存在,RAIM算法抑制粗差的性能受到了影響,定位精度有所下降。當兩個粗差出現時, RAIM算法的定位精度明顯下降,傳統M估計的Z軸方向方差明顯增大。與前兩種算法相比,無論是在單個還是多個粗差情況下改進后的M估計表現出更強的抗差性。

表2 定位精度的比較Tab．2 The comparison of the positioning RMS m

為了進一步驗證初值誤差對抗差性能的影響,將初值相對于真值的誤差進行了統計。表3給出了M估計算法改進前和改進后初值的三維誤差。從表3中可以看出,傳統M估計的初值容易受到粗差影響,隨著粗差個數增多,初值偏差顯著增大;而改進后的M估計初值具有穩(wěn)定的粗差抑制性能,受粗差數量變化影響較小。

表3 M估計算法改進前后的初值誤差比較Tab．3 Comparison of the initial errors m

圖1和圖2分別給出了傳統M估計和改進后M估計兩種算法在不同粗差個數下定位誤差隨時間變化的趨勢,其中定位精度用球誤差概率來衡量,圖示中的SEP95表示95%的球誤差概率,即以真值為圓心,SEP95為半徑的圓球,估算的三維位置誤差落在球內的概率是95%。在圖1中,傳統M估計的定位精度隨著粗差個數的增加而急劇惡化。相比之下,圖2中定位結果具有很強的抗差性。同時,在圖2(c)中出現了一段時間內定位精度較差的情況,這是因為該時間段歷元數據中可用衛(wèi)星數小于12,衛(wèi)星數的減少使得改進的M估計算法無法有效抑制3個粗差。

結合以上數據處理結果,分析總結如下:

(1)無粗差時,3種算法的定位精度相差不大,由于M估計是次優(yōu)估計,所以傳統的和改進后的M估計的位置誤差都略大于帶RAIM功能LS算法。兩種M估計均以精度的略微下降來換取抗差性能的提高。

(2)當有1顆衛(wèi)星含粗差時,傳統的和改進后的M估計的定位精度基本不受粗差影響。當粗差污染嚴重增至兩顆衛(wèi)星時,帶RAIM功能的定位精度變差;傳統M估計在少部分歷元上抑制了粗差,但是產生100 m以上偏差的概率明顯增加,定位精度無法滿足需要;改進后的M估計有效地抵制了兩個粗差,最大三維定位偏差仍在10 m以內,且精度與0和1個粗差時基本一致。

圖1 傳統M估計Fig．1 The traditional M-estimate

(3)隨著可見衛(wèi)星數的增多,改進后的M估計抗差性能增強。在含3個粗差的圖2(c)中,當衛(wèi)星數為12時,三維定位偏差在5 m以下;當衛(wèi)星數較少為11時,定位誤差達20 m左右,當衛(wèi)星數較少到10顆時,定位誤差將近300 m。同時該結果也初步驗證了式(13)的準確性:從式(13)可知當衛(wèi)星數大于等于12顆時,改進后的M估計才能有效抑制3個粗差。

4 結 論

在高斯噪聲假設下,基于LS的定位解算是近似最優(yōu)的定位算法,但是容易受到粗差的影響。傳統RAIM算法和M估計無法有效抵制多個粗差,改進后的抗差定位算法提高了M估計的抗差性能,實現了多個粗差的抑制,并且給出了可見衛(wèi)星數與所能抑制粗差數之間的關系。GPS數據處理結果驗證了該方法的有效性。如何檢驗粗差個數未超出上限以及提供類似RAIM的保護水平值得進一步深入研究。

圖2 改進后的M估計Fig．2 The modified M-estimate

[1] WANG J L,WANG J．Mitigating the Effect of Multiple Outliers on GNSS Navigation with M-Estimation Schemes [C]∥International Global Navigation Satellite Systems Society Symposium．Sydney:Menary Pty Limited,2007: 1-9．

[2] ANDERSEN R．Modern Methods for Robust Regression [M]．London:Corwin Press,2008:48-70．

[3] HUBER P J,RONCHETTI E M．Robust Statistics[M]．2nd ed,Hoboken:John Wiley&Sons Inc,2009:125-144．[4] KHODABANDEH A,AMIRI-SIMKOOEI A R．GPS Position Time-series Analysis Based on Asymptotic Normality of M-estimation[J]．Journal of Geodesy,2012,86(1): 15-33．

[5] YANG Y．Robust Estimation for Dependent Observations [J]．Manuscripta Geodaetica,1994,19(1):10-17．

[6] YANG Y,SONG L,XU T．Robust Estimator for Correlated Observations Based on Bifactor Equivalent Weights[J]．Journal of Geodesy,2002,76(6-7):353-358．

[7] WIESER A,FRITZ B K．Short Static GPS Sessions:Robust Estimation Results[J]．GPS Solutions,2002,5(3): 70-79．

[8] CHANG X W,GUO Y．Huber＇s M-estimation in Relative GPS Positioning:Computational Aspects[J]．Journal of Geodesy,2005,79(6-7):351-362．

[9] YANG Yuanxi．Adaptively Robust Least Squares Estimation [J]．Acta Geodaetica et Cartographica Sinica,1996,25 (3):206-211．(楊元喜．自適應抗差最小二乘估計[J]．測繪學報,1996,25(3):206-211．)

[10] KOCH K R,YANG Yuanxi．Robust Kalman Filter for Rank Deficient Observation Models[J]．Journal of Geodesy, 1998,72(7-8):436-441．

[11] YANG Y,HE H,XU G．Adaptively Robust Filtering for Kinematic Geodetic Positioning[J]．Journal of Geodesy, 2001,75(2-3):109-116．

[12] CROUX C,DHAENE G．Robust Standard Errors for Robust Estimators[J]．CES-Discussion Paper Series(DPS), 2004,16(3):1-20．

[13] TONG H,ZH ANG G,OU G．Iterative Reweighted Recursive Least Squares for Robust Positioning[J]．Electronics Letters,2012,48(13):789-791．

[14] CHEN Mingjian,ZHOU Fengqi,HAO Jinming．Improved Robust Filtering Based on Residuals in Real Time Single Point Positioning[J]．Journal of Chinese Inertial Technology, 2009,17(6):692-696．(陳明劍,周鳳岐,郝金明．在實時單點定位中基于殘差信息的改進抗差濾波方法[J]．中國慣性技術學報,2009,17(6):692-696．)

[15] LI Deren,YUAN Xiuxiao．Error Processing and Reliability Theory[M]．Wuhan:Wuhan University Press,2002．(李德仁,袁修孝．誤差處理與可靠性理論[M]．武漢:武漢大學出版社,2002．)

[16] KNIGHT N L,WANG J．A Comparison of Outlier Detection Procedures and Robust Estimation Methods in GPS Positioning[J]．The Journal of Navigation,2009,62(4): 699-709．

[17] SALIBIAN M,YOH AI V J．A Fast Algorithm for SRegression Estimates[J]．Journal of Computational and Graphical Statistics,2006,15(2):414-427．

[18] BERRENDERO J R,MENDES B M．On the Maximum Bias Functions of MM-estimates and Constrained M-estimates of Regression[J]．The Annals of Statistics, 2007,35(1):13-40．

[19] MARONNA R A．Robust Statistics Theory and Method [M]．Singapore:John Wiley&Sons,2006:58-62．

[20] DO Juyong．Stanford GPS/GNSS Matlab Platform[EB/ OL]．[2013-01-08]．http:∥waas．stanford．edu/sgmp/ sgmp．html．

(責任編輯:宋啟凡)

Robust Positioning Algorithm with Modified M-estimation for Multiple Outliers

TONG Haibo,ZHANG Guozhu
College of Electronic Science and Engineering,National University of Defense Technology,Changsha 410073,China

The possibility of multiple outliers should not be neglected in measurements for the increment number of the navigation satellites,and the RAIM based on the single-outlier hypothesis cannot provide effective restraining against multiple outliers．The robust estimation has gained attention widely．A robust positioning algorithm is presented for resisting multiple outliers which make the traditional M-estimator ineffective．The robustness can be improved by the modified positioning algorithm．To handle the biased convergence problem,the robust estimate of the initial values is realized by modifying the S-estimation with available satellite number in real time．Positioning results of the actual GPS measurements show that the proposed method resists multiple outliers effectively．

global navigation satellite system;positioning;robust estimation;fault detection;M-estimation;S-estimation

TONG Haibo(1984—),male,PhD candidate,majors in GNSS receiver autonomous integrity monitoring(RAIM),robust positioning algorithms, weak signal acquisition and tracking technology．

TN927．23

A

1001-1595(2014)04-0366-06

2013-03-13

仝海波(1984—),男,博士生,研究方向為衛(wèi)星導航接收機自主完好性監(jiān)測、抗差定位、弱信號的捕獲與跟蹤技術。

E-mail:hbo．tong＠gmail．com

TONG Haibo,ZHANG Guozhu．Robust Positioning Algorithm with Modified M-estimation for Multiple Outliers[J]．Acta Geodaetica et Cartographica Sinica,2014,43(4):366-371．(仝海波,張國柱．改進M估計的抗多個粗差定位解算方法[J]．測繪學報,2014, 43(4):366-371．)

10．13485/j．cnki．11-2089．2014．0055

中國第二代衛(wèi)星導航系統重大專項

修回日期:2013-12-13