大连理工大学时间序列作业.docx
- 文档编号:3330154
- 上传时间:2023-05-05
- 格式:DOCX
- 页数:19
- 大小:515.52KB
大连理工大学时间序列作业.docx
《大连理工大学时间序列作业.docx》由会员分享,可在线阅读,更多相关《大连理工大学时间序列作业.docx(19页珍藏版)》请在冰点文库上搜索。
大连理工大学时间序列作业
时间序列分析实验
一、时序数据表(去均值)
此处只列出部分时序数据,如下:
-0.08959
-0.01775
0.0625
0.08103
-0.02341
0.02499
-0.07067
-0.07357
0.01754
0.05046
0.02298
-0.06761
0.14238
-0.00749
0.01607
-0.05038
-0.02956
0.02935
0.05161
0.03473
-0.17435
0.01329
-0.0503
0.02731
0.04063
-0.10376
-0.06801
0.01264
-0.11361
0.01586
0.06843
-0.03552
0.03008
0.00465
-0.04569
0.02872
-0.10188
0.11759
-0.10946
0.00083
0.0836
-0.03002
-0.05888
0.01277
0.02264
0.02442
0.05563
0.01102
-0.04836
-0.01016
-0.06611
-0.00932
-0.00739
0.01196
-0.09943
0.04133
0.02473
0.08021
0.04504
-0.0499
-0.06432
-0.08612
0.01965
-0.02935
0.01568
0.05521
0.00598
0.01475
0.068
0.00074
0.07667
0.07129
0.00556
0.01817
-0.08356
0.05343
-0.03279
0.01009
-0.01803
0.00839
0.00781
0.00352
-0.08454
0.16242
0.04453
-0.06626
-0.00061
0.07302
-0.06301
-0.06835
0.03202
0.03051
0.00245
0.05025
0.0347
0.03423
0.0537
-0.06215
-0.06907
0.0549
-0.10403
0.09124
0.06894
0.12646
-0.0923
0.0188
-0.07684
0.01556
0.06613
-0.04629
0.03165
0.03763
-0.04739
-0.0387
-0.02969
-0.02256
0.01459
-0.02837
-0.03985
0.05124
-0.07758
0.01995
-0.08152
-0.02503
0.08058
-0.10023
-0.03022
0.02673
-0.00058
0.00804
0.02834
0.07375
-0.12752
-0.03832
-0.02946
0.02818
-0.03897
0.00301
0.01513
-0.04519
0.01747
-0.02511
0.03847
-0.03959
-0.03026
-0.02406
-0.04652
0.03301
0.01706
-0.02618
-0.03203
0.04173
0.10232
-0.07372
-0.01493
0.01614
0.04232
0.00101
0.00996
0.03327
-0.11276
-0.01314
0.16323
0.00643
-0.01035
0.01119
0.01544
0.07149
-0.00377
0.0476
0.15212
0.03364
-0.0602
-0.03387
0.00612
-0.05429
-0.04
-0.01769
-0.06824
0.0281
-0.02903
0.04163
-0.01799
0.02477
-0.03745
-0.00702
0.0415
-0.04943
0.02159
0.03184
0.12661
-0.07671
0.05233
-0.07138
-0.07771
0.08634
-0.09408
0.04261
0.03566
0.03936
-0.06293
-0.01411
-0.00687
-0.09608
-0.08268
-0.02495
-0.01931
-0.00072
0.1062
-0.00436
0.00466
-0.09427
-0.10789
0.02254
0.0492
-0.03173
0.05783
0.01489
0.10924
-0.1076
0.02668
-0.01664
-0.1014
0.01382
0.03352
-0.01419
0.00827
0.09193
-0.0703
0.09211
-0.01847
0.03345
-0.00607
0.03516
0.00086
-0.03282
0.08057
0.00542
0.06062
-0.09209
0.01551
-0.0785
0.04286
0.03583
0.0619
0.02313
0.04226
0.07538
0.03936
-0.08119
0.00931
0.04798
0.03315
-0.00165
-0.00221
-0.0479
-0.07397
0.03803
-0.12415
0.03993
0.02224
0.00792
0.00632
0.06704
-0.0104
-0.02566
0.02387
0.04563
0.05906
-0.08543
0.08801
0.04947
0.00385
0.00583
0.04095
0.10594
-0.02018
0.01862
-0.00574
0.02496
二、时序图形(去均值)
Matlab程序如下:
%去均值
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
clc
closeallhidden
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%循环读入数据
s=1;
fid11=fopen(['D:
\matlab2009\work\时间序列数据\Time(Signal1)-Input0',num2str(s),'.txt'],'r');
i=1;
while~feof(fid11)
tline1=fgetl(fid11);
ifdouble(tline1
(1))>=48&&double(tline1
(1)<=57)
tline11=str2num(tline1);
t(i)=tline11
(2);
a(i)=tline11(3);
i=i+1;
continue
end
end
fclose(fid11);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
NN=length(aa);
k=1;
forj=1:
NN
ifrem(j,4)==0
t(k)=tt(j);
a(k)=aa(j);
k=k+1;
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
N=length(a);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%去均值
d=mean(a,2);
x=a-d;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%绘制去均值和未去均值时的振动信号图
figure
(1)
plot(t,a);
xlabel('时间t');
ylabel('加速度a');
title('未去均值的振动信号');
gridon
figure
(2)
plot(t,x);
xlabel('时间t');
ylabel('加速度a');
title('去均值的振动信号');
gridon
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%以写的方式建立一个新文件,存储原始的数据
fid1=fopen(['yuanshixushujubiao',num2str(s),'.txt'],'w');
%输出原始的数据
fork=1:
N;
%每行输出1个实型的原始数据
fprintf(fid1,'%10.5f\n',x(1,k));
end
fclose(fid1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%以写的方式建立一个新文件,存储去均值后的数据
fid2=fopen(['qujunzhishixushujubiao',num2str(s),'.txt'],'w');
%输出去均值后的数据
fork=1:
N;
%每行输出1个实型数据,去均值振动信号的加速度
fprintf(fid2,'%10.5f\n',a(k));
end
fclose(fid2);
利用matlab读出数据,并显示出其原始时序图和去均值后的时序图,如下:
未去均值的原始数据的时序图
去均值后的时序图
三、用最小二乘法建立AR(n)模型,并利用BIC准则对模型进行适应性检查
最小二乘法公式为:
即
BIC准则公式为:
BIC(n)=Nln
,使得BIC(n)的值为最小即可。
Matlab程序如下:
%最小二乘估计法估计AR模型参数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
clc
closeallhidden
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%读入去均值后的振动信号
s=1;
fid=fopen(['D:
\matlab2009\work\qujunzhishixushujubiao',num2str(s),'.txt'],'r');
h=fscanf(fid,'%f');
fclose(fid);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%估计AR模型参数,建立AR模型
N=length(h);
forn=1:
20
y=h(n+1:
N);
%y=h;%构造矩阵y
%y(1:
n-1)=[];
x=zeros(N-n,n);
fori=1:
N-n%构造系数矩阵x
forj=1:
n
x(i,j)=h(n+i-j);
end
end
b=inv(x'*x)*x'*y;%估计参数φ
c{n,:
}=b;%将每次估计的参数φ存在二维c中
%A=y-x*b;%残差矩阵
%af=A'*A/(N-n);%残差的方差
suma=0;%估计残差a的方差
fort=n+1:
N
sumbh=0;
fori=1:
n
bh=b(i)*h(t-i);
sumbh=sumbh+bh;
end
a=h(t)-sumbh;
suma=suma+a.^2;
end
af=suma/(N-n);
d(n)=af;%将每次估计的方差存在数组d中
%模型的适用性检验
fpe(n)=(N+n)*af/(N-n);%FPE准则
aic(n)=N*log(af)+2*n;%AIC准则
bic(n)=N*log(af)+n*log(N);%BIC准则
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%绘制a的方差和BIC准则曲线
figure
(1)
plot(d);
xlabel('n的值');
ylabel('残差的方差');
title('估计的a的方差的曲线');
gridon
figure
(2)
plot(bic);
xlabel('n的值');
ylabel('BIC值');
title('BIC准则的曲线');
gridon
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%利用BIC准则判断AR(n)模型的阶数,及模型下的各参数,并显示结果
min=bic
(1);
forn=1:
20
ifmin>bic(n)
min=bic(n);
m=n;
end
end
fprintf('AR模型的阶数为%d\n',m);
fprintf('参数φ为\n');
fprintf('%10.5f\n',c{m,1});
fprintf('残差的方差为%10.5f',d(m));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%以写的方式建立一个新文件,存储AR模型的参数
fid=fopen(['AR模型参数',num2str(s),'.txt'],'w');
%输出AR模型的参数
fprintf(fid,'%d\n',m);
fprintf(fid,'%10.5f\n',d(m));
fprintf(fid,'%10.5f\n',c{m,1});
fclose(fid);
运行matlab程序,可得如下结果:
残差的方差曲线
BIC准则曲线
AR模型的阶数为15
参数φ为
[0.235180.05443-0.030540.18497-0.08442-0.196210.001950.02479
-0.02883-0.05662-0.117060.03522-0.105700.087100.21186],
残差的方差为0.00254
也就是说n=15最合适,此时参数φ和残差的方差的值如上所示。
四、给出所建的AR(n)模型
AR(n)模型中的各参数分别为:
n=15,
φ=[0.235180.05443-0.030540.18497-0.08442-0.196210.001950.02479
-0.02883-0.05662-0.117060.03522-0.105700.087100.21186],
。
将以上各参数代入就可求出AR模型。
五、进行AR谱分析
AR谱的公式为:
。
Matlab程序如下:
%进行AR谱分析
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
clc
closeallhidden
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%AR(n)模型的初始参数
s=1;
fid=fopen(['D:
\matlab2009\work\AR模型参数',num2str(s),'.txt'],'r');
h=fscanf(fid,'%f');
fclose(fid);
n=h
(1);
af=h
(2);
b=h(3:
(n+2));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%计算AR谱值
drt=2.5625e-004;
fs=6400;%采样频率
f=0:
6.25:
(fs-6.25);
w=2*pi*f;
sumf=0;
forj=1:
n
sumf=sumf+b(j)*exp(-w*j*drt*i);
end
sarws=af./((abs(1-sumf)).^2);
sarw=2*sarws;
N=length(sarw);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%绘制AR模型的谱图
figure;
plot(f(1:
N/2),sarw(1:
N/2));
xlabel('频率(Hz)');
ylabel('功率谱值');
title('AR谱图');
gridon;
运行matlab程序,可得AR谱图如下:
AR模型的AR谱图
六、进行周期图谱(FFT)分析
AR模型的周期图谱的公式为:
Matlab程序如下:
%进行周期图谱(FFT)分析
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
clc
closeallhidden
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%读入原始数据
s=1;
fid=fopen(['yuanshixushujubiao',num2str(s),'.txt'],'r');
h=fscanf(fid,'%f');
fclose(fid);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%计算周期图谱
N=length(h);
fs=6400;%采样频率
spews=fft(h,N);%原信号的FFT变换
spew=2*abs(spews).^2/(2*pi*N);
f=0:
fs/N:
(fs-fs/N);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%绘制周期图谱
figure;
plot(f(1:
N/2),spew(1:
N/2));
xlabel('频率(Hz)');
ylabel('周期图谱值');
title('原始数据的周期图谱');
gridon;
运行matlab程序,得AR模型的周期图谱如下:
AR模型的周期图谱
比较:
从AR谱和周期图谱两图中可以看出,AR谱比周期图谱平滑,且无毛刺,谱峰突出。
七、程序流程图
八、计算Green函数
利用递推算式计算Green函数:
Matlab程序如下:
%计算Green函数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
clc
closeallhidden
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%AR(n)模型的初始参数
s=1;
fid=fopen(['D:
\matlab2009\work\AR模型参数',num2str(s),'.txt'],'r');
h=fscanf(fid,'%f');
fclose(fid);
n=h
(1);
b=h(3:
(n+2));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%利用递推法算式计算Green函数
G0=1;
G1=b
(1)*G0;
G(1:
40)=0;
G
(1)=G0;
G
(2)=G1;
forj=3:
40
fori=1:
j-1
ifj<=n
G(j)=G(j)+b(i)*G(j-i);
end
end
ifj>n
fori=1:
n
G(j)=G(j)+b(i)*G(j-i);
end
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%绘制Green函数
stem(G);
xlabel('n的值');ylabel('G的值');
title('Green函数');
gridon
运行matlab程序,可得Green函数如下:
AR模型的Green函数
Green函数为G=
1.00000.23520.10970.00810.1857-0.0002-0.1859
-0.10140.0034-0.0391-0.1374-0.17050.0051-0.1080
0.00720.18550.16820.09300.05520.12480.0553
-0.0423-0.03990.0029-0.0022-0.0982-0.0943-0.0411
-0.0569-0.04430.00760.06010.04750.03610.0601
0.05360.0207-0.00070.01520.0107
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 大连理工大学 时间 序列 作业