前言
本欄前兩節(jié)經(jīng)典譜估計(jì)中提到:經(jīng)典譜估計(jì)下,方差和分辨率是一對(duì)矛盾。這是因?yàn)榻?jīng)典譜估計(jì)將數(shù)據(jù)進(jìn)行了加窗,自相關(guān)法還對(duì)自相關(guān)進(jìn)行了加窗(二次加窗),這就讓我們想到把原始數(shù)據(jù)藏在一個(gè)系統(tǒng)H(Z)中,讓這個(gè)系統(tǒng)包含這組數(shù)據(jù)的特性,這樣一來(lái),系統(tǒng)中的系數(shù)就可以表示系統(tǒng)反映的數(shù)據(jù)。這就是現(xiàn)代功率譜密度估計(jì)-參數(shù)模型法的思想。按照書(shū)本的就是先根據(jù)數(shù)據(jù)的自相關(guān)函數(shù)r(m)求出H(Z)系數(shù),再通過(guò)H(Z)進(jìn)行譜估計(jì)。
參數(shù)模型法有AR,MA,ARMA模型,其性質(zhì)為:
AR | MA | ARMA | |
---|---|---|---|
H(Z) | ![]() |
![]() |
![]() |
線性/非線性 | 線性 | 非線性 | 非線性 |
反映頻譜特性 | 峰值 | 谷值 | 兼顧 |
一、概率梳理
由于AR模型是線性方程,也可以等效預(yù)測(cè)模型,所以AR模型遠(yuǎn)比另外兩種實(shí)用,所以本文只梳理和仿真了AR模型。
首先模型中令輸入是白噪聲,需要求解H(Z)的系數(shù)也就是ak,k=1,2…p。也就是要通過(guò)數(shù)據(jù)的自相關(guān)與ak的關(guān)系進(jìn)行求解,也就是需要通過(guò)正則方程(Yule-Walker方程),正則方差的推導(dǎo)過(guò)程如下:
其中,正則方差可以用Lecinson-Durbin遞推快速算法計(jì)算,也就是自相關(guān)法的方法。后面還會(huì)說(shuō)其他的方法以及比較。
并且在預(yù)測(cè)模型中,二者系統(tǒng)的系數(shù)是相等的,其最小誤差等效于AR模型輸入端白噪聲的方差。其關(guān)系如下:
二、AR模型的幾種方法
比較常見(jiàn)的有剛剛提及的自相關(guān)法,還有burg法和改進(jìn)后的協(xié)方差法,它們之間的區(qū)別如下:
自相關(guān)法 | burg法 | 改進(jìn)后的協(xié)方差法 | |
---|---|---|---|
預(yù)測(cè)方式 | 前向預(yù)測(cè) | 前后向預(yù)測(cè) | 前向預(yù)測(cè) |
加窗 | 前后加窗 | 前后都不加窗 | 前后都不加窗 |
Levinson遞推算法 | 可以使用 | 可以使用 | 不可使用 |
除此之外,還有常會(huì)用到的最大熵譜估計(jì),由于數(shù)據(jù)可能相對(duì)于原始數(shù)據(jù)還是有截短。之前的經(jīng)典譜估計(jì)是將兩邊直接為零,而這里是將兩邊加上最隨機(jī)的數(shù),也就是最大熵的準(zhǔn)則。
三、AR模型的方法與具體仿真
為了和經(jīng)典功率譜估計(jì)比較,采用的原數(shù)據(jù)和前兩節(jié)經(jīng)典功率譜估計(jì)一樣的,仿真采取的階數(shù)均為50文章來(lái)源:http://www.zghlxwxcb.cn/news/detail-465238.html
%%現(xiàn)代功率譜估計(jì)的一些方法
clear all;
clc;close all;clear;
N=200; Nfft=20000; Fs=1; n=0:N-1;
x=cos(0.3*pi*n)+cos(0.32*pi*n)+0.1*randn(size(n));
fn=-0.5:2/N:0.5;
figure;
%% burg法
% 使用 Burg 算法得到功率譜估計(jì);
xpsd=pburg(x,50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(221);
plot(fn,fftshift(xpsd));grid on;title('burg法');
%% 協(xié)方差法
xpsd=pcov(x,50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(222);
plot(fn,fftshift(xpsd));grid on;title('協(xié)方差法');
%% 協(xié)方差的改進(jìn)法
xpsd=pmcov(x,50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(223);
plot(fn,fftshift(xpsd));grid on;title('改進(jìn)的協(xié)方差法');
%% 最大熵法
xpsd=pmem(x,50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(224);
plot(fn,fftshift(xpsd));grid on;title('最大熵法');
%% 自相關(guān)
figure;
% 使用自相關(guān)矩陣分解的 MUSIC 算法得到功率譜估計(jì);
xpsd=pmusic(x',50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(311);
plot(fn,fftshift(xpsd));grid on;title('自相關(guān)矩陣分解的MUSIC算法');
% 使用自相關(guān)矩陣分解的特征向量算法得到功率譜估計(jì);
[xpsd,F,V,E]=peig(x',50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(312);
plot(fn,fftshift(xpsd));grid on;title('自相關(guān)矩陣分解的特征向量算法');
% 使用自相關(guān)法得到功率譜估計(jì);
xpsd=pyulear(x,50,N);
pmax=max(xpsd);
xpsd=xpsd/pmax;
xpsd=10*log10(xpsd+0.000001);
subplot(313);
plot(fn,fftshift(xpsd));grid on;title('自相關(guān)法');
文章來(lái)源地址http://www.zghlxwxcb.cn/news/detail-465238.html
到了這里,關(guān)于現(xiàn)代信號(hào)處理-現(xiàn)代功率譜密度估計(jì)AR模型的文章就介紹完了。如果您還想了解更多內(nèi)容,請(qǐng)?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章,希望大家以后多多支持TOY模板網(wǎng)!