劉勝,牛鴻敏,張?zhí)m勇,郭曉杰
(哈爾濱工程大學(xué) 自動(dòng)化學(xué)院,黑龍江 哈爾濱 150001)
在工程設(shè)計(jì)中,非線性濾波問題普遍存在于信號(hào)處理、目標(biāo)跟蹤等領(lǐng)域,對(duì)于非線性系統(tǒng)而言,要得到精確的最優(yōu)解較為困難,需要精確已知系統(tǒng)狀態(tài)的后驗(yàn)狀態(tài)分布。傳統(tǒng)的線性濾波方法[1]不能滿足非線性系統(tǒng)的濾波要求,因此需要對(duì)濾波方法進(jìn)行改進(jìn)以滿足非線性系統(tǒng)濾波要求,提高濾波精度。一般常用的非線性濾波方法包括擴(kuò)展卡爾曼濾波,無跡卡爾曼濾波和粒子濾波等。
擴(kuò)展卡爾曼濾波結(jié)構(gòu)簡單,易于執(zhí)行,但是該算法采用一階泰勒級(jí)數(shù)對(duì)非線性方程進(jìn)行線性局部近似,因此對(duì)于階段誤差較為敏感,對(duì)于高階系統(tǒng)進(jìn)行線性近似會(huì)帶來較大的誤差[2-4]。無跡卡爾曼濾波采用UT變換,以確定性采樣策略逼近非線性系統(tǒng)的狀態(tài)后驗(yàn)分布,雖然能夠減少非線性產(chǎn)生的誤差影響[5-6],但其容易受到初始觀測(cè)誤差和參數(shù)變化的影響,當(dāng)系統(tǒng)維數(shù)較高時(shí),可能導(dǎo)致濾波發(fā)散等現(xiàn)象。粒子濾波以不確定采樣方法逼近狀態(tài)的條件概率密度函數(shù),能夠處理任意非線性和非高斯的估計(jì)問題,但粒子濾波實(shí)現(xiàn)過程中計(jì)算量大,不適用于反應(yīng)速度快的反饋控制系統(tǒng);并且隨著迭代次數(shù)的增加,容易產(chǎn)生粒子退化和局部最優(yōu)的問題。近年來,容積非線性卡爾曼濾波的發(fā)展受到廣泛的關(guān)注[7],為非線性系統(tǒng)的狀態(tài)估計(jì)提供新的實(shí)現(xiàn)方式。容積卡爾曼濾波使用容積數(shù)值積分原則計(jì)算非線性變換后的隨機(jī)變量的均值和協(xié)方差,采用一組等權(quán)值的容積點(diǎn)對(duì)最優(yōu)狀態(tài)的后驗(yàn)分布進(jìn)行逼近,相比于其他非線性卡爾曼濾波具有更高的非線性逼近性能和估計(jì)精度,并且其運(yùn)行速度快。
標(biāo)準(zhǔn)容積卡爾曼濾波 (cubature Kalman filter, CKF)算法設(shè)計(jì)需要已知精確的系統(tǒng)模型,然而實(shí)際情況下,噪聲為非零均值、非高斯噪聲,系統(tǒng)往往受到外界環(huán)境干擾以及計(jì)算機(jī)計(jì)算字長限制等,不可避免地存在非高斯非零均值的噪聲以及系統(tǒng)參數(shù)攝動(dòng)的情況,影響濾波的穩(wěn)定性和精度。因此,近年來一些學(xué)者提出了平方根容積卡爾曼濾波方法[8-9]來避免系統(tǒng)發(fā)散,然而其計(jì)算復(fù)雜,濾波速度受到限制。進(jìn)而,為了提高濾波系統(tǒng)的魯棒性,H∞濾波方法被提出[10],H∞濾波基于魯棒理論設(shè)計(jì),對(duì)于系統(tǒng)的攝動(dòng)和未知特性具有一定的魯棒性,近年來一些學(xué)者結(jié)合擴(kuò)展卡爾濾波、無跡卡爾曼濾波應(yīng)用到電機(jī)以及電力系統(tǒng)狀態(tài)估計(jì)等,取得了較好的效果[11-15]。
此外,在實(shí)際工程問題中系統(tǒng)的先驗(yàn)噪聲統(tǒng)計(jì)特性往往不能準(zhǔn)確獲得,如果直接應(yīng)用到CKF算法中,將引起較大的誤差。為了對(duì)未知噪聲特性進(jìn)行估計(jì),噪聲統(tǒng)計(jì)估值器和基于極大后驗(yàn)估計(jì)的方法得到應(yīng)用[16-18],然而其遞推迭代過程計(jì)算量大,并且估值器對(duì)噪聲方差的初始值敏感。模糊自適應(yīng)通過設(shè)計(jì)模糊規(guī)則,根據(jù)估計(jì)誤差自適應(yīng)噪聲統(tǒng)計(jì)特性,其對(duì)噪聲方差的初值不敏感,并且在存在較大估計(jì)誤差的情況下,也能實(shí)現(xiàn)較好的估計(jì)效果[19-20]。
因此,本文提出一種模糊自適應(yīng)H∞容積卡爾曼濾波方法(H∞Kalman filtering, HCKF),對(duì)系統(tǒng)噪聲和觀測(cè)噪聲進(jìn)行實(shí)時(shí)估計(jì),使得濾波算法對(duì)噪聲特性具有一定的自適應(yīng)能力,并且對(duì)于系統(tǒng)的模型攝動(dòng)和不確定噪聲具有魯棒性,提高濾波的收斂速度和穩(wěn)定性。
考慮如下非線性動(dòng)態(tài)系統(tǒng)狀態(tài)方程:
(1)
式中:xk∈Rn為系統(tǒng)的狀態(tài)估計(jì)向量;zk∈Rl為系統(tǒng)的觀測(cè)輸出向量;uk∈Rp為控制輸入向量,假設(shè)過程噪聲wk-1和量測(cè)噪聲vk相互獨(dú)立,且均值和方差滿足:
(2)
容積卡爾曼濾波算法將非線性濾波歸結(jié)為求解非線性函數(shù)與高斯概率密度乘積的積分問題,并采用一組等權(quán)值的2n個(gè)容積點(diǎn)對(duì)狀態(tài)后驗(yàn)概率密度進(jìn)行逼近。使用三階容積原則獲得的基本容積點(diǎn)和所對(duì)應(yīng)的權(quán)值:
(3)
(4)
實(shí)際生產(chǎn)應(yīng)用中,非線性系統(tǒng)存在參數(shù)攝動(dòng)、未知輸入、觀測(cè)噪聲和系統(tǒng)噪聲為非零均值非高斯等情況,從而影響容積卡爾曼濾波性能和穩(wěn)定性。與傳統(tǒng)的卡爾曼濾波最小化誤差均方差不同,H∞卡爾曼濾波是基于最小化最大估計(jì)誤差協(xié)方差進(jìn)行設(shè)計(jì)的。因此,采用H∞卡爾曼濾波與容積卡爾曼濾波相結(jié)合,提高系統(tǒng)對(duì)不確定擾動(dòng)、噪聲的魯棒性和估計(jì)精度。H∞容積卡爾曼濾波設(shè)計(jì)分為以下4個(gè)步驟。
1)時(shí)間更新。
進(jìn)行容積點(diǎn)的計(jì)算,通過非線性函數(shù)進(jìn)行容積點(diǎn)傳播,以及計(jì)算預(yù)測(cè)狀態(tài)向量和系統(tǒng)協(xié)方差矩陣。
初始化系統(tǒng)狀態(tài)及誤差協(xié)方差,系統(tǒng)噪聲均值和方差,測(cè)量噪聲均值和方差:
(5)
(6)
設(shè)期望獲取的觀測(cè)信息為:
yk=Ckxk
(7)
若Ck=I,則實(shí)際獲取的觀測(cè)量為全部狀態(tài)輸出量。根據(jù)式(4)產(chǎn)生基本容積點(diǎn),采用Cholesky方法對(duì)協(xié)方差進(jìn)行分解:
Sk-1|k-1=Chol(Pk-1|k-1)
(8)
計(jì)算容積點(diǎn):
(9)
計(jì)算通過狀態(tài)函數(shù)傳播的容積點(diǎn):
(10)
計(jì)算狀態(tài)估計(jì)值:
(11)
計(jì)算估計(jì)誤差協(xié)方差矩陣:
(12)
2)測(cè)量更新。
對(duì)容積點(diǎn)進(jìn)行更新,通過量測(cè)方程傳播容積點(diǎn),預(yù)測(cè)觀測(cè)向量及協(xié)方差。
首先對(duì)誤差協(xié)方差進(jìn)行分解:
(13)
計(jì)算更新容積點(diǎn):
(14)
計(jì)算通過量測(cè)方程傳播的容積點(diǎn):
(15)
進(jìn)行量測(cè)預(yù)測(cè):
i=1,2,…,L
(16)
計(jì)算新息的方差矩陣:
(17)
計(jì)算估計(jì)協(xié)方差矩陣:
(18)
3)為了便于進(jìn)行H∞濾波器設(shè)計(jì),將容積卡爾曼濾波通過統(tǒng)計(jì)線性化等效為線性回歸形式。通過將統(tǒng)計(jì)線性化應(yīng)用于狀態(tài)方程和測(cè)量方程,得到線性化近似形式。
證明:將非線性函數(shù)g=σ(y)在L個(gè)容積點(diǎn)(χi,ζi)處傳遞:
ζi=σ(χi),i=1,2,…,L
(19)
設(shè)非線性函數(shù)通過統(tǒng)計(jì)線性化近似為:
g=Ay+B+μ
(20)
(21)
μi=ζi-(Ayi+B)
(22)
對(duì)目標(biāo)函數(shù)分別對(duì)A、B求導(dǎo),得:
(23)
得到:
(24)
(25)
從而計(jì)算得到估計(jì)誤差方差陣為:
(26)
(27)
現(xiàn)通過求取期望值以及統(tǒng)計(jì)線性模型近似方法求得模型的后驗(yàn)統(tǒng)計(jì)規(guī)律,即:
(28)
(29)
因此可以得出通過統(tǒng)計(jì)線性化推導(dǎo)出的形式與通過容積點(diǎn)傳播得到的均值和誤差方差相同。
(30)
Fk=(PxX′,k|k-1)T(Pxx,k-1|k-1)-1
(31)
(32)
式中δk-1為統(tǒng)計(jì)線性近似誤差。
同理對(duì)觀測(cè)方程進(jìn)行統(tǒng)計(jì)等效線性化得到統(tǒng)計(jì)回歸矩陣:
(33)
(34)
式中εk為統(tǒng)計(jì)線性近似誤差,用于補(bǔ)償系統(tǒng)的非線性項(xiàng):
(35)
4)H∞容積卡爾曼濾波設(shè)計(jì)。
設(shè)計(jì)一個(gè)濾波器對(duì)濾波器估計(jì)誤差進(jìn)行限制,保證估計(jì)誤差有界,在任意wk、vk、x0的情況下,使得觀測(cè)信息誤差最小化,盡量減少估計(jì)誤差上界。設(shè)計(jì)濾波器滿足:
(36)
(37)
則:
(38)
將式(36)轉(zhuǎn)化為不定式形式:
(39)
進(jìn)行狀態(tài)更新:
(40)
(41)
Pxx,k|k-1
(42)
(43)
為了表示誤差協(xié)方差和調(diào)整參數(shù)之間的關(guān)系,應(yīng)用矩陣求逆定理將式(42)求逆得到:
(44)
β為調(diào)節(jié)因子,權(quán)衡H∞濾波和最小均方根性能。
為了保證誤差方差的正定性,滿足:
(45)
即
(46)
式中eig(·)代表均值的特征值。
噪聲協(xié)方差矩陣對(duì)算法的穩(wěn)定性有著重要影響。為了保證濾波算法的穩(wěn)定性,許多學(xué)者通過在誤差方差和新息方差矩陣中引入附加的正定矩陣[16]:
Qk-1+ΔQk-1
(47)
通過增大Qk來增強(qiáng)穩(wěn)定性。然而其會(huì)增大誤差方差陣Pxx,k|k-1,同時(shí)增大增益矩陣Kk,影響估計(jì)精度。
標(biāo)準(zhǔn)的容積卡爾曼濾波需要假設(shè)噪聲的均值為零,且精確已知噪聲統(tǒng)計(jì)特性。一般情況下過程噪聲和觀測(cè)噪聲的統(tǒng)計(jì)特性不易準(zhǔn)確獲得。較大的Qk-1,可以保證系統(tǒng)的穩(wěn)定性,但估計(jì)精度下降;而Qk-1較小,則會(huì)是濾波器不穩(wěn)定,從而降低濾波性能。因此采用模糊規(guī)則對(duì)過程噪聲和觀測(cè)噪聲進(jìn)行估計(jì),通過增強(qiáng)對(duì)噪聲的自適應(yīng)力來提高濾波精度和收斂速度。
設(shè)過程噪聲方差和測(cè)量噪聲方差的模糊自適應(yīng)律為:
(48)
(49)
(50)
(51)
則自適應(yīng)H∞容積卡爾曼濾波系統(tǒng)結(jié)構(gòu)圖1。
圖1 模糊自適應(yīng)HCKF系統(tǒng)結(jié)構(gòu)Fig.1 The structure of the fuzzy adaptive HCKF system
為了驗(yàn)證模糊自適應(yīng)HCKF算法的有效性,選用以下一般非線性系統(tǒng)進(jìn)行仿真驗(yàn)證,假設(shè)系統(tǒng)的過程噪聲和觀測(cè)噪聲不準(zhǔn)確,采用標(biāo)準(zhǔn)的CKF和模糊自適應(yīng)HCKF的算法對(duì)系統(tǒng)狀態(tài)進(jìn)行估計(jì)比較,并驗(yàn)證統(tǒng)計(jì)線性化的有效性,從而具有一定的普遍性和適用性。
系統(tǒng)的狀態(tài)方程為:
x(k)=0.5x(k-1)+8cos(1.2(k-1))+
w(k)+2.5x(k-1)/(1+x2(k-1))
(52)
觀測(cè)方程為:
z(k)=x2(k)/20+v(k)
(53)
首先對(duì)統(tǒng)計(jì)線性化算法(30)~(35)的有效性進(jìn)行驗(yàn)證,系統(tǒng)的噪聲統(tǒng)計(jì)特性為q=0,r=0,Q=0.06,R=0.01,P0=10。仿真結(jié)果如圖2、圖3所示。
從圖2和圖3的仿真結(jié)果可以得出,通過統(tǒng)計(jì)線性化方法對(duì)非線性系統(tǒng)的近似可以達(dá)到較好的效果,近似值與真實(shí)值基本相符,誤差較小,因此可以通過統(tǒng)計(jì)線性化方法來近似非線性系統(tǒng),從而進(jìn)一步進(jìn)行濾波器的設(shè)計(jì)研究和分析。
圖3 統(tǒng)計(jì)線性化近似誤差Fig.3 Error of approximation
圖2 統(tǒng)計(jì)線性化近似Fig.2 Statistical linearization approximation
為了驗(yàn)證模糊自適應(yīng)HCKF算法的有效性和魯棒性,分別采用標(biāo)準(zhǔn)的容積卡爾曼濾波及模糊自適應(yīng)HCKF濾波方法在以下2種情況下進(jìn)行仿真實(shí)驗(yàn)。
Case1:q=0,r=0,Q=0.06,R=0.01。
Case2:q=0.4,r=0.2,Q=0.06,R=0.01。
Case3:q=0.4,r=0,Q=0.06,R=0.01。
Case4:q=0,r=0.2,Q=0.06,R=0.01。
實(shí)際的噪聲協(xié)方差為Qe=0.05,Ra=0.01,采用模糊自適應(yīng)算法對(duì)噪聲協(xié)方差進(jìn)行估計(jì)。由于篇幅限制,本文僅給出Case1和Case2的結(jié)果仿真圖。Case1,Case2,Case3和Case4的仿真誤差統(tǒng)計(jì)結(jié)果在表中給出,以進(jìn)行結(jié)果對(duì)比,得出相關(guān)結(jié)論。
在Case 1 情況下對(duì)系統(tǒng)的狀態(tài)進(jìn)行估計(jì)比較,在CKF和模糊自適應(yīng)HCKF 2種算法下的估計(jì)結(jié)果和估計(jì)誤差如圖4、圖5所示。
圖4 Case 1情況下不同濾波方法估計(jì)Fig.4 Estimation of different filtering method of Case 1
圖5 Case 1情況下不同濾波方法的估計(jì)誤差Fig.5 The estimation error of different filtering method of Case 1
在Case 2 情況下對(duì)系統(tǒng)的狀態(tài)進(jìn)行估計(jì),在CKF和模糊自適應(yīng)HCKF這2種算法下的估計(jì)結(jié)果和估計(jì)誤差如圖6及圖7所示。
圖6 Case 2情況下不同濾波方法估計(jì)Fig.6 Estimation of different filtering method of Case 2
圖7 Case 2情況下不同濾波方法的估計(jì)誤差Fig.7 The estimation error of different filtering method of Case 2
從仿真結(jié)果圖3~圖6研究中可以發(fā)現(xiàn)當(dāng)系統(tǒng)噪聲和觀測(cè)噪聲均值為零時(shí),模糊自適應(yīng)HCKF和CKF 算法的估計(jì)效果相當(dāng),估計(jì)誤差均值和均方差相差較少,都能較好的估計(jì)系統(tǒng)的狀態(tài)值;當(dāng)系統(tǒng)噪聲和觀測(cè)噪聲均值不為零時(shí),模糊自適應(yīng)HCKF 的估計(jì)效果優(yōu)于CKF。
為了便于更精確的比較2種算法的估計(jì)結(jié)果,本文分別對(duì)2種算法的誤差均值和誤差均方差值進(jìn)行計(jì)算和比較,統(tǒng)計(jì)結(jié)果如表1所示。
(54)
(55)
由表1和表2可知,在分別存在系統(tǒng)噪聲、觀測(cè)噪聲以及同時(shí)存在系統(tǒng)噪聲和觀測(cè)噪聲的情況下,HCKF誤差均值和均方根值小于標(biāo)準(zhǔn)容積卡爾曼,因而提出的HCKF濾波方法能夠有效的減少噪聲產(chǎn)生的影響,在存在非零均值噪聲情況下具有一定的魯棒性和較高的估計(jì)精度,并且在噪聲統(tǒng)計(jì)特性不準(zhǔn)確的情況下,模糊自適應(yīng)可以有效的估計(jì)不準(zhǔn)確的噪聲協(xié)方差,提高估計(jì)效率和收斂速度,防止濾波發(fā)散。
表1 Case 1和Case 2的估計(jì)誤差及均方差Table 1 The AVE and RMSE of Case 1and Case 2
表2 Case 3和Case 4的估計(jì)誤差及均方差Table 2 The AVE and RMSE of Case 3 and Case 4
1)在濾波過程中,在存在噪聲統(tǒng)計(jì)特性不準(zhǔn)確及非零均值噪聲協(xié)方差的情況下,與標(biāo)準(zhǔn)的容積卡爾曼濾波相比,本文提出模糊自適應(yīng)HCKF能夠更好的估計(jì)非線性系統(tǒng)的狀態(tài)值,具有較高的估計(jì)精度和魯棒性。
2)提出的模糊自適應(yīng)規(guī)則能夠較好的估計(jì)不準(zhǔn)確的噪聲協(xié)方差,從而避免濾波發(fā)散,提高收斂速度。
3)本文提出的濾波方法可以較好地推廣到任意非線性系統(tǒng),應(yīng)用廣泛,具有一定的工程應(yīng)用價(jià)值。