戴 卿，許輝熙，馬云峰

(1.西南交通大學(xué) 地球科學(xué)與環(huán)境工程學(xué)院，成都 610031；2.四川建筑職業(yè)技術(shù)學(xué)院 博士后創(chuàng)新實踐基地，四川 德陽 618000；3. 重慶水利電力職業(yè)技術(shù)學(xué)院 建筑工程系，重慶 402100)

0 引言

以GNSS(Global Navigation Satellite System，全球?qū)Ш叫l(wèi)星系統(tǒng))/SINS(Strapdown Inertial Navigation System，捷聯(lián)慣導(dǎo)系統(tǒng))緊耦合系統(tǒng)為研究熱點的現(xiàn)代導(dǎo)航技術(shù)，因其具有靈活性大和動態(tài)導(dǎo)航定位精度高的特點，近年來受到越來越多測量學(xué)家的關(guān)注。然而導(dǎo)航系統(tǒng)在實際工作中受衛(wèi)星可見性和外界干擾等因素影響，量測噪聲隨機模型的統(tǒng)計特性隨觀測條件會發(fā)生改變，如不加以考慮就會引起系統(tǒng)穩(wěn)定性和精度的降低[1]。為此，利用函數(shù)模型和隨機模型作為補償手段的自適應(yīng)濾波理論相繼被提出[2]。Sage-Husa濾波能夠?qū)崟r獲取當前歷元量測噪聲的統(tǒng)計特性，但該方法要求各歷元量測信息同類、同維、同分布，且當運動載體產(chǎn)生大的擾動時，會影響估計效果[3]。漸消濾波通過引入一個漸消因子使濾波算法符合最優(yōu)估計條件，但這種方法只能用于時變過程噪聲，不能處理量測噪聲[4]。抗差自適應(yīng)濾波可用于量測噪聲發(fā)生擾動或異常的情況，通過調(diào)整自適應(yīng)因子和等價權(quán)矩陣因子達到穩(wěn)定狀態(tài)估值的目的，但該算法需要有冗余的量測值[5]。為更好地解決濾波中量測噪聲協(xié)方差未知或不準確的問題，有學(xué)者提出包含量測噪聲的最優(yōu)Bayesian濾波，但該算法在實際應(yīng)用中卻面臨被積函數(shù)復(fù)雜、變量維數(shù)高的問題，且噪聲方差陣的引入使求解更加困難，從而影響其推廣使用[6]。研究人員基于機器學(xué)習(xí)的思想提出了變分Bayesian學(xué)習(xí)(Variational Bayesian, VB)，即通過近似的方法求得次優(yōu)解，以較少的運算量進行較高精度的參數(shù)估計，解決Bayesian估計在實際應(yīng)用中被積函數(shù)計算復(fù)雜的問題，具有收斂速度快的特點，適應(yīng)于高維觀測環(huán)境，可有效提升估計速度[7]。變分Bayesian學(xué)習(xí)在隨機系統(tǒng)參數(shù)估計、盲源分離、語音識別和獨立成分分析等領(lǐng)域取得了一定的應(yīng)用和進展[8]。

受上述研究工作啟發(fā)，本文考慮GNSS/SINS緊耦合組合導(dǎo)航系統(tǒng)中衛(wèi)星量測噪聲的時變特性，對濾波隨機模型做進一步精化處理，在變分Bayesian理論框架下提出了一種基于逆Gamma分布的變分自適應(yīng)濾波算法，利用變分Bayesian理論來解決量測噪聲協(xié)方差未知或不準確時的參數(shù)自適應(yīng)估計問題，并通過數(shù)據(jù)測試驗證了該算法的魯棒性、實時性和準確性，為其今后在組合導(dǎo)航系統(tǒng)及其擴展應(yīng)用提供了一定的理論支持。

1 問題描述

最優(yōu)Bayesian濾波框架推導(dǎo)的卡爾曼濾波，需要以系統(tǒng)模型準確已知為前提條件。當過程噪聲已知，量測噪聲未知時，可導(dǎo)出時刻包含量測噪聲的最優(yōu)Baysian濾波形式為

p(xk,Rk|xk-1,Rk-1)=p(xk|xk-1)p(Rk|Rk-1) (1)

式中：p(·)代表似然概率；x為狀態(tài)量；R為量測噪聲方差。設(shè)量測量Yk-1={y1,y2，…,yk-1}，由此得到時間更新方程為

量測更新方程分別為

式(1)～(3)構(gòu)成了量測噪聲方差未知的最優(yōu)Bayesian濾波過程。分析發(fā)現(xiàn)上述計算中存在多維數(shù)值積分的問題，且量測噪聲方差的引入使得式(1)～(3)難以獲得解析解。

2 基于逆Gamma分布的變分自適應(yīng)濾波

2.1 算法原理

當量測噪聲難以準確描述時，可用逆Gamma分布(Inverse Gamma, IG)來模擬協(xié)方差，并進行近似求解[9]。為此，將式(2)中的量測噪聲通過逆Gamma分布近似

由于式(3)中p(yk|xk,Rk)噪聲方差與量測噪聲具有相關(guān)性，難以獲得解析形式，故將其改寫為

p(xk,Rk|Yk)≈Qx(xk)QR(Rk) (5)

并用KL(Kullback-Leibler)散度進行度量，得

將上式分別對Qx(xk)和QR(Rk)求最小極值

進一步求得量測更新方程

假設(shè)存在線性化模型

xk=Fk-1xk-1+wk-1(10)

yk=Hkxk+vk(11)

式中：Fk-1為狀態(tài)轉(zhuǎn)移陣；wk-1為過程噪聲；Hk-1為量測矩陣；vk-1為量測噪聲。

式中：i=1,2,…,l，l為量測方差維數(shù)；p∈(0,1]為預(yù)測加權(quán)系數(shù)，p=1時，噪聲方差為固定值，p越小方差參數(shù)和α與β前一時刻值的關(guān)系越弱。

上式中νi難以從量測量中求得，故根據(jù)式(11)得到k時刻的量測噪聲精確表達式為

vk=yk-Hkxk(14)

式(14)說明vk與xk有關(guān)，但xk只能由卡爾曼濾波的中間變量近似得到。

綜上，得到基于逆Gamma分布的自適應(yīng)卡爾曼濾波遞推過程，即

時間更新階段

參數(shù)迭代更新

量測更新階

式中：為量測噪聲方差陣的估值；αk=[αk,1,αk,2,…,αk,l]T；βk=[βk,1,βk,2,…,βk,l]T；ρ=[ρ1,ρ2,…,ρl]T。計算迭代次數(shù)≤3次。

2.2 性能分析

借鑒文獻[5]中基于殘差的自適應(yīng)估計(Residual Adaptive Estimation, RAE)算法，量測噪聲方差估計采用下式表述

式中，N=k-i0+1為殘差序列總數(shù)。與本文基于逆Gamma分布的算法相比較，兩種算法的自適應(yīng)調(diào)整均依賴于殘差序列和量測協(xié)方差陣，但本文討論的算法避免了RAE算法中滑動窗中新息序列存儲及求和運算，其k時刻估值僅與k-1時刻有關(guān)，易于迭代實現(xiàn)。另一方面，基于逆Gamma分布的算法在噪聲方差估計中只求解對角元素，并利用ρ實現(xiàn)對獨立變化噪聲的跟蹤調(diào)節(jié)，較RAE算法更為靈活。因此，兩種算法在估計精度上接近，但在實用性上本文所研究的優(yōu)化算法更有優(yōu)勢。

3 數(shù)據(jù)測試

為驗證算法有效性，本文采用MATLAB R2014a 編寫的GNSS/SINS緊耦合仿真平臺進行數(shù)據(jù)測試。設(shè)置部分仿真參數(shù)如下：運動載體初始位置經(jīng)度為125°，緯度為50°，高度為300 m；初始失準角為5°5′20″；陀螺儀零偏為20°/h，隨機游走為0.07°/(h·Hz1/2)；加速度計零偏為50 mg，隨機游走為6 mg/(h·Hz1/2)；機動時間為1000 s；仿真可見衛(wèi)星為5顆。借鑒文獻[5]的思路，雖然實際觀測中每顆衛(wèi)星受到的干擾程度不同，但受制于相同無線電信號下量測噪聲方差變化趨勢相同，因此為考察算法對衛(wèi)星量測噪聲方差變化的跟蹤能力，對PRN-8號衛(wèi)星偽距的量測噪聲在600 s至700 s設(shè)置了協(xié)方差劇烈變動，如圖1所示。采用如下兩種不同算法方案進行性能驗證：

方案1：基于殘差的自適應(yīng)估計濾波算法；

方案2：基于逆Gamma分布的變分自適應(yīng)濾波算法。

圖1～圖4所示為數(shù)據(jù)測試結(jié)果。當量測噪聲協(xié)方差估計不準確或噪聲方差發(fā)生突變時，圖1中方案1和方案2的自適應(yīng)估計策略均能較好的跟蹤方差變化。圖2～圖4所示兩方案濾波結(jié)果仍能保持收斂狀態(tài)。由此說明，方案1和方案2能有效辨識噪聲變化并加以利用，維持濾波的穩(wěn)定性和魯棒性。

圖1 PRN-8偽距噪聲方差估計

圖2 姿態(tài)誤差曲線

圖3 速度誤差曲線

圖4 位置誤差曲線

在實驗條件不變的情況下，進行了50次蒙特卡羅仿真測試，用姿態(tài)、速度和位置的均方根誤差(Root Mean Square Error，RMSE)來衡量精度，用單位歷元耗時衡量計算效率，并引入擴展卡爾曼濾波(Extended Kalman Filter, EKF)作為對比項，結(jié)果如表1所示。EKF算法由于對量測噪聲建模不準確導(dǎo)致算法估計精度急劇下降，而方案1和方案2估計精度優(yōu)勢明顯。方案1計算負擔(dān)偏大，這是因為滑動窗中新息序列存儲及求和運算占用了計算資源，方案2在噪聲方差估計中只求解對角元素，通過可實現(xiàn)對獨立變化的噪聲進行跟蹤調(diào)節(jié)，計算量相對減少。

表1 不同算法估計性能比較(RMSE)

為了進一步討論算法性能，引入另一評價指標即新息加權(quán)指數(shù)

當值接χ2近于1，濾波對噪聲估計準確，濾波趨向最優(yōu)。χ2不僅能反應(yīng)濾波對噪聲的估計精度，還能體現(xiàn)對噪聲跟蹤能力的好壞。圖5中分別設(shè)置方案1中的ρ值和方案2中的N值，發(fā)現(xiàn)ρ和N的取值與χ2參數(shù)值相關(guān)，不合適的取值會直接影響自適應(yīng)濾波的跟蹤性能，甚至引起濾波發(fā)散。

綜上所述，時變噪聲環(huán)境下的GNSS/SINS緊耦合導(dǎo)航系統(tǒng)中，本文討論的優(yōu)化算法能夠保證濾波穩(wěn)定性和精度，且更適合工程計算應(yīng)用。

4 結(jié)論

本文討論了一種基于逆Gamma分布的變分自適應(yīng)濾波算法，該算法顧及量測噪聲的時變性對獨立變化的量測噪聲進行跟蹤調(diào)節(jié)，精化了系統(tǒng)隨機模型，提高了算法的魯棒性。數(shù)據(jù)測試表明該方法能有效實現(xiàn)噪聲量測方差不準確或突變時的參數(shù)自適應(yīng)估計，擁有較高的濾波精度，且計算負擔(dān)小易于編程實現(xiàn)。將該算法用于含時變觀測噪聲的緊組合導(dǎo)航系統(tǒng)中，可為導(dǎo)航工程研究提供了一定的參考價值。