非高斯環境下船變測量最大熵卡爾曼濾波方法

龍子旋, 周 琪, 彭俠夫, 張霄力,*

(1.廈門大學航空航天學院, 福建 廈門 361101;2.中國航空工業集團公司西安飛行自動控制研究所, 陜西 西安 710076)

0 引 言

慣性量匹配方法因其精度高、實時性好而在船體變形測量中得到廣泛應用,其本質就是建立船變測量過程的卡爾曼濾波方程。目前主流的慣性匹配算法包括角速度匹配法[1-3]、姿態陣匹配法[4-5]、姿態角匹配法[6]和姿態四元數匹配法[7-8]。

文獻[9]指出,當噪聲建模失準或測量系統受到未知噪聲干擾時,傳統卡爾曼濾波或無法準確跟蹤角形變。在船體變形測量中,為簡化建模過程,往往將量測噪聲理想化為模型準確已知的高斯白噪聲,對模型失準系統受到有色噪聲干擾下的濾波研究較少。文獻[10]將系統噪聲建模為白噪聲驅動的有色噪聲,然后修改了卡爾曼濾波器的迭代方程來完成非高斯濾波。文獻[11]提出了一種強跟蹤的最大互相關熵(簡稱為最大熵)卡爾曼濾波,即將自適應因子引入傳統的最大熵濾波器改進姿態陣匹配法,減小模型失配誤差。

針對噪聲統計特性建模失準問題,學者們提出了一種自適應卡爾曼濾波,其核心思想是,利用量測數據自動進行相關迭代矩陣的在線估計和修正。文獻[12-13]介紹了一種Sage-Husa自適應卡爾曼濾波,該濾波器引入單一的遺忘因子來減小歷史數據的影響,直接估計噪聲的一階和二階統計特性。此后,單自適應漸消因子的概念被提出,其本質是強迫新息協方差的理論值與估計值相等以穩定濾波過程[14]。文獻[15-18]引入多通道自適應漸消因子對狀態一步預測誤差矩陣進行實時更新調整以實現對多通道時變噪聲的自適應。此外,文獻[19-20]還引入了χ2檢驗來提高自適應漸消因子引入時機,抑制可能出現的狀態跳變,成功應用于捷聯慣性導航系統(strapdown inertial navigation system, SINS)初始對準問題,并推廣到非線性系統中[19]。

針對非高斯白噪聲的濾波問題,在信息論(information theoretical learning, ITL)中,有關學者提出了另一種最優理論,即最大熵準則(maximum correntropy criterion, MCC)[21-22]。其核心思想是,在卡爾曼濾波迭代過程中引入MCC構造核函數,以此作為卡爾曼濾波器的目標函數,重新求解濾波器的增益矩陣,充分利用量測值的高階信息,這便是用來處理非高斯噪聲的最大相關熵卡爾曼濾波器。文獻[23]在進行增益矩陣的求解時,近似認為狀態值與一步預測值相等,得到基于MCC的卡爾曼濾波迭代(MCC Kalman filter, MCKF)算法,由于其近似處理較多,故而其算法簡單但濾波器精度較差。文獻[24]通過改進近似處理過程并融合加權最小二乘(weighted least squares, WLS)理論修改高斯核函數得到性能更優的非高斯濾波器,其濾波收斂性更好、精度更高、穩定性更好。此外,基于擴展卡爾曼濾波器(extend Kalman filter, EKF)和無跡卡爾曼濾波器(unscented Kalman filter, UKF)的最大熵濾波器也相繼被提出以解決工程遇到的實際問題[25-26],將MCC推廣到非線性系統,進一步說明了MCC在非線性環境中同樣適用。

受相關論文啟發,本文介紹一種基于權重矩陣的MCCKF算法解決船體變形角估計中遇到的實際工程問題。

首先，介紹了最大熵概念,推導了傳統MCKF和改進MCCKF迭代公式。然后,引入兩種主流自適應濾波器:多重自適應漸消因子濾波器(簡稱為MAKF)和基于χ2檢驗的多重自適應漸消因子濾波器(簡稱為TAKF),并與MCKF一同作為對比算法，基于船體變形角估計應用場景,在兩種典型非高斯噪聲和兩種船體靜態變形建模情形下,針對4種濾波器開展多組對比實驗。此外,還探討了MAKF、TAKF實際使用受限問題。最后，通過實驗結果分析,驗證了本文介紹的濾波算法在收斂性、變形估計精度、實用性等綜合性能上相對于主流自適應和最大熵濾波器具備一定的優越性。

1 MCKF

1.1 傳統MCKF

MCC準則可充分利用信息的高階統計特性來衡量兩個隨機變量之間的相關性[27-28],在非高斯環境中的濾波有廣泛的應用[29-30]。最大相關熵濾波器源于對經典卡爾曼濾波迭代公式中的代價函數(又稱目標函數)進行替換繼而推導得到的。首先,采用如下形式的高斯核函數替代目標函數:

(1)

式(1)表明新的高斯核函數由兩個核函數組成,其中,

(2)

類似于標準卡爾曼濾波器的最小均方誤差(minimum square error, MSE)理論求解最優解的方式,求解式(2)的最大熵,令

(3)

得

(4)

即

(5)

為求狀態變量更新公式,通常在較短濾波更新周期內認為

(6)

則,式(5)中可有以下相關近似項:

(7)

式中:I為n維的單位矩陣。

聯立式(5)與式(7)得到狀態量更新公式:

(8)

式中：增益矩陣可表示為

(9)

(10)

(11)

(12)

最后求得改進后的MCKF狀態量更新公式:

(13)

式(13)還可以得到狀態預測的增益矩陣:

(14)

用式(13)和式(14)替換掉經典卡爾曼濾波狀態更新公式便可得到MCKF的迭代算法。算法流程如下所示。

算法1 MCKF算法濾波過程迭代公式初始條件x0=x0,P0=P0一步預測^xk|k-1=Φk|k-1^xk-1Pk|k-1=Φk|k-1Pk-1ΦTk|k-1+Γk-1Qk-1ΓTk-1濾波增益 Kk=G1(zk-Hk^xk|k-1)HTkI+G1(zk-Hk^xk|k-1)HTkHk濾波更新^xk=^xk|k-1+ KkεkPk=(I- KkHk)Pk|k-1

1.2 改進型MCCKF

(15)

類似式(3)方式求其最優解,得

(16)

將式(16)改寫為

(17)

其中,

(18)

類似式(11)的處理,得

(19)

最后求得MCCKF狀態量更新公式:

(20)

類似式(14),得到濾波過程的增益矩陣為

(21)

同樣, MCCKF迭代算法流程如下所示。

算法2 MCCKF算法濾波過程迭代公式初始條件x0=x0,P0=P0一步預測^xk|k-1=Φk|k-1^xk-1Pk|k-1=Φk|k-1Pk-1ΦTk|k-1+Γk-1Qk-1ΓTk-1濾波增益^Kk=GkHTkR-1kP-1k|k-1+GkHTkR-1kHk濾波更新^xk=^xk|k-1+^KkεkPk=(I-^KkHk)Pk|k-1

2 實驗與分析

目前角速度匹配法在慣性匹配方式體系中應用較廣、技術成熟,因此以下變形濾波驗證實驗亦基于此。另外,因試驗無法提供基準及真實數據,以下驗證實驗均以計算機仿真平臺為主。

2.1 角速度匹配法和實驗參數設置

船體變形角包括靜態變形角和動態變形角。靜態變形視為“準靜態”模型,用隨機游走模型近似代替其長周期的緩變過程:

(22)

動態變形θi視為白噪聲信號通過二階濾波器得到,其模型如下:

(23)

式中：i=x,y,z。光纖陀螺器件漂移模型分為兩類,即常值漂移εc和偏置穩定性漂移εr,前者視為常量,后者視作一階Markov過程:

(24)

式中：j=1,2。式(23)和式(24)相關參數含義如表1所示。

表1 船體變形模型關鍵參數

綜合以上變形模型、陀螺漂移模型,將角速度匹配方法的狀態方程組表示如下:

(25)

另外,量測方程為

(26)

(27)

其中,各個矩陣的含義由于篇幅有限這里不作展開。實驗給定的船變參數如表1所示。對搖擺譜的選擇需要符合船舶航行的海浪實況,實測數據表明搖擺激勵越小,觀測量越小,變形結果耦合的地速就越大,精度越低。表1中船體靜、動態變形角示意圖如圖1所示。

圖1 船體變形角示意圖

2.2 實驗與分析

本節通過兩種典型的非高斯噪聲來對比驗證MCKF、MAKF、TAKF、MCCKF等4種濾波器的船體變形角估計性能:沖擊噪聲和高斯混合噪聲,又稱強拖尾噪聲。此外,為確保角度估計長時間的收斂穩定性,整個濾波時間較長;且為避免角度估計結果的偶然性,所有實驗的RMSE值均采用50次蒙特卡羅獨立重復實驗計算其均值,并統一從濾波穩定之后開始計算(5 min)。

2.2.1 實驗一: 沖擊噪聲

假設量測量Z每個通道均受到如下非高斯分布的劇烈沖擊噪聲(單位:rad/s):

Δvshot=5×10-4,t≥100 s

(28)

該沖擊信號次數設置為1 000次,且隨機施加給系統。單通道噪聲示意圖如圖2所示,4種濾波器各軸的變形角跟蹤曲線和誤差曲線分別如圖3和圖4所示。

圖2 沖擊噪聲

圖3 變形角跟蹤曲線(沖擊)

圖4 變形角估計誤差(沖擊)

各軸變形角估計誤差的RMSE值如表2所示。

表2 船體變形估計RMSE(沖擊)

從圖3和圖4可以看出,在系統受到隨機的沖擊噪聲下,傳統的MCKF濾波效果最差,在規定時間內甚至無法收斂,這是因為長時間的濾波誤差積累導致發散;MAKF可以收斂,但由于自適應因子的引入時機可能有誤,導致個別通道(橫扭角)狀態估計誤差跳動范圍大,穩定性較差;TAKF解決了這一問題,但濾波開始3 min內角度跟蹤誤差曲線出現巨大跳變;而MCCKF在穩定性上顯著優于前者,但快速收斂之后抗隨機噪聲干擾的能力要稍弱于自適應濾波器。根據表2顯示,TAKF在精度上具有微弱優勢,穩定后的RMSE值要略大于前者0.4″,這在更為強調穩定性的工程上幾乎忽略不計。可以認為,在收斂精度滿足要求時,可以選擇MCCKF以期實現更穩定、更快的變形測量。

2.2.2 實驗二: 混合噪聲

高斯混合噪聲,即將兩種非零均值、非高斯分布的噪聲混合,同樣施加給量測信息(均值單位:rad/s;方差單位:(rad/s)2),表達式為

(29)

式中：N(·，·)為正態分布。為方便起見,其中相關參數設置下:

(30)

(31)

不同于實驗一隨機施加的沖擊噪聲,混合噪聲存在于整個實驗過程,其單通道示意圖如圖5所示,各濾波器的變形角跟蹤曲線如圖6所示。

圖5 混合噪聲

圖6 變形角跟蹤曲線(混合)

相應的誤差收斂曲線也給出,如圖7所示。

圖7 變形角估計誤差(混合)

收斂曲線表明自適應濾波器和MCCKF均可快速達到穩定狀態,但明顯可見,MCCKF誤差更小,波動閾值要小于MAKF和TAKF。

各軸變形估計的RMSE結果如表3所示。

表3 船體變形估計RMSE(混合)

混合噪聲擾動下,與實驗一不同的是,此時MCCKF的收斂精度要高于TAKF,這是因為混合噪聲在整個濾波過程均存在,自適應濾波器在這種情況下不具備自適應抗隨機擾動優勢,而MAKF不存在判斷失誤的情況因此并未出現誤差跳變的情形。從實驗結果可以看出,MCCKF在收斂精度和穩定性上也要比MAKF、TAKF優1″左右。此外,當逐步加大噪聲幅值,MCCKF的RMSE值要比TAKF小至數個角秒,表明強非高斯環境下,MCCKF具備更明顯的優勢。

2.2.3 實驗三: 常值靜態變形角

此外,針對船體靜態變形角建模為常值的情形,進一步驗證MCCKF濾波算法的魯棒性。靜態變形角大小均設置為0.1°,以沖擊噪聲為例,靜態變形角示意圖、變形角跟蹤曲線分別如圖8和圖9所示。

圖8 常值靜態變形角

圖9 變形角跟蹤曲線(常值靜態變形角+沖擊噪聲)

相應的誤差收斂曲線、變形角估計誤差的RMSE值分別如圖10和表4所示。

圖10 變形角估計誤差(常值靜態變形角+沖擊噪聲)

表4 船體變形估計RMSE

當靜態變形角建模為常值情形時,實驗結果再次驗證了實驗一的結論,但MCKF有收斂趨勢,延長濾波時間或可收斂。實驗表明本文介紹的MCCKF收斂性能始終優于MCKF,且對不同的變形測量模型也具有較好的魯棒性。

此外,當加大實驗三的沖擊噪聲頻率次數,MCCKF的RMSE值與TAKF差距越來越小,再次驗證MCCKF相對于自適應濾波器,對強非高斯環境濾波更加具備優勢。

最后,從算法本身來講,基于權重矩陣的MCCKF相較于兩種自適應濾波器MAKF、TAKF,迭代算法簡潔和直觀,無需設置判斷語句,更方便編程。并且,后者涉及自適應漸消因子縮放系數的選取,選取結果直接影響自適應濾波性能。考慮到最大熵濾波器對隨機外擾、強非高斯均具備不亞于自適應濾波器的性能,本文介紹的MCCKF通用性也更強。

3 結 論

本文介紹了一種基于最大熵準則和權重矩陣核函數的改進卡爾曼濾波器,并成功應用于非高斯白噪聲環境下的船體變形測量。重點介紹了基于最大熵的卡爾曼濾波器算法迭代公式,并對比兩種典型的自適應濾波器和傳統最大熵濾波器來驗證MCCKF濾波性能。各項實驗的結果證實,新算法在濾波精度和收斂速度上滿足相關收斂指標前提下,還可更好地解決噪聲模型失配和非高斯噪聲環境下的變形角估計問題,且在不同實驗條件下仍保持良好的魯棒性。此外,算法不僅可直接應用于復雜干擾、噪聲建模失準環境下的光纖陀螺組件船體變形測量問題,對于其他用到卡爾曼濾波算法的實際工程問題也具有一定的指導意義,如運載體跟蹤、導航和定位,組合導航的組合濾波算法等。