蔡 奇 劉 寧 戴吾蛟 楊德地

1 中南大學地球科學與信息物理學院，長沙市麓山南路932號，410083

可靠的Kalman濾波算法是基于可靠的函數(shù)模型和隨機模型［1］，當系統(tǒng)狀態(tài)模型不準確或觀測誤差中含有粗差時，Kalman 濾波算法性能將大為降低甚至導致濾波發(fā)散。這種情況下，一般將抗差估計與Kalman 濾波相結(jié)合，建立抗差Kalman濾波算法進行濾波估計［2］。

抗差估計的設計一般通過建立等價權(quán)函數(shù)來實現(xiàn)［3］?，F(xiàn)有的等價權(quán)主要有Turkey、Huber、Hampel、IGGⅠ、IGGⅢ等［4－5］，其共同特點為分段函數(shù)，涉及分段降權(quán)臨界值的選取問題，而臨界值的選取一般根據(jù)誤差檢驗的置信水平和經(jīng)驗來確定［6］。為克服憑經(jīng)驗確定統(tǒng)計量臨界值所帶來的抗差效果風險，歐吉坤［6］提出一種將抗差初值、綜合抗差以及優(yōu)化抗差效率結(jié)合起來的三步抗差方案，楊元喜等［7］則利用學生化殘差統(tǒng)計量構(gòu)建一種臨界值可變的抗差估計等價權(quán)函數(shù)。而對于不同的抗差環(huán)境，除臨界值可變外，等價權(quán)的函數(shù)形態(tài)也應不同。為此，本文從粗差檢驗統(tǒng)計量的分布規(guī)律出發(fā)，提出一種等價權(quán)函數(shù)為連續(xù)函數(shù)且函數(shù)形態(tài)可變（continuous variable，簡稱CV）的抗差估計方法，其抗差性能優(yōu)于傳統(tǒng)模型。

1 連續(xù)可變等價權(quán)函數(shù)設計與分析

1.1 連續(xù)可變等價權(quán)函數(shù)的設計

抗差Kalman濾波算法的關(guān)鍵是抗差等價權(quán)函數(shù)的設計以及定位參數(shù)、尺度參數(shù)的確定［8］。

等價權(quán)函數(shù)的函數(shù)形態(tài)應含有3個階段——高端平緩區(qū)、低端平緩區(qū)以及介于兩者之間的下降收斂區(qū)，即“反S型”曲線。據(jù)此，設計連續(xù)可變等價權(quán)函數(shù)模型：

其函數(shù)形態(tài)為“反S型”，A、B、C、D為待定參數(shù)，S為形態(tài)參數(shù)。

假定該函數(shù)有3 個特征點，當x＝k0時，f（x）＝1，使函數(shù)從1 開始下降；當x＝a時，f（x）＝1／2，用來控制其收斂速度；當x＝k1時，f（x）＝0，使函數(shù)收斂為0。并且設a∝S，由此可以聯(lián)立以下3個方程：

解得：

由標準正態(tài)分布的性質(zhì)可知，k0＝0，k1＝∞。將k0＝0與式（3）代入式（1），得：

又當k1＝∞時，

在條件一定的情況下，S視為常數(shù)，那么當x增大到一定程度時，a近似于x的線性函數(shù)，使得函數(shù)形態(tài)向右移，但同時可能會出現(xiàn)x增大時，f（x）函數(shù)尾部向0收斂過慢。

為解決函數(shù)尾部收斂過慢問題，現(xiàn)對函數(shù)進行一定的修改。假設系統(tǒng)的模型誤差為υm，當a＞υm時，因為已經(jīng)超出模型誤差改正的范圍，權(quán)函數(shù)應立即收斂為0，而由于上述情況的發(fā)生，有可能導致權(quán)函數(shù)并未收斂。所以，現(xiàn)構(gòu)造f（x）k函數(shù)，通過k值來控制收斂。

由于0＜f（x）＜1，當x＞υm時，應使k＞1來控制f（x）k值收斂下降；當x＜υm時，應使0＜k＜1來控制f（x）k值不收斂。濾波過程中υm未知，在實際數(shù)據(jù)中隨著υm模型誤差增大，導致粗差檢驗統(tǒng)計量右移，因此抗差過程中，需要將抗差降權(quán)區(qū)域右移。由于f（x）中S形態(tài)參數(shù)越大，對應的函數(shù)形態(tài)最初高端平緩區(qū)越長，即抗差區(qū)域右移，因此υm與S應具有正比關(guān)系。不妨設υm∝S2，取：

綜上所述，最終所得的函數(shù)為：

其中，a＝（x2＋S2）／（xS）。

將該函數(shù)模型應用到等價權(quán)函數(shù)中，得CV抗差方案模型：

1.2 等價權(quán)函數(shù)模型分析

1）該抗差模型是根據(jù)實際情況下的標準化殘差的分布以及形態(tài)系數(shù)來確定函數(shù)模型的形態(tài)，以此設定一個較好的函數(shù)形態(tài)來處理數(shù)據(jù)。通過調(diào)整函數(shù)形態(tài)系數(shù)S，可以靈活地處理不同情況下的數(shù)據(jù)。

2）從圖1可知，當S＝0.8時CV 函數(shù)模型與Turkey函數(shù)模型形態(tài)近似一致。從圖2可知，當S＝2時CV 函數(shù)模型與Huber函數(shù)模型形態(tài)近似一致。當0＜S＜1 時，適用于線性模型；當S＞1時，適用于具有系統(tǒng)誤差的非線性模型。隨著S的增大，CV函數(shù)模型可以由Turkey 到Huber轉(zhuǎn)換。假設S是時間t的函數(shù)，可以根據(jù)不同時刻的具體環(huán)境動態(tài)改變CV 函數(shù)形態(tài)，從而達到更好的抗差效果。

圖1 抗差方案權(quán)因子函數(shù)模型（CV形態(tài)參數(shù)S＝0.8）Fig.1 The model plot of robust weight function（S＝0.8）

圖2 抗差方案權(quán)因子函數(shù)模型（CV形態(tài)參數(shù)S＝2）Fig.2 The model plot of robust weight function（S＝2）

2 實驗分析

2.1 模擬實驗數(shù)據(jù)

為比較分析基于連續(xù)可變等價權(quán)函數(shù)的抗差Kalman濾波算法的濾波處理效果，先利用兩組模擬數(shù)據(jù)，采用以下5種方案進行對比分析：IGGⅢ法、IGGⅠ法、Turkey法、Huber法、連續(xù)可變等價權(quán)函數(shù)法。

1）第一組模擬數(shù)據(jù)。首先利用式（10）非線性運動模型計算得到5 000個觀測真值序列，然后在所有觀測值上加入中誤差為5的隨機誤差，并隨機選擇50個觀測值加入大粗差，作為觀測值：

2）第二組模擬數(shù)據(jù)。利用式（11）模型計算得到5 000個觀測真值序列，然后在所有觀測值上加入中誤差為0.5的隨機誤差，并隨機選擇50個觀測值加入大粗差，作為觀測值：

2.2 真實實驗數(shù)據(jù)

這里引入五強溪大壩引張線變形觀測數(shù)據(jù)——22個控制點的2 978個時刻的時序數(shù)據(jù)?？紤]到數(shù)據(jù)處理的時間效率，選擇1～5號控制點200個時刻的數(shù)據(jù)進行分析。在實際處理中，選擇一段數(shù)據(jù)比較平穩(wěn)的時段，人為在50、100、150時刻加入粗差污染。

2.3 實驗成果與分析

標準Kalman與CV 濾波殘差比較見圖3，非線性模擬數(shù)據(jù)處理（CV 形態(tài)參數(shù)為5）見圖4，線性模擬數(shù)據(jù)處理（CV 形態(tài)參數(shù)為0.5）見圖5。

1）由圖3可知，標準Kalman濾波抗差效果較差，當觀測信息存在較大粗差時，標準Kalman濾波值不可靠。

2）從表1可知，對于非線性數(shù)據(jù)，濾波時存在系統(tǒng)誤差，由圖4以及不同抗差方法RMS比較可以看出，在較小系統(tǒng)誤差情況下，不宜將觀測值的權(quán)迅速降低。IGGⅠ、Huber和CV（S＝5）在處理非線性數(shù)據(jù)上具有較好的效果，而Turkey 的處理效果不佳。從圖4（c）和表2可知，由于Turkey 整體降權(quán)導致過抗差現(xiàn)象，因此抗差濾波發(fā)散。

3）對于線性數(shù)據(jù)，根據(jù)圖5可知，Turkey和CV（S＝0.5）數(shù)據(jù)處理效果更好。由于模型建立準確，因此整體降權(quán)下的濾波效果更好。

4）本文提出的抗差方法可以根據(jù)實際，靈活地對函數(shù)形態(tài)進行調(diào)整，從而更好地適用于Turkey和Huber的濾波環(huán)境。而且不用過多地嘗試閾值的設置，在抗差濾波使用上較為方便，抗差濾波的結(jié)果穩(wěn)健性較好。

圖3 標準Kalman與CV 濾波殘差比較Fig.3 The residual plot of Kalman filtering ＆CV filtering

圖4 非線性模擬數(shù)據(jù)處理（CV 形態(tài)參數(shù)為5）Fig.4 Non－linear analog data calculation（S＝5）

圖5 線性模擬數(shù)據(jù)處理（CV 形態(tài)參數(shù)為0.5）Fig.5 Linera analog data calculation（S＝0.5）

表1 非線性數(shù)據(jù)與線性數(shù)據(jù)時各類模型之間的RMS值Tab.1 The RMS of linear and non－linear data for different models

表2 真實數(shù)據(jù)時各個濾波方案RMS值比較Tab.2 The comparison for RMS of different schemes under true data

3 結(jié) 語

針對目前抗差估計在Kalman濾波中的不同算法，總結(jié)了各個抗差方案的適用范圍，并通過分析抗差處理中存在的問題，提出一種基于連續(xù)可變等價權(quán)函數(shù)的抗差方案——CV 法。通過兩組模擬實驗和真實數(shù)據(jù)處理表明：

1）Turkey適用于線性數(shù)據(jù)，Huber適用于非線性數(shù) 據(jù)，而IGG Ⅰ、IGG Ⅲ和CV 法 抗 差 適 用范圍更廣；

2）從靈活度上來看，CV 法可以根據(jù)系統(tǒng)實際情況及時調(diào)整函數(shù)形態(tài)，從而進行較為可靠的抗差濾波。

［1］付夢印，鄧志紅，張繼偉.動態(tài)Kalman濾波理論及其在導航系統(tǒng)中的應用［M］.北京：科學出版社，2003（Fu Mengyin，Deng Zhihong，Zhang Jiwei.Thesis on Dynamic Kalman Filtering and Its Application in Navigation System［M］.Beijing：Science Press，2003）

［2］宋迎春.動態(tài)定位中的卡爾曼濾波研究［D］.長沙：中南大學，2006（Song Yingchun.Research on Kalman Filtering in Dynamic Positioning［D］.Changsha：Central South University，2006）

［3］李德仁，袁修孝.誤差處理與可靠性理論［M］.武漢：武漢大學出版社，2002（Li Deren，Yuan Xiuxiao.Error Processing and Reliability Theory［M］.Wuhan：Wuhan University Press，2002）

［4］楊元喜.抗差估計理論及其應用［M］.北京：測繪出版社，1993（Yang Yuanxi.Robust Estimation Principle and Its Application［M］.Beijing：Surveying and Mapping Press，1993）

［5］周江文.經(jīng)典誤差理論與抗差估計［J］.測繪學報，1989，18（2）：116－118（Zhou Jiangwen.Classical Error Theory and Robust Estimation［J］.Acta Geodaetica et Cartographica Sinica，1989，18（2）：116－118）

［6］歐吉坤.一種三步抗差方案的設計［J］.測繪學報，1996，25（3）：173－179（Ou Jikun.A Design on Three－Step Robust Estimation Project［J］.Acta Geodaetica et Cartographica Sinica，1996，25（3）：173－179）

［7］楊元喜，吳富梅.臨界值可變的抗差估計等價權(quán)函數(shù)［J］.測繪科學技術(shù) 學報，2006，23（5）：317－320（Yang Yuanxi，Wu Fumei.Robust Estimation Equivalent Weight Function of Variable Critical Value［J］.Journal of Geomatics Science and Technology，2006，23（5）：317－320）

［8］楊玲，沈云中，樓立志.基于中位參數(shù)初值的等價權(quán)抗差估計方法［J］.測繪學報，2011，40（1）：28－32（Yang Ling，Shen Yunzhong，Lou Lizhi.Method of Equivalent Weight Robust Estimation Based on the Initial Value of Neutral Position Parameters［J］.Acta Geodaetica et Cartographica Sinica，2011，40（1）：28－32）