AR预测模型

leadersheng

AR预测模型
% Akaike信息检验准则率定模型阶数p
% 2013/6/14
% mingbo1501

%%
clear,clc,close all;

p=12;%模型阶数
num=500;%检验个数

%----------------生成样本--------------------%
data0=load('runoff.txt');
[m,n]=size(data0);
data=zeros(m-p,p+1);
for i=1:m-p
    data(i,:)=data0(i:i+p);%前p列输入，第p+1列输出
end
[m1,n1]=size(data);
X1=data(1:m1-num,1:p);X1=[ones(m1-num,1),X1];
Y1=data(1:m1-num,p+1);%训练样本
X2=data(m1-num+1:end,1:p);X2=[ones(num,1),X2];
Y2=data(m1-num+1:end,p+1);%检验样本

%-----------------估计参数--------------------%
fai=(X1'*X1)\X1'*Y1;%最小二乘估计
yuce1=X1*fai;
err1=Y1-X1*fai;%残差
sigma=var(err1);%残差的方差

figure
plot(1:m1-num,yuce1,'b',1:m1-num,Y1,'r');
legend('预测值','实测值');
rate1=length(find((yuce1-Y1)./Y1*100<10))/(m1-num)*100;%合格率
DC1=1-sum((Y1-yuce1).^2)./sum((Y1-mean(Y1)).^2);%确定性系数

%-----------------预报精度评定--------------------%
Y_rand=zeros(1,num);
for i=1:num
    Y_rand(i)=normrnd(0,sqrt(sigma));%随机项
end
yuce2=X2*fai;
figure
plot(1:num,yuce2,'b',1:num,Y2,'r');
legend('预测值','实测值');
rate2=length(find((yuce2-Y2)./Y2*100<10))/num*100;%合格率
DC2=1-sum((Y2-yuce2).^2)./sum((Y2-mean(Y2)).^2);%确定性系数

sciren

预测啥的  。。。。能用吗