function[yucezhi]=shijianxulie(data,cishu);[hanglie]=size(data);yucezhi=zeros(hang,lie);z=zeros(hang,cishu);fornt=1:hangx=data(nt,:);xx=x;x=[0diff(x)];[NnNNn]=size(x);%%%%%%用AR模型forn=1:sqrt(Nn)x2=x((n+1):Nn);forj=1:(Nn-n)forji=1:nx1(j,ji)=x(n-ji+j);endenddd=pinv((x1'*x1))*x1'*x2';p=dd;%%%%%%%%p为行向量forj=n+1:Nnforii=1:nax(ii)=p(ii)*x(j-ii);%%%%拟合出的xtendct(j)=x(j)-sum(ax);%%%%%%%%白噪声序列endt=[];t=ct(n+1:Nn).^2;b(n)=sum(t)/(Nn-n);AIC(n)=log(b(n))+2*n/Nn;%%%%%%%%AIC准则x1=[];x2=[];end%%%%%%%%%%%%%%%%%求AIC的最小值pp=AIC(1);nn=1;fori=1:sqrt(Nn)ifppAIC(i)pp=AIC(i);nn=i;endend%%%%%%%%%%%%%%%AIC中最小值为pp,回归节数为nn%%%%%%%%%求白噪声序列%%%%%%forj=nn+1:Nnforii=1:nnax(ii)=p(ii)*x(j-ii);endct(j)=x(j)-sum(ax);%%%%%%%%%白噪声序列%%%%%%endat0(nn+1:Nn)=ct(nn+1:Nn);%%%%%%%%%%%%%估计参数aabbjies=12;ifjies=Nnjies=Nn-1;endforp=1:jiesforq=1:jieshh=[p,q,nn];L=max(hh);X=[];e1=[];YY=[];Y=[];m=[];n=[];aabb=[];pa=[];pb=[];forj=1:Nn-Lforkm=1:pX(j,km)=x(L+j-km);endendforj=1:Nn-Lforjm=1:qe1(j,jm)=at0(L-jm+j);endendfori=1:Nn-LYY(i)=x(L+i);endY=YY';XX=[Xe1];aabb=pinv(XX)*Y;pa=aabb(1:p,1);%%%%%%%%%%%pa为自回归系数pb=aabb(p+1:p+q,1);%%%%%%%%%%%pb为移动平均系数%%%%%%%%%%%%计算AICB(p,q)forj=L+1:Nnforii=1:pax(ii)=pa(ii)*x(j-ii);endforkk=1:qbx(kk)=pb(kk)*x(j-kk);endabtt(j)=x(j)-sum(ax)-sum(bx);%%%%%%白噪声序列%%%%%%endt=[];t=abtt(L+1:Nn).^2;ab(p,q)=sum(t)/(Nn-L);AICB(p,q)=log(ab(p,q))+2*(p+q)/Nn;endend%%%%%%%%%%%%%%%%求AICB(p,q)的最小值pppp=AICB(1,1);nna=1;nnb=1;forik=1:jiesforjk=1:jiesifppppAICB(ik,jk)pppp=AICB(ik,jk);nna=ik;nnb=jk;endendend%%%%%%%%%%%%%%%%确定模型为ARMA(p,q)%%%%%%%%%%求自回归系数和移动平均系数X=[];e1=[];p=nna;q=nnb;hh=[p,q];L=max(hh);forj=1:Nn-Lforkm=1:pX(j,km)=x(L+j-km);endendforj=1:Nn-Lforjm=1:qe1(j,jm)=at0(L-jm+j);endendfori=1:Nn-LYY(i)=x(L+i);endY=YY';XX=[Xe1];aabb=pinv(XX)*Y;abc=aabb;pma=[];pmb=[];pma=abc(1:p,1)';%%%%%%%%%自回归系数(行向量pmb=abc(p+1:p+q,1)';%%%%%%%%%%移动平均系数(行向量%%%%%%%%%%%模型ARMA(p,q)%%%%%%自回归系数pma%%%%%%%移动平均系数pmb%%%%%signa=sum(Y-X*pma'-e1*pmb').^2/(Nn-L)%%%%%%%%%%%%%%%预测%%%%%%%%%%%%%%%yyuce=[x(1:L)';X*pma'+e1*pmb'];ee=(Y-X*pma'-e1*pmb');ee=[zeros(L,1);ee;zeros(cishu,1)]';%%%%%%%白噪声(L+1:Nn)(行向量c=zeros(1,p);cc=zeros(1,q);fori=1:pc(i)=x(Nn-i+1);endfori=1:qcc(i)=ee(Nn-i+1);endz1=c*pma'+cc*pmb';z(nt,1)=z1;fork=2:cishuifp==1z2=z1;z1=z2*pma'+cc*pmb';endifp1fori=1:p-1c(p-i+1)=c(p-i);endc(1)=z1;fori=1:qcc(i)=ee(Nn-i+k);endz1=c*pma'+cc*pmb';endz(nt,k)=z1;endz(nt,1)=z(nt,1)+xx(Nn);fork=2:cishuz(nt,k)=z(nt,k)+z(nt,k-1);endyyuce=(yyuce)';yuce(1)=xx(1);fori=1:Nn-1yuce(i+1)=yyuce(i)+xx(i+1);endyucezhi(nt,:)=yuce;endyucezhi=[yucezhiz];end