基于ARMA模型的功率譜估計(jì)

1、LOGO 賈鑫賈鑫2012.12.08 目錄目錄一、一、ARMAARMA過程基本理論過程基本理論二、平穩(wěn)二、平穩(wěn)ARMAARMA過程功率譜過程功率譜三、平穩(wěn)三、平穩(wěn)ARMAARMA過程譜估計(jì)過程譜估計(jì)四、四、ARAR模型辨識(shí)模型辨識(shí)五、算例五、算例目錄目錄一、一、ARMAARMA過程基本理論過程基本理論二、平穩(wěn)二、平穩(wěn)ARMAARMA過程功率譜過程功率譜三、平穩(wěn)三、平穩(wěn)ARMAARMA過程譜估計(jì)過程譜估計(jì)四、四、ARAR模型辨識(shí)模型辨識(shí)五、算例五、算例1、相當(dāng)多的平穩(wěn)隨機(jī)過程都可以通過用白噪聲激勵(lì)一線性時(shí)不變系統(tǒng)來產(chǎn)生。2、線性系統(tǒng)可以用線性差分方程（ARMA模型）進(jìn)行描述。3、任何一個(gè)有理式

2、的功率譜密度都可以用一個(gè)ARMA隨機(jī)過程的功率譜來精確逼近。將廣義的平穩(wěn)過程x(n)表示成一個(gè)輸入序列u(n)(白噪聲)激勵(lì)線性系統(tǒng)H(z)(ARMA模型)的輸出由H(z)的輸出功率譜來估計(jì)x(n)的功率譜ARMA過程定義過程定義( )( ) ie nx nh 利用已知的x(n)來估計(jì)H(z)的參數(shù)將廣義的平穩(wěn)過程x(n)表示成一個(gè)輸入序列u(n)(白噪聲)激勵(lì)線性系統(tǒng)H(z)(ARMA模型)的輸出由H(z)的輸出功率譜來估計(jì)x(n)的功率譜ARMA過程定義過程定義 離散隨機(jī)過程離散隨機(jī)過程 服從線性差分方程：服從線性差分方程： 為離散白噪聲，則稱為離散白噪聲，則稱 為為ARMA過程。過程。自

3、回歸自回歸 (autoregressive)滑動(dòng)平均滑動(dòng)平均(moving average)過程過程 ( )e n ( )x n ( )x n11( )(1)()( )(1)()pqx na x na x npe nbe nb e nq11( )()( )()pqijijx na x n ie nb e njARAR階數(shù)階數(shù)ARAR參數(shù)參數(shù)MAMA階數(shù)階數(shù)MAMA參數(shù)參數(shù)ARMA過程定義過程定義2( ) (0,)e nN( )()jz x nx nj后向移位算子：11( )1ppA za za z 其中：00()()pqijija x nib e nj( ) ( )( ) ( )A z x n

4、B z e n11( )1qqB zb zb z ARMA過程定義過程定義( )( )( )nknkx ne k he nhARMA模型描述的線性時(shí)不變（模型描述的線性時(shí)不變（LTI）系統(tǒng)系統(tǒng)傳遞函數(shù)：傳遞函數(shù)：( )( ) ie nx nh ( )( )( )iiiB zH zh zA zARMA過程定義過程定義沖擊響應(yīng)沖擊響應(yīng)系數(shù)系數(shù)滿足滿足ARMAARMA模型的條件：模型的條件：(1)(1)沖激響應(yīng)系數(shù)必須絕對(duì)可求和：沖激響應(yīng)系數(shù)必須絕對(duì)可求和： ( (系統(tǒng)穩(wěn)定系統(tǒng)穩(wěn)定) )(2)(2)A(z)和和B(z)無公共因子無公共因子(p,q唯一唯一)(3)(3)系統(tǒng)是物理可實(shí)現(xiàn)的系統(tǒng)是物理可實(shí)

5、現(xiàn)的( (因果系統(tǒng)因果系統(tǒng)) )u極點(diǎn)的作用：決定系統(tǒng)的穩(wěn)定性和因果性u(píng)因果性：稱x(n)是e(n)的因果函數(shù)，若即因果系統(tǒng)要求極點(diǎn)在單位圓以內(nèi)，A(z)的根|z|1kkh ( )( )( )B zH zA z零點(diǎn)部分零點(diǎn)部分極點(diǎn)部分極點(diǎn)部分0 ( )()iiiihx nh e ni ARMA過程性質(zhì)過程性質(zhì)u零點(diǎn)的作用：決定系統(tǒng)的可逆性，即 u可逆性：稱e(n)是x(n)的可逆函數(shù)，若 (1)存在序列 ，并滿足 (2) 可逆系統(tǒng)的穩(wěn)定 可逆性條件11( )( )( )( )A zHzH zB z0( )()iie nx ni iii ARMA過程性質(zhì)過程性質(zhì)( )1A z 11( )( )1

6、( )iqiqiB zH zh zb zb zA z ( ) ( )( ) ( )A z x nB z e n特例一：MA過程( )( ) ( )x nB z e n 1,ihiq抽頭有限沖激響應(yīng)(FIR)系統(tǒng)ARMA過程特例過程特例MA過程過程滑動(dòng)平均滑動(dòng)平均1( )( )H zA z2( )1, ( )WN(0,)eB ze n特例二：AR過程中含有中含有 的無數(shù)多項(xiàng)的無數(shù)多項(xiàng)1z無限沖激響應(yīng)(IIR)系統(tǒng)ARMA過程特例過程特例AR過程過程自回歸自回歸( ) ( )( ) ( )A z x nB z e nARMA過程的過程的Wold分解定理分解定理vWold分解定理：分解定理： 任何一

7、個(gè)具有有限方差的ARMA或MA過程，可以表示成唯一的、階數(shù)有可能無窮大的AR過程； 同樣，任何一個(gè)ARMA或AR過程也可以表示成一個(gè)階數(shù)可能無窮大的MA過程。目錄目錄一、一、ARMAARMA過程基本理論過程基本理論二、平穩(wěn)二、平穩(wěn)ARMAARMA過程功率譜過程功率譜三、平穩(wěn)三、平穩(wěn)ARMAARMA過程譜估計(jì)過程譜估計(jì)四、四、ARAR模型辨識(shí)模型辨識(shí)五、算例五、算例1*11*1()1( )()1( )ppqqA za za zA zB zb zb zBz 21221( )( ) ()( )( )( ) ()xjwjwz ez eB zB z B zPA zA z A z 則功率譜則功率譜 其中其

8、中( ) ( )( ) ( )A z x nB z e n2( )(0,)e nNARMA過程功率譜定義過程功率譜定義ARMA過程功率譜定義過程功率譜定義v證明證明 設(shè) 是零均值離散時(shí)間平穩(wěn)過程，取ARMA過程 則：對(duì)上式兩邊取數(shù)學(xué)期望 計(jì)算自相關(guān)函數(shù)ARMA過程功率譜定義過程功率譜定義 由上式計(jì)算功率譜密度函數(shù) 取 為白噪聲 ，則有 （白噪聲功率譜密度為常數(shù)， ） 固有目錄目錄一、一、ARMAARMA過程基本理論過程基本理論二、平穩(wěn)二、平穩(wěn)ARMAARMA過程功率譜過程功率譜三、平穩(wěn)三、平穩(wěn)ARMAARMA過程譜估計(jì)過程譜估計(jì)四、四、ARAR模型辨識(shí)模型辨識(shí)五、算例五、算例vWold定理表明

9、：定理表明： 一個(gè)ARMA模型可以用一個(gè)階數(shù)足夠大的AR模型來近似。 相比于ARMA模型不僅需要確定AR階數(shù)和MA階數(shù)，還需要估計(jì)AR參數(shù)和MA參數(shù)（MA參數(shù)估計(jì)必須求解非線性方程組），AR模型相對(duì)簡(jiǎn)單，故工程上常用AR模型作近似。vARMA功率譜的線性估計(jì)方法的基本思路都是首功率譜的線性估計(jì)方法的基本思路都是首先解線性方程估計(jì)出先解線性方程估計(jì)出AR參數(shù)，再通過一定的方參數(shù)，再通過一定的方法，將功率譜表達(dá)式轉(zhuǎn)換成只需要法，將功率譜表達(dá)式轉(zhuǎn)換成只需要AR參數(shù)，而參數(shù)，而不需要不需要MA具體參數(shù)值的計(jì)算表達(dá)式。具體參數(shù)值的計(jì)算表達(dá)式。估計(jì)方法估計(jì)方法21221( )( ) ()( )( )( )

10、 ()xjwjwz ez eB zB z B zPA zA z A z估計(jì)方法估計(jì)方法AR過程的過程的實(shí)現(xiàn)方法實(shí)現(xiàn)方法ARMA過程過程的實(shí)現(xiàn)方法的實(shí)現(xiàn)方法定階定階p&q估計(jì)估計(jì)AR、MA參數(shù)參數(shù)功率譜功率譜計(jì)算計(jì)算將ARMA功率譜密度分解為兩部分之和：線性化方法一線性化方法一: Cadzow譜估計(jì)子譜估計(jì)子其中，?。毫硪环矫婀β首V可做如下類似分解：其中，取：線性化方法一線性化方法一: Cadzow譜估計(jì)子譜估計(jì)子可以得到：從而可以計(jì)算ARMA模型的功率譜：線性化方法二線性化方法二: Kaveh譜估計(jì)子譜估計(jì)子將ARMA功率譜密度公式作如下變形：為了保證上式中第二個(gè)等號(hào)相等，有：可以看出

11、， 具有對(duì)稱性，即：從上式中第三個(gè)等式，有：線性化方法二線性化方法二: Kaveh譜估計(jì)子譜估計(jì)子比較上式兩邊同冪次項(xiàng)的系數(shù)，可以得到：從而可以計(jì)算ARMA模型的功率譜：目錄目錄一、一、ARMAARMA過程基本理論過程基本理論二、平穩(wěn)二、平穩(wěn)ARMAARMA過程功率譜過程功率譜三、平穩(wěn)三、平穩(wěn)ARMAARMA過程譜估計(jì)過程譜估計(jì)四、四、ARAR模型辨識(shí)模型辨識(shí)五、算例五、算例AR模型階數(shù)確定模型階數(shù)確定vFPE(Final Prediction Error)準(zhǔn)則函數(shù)準(zhǔn)則函數(shù)vAIC(An Information Criterion)準(zhǔn)則函數(shù)準(zhǔn)則函數(shù)vMDL(Minimum Descriptio

12、n Length)準(zhǔn)則函數(shù)準(zhǔn)則函數(shù)v在各自準(zhǔn)則取得最小值時(shí)的模型為適用模型在各自準(zhǔn)則取得最小值時(shí)的模型為適用模型 為AR模型階數(shù)， 為激勵(lì)方差， 為樣本點(diǎn)數(shù)。赤池，日本，1969赤池，日本，1974Rissanen，芬蘭，1983AR模型參數(shù)估計(jì)模型參數(shù)估計(jì)ARMA過程可以表示為：其自相關(guān)函數(shù)為：由白噪聲，有：因此，可得：AR模型參數(shù)估計(jì)模型參數(shù)估計(jì)由ARMA過程的定義式，有：從而可以得到下式：注意，對(duì)于一個(gè)ARMA過程而言，其MA參數(shù)在q階以上為零，即有：ARMA過程的自相關(guān)函數(shù)可總結(jié)為如下結(jié)構(gòu):式中，p和q分別是AR和MA的階數(shù)，ai和bj分別是AR參數(shù)和MA參數(shù)，r(k)是輸入信號(hào)的自相

13、關(guān)函數(shù)，h是ARMA模型的參數(shù)，當(dāng)h下標(biāo)小于0時(shí)，h均取零。該式是很多AR(MA)過程確定AR系數(shù)估計(jì)器的基礎(chǔ)。AR模型參數(shù)估計(jì)模型參數(shù)估計(jì)0001)0()2() 1()()2()0() 1 ()2() 1() 1 ()0() 1 ()()2() 1 ()0(221pxxxxxxxxxxxxxxxxaaarprprprprrrrprrrrprrrr解上述方程，就可以求出功率譜計(jì)算公式中的所需參數(shù)，進(jìn)而求出功率譜。對(duì)于該方程，可以采用直接解法，也可以采用Levinson-Durbin或Delsarte-Genin等階遞推算法來減小計(jì)算量。AR模型參數(shù)估計(jì)模型參數(shù)估計(jì)-Yule-Walker方法方

14、法AR模型參數(shù)估計(jì)模型參數(shù)估計(jì)-最小二乘方法最小二乘方法AR模型參數(shù)估計(jì)模型參數(shù)估計(jì)-最小二乘方法最小二乘方法取目標(biāo)函數(shù)：求解方程組： 可得： 令：即可確定AR模型參數(shù)。引申：當(dāng)同時(shí)考慮A和b二者的誤差或擾動(dòng)時(shí)，可獲得AR參數(shù)估計(jì)的總體最小二乘法。目錄目錄一、一、ARMAARMA過程基本理論過程基本理論二、平穩(wěn)二、平穩(wěn)ARMAARMA過程功率譜過程功率譜三、平穩(wěn)三、平穩(wěn)ARMAARMA過程譜估計(jì)過程譜估計(jì)四、四、ARAR模型辨識(shí)模型辨識(shí)五、算例五、算例算例算例1、利用、利用matlab自帶的計(jì)算函數(shù)，實(shí)現(xiàn)了對(duì)信號(hào)的自帶的計(jì)算函數(shù)，實(shí)現(xiàn)了對(duì)信號(hào)的 AR功率譜估計(jì)。功率譜估計(jì)。2、利用、利用mat

15、lab自帶的自帶的AR模型參數(shù)計(jì)算函數(shù)，結(jié)模型參數(shù)計(jì)算函數(shù)，結(jié) 合合YuleWalker方程，實(shí)現(xiàn)了對(duì)信號(hào)的功率方程，實(shí)現(xiàn)了對(duì)信號(hào)的功率譜譜 估計(jì)。估計(jì)。3、利用、利用YuleWalker方法，首先對(duì)方法，首先對(duì)ARMA模型的模型的 AR參數(shù)進(jìn)行計(jì)算，并利用參數(shù)進(jìn)行計(jì)算，并利用LevinsonDurbin 算法實(shí)現(xiàn)對(duì)算法實(shí)現(xiàn)對(duì)MA參數(shù)的估計(jì)，完成對(duì)信號(hào)的功參數(shù)的估計(jì)，完成對(duì)信號(hào)的功 率譜估計(jì)。率譜估計(jì)。4、利用、利用YuleWalker方法，首先對(duì)方法，首先對(duì)ARMA模型的模型的 AR參數(shù)進(jìn)行計(jì)算，并利用參數(shù)進(jìn)行計(jì)算，并利用Kaveh譜估計(jì)子算譜估計(jì)子算 法，實(shí)現(xiàn)對(duì)信號(hào)的功率譜估計(jì)。法，實(shí)現(xiàn)對(duì)

16、信號(hào)的功率譜估計(jì)。算例算例估計(jì)信號(hào)如下：估計(jì)信號(hào)如下：x = cos(2*pi*40*n)+3*cos(2*pi*90*n)+randn(size(n);采樣頻率：采樣頻率：Fs=1024傅里葉變換數(shù)：傅里葉變換數(shù)：nfft = 512;算例算例1、 AR matlab各種自帶函數(shù)各種自帶函數(shù)算例算例1 AR matlab自帶函數(shù)自帶函數(shù)v 可以看出利用可以看出利用matlab的自帶函數(shù)，各種估計(jì)方法所得的的自帶函數(shù)，各種估計(jì)方法所得的結(jié)果非常接近結(jié)果非常接近v 接下來我們就考察接下來我們就考察AR參數(shù)的不同對(duì)估計(jì)結(jié)果的影響了，參數(shù)的不同對(duì)估計(jì)結(jié)果的影響了，有了上一結(jié)論，我們可以只采用一種算法來進(jìn)行不同有了上一結(jié)論，我們可以只采用一種算法來進(jìn)行不同AR階數(shù)的