典型时间序列模型分析.doc

实验1 典型时间序列模型分析 1、实验目的 熟悉三种典型的时间序列模型：AR 模型，MA 模型与ARMA 模型，学会运用Matlab工具对对上述三种模型进行统计特性分析,通过对2 阶模型的仿真分析，探讨几种模型的适用范围，并且通过实验分析理论分析与实验结果之间的差异。 实验原理 AR 模型分析： 设有 AR(2)模型， X(n)=-0.3X(n-1)-0.5X(n-2)+W(n) 其中：W(n)是零均值正态白噪声，方差为4。 （1）用MATLAB 模拟产生X(n)的500 观测点的样本函数，并绘出波形 （2）用产生的500 个观测点估计X(n)的均值和方差 （3）画出理论的功率谱 （4）估计X(n)的相关函数和功率谱 【分析】给定二阶的AR 过程，可以用递推公式得出最终的输出序列。或者按照一个白噪声 通过线性系统的方式得到，这个系统的传递函数为： 这是一个全极点的滤波器，具有无限长的冲激响应。 对于功率谱，可以这样得到， 可以看出，完全由两个极点位置决定。 对于 AR 模型的自相关函数，有下面的公式： 这称为 Yule-Walker 方程，当相关长度大于p 时，由递推式求出： 这样，就可以求出理论的 AR 模型的自相关序列。 产生样本函数，并画出波形 题目中的 AR 过程相当于一个零均值正态白噪声通过线性系统后的输出，可以按照上面的方法进行描述。 clear all; b=[1]; a=[1 0.3 0.5]; % 由描述的差分方程，得到系统传递函数 h=impz(b,a,20); % 得到系统的单位冲激函数，在20 点处已经可以认为值是0 randn(state,0); w=normrnd(0,2,1,500); % 产生题设的白噪声随机序列，标准差为2 x=filter(b,a,w); % 通过线形系统，得到输出就是题目中要求的2 阶AR 过程 plot(x,r); ylabel(x(n)); title(邹先雄——产生的AR 随机序列); grid on; 得到的输出序列波形为： 估计均值和方差 可以首先计算出理论输出的均值和方差，得到，对于方差可以先求出理论自相 关输出，然后取零点的值。 并且，，带入有 在最大值处输出的功率，也就是方差，为 对实际数据进行估计，均值为mean(x)=-0.0703，而方差为var(x)=5.2795，两者合理论值吻合得比较好。 程序及运行结果图如下，其中y_mean表示均值，y_var表示方差。 画出理论的功率谱密度曲线 理论的功率谱为， 用下面的语句产生： delta=2*pi/1000; w_min=-pi; w_max=pi; Fs=1000; w=w_min:delta:w_max; % 得到数字域上的频率取样点，范围是[-pi,pi] Gx=4*(abs(1./(1+0.3*exp(-i*w)+0.5*exp(-2*i*w))).^2); % 计算出理论值 Gx=Gx/max(Gx); % 归一化处理 f=w*Fs/(2*pi); % 转化到模拟域上的频率 plot(f,Gx); title(邹先雄——理论功率谱密度曲线); grid on; 得到的图形为： 可以看出，这个系统是带通系统。 估计自相关函数和功率谱密度 用实际数据估计自相关函数和功率谱的方法前面已经讨论过，在这里仅给出最后的仿真图形。 Mlag=20; % 定义最大自相关长度 Rx=xcorr(x,Mlag,coeff); m=-Mlag:Mlag; stem(m,Rx,r.); title(邹先雄——自相关函数); 最终的值为 可以看出，它和上面的理论输出值吻合程度很好。实际的功率谱密度可以用类似于上面的方法进行估计， window=hamming(20); % 采用hanmming 窗，长度为20 noverlap=10; % 重叠的点数 Nfft=512; % 做FFT 的点数 Fs=1000; % 采样频率，为1000Hz b=[1]; a=[1 0.3 0.5]; % 由描述的差分方程，得到系统传递函数 h=impz(b,a,20); % 得到系统的单位冲激函数，在20 点处已经可以认为值是0 randn(sta