《MATLAB中AR模型功率谱估计中AR阶次估计的实现》由会员分享,可在线阅读,更多相关《MATLAB中AR模型功率谱估计中AR阶次估计的实现(8页珍藏版)》请在金锄头文库上搜索。
1、MATLAB中AR模型功率谱估计中AR阶次估计的实现(最近看了几个关于功率谱的问题,有关AR模型的谱估计,在此分享一 下,希望大家不吝指正)(声明:本文内容摘自我的毕业论文心率变异信号的预处理及功率 谱估计)(按:AR模型功率谱估计是对非平稳随机信号功率谱估计的常用方法, 但是其模型阶次的估计,除了 HOSA工具箱里的arorder函数外,没有现成的 函数可用,arorder函数是基于矩阵SVD分解的阶次估计方法,为了比较各 种阶次估计方法的区别,下面的函数使用了FPE, AIC, MDL, CAT 集 中准则一并估计,并采用试验方法确定那一个阶次更好。)以上省 略假设原始数据序列为x,那么n
2、阶参数使用最小二乘估计在MATLAB中实 现如下:复制内容到剪贴板代码:Y = x;Y(1:n) = ;m = N-n;X = ;% 构造系数矩阵for i = 1:mfor j = 1:nX(i,j) = xt(n+i-j);endendbeta = inv(X*X)*X*Y;beta即为用最小二乘法估计出的模型参数。此外,还有估计AR模型参数的Yule-Walker方程法、基于线性预测理论 的Burg算法和修正的协方差算法等26。相应的参数估计方法在MATLAB中 都有现成的函数,比如aryule、arburg以及arcov等。4.3.3 AR模型阶次的选择及实验设计文献26冲介绍了五种不
3、同的AR模型定阶准则,分别为矩阵奇异值分 解(Singular Value Decomposition, SVD)定阶法、最小预测定误差阶准 则(Final Prediction Error Criterionn, FPE)、AIC定阶准则(Akaikas Information theoretic Criterion, AIC)、MDL 定阶准则以及 CAT 定阶准 则。文献28中还介绍了一种BIC定阶准则。SVD方法是对Yule-Walker方 程中的自相关矩阵进行SVD分解来实现的,在MATLAB工具箱中arorder函数 就是使用的该算法。其他五种算法的基本思想都是建立目标函数,阶次估
4、计 的标准是使目标函数最小化。以上定阶准则在MATLAB中也可以方便的实现,下面是本文实现FPE、AIC、 MDL、CAT定阶准则的程序(部分):复制内容到剪贴板代码:for m = 1:N-1%判断是否达到所选定阶准则的要求if strcmp(criterion,FPE) objectfun(m+1) = (N+(m+1)/(N-(m+1)*E(m+1);elseif strcmp(criterion,AIC) objectfun(m+1) = N*log(E(m+1)+2*(m+1);elseif strcmp(criterion,MDL) objectfun(m+1) = N*log(E
5、(m+1)+(m+1)*log(N);elseif strcmp(criterion,CAT)for index = 1:m+1temp = temp+(N-index)/(N*E(index);end objectfun(m+1) = 1/N*temp-(N-(m+1)/(N*E(m+1);endif objectfun(m+1) = objectfun(m) orderpredict = m; break;endendorderpredict 变量即为使用相应准则预测的 AR 模型阶次。(注:以上代码为结合 MATLAB 工具箱函数 pburg,arburg 两个功率谱估计函数增加而得,修
6、改后的Pburg等函数会在附件中示意,名为 pburgwithcriterion)e-rrorargchkt Lu g nafgini iif Tiaxgili = 1crtt eiLcn = J FPE-1 ;endtiXjIix = siae:n)J t:k mH 矩(n);X = X : 1 .1ST昱占达拥臓定断淮山樓卞 if strciip CELteriarb FPE)ab jtctfiuitnrl-) = (ff+tiH-JJJ/CK-tnL+LJJE-tid-);1 sei st rcnp t.crzt ec ion.,AlC.)ob jectfuriCiri-) = N:E
7、1 irides);下面本文使用3.2.1实验设计的输出结果即20例经预处理的HRV信号序 列作为实验对象,分别使用FPE、AIC、MAL和CAT定阶准则预测AR模型阶次,图(见下页)为其中一例典型信号使用不同预测准则其目标函数随阶次的变 化情况。从图中可以看出,使用FPE、AIC以及MDL定阶准则所预测的AR模 型阶次大概位于10附近,即阶次10左右会使相应的目标函数最小化,符合 定阶准则的要求,使用CAT定阶准则预测的阶次较小,在510之间。图(见 下页)为另一例信号的阶次估计情况,从中也可以得到同样的结论。搜索更多相关主题的帖子:MATLAB AR模型阶次估计本帖最近评分记录mengze
8、 财富-1 2009-9-8 10:43mooni财富+12鼓励原创内容2009-8-28 20:592#大中小发表于2009-8-28 20:54 只看该作者本文分别选择6、8、9、10、11、12、16阶AR模型进行上述20例信号 序列进行频谱分析,图(见下页)为其中一例信号使用不同阶次AR模型谱 估计的效果,从图中可以看出,使用6阶AR模型已经可以获得比较光滑的 功率谱估计,但是频率分辨率不高,谱峰不容易分辨,随着模型阶次的增大, 频率分辨率逐渐增加,到16阶AR模型时,虽然频率分辨率较大了,但是随 之牺牲的是频谱的方差特性,因此综合权衡谱图的方差特性以及频率分辨率 的要求,本文选择11
9、阶AR模型对HRV信号进行谱估计。關25201510MatlabififeiLoveMatlabnKU省 略4.5.2 HRV频谱分析实验结果及分析图为从20例信号序列中抽取的一组典型信号用各种算法绘制的功率谱 图形。从图中可以看到,傅里叶变换方法得到的频谱(图,见下页)频率分 辨率较高,但是方差特性很差,Lomb周期图谱(图,见下页)也存在方差特 性差的问题;周期图法中不管是Welch方法还是多窗估计法(图),见下页, 它们对方差特性进行了改善,但是牺牲了频率分辨率,估计效果也不是很好; EMD分解方法可以获得不同频段的频域波形(图),并且频率分辨率高,并 且每个IMF分量频谱可以使用各种算
10、法进行,在进行频域参数的分析时可以 方便的获得各频段的功率,但是整体的目测效果不好,不容易对HRV频谱获 得直观的理解;相比之下,AR模型获得的功率谱图形(图),在频率分辨率 和方差特性上得到折中,功率谱图光滑,不管是分析还是目测效果都比较好, 是HRV谱估计的较好算法。252010Yule-WalkedBug笄、去 改进的樹方羞算注1050 10D.20.25030 35实际频率(Hz)MatlabifiiLaYeMl(上面是使用11阶AR模型对HRV信号的功率谱估计结果。)TOP3#大中小发表于2009-8-28 20:55只看该作者(下面的内容不是毕业设计中的了,是根据本贴的主题做的实验
11、)下面使用自编函数pburgwithcriterion函数预测AR模型系数) 复制内容到剪贴板代码:%功能:验证自编函数arburgwithcriterion估计AR模型阶次的准确性能%方法:使用已知的AR模型的输出,然后使用该函数预测AR参数进行比较% lskypaorig = ;for k = 2:10atemp = xlsread(,A num2str(k);%xls 文件第 1 列从第 2 行开始存储原始AR模型系数,后面进行估计if isempty(atemp)break;endaorig = aorig atemp;endx=randn(1,512);y=filter(1,aori
12、g,x);apredict,e,k,order_predict = arburgwithcriterion(y,FPE); xlswrite(,apredict,B2:B num2str(length(apredict)+1);% 图形显示figure (1)arorder = length(aorig);subplot(211);pburg(y,arorder);subplot(212);pburgwithcriterion(y);三KB)有些地方不是很完善,希望大家予以指正,上面的几个函数为根据 MATLAB信号处理工具箱中的几个函数改写,增加了需要的东西,如果有不 当的地方,还希望大家实验之中改正。26 祈才君数字信号处理技术的算法分析与应用M.北京:机械工业出版 社,2005:196-197.27 赵静,刘琦.Welch法谱估计和参数模型谱估计的MATLAB分析J.水利 电力机械,2006,28(4):48-50.28 陈国强,赵俊伟,黄俊杰,等基于Matlab的AR模型参数估计J.工具技 术,2005,39(4):39-40.