ar_smart_ensemble.m

第一次上载程序

clear;
%%%AR
% load   Eunitefhgz.mat    %%  Eunitef competition-data 2001
load  Eunitefhgz97.mat
e=0.001;                         %矩阵的系数，一般是小于1的数，但也不宜过小，可能导致计算过程发散
NN=30;                           %采样次数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%fhysh=fhysh/100;
for h=1:1:48
    for i=1:1:250
      y(h,i)=fhysh((i-1)*48+h);  % 按每天48小时组成矩阵
    end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
yy=y;                           
noise=randn(1,20*250);           %产生随机数
for n=1:1:20                     %循环20次取平均


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for d=1:1:48                     %每天48小时


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%    
 nnn=noise((n-1)*250+1:n*250);   %分20段
 std2=std(nnn);                  %取标准差    
 yyy=diff(yy(d,:));              %把负荷序列做差分
 
 std1=std(yyy);                  %求差分完的负荷序列标准差
 std3=std1*0.02;                 %信噪比
 beishu=std2/std3;               
 nnn=nnn/beishu;                 %相应信噪比的噪声
 y(d,:)=y(d,:)+nnn;   
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%    
  p=5;                           %AR模型的阶数
  start=16;                      %序列的起始点
  Y=0;                           %初始化成0，为了下一次循环进入
  canshu=0;                      %初始化
  P=0;                           %初始化
  P=1/e*eye(p,p);                %初始P矩阵                
  
  for i=1:1:p                    %初始化 canshu(:,1);
    canshu(i,1)=0;
  end
  
  for N=1:1:NN                   % 采样次数
    
    for i=1:1:p
      Y(i,N+1)=y(d,start+(N-1)-i+1);%%% 顺序是从前向后,最后一个数是114，滚动的取值，每循环一次N，多一个
    end
    
    %%%%%%%%%%%%%%%%% ＃＃求最佳参数向量fai＃＃ %%%%%%%%%%%%%%%%%%%%%%%%%%%
    ZZ=1/(1+Y(:,N+1)'*P(:,:,N)*Y(:,N+1));%% 实数
    P(:,:,N+1)=P(:,:,N)-ZZ*P(:,:,N)*Y(:,N+1)*Y(:,N+1)'*P(:,:,N);  %
    canshu(:,N+1)=canshu(:,N)+P(:,:,N+1)*Y(:,N+1)*(y(d,start+(N-1)+p)-Y(:,N+1)'*canshu(:,N));%%Yt
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%   
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    for i=1:1:p  
        fai(p,N+1,i)=canshu(i,N+1);%% 记录一下，用于以后预报
    end
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

  end%% NN
  
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ＃＃ 预报 ＃＃ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    CC=0;
    DD=0;
    EE=0;

    for i=1:1:p
      CC=CC+fai(p,NN+1,i)*y(d,NN+start-i);  %%接近预报那天的前P个数来预报，不管优化的P是多大，都是预报最大P的后面的三天，而fai是倒着取的，假如fai比P小，就取后面几个数作为fai
    end
    yb(d,NN+start)=CC; 
    for i=2:1:p
      DD=DD+fai(p,NN+1,i)*y(d,NN+start+1-i);
    end
    yb(d,NN+start+1)=DD+fai(p,NN+1,1)*yb(d,NN+start);
    for i=3:1:p
      EE=EE+fai(p,NN+1,i)*y(d,NN+start+2-i);
    end
    yb(d,NN+start+2)=EE+fai(p,NN+1,1)*yb(d,NN+start)+fai(p,NN+1,2)*yb(d,NN+start+1);
   
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 
end%% d



%%%%%%%%%%%%%%%%%%%%%%%% ＃＃计算误差＃＃ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%   
for i=1:1:48
     wucha1(i)=(y(i,NN+start)-yb(i,NN+start))/y(i,NN+start)*100;
     wucha11(n,i)=wucha1(i);
end
for i=1:1:48
     wucha2(i)=(y(i,NN+start+1)-yb(i,NN+start+1))/y(i,NN+start+1)*100;
     wucha22(n,i)=wucha2(i);
end
for i=1:1:48
     wucha3(i)=(y(i,NN+start+2)-yb(i,NN+start+2))/y(i,NN+start+2)*100;
     wucha33(n,i)=wucha3(i);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


% wucha=[wucha1 wucha2 wucha3];

% i=1:1:48;plot(i,wucha1(i),i,wucha2(i),i,wucha3(i),'-');grid on;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end %n


%%%%%%%%%%%%%%%%%%%%%% ＃＃计算平均误差＃＃ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fff=0;
for i=1:1:48
  for n=1:1:20
     fff=fff+wucha11(n,i);
  end
     wucha111(i)=fff/20;
     fff=0;
end

for i=1:1:48
  for n=1:1:20
     fff=fff+wucha22(n,i);
  end
     wucha222(i)=fff/20;
     fff=0;
end

for i=1:1:48
  for n=1:1:20
     fff=fff+wucha33(n,i);
  end
     wucha333(i)=fff/20;
     fff=0;
end

i=1:1:48;plot(i,wucha111(i),'.-',i,wucha222(i),'.-',i,wucha333(i),'.-');grid on;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%