蔡 奇 劉 寧 戴吾蛟 楊德地

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

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

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

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

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

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

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

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

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

解得：

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

又當(dāng)k1＝∞時，

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

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

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

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

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

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

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

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

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

圖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種方案進(jìn)行對比分析：IGGⅢ法、IGGⅠ法、Turkey法、Huber法、連續(xù)可變等價權(quán)函數(shù)法。

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

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

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

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

2.3 實驗成果與分析

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

1）由圖3可知，標(biāo)準(zhǔn)Kalman濾波抗差效果較差，當(dāng)觀測信息存在較大粗差時，標(biāo)準(zhǔn)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)導(dǎo)致過抗差現(xiàn)象，因此抗差濾波發(fā)散。

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

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

圖3 標(biāo)準(zhǔn)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)，從而進(jìn)行較為可靠的抗差濾波。

［1］付夢印，鄧志紅，張繼偉.動態(tài)Kalman濾波理論及其在導(dǎo)航系統(tǒng)中的應(yīng)用［M］.北京：科學(xué)出版社，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］.長沙：中南大學(xué)，2006（Song Yingchun.Research on Kalman Filtering in Dynamic Positioning［D］.Changsha：Central South University，2006）

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

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

［5］周江文.經(jīng)典誤差理論與抗差估計［J］.測繪學(xué)報，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］歐吉坤.一種三步抗差方案的設(shè)計［J］.測繪學(xué)報，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］.測繪科學(xué)技術(shù) 學(xué)報，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］.測繪學(xué)報，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）