王哲,李國輝,趙善祿
(空軍航空大學(xué)軍事仿真研究所,長春130022)
經(jīng)典譜估計雖然計算效率高,但是由于窗函數(shù)的局限性,在計算過程中不可避免地出現(xiàn)頻譜泄露、譜峰轉(zhuǎn)移、分辨率低等缺點(diǎn)。針對經(jīng)典譜估計的局限,現(xiàn)代譜估計方法應(yīng)運(yùn)而生,其中最大熵譜估計法是一種常用的現(xiàn)代譜估計方法。
1967年,Burg最早將“熵”這一概念從信息論中引入到譜估計中,提出最大熵譜估計[1]。最大熵譜估計具有分辨率高,使用范圍廣的特點(diǎn)。
“熵”是對于隨機(jī)性進(jìn)行度量的量,熵值越大,則隨機(jī)性越大。假設(shè)p(x)為連續(xù)隨機(jī)變量的概率函數(shù)密度,則“熵”可以定義為:

如果平穩(wěn)時間序列{xt,t∈Z}是正態(tài)隨即過程,其功率譜為S(ω),那么功率譜熵為:

與熵的原始定義相似,I(S)表示隨即序列不確定性的度量,I(S)值越大,表明隨機(jī)性越強(qiáng)。根據(jù)譜估計知識可知,樣本序列自相關(guān)函數(shù)和功率譜密度構(gòu)成傅里葉變換對,在已知相關(guān)函數(shù)r(m),m=0,1,2,…,p情況下,對m>p的自相關(guān)函數(shù)進(jìn)行外推,當(dāng)功率譜S(ω)滿足如下約束條件:

當(dāng)在上數(shù)條件下使得H(S(ω))最大時,這樣的過程稱為最大熵譜估計,估計出的S?()ω為最大熵譜。可以證明,正態(tài)隨機(jī)過程的最大熵譜估計與AR模型譜估計是等效的[2]。AR模型的功率譜為:要想求信號功率譜密度pn(ω),需要求出上式中的參數(shù)φ和σ。可以通過Yule-Walker法或Burg法求得相關(guān)參數(shù),具體過程見參考文獻(xiàn)[3]。當(dāng)確定相關(guān)參數(shù)后,就可通過4式計算功率譜。

對最大熵譜估計進(jìn)行一致性檢驗(yàn),有常規(guī)的非參數(shù)檢驗(yàn)、窗譜估計一致性檢驗(yàn)和最大熵譜估計一致性檢驗(yàn)等多種方法[4]。本文采用最大熵譜估計一致性檢驗(yàn)法與模型相結(jié)合進(jìn)行檢驗(yàn)。
最大熵譜估計的統(tǒng)計結(jié)果是有限的近似結(jié)果,是漸進(jìn)無偏與漸進(jìn)正態(tài)的[5]。當(dāng)采樣數(shù)目N與AR模型的階數(shù)M足夠大時,最大熵譜的均值與方差為:

對 S?(ω )作對數(shù)變換得到logS?(ω ),由于 S?(ω)是正態(tài)的,所以 logS?(ω )仍為正態(tài)的。logS?(ω)的均值和方差為:

該對數(shù)譜 logS?(ω)的100(1 -α)%置信區(qū)間近似為

對于兩個獨(dú)立時間序列的功率譜密度均漸進(jìn)正態(tài)情況下,第i個帶寬中對數(shù)的抽樣分布近似為:

假設(shè)兩個時間序列具有相同功率譜密度函數(shù),即S1(ω1)=S2(ω2),根據(jù)式8可以得出:

該假設(shè)檢驗(yàn)的接受域是[ ]-Zα/2≤D≤Zα/2,α 為該檢驗(yàn)的顯著水平,N1、N2是兩次試驗(yàn)的樣本容量,M1、M2是AR模型的階次。
在檢驗(yàn)過程中,原則上對每個頻率點(diǎn)都要進(jìn)行檢驗(yàn),如果每個頻率點(diǎn)都接受該檢驗(yàn),則說明兩序列相容。檢驗(yàn)的公式為:

但是,與窗譜估計一致性檢驗(yàn)類似,在實(shí)際驗(yàn)證過程中,由于頻率點(diǎn)過多并且有的點(diǎn)功率譜較弱,所以工程中常常關(guān)注于頻率譜較強(qiáng)的點(diǎn)。
由前面章節(jié)的知識可知,對于正態(tài)隨機(jī)過程,可以通過建立與之等效的AR模型進(jìn)行最大熵譜估計。本節(jié)仍采用飛機(jī)水平加速的縱向過載nx隨飛行速度ViV變化的仿真數(shù)據(jù)與試飛數(shù)據(jù)進(jìn)行一致性檢驗(yàn)。故需先建立AR模型。
在對仿真數(shù)據(jù)與試飛數(shù)據(jù)進(jìn)行最大熵譜估計前,首先需要建立AR模型。由于經(jīng)過數(shù)據(jù)預(yù)處理,仿真數(shù)據(jù)與試飛數(shù)據(jù)均可視為正態(tài)隨機(jī)過程,可以參考文獻(xiàn)[6]中時序建模的方法建立AR模型。對于模型階數(shù)的確定,目前有多種準(zhǔn)則與方法,主要有赤池FPE準(zhǔn)則、AIC準(zhǔn)則和Parzen的自回歸傳遞函數(shù)準(zhǔn)則[7]。本文對于模型階數(shù)的確定,采用工程上常用的經(jīng)驗(yàn)公式由于采樣數(shù)目為N1,N2=3000,則M≈18,采樣頻率為
階數(shù)為 P,模型參數(shù)為 φ1,φ2,…,φP的 AR(p)可以表示為

其中,{εt}是白噪聲序列,滿足 EXsεt=0 對一切s 采用Marple算法解yule-walker方程,利用MATLAB軟件計算出相關(guān)參數(shù),得到仿真數(shù)據(jù)的AR模型為: Xt=1.5761Xt-1-0.7562Xt-2-0.5123Xt-3-0.3165Xt-4+0.2119Xt-5+0.0352Xt-6-0.0158Xt-7-0.0264Xt-8+0.0322Xt-9+0.0286Xt-10-0.0346Xt-11+0.0153Xt-12-0.0234Xt-13+0.0315Xt-14+0.0446Xt-15-0.0061Xt-16+0.0019Xt-17+0.0211Xt-18 試飛數(shù)據(jù)的AR模型為: Xt=1.3213Xt-1-0.8613Xt-2-0.7449Xt-3-0.4881Xt-4-0.3591Xt-5-0.1011Xt-6-0.0695Xt-7+0.0987Xt-8+0.0413Xt-9-0.0363Xt-10-0.0303Xt-11-0.0134Xt-12+0.0006Xt-13+0.0201Xt-14+0.0362Xt-15-0.0061Xt-16+0.0109Xt-17+0.0211Xt-18 根據(jù)4式對仿真數(shù)據(jù)與試飛數(shù)據(jù)進(jìn)行最大熵譜估計,得到如圖1所示功率譜對比曲線。 圖1nx隨Vi變化最大熵譜估計對比曲線 根據(jù)歸一化處理后的功率譜密度,確定主要頻率區(qū)間為[0,0.025],因此可以對該區(qū)間內(nèi)的頻率點(diǎn)進(jìn)行一致性檢驗(yàn)。要驗(yàn)證nx隨Vi變化的試飛數(shù)據(jù)功率譜S1(ω1)和仿真數(shù)據(jù)功率譜S2(ω2)的相容性,需要檢驗(yàn)S1(ω1)=S2(ω2),構(gòu)造統(tǒng)計量D: 其中,M1,M2為AR模型階次,N1,N2為采樣點(diǎn)數(shù)量,α為顯著水平,統(tǒng)計量D對于該假設(shè)檢驗(yàn)的接受域?yàn)閇 ] -Zα/2,Zα/2,按公式 11進(jìn)行檢驗(yàn): 當(dāng)顯著水平 α取 0.05時,查表得 Zα/2=1.96,N1=N2=3000,M1=M2=18,計算得B≈0.30。因此,只需要驗(yàn)證在主要頻率區(qū)間內(nèi)每個頻率點(diǎn)處A的值小于0.30,就可以說明試飛數(shù)據(jù)與仿真數(shù)據(jù)的功率譜密度具有一致性。圖2為A在[0,0.025]區(qū)間內(nèi)的曲線圖。 圖2nx隨Vi變化最大熵譜估計驗(yàn)證曲線 由圖2可知,在主要頻率區(qū)間內(nèi),有A≤B,所以在置信水平為95%情況下,S1(ω1)=S2(ω2)成立,即 nx隨Vi變化的試飛數(shù)據(jù)與仿真數(shù)據(jù)具有一致性。 本文以飛機(jī)在某高度平飛加速過程中過載隨速度的變化為例,采用最大熵譜估計法對仿真數(shù)據(jù)和試飛數(shù)據(jù)進(jìn)行一致性檢驗(yàn),取譜峰附近[0,0.025dB]作為主要頻率區(qū)間,檢驗(yàn)結(jié)果顯示,功率譜的主要頻率點(diǎn)在95%的置信水平下,nx隨Vi變化的仿真數(shù)據(jù)與試飛數(shù)據(jù)在統(tǒng)計意義下一致,該方法適用于其他動態(tài)飛行性能的一致性檢驗(yàn)。 [1]李鵬波.時間序列樣本的總體一致性檢驗(yàn)——頻域方法[J].飛行器測控學(xué)報,1999,18(4). [2]曾鳴,李雪青,謝保川,等.基于飛參數(shù)據(jù)的飛行仿真模型驗(yàn)證[J].指揮控制與仿真,2011,6(12). [3]徐慧娟.自回歸AR模型的整體最小二乘分析研究[D].上海:東華理工大學(xué),2012. [4]李鶴.基于試飛數(shù)據(jù)的模型飛行模擬器飛行性能驗(yàn)證研究[D].空軍工程大學(xué),2009. [5]李鵬波,高霞.應(yīng)用最大熵譜估計進(jìn)行導(dǎo)彈系統(tǒng)的仿真模型驗(yàn)證[J].國防科技大學(xué)學(xué)報,1999,21(2). [6]李鶴,呂巖,李國輝.飛行仿真動態(tài)數(shù)據(jù)驗(yàn)證的時序建模方法[J].兵工自動化,2008,27(12). [7]王建華,符文星,董敏周.最大墑譜估計在空空導(dǎo)彈仿真模型驗(yàn)證中的應(yīng)用[J].彈箭與制導(dǎo)學(xué)報,2005,25(4).




3 結(jié)語