ar_year(old).m

用matlab做的多元线形回归预测的小程序

function [FlowReal,FlowFore,Error,p]=AR_year(year,num)
% 时间序列AR(p)    模型的日径流预报
% year            待预报日期
% num             历史数据样本容量
if  year<1984|year>2000
    input('时间选择错误，历史数据范围1984年~2000年');
    return
end
if (num<10|num>(year-1984-1))
    input('样本容量溢出，请更改样本容量值！');
    return
end
load Qx
Data=[];
for i=1:17    
    Data=[Data;[i,sum(Qx(find(Qx(:,2)==(1983+i)),5))]];
end
N=num; %数据个数
se2=year-1984-1;
se1=se2-num+1;
xt=Data(se1:se2,2);
OutCompare=[];fai=[];
for n=1:5
  Y=xt;
  Y(1:n)=[];
  m=N-n;
  X=[]; %构造系数矩阵x
  for i=1:m
    for j=1:n
      X(i,j)=xt(n+i-j);
    end
  end
  f=inv(X'*X)*(X'*Y);
%方差
  Delta(n)=(Y-X*f)'*(Y-X*f)/(N-n);
%FPE准则
  criterion(n,1)=(N+n)*Delta(n)/(N-n);
%AIC准则
%criterion(n,2)=N*log(Delta(n))+2*n;
%BIC准则
%criterion(n,3)=N*log(Delta(n))+n*log(log(N));
  fa=zeros(5,1);
  fa(1:n)=f;
  fai=[fai,fa];
  OutCompare=[OutCompare;[n,Delta(n),criterion(n,1)]];
end
p=OutCompare(find(OutCompare(:,3)==min(OutCompare(:,3))),1);
fai=fai(1:p,p) %AR（p）模型的系数fai
DataInput=xt(num-p+1:num);
FlowFore=DataInput'*fai
FlowReal=sum(Qx(find(Qx(:,2)==year),5))
Error=abs(FlowFore-FlowReal)/FlowReal