模型识别——确定线性模型的类别和阶数

ling0

1. function [p,Q]=mxsb(Z)
   
2. %《随机过程》(汪荣鑫编)P132.
   
3. %Z为时间序列的向量,p为自相关系数向量，Q为偏相关系数向量
   
4. %蓝色线是自相关系数向量p的图像；红色线是偏相关系数向量Q的图像
   
5. %书上例题1中的Z=[15600,8960,10400,10600,10800,9880,9850,10900,8810,9960,12200,7510,8640,6380,6810,8820,14400,7440,7240,6430,11000,7340,9260,5290,9130,7480,6980,9650,7260,8750,9900,7310,9040,7310,8850,7840,10700,6190,9610,7580,9990,6150,8250,6030,8980,6180,9630,9490,2340,11100,5090,10900,6490,12600,6640,7430,6760,10000,9300];
   
6. W=Z-mean(Z);
   
7. n=length(Z);K=15;
   
8. r=ones(1,K);p=ones(1,K);
   
9. r0=W(1:n)*W(1:n)'/n;
   
10. for k=1:K
    
11.     r(k)=W(1:n-k)*W(1+k:n)'/n;
    
12.     p(k)=r(k)/r0;
    
13. end
    
14. q=ones(K,K);Q=ones(1,K);
    
15. q(1,1)=p(1);Q(1)=q(1,1);
    
16. for k=1:K-1
    
17.     q(k+1,k+1)=(p(k+1)-p(k:-1:1)*q(k,1:k)')/(1-p(1:k)*q(k,1:k)');
    
18.     for j=1:k
    
19.         q(k+1,j)=q(k,j)-q(k+1,k+1)*q(k,k-j+1);         
    
20.     end
    
21.     Q(k+1)=q(k+1,k+1);
    
22. end
    
23. plot(0:K,[1,p])
    
24. hold on
    
25. plot(0:K,[1,Q],'r:p')
    
26. legend('p(k)','Q(k)')
    
27. title('模型识别——确定线性模型的类别和阶数')
    
28. xlabel('k')
    
29. ylabel('p(k),Q(k)')
    
30. grid on   

复制代码