Yule-Walker方程

实验报告课程名称：学院：专业：年级：班级：姓名：学号：指导教师：实验五Yule-Walker方程（一）实验目的学习求解Yule-Walker方程，建立随机信号的AR模型。（二）实验原理随机信号可以看作是由当前激励白噪声w(n)以及若干次以往信号x(n-k)的线性组合产生，即所谓自回归模型（AR模型）pkkkpkkzazHknxanwnx1111模型参数满足Yule-Walker方程01mkmRamRpkxxkxx矩阵形式pRRRaaaRpRpRpRRRpRRRp2102120111021求解Yule-Walker方程，就可以得到AR模型系数当模型阶次较大时，直接用矩阵运算求解的计算量大，不利于实时运算。利用系数矩阵的特性，人们提出了如L-D算法等快速算法。（三）实验内容和步骤编写求解Yule-Walker方程的程序，并对实际生理信号（例如脑电）建立AR模型。对同一数据，使用matlab信号处理工具箱自带函数aryule计算相同阶数AR模型系数，检验程序是否正确。用伪随机序列（白噪声）驱动AR模型，观察输出是否与真实脑电信号相似，对比真实信号与仿真信号的功率谱。流程图：实验流程图如下：输入信号计算x信号的自相关函数和Toeplitz型矩阵L-D算法，（aryule函数或用arburg）Yule-Woker方程参数估计向量参数估计向量a两种不同类型方法的参数估计值比较向量da判断自编函数是否正确利用L-D算法或者Burg算法产生的参数变量进行脑电信号、心电信号处理输出信号最小均方误差E计算最小均方误差选定AR模型的阶数p源程序：clear;clc;loadeegdata;x=eegdata(1:1024);%长度可以任意选择，但信号越长计算量越大loadecgdata;x=ecgdata(1:1024);%长度可以任意选择，但信号越长计算量越大p=12;%尝试改变模型阶数，观察效果Rxx=xcorr(x,'biased');Rtemp=zeros(1,p);Rl=zeros(p,1);fork=1:length(Rtemp)Rtemp(k)=Rxx(length(x)-1+k);Rl(k)=Rxx(length(x)+k);endRs=toeplitz(Rtemp);%生成自相关系数矩阵(Toeplitz型)A=-inv(Rs)*Rl;%AR模型系数估计Sw=[Rtemp(1),Rl']*[1;A];%白噪声方差估计%采用malab自带函数估计模型系数[a,E]=aryule(x,p);%a--系数，E--预测误差，k--反射系数da=a(2:end)-A'%自编程序求解是否正确？Stem(da);title(‘参数估计偏差’)w=randn(size(x));x2=filter(1,a,w);%仿真数据figure;subplot(3,1,1);plot(x);title('真实数据');subplot(3,1,2);plot(x2);title('仿真数据');error=mean((x2-x).^2);Rxx2=xcorr(x2,'biased');Px=abs(fft(Rxx));Px2=abs(fft(Rxx2));figure;subplot(2,1,1);plot(-1023:1023,Px);title('真实信号功率谱');subplot(2,1,2);plot(-1023:1023,Px2);title('仿真信号功率谱');%%绘制Sw、E随p变化散点图，以下是统计得到的结果p=[8910111213141516];Sw=[1.05271.03011.01500.99610.99430.99420.98960.98770.9872];E=[3.87413.43373.32003.73143.36833.48093.58263.62533.4176];subplot(2,1,1);stem(p,Sw);xlabel('p');ylabel('Sw');title('Sw随p变化散点图');subplot(2,1,2);plot(p,E,'-o');xlabel('p');ylabel('E');title('E随p变化散点图');实验结果：1、p=12时，比较心电信号与脑电信号的da值，并作出散点图。024681012-2.5-2-1.5-1-0.500.51x10-16参数估计偏差024681012-1.5-1-0.500.511.522.5x10-13参数估计偏差脑电信号心电信号结果分析：由以上两个图可以看出脑电信号的参数估计偏差比心电信号的参数估计偏差要小许多，脑电信号的自编程序跟MATLAB信号处理工具箱自带函数aryule的处理更接近。由此可知，L-D算法与自编程序相比较，自编程序对估计的参数比较精确一些。2、p=12时，比较心电信号与脑电信号的真实数据与仿真数据020040060080010001200-10-50510真实数据020040060080010001200-50510仿真数据020040060080010001200-1-0.500.51真实数据020040060080010001200-2002040仿真数据脑电信号心电信号结果分析：由上面两幅图可以看出，对心电信号与脑电信号建模后，心电信号的真实数据和仿真数据相差很大，该模型产生的信号更能真实的反映脑电信号的特征。不同噪声的激励得到的信号时不同的。3、p=12时，比较心电信号与脑电信号的功率谱，-1500-1000-500050010001500050100150真实数据功率谱-1500-1000-5000500100015000204060仿真数据功率谱-1500-1000-500050010001500050100150真实数据功率谱-1500-1000-5000500100015000123x104仿真数据功率谱脑电信号心电信号结果分析：在比较脑电信号与心电信号的图之后不难发现，对信号进行功率谱估计时，脑电信号的功率谱的纵坐标相差不大，但是相对而说心电信号的相差是极大的。4、L-D算法与Burg算法比较4.1、脑电信号L-D算法与Burg算法比较020040060080010001200-10-50510真实数据020040060080010001200-10-505仿真数据020040060080010001200-10-50510真实数据020040060080010001200-4-2024仿真数据L-D算法、error=3.8118Burg算法、error=4.0522024681012-2.5-2-1.5-1-0.500.51x10-16参数估计偏差024681012-1.5-1-0.500.51x10-3参数估计偏差Da值L-D-0.2220-0.09710.0312-0.0555-0.05550.06940.0833-0.0104-0.06250.0347-0.1110-0.0833Burg0.0007-0.0000-0.00020.0001-0.00100.00060.0009-0.0004-0.0011-0.00100.0006-0.0008结果分析：当AR模型的阶数p为一个确定值时，对于脑电信号，L-D算法的均方误差比起arburg算法的均方误差要小，但是它的噪声方差估计不变，4.2、心电信号L-D算法与Burg算法比较020040060080010001200-1-0.500.51真实数据020040060080010001200-60-40-20020仿真数据020040060080010001200-1-0.500.51真实数据020040060080010001200-50050仿真数据L-D算法、error=162.2773Burg算法、error=190.7024024681012-1.5-1-0.500.511.522.5x10-13参数估计偏差024681012-0.3-0.2-0.100.10.20.30.40.50.6参数估计偏差Da值L-D0.0235-0.05720.0490-0.14030.01180.00830.0136-0.01150.03100.2091-0.1227-0.0116Burg-0.28380.5068-0.0974-0.27940.10510.06590.0240-0.07920.03140.0121-0.02210.0149结果分析：由数据可知，心电信号的AR模型的阶数与均方误差和噪声方差估计之间没有必然的联系。与L-D算法相比，Burg算法的AR模型的均方误差比较小，这是因为Burg算法具有较高的计算效率，该算法从总体上选择所有的模型参数达到总的均方误差最小4.3、改变p值，比较L-D算法与Burg算法（脑电信号）为观察方便我做了下面的表格P值8910111213141516Burg算法Sw1.05271.03011.01500.99610.99430.99420.98960.98770.9872E3.87413.43373.32003.73143.36833.48093.58263.62533.4176L-D算法Sw1.05271.03011.01500.99610.99430.99420.98960.98770.9872E3.46873.20103.85303.38883.28743.39453.93154.15733.3942891011121314151600.511.5pSwSw随p变化散点图89101112131415163.23.43.63.84pEE随p变化散点图891011121314151600.511.5pSwSw随p变化散点图891011121314151633.544.5pEE随p变化散点图Burg算法L-D算法实验分析：观察图和表，不难发现，先不管是Burg算法还是L-D算法在改变了p的大小之后，噪声方差估计Sw和均方误差error都有在改变，基于我们取的值来看，噪声方差估计呈现一个先减小后增大的趋势，但是改变的很微弱，不明显。对于均方误差显得杂乱。不管是Burg算法还是L-D算法，噪声方差估计在p值不变的情况下是相同的，只是均方误差在一直改变，而白噪声方差估计值则随阶数的增大而减小，可发现当阶数为12时，白噪声方差估计最接近1，所以，p取12最为合理。四、思考题对ECG、EEG建模后，该模型产生的信号能否反映ECG、EEG信号的特征？答：通过上述各实验，可得出的结论为：心电信号模型产生的信号不能精确的反映ECG信号的特征，因为ECG信号有一定的规律性，而AR模型不能用在信号具有确定性的场合，对于脑电信号信号，因为EEG信号没有规律性，脑电信号可以由均值为0，方差为1的白噪声激励而成，所以它是不确定信号，因而在AR模型中得到它的误差比较小，所以AR模型基本能反应EEG信号的特征。五、实验总结最后一次实验了，还是感觉没有整体的学好这门课，也不知道下次还有没有机会再接触了。