wxyhome353.m

程序包含的设置包括：列数、样本总个数、建模样本数、预报因子数、预报对象数、学习因子、 动量因子、训练次数、总体误差、隐节点数。

%'因子数n1,预报对象数 n2,试报年数n,学习样本数 m,，mm为总样本数，学习因子a1，动量因子 a2,总体误差ee1
%个体误差ee2,迭代次数 num,隐层节点数 yj，
close all
clear  all
clc
[file,path]=uigetfile('*.*','打开数据文件');
prompt={'包括Y的列数','样本总个数','建模样本数','预报因子数','预报对象数','学习因子','动量因子','训练次数','总体误差','隐节点数'};
def={'4','45','40','3','1','0.9','0.7','5','0.001','3'};
dlgTitle='参数输入';
lineNo=1;
answer=inputdlg(prompt,dlgTitle,lineNo,def);
answer=char(answer);
h=str2num(answer(1,:));
mm=str2num(answer(2,:));
m=str2num(answer(3,:));
n1=str2num(answer(4,:));
n2=str2num(answer(5,:));
a1=str2num(answer(6,:));
a2=str2num(answer(7,:));
num=str2num(answer(8,:));
ee1=str2num(answer(9,:));
yj=str2num(answer(10,:));
filename1=strcat(path,file);
fid1=fopen(filename1,'r');
a=fscanf(fid1,'%f',h*mm);
fclose(fid1);
a=reshape(a,h,mm);
a=a';
p=a;
[file,path]=uiputfile('*.*','保存');
filename2=strcat(path,file);
fid2=fopen(filename2,'w');
n=mm-m;
hdata=p;
hda(3,(h))=0;
fdata(mm,(h))=0;                                                                                                                                                                                                                                                                                                                                                                           
s1=min(hdata(1:m,:));
s2=max(hdata(1:m,:));
hda(1,1:h)=mean(hdata(1:m,1:(h)));
mean(hdata(1:m,1:(h)));
hda(2,1:h)=(9*s2-s1)./8;
hda(3,1:h)=(9*s1-s2)./8;
for i=1:mm
    for j=1:(n1+n2)
        if hdata(i,j)~=-999
        fdata(i,j)=(hdata(i,j)-hda(3,j))/(hda(2,j)-hda(3,j));
        else
            fdata(i,j)=(hda(1,j)-hda(3,j))/(hda(2,j)-hda(3,j));
        end
    end
end
%-----------------------------------------fdata 中放的是标准化的资料，gdata 中放的是预报量资料
 gdata=fdata(1:mm,(n1+1):h)  
 fdata=fdata(:,1:n1) 
 ggdata=p(1:mm,(n1+1):h)
 %----------------------------------------------------------  
v=randn(n1,yj);
vc=randn(yj);
w=randn(yj,n2);
rc=randn(n2);
f=fdata;
g=gdata;
m_data(m+n,n2)=0;
n_dd=0;
ee=0.1;
eee=0
%计算
n_dd=1
bar=waitbar(0,'please wait....');
for n_dd=1:num
 waitbar(n_dd/num,bar);
        %给连接权付初值0
        v1=zeros(n1,yj);
        vc1=zeros(yj);
        w1=zeros(yj,n2);
        rc1=zeros(n2);
                                              
 for i=1:m          
    % '计算隐含层的激活值  
    for j=1:yj                             
       xx=0 ;                               
        for k=1:n1
          xx=xx+f(i,k)*v(k,j);  
        end
        xx=xx+vc(j);
       b(j)=1/(1+exp(-xx));
    end                                      
    % 计算输出层单元的激活值                                       
       for j=1:n2
       xx=0;
       for k=1:yj
       xx=xx+b(k)*w(k,j);
       end
       xx=xx+rc(j);
       c(j)=1/(1+exp(-xx));
       mdata(i,j)=c(j);
       end
    % 计算输出层单元的一般化误差
    for j=1:n2
        d(j)=c(j)*(1-c(j))*(g(i,j)-c(j));
    end
     % 计算隐含层对于每个dj误差
     for j=1:yj
         xx=0;
         for k=1:n2
             xx=xx+w(j,k)*d(k);
             
         end
         e(j)=b(j)*(1-b(j))*xx;
     end
      %调整隐含层单元 到输出层单元的连接权
      for j=1:yj
          for k=1:n2
              w1(j,k)=w1(j,k)+a1*b(j)*d(k);
          end
      end
       %调整 输出层单元的阈值
       for j=1:n2
         rc1(j)=rc1(j)+a1*d(j);
       end
      % 调整输入层单元 到隐含层单元的连接权
      for j=1:n1
          for k=1:yj
              v1(j,k)=v1(j,k)+a2*f(i,j)*e(k);
          end
      end
       % 计算隐含层单元的阈值
       for j=1:yj
           vc1(j)=vc1(j)+a2*e(j);
       end
 end
 %-------------------------------------  
   labuda=0.15;
   %调整连接权
   for j=1:n1
       for k=1:yj
           v(j,k)=v(j,k)+labuda*v1(j,k);
       end
   end
   for j=1:yj
       vc(j)=vc(j)+labuda*vc1(j);
   end
   for j=1:yj
       for k=1:n2
           w(j,k)=w(j,k)+labuda*w1(j,k);
       end
   end
   for j=1:n2
       rc(j)=rc(j)+labuda*rc1(j);
   end
   %计算总体误差
   
   ee=0;
   for i=1:m
       for j=1:n2
           ee=ee+(mdata(i,j)-g(i,j))^2/2;
       end
   end
   ee=ee/m;
   eee(n_dd)=ee;
   
if ee<ee1
   break
end
end
close(bar);
ee;
plot(eee);
fprintf(fid2,'包括预报量的列数：%1d,样本总数：%2d,建模数：%1d,预报因子数：%1d,预报量数：%1d,隐节点数：%1d\n',h,mm,m,n1,n2,yj);
fprintf(fid2,'训练次数：%d:,动量因子：%4.1f,学习因子:%4.1f,总体误差：%f\n',num,a1,a2,ee1);
fprintf(fid2,'收敛误差为：%8.6f\n\n',ee) ;
sprintf('收敛误差为：%12.6f',ee) 
n_dd;
%预报
for i=1:mm
    %计算隐含层新的激活值
   for j=1:yj
       xx=0;
       for k=1:n1
           xx=xx+f(i,k)*v(k,j);
       end
       xx=xx+vc(j);
       b(j)=1/(1+exp(-xx));
   end
   %
     for j=1:n2
       xx=0;
       for k=1:yj
           xx=xx+b(k)*w(k,j);
       end
       xx=xx+rc(j);
       c(j)=1/(1+exp(-xx));
      mdata(i,j)=c(j);
     end
 end

mdata;
  % transfering the original value
for i=1:mm
    for j=1:n2
        if mdata(i,j)~=-999
          mdata(i,j)=mdata(i,j)*(hda(2,h)-hda(3,h))+hda(3,h)  ;
    
      elseif mdata(i,j)==-999
            
            mdata(i,j)=(hda(1,j)-hda(3,h))/(hda(2,h)-hda(3,h));
        end
    end
end
mdata
fprintf(fid2,'predictand  fitting\n\n');
  fprintf(fid2,'              实况        预测       误差       相对误差\n');
for i=1:m
    for j=1:n2
        xd(i,j)=abs(ggdata(i,j)-mdata(i,j))/ggdata(i,j)*100;
        fprintf(fid2,'        %12.4f%12.4f%12.4f%12.4f%%',ggdata(i,j),mdata(i,j),(ggdata(i,j)-mdata(i,j)),abs(ggdata(i,j)-mdata(i,j))/ggdata(i,j)*100);
    end
    fprintf(fid2,'\n');
end
jd1=mean(abs(ggdata(1:m,:)-mdata(1:m,:)));
fprintf(fid2,'average %12.4f%12.4f%12.4f%12.4f%%\n\n',mean(ggdata(1:m,:)),mean(mdata(1:m,:)),jd1,mean(xd(1:m,:)));
    fprintf(fid2,'\npredictand  prediction\n\n');
for i=m+1:mm
    for j=1:n2
        xd(i,j)=abs(ggdata(i,j)-mdata(i,j))/ggdata(i,j)*100;
        fprintf(fid2,'        %12.4f%12.4f%12.4f%12.4f%%',ggdata(i,j),mdata(i,j),(ggdata(i,j)-mdata(i,j)),abs(ggdata(i,j)-mdata(i,j))/ggdata(i,j)*100);
    end
    fprintf(fid2,'\n');
end
    fprintf(fid2,'\n\n');    

jd2=mean(abs(ggdata((m+1):mm,:)-mdata((m+1):mm,:)));
       
       fprintf(fid2,'\n\n  average%12.4f%12.4f%12.4f%12.4f%%\n\n',mean(ggdata(1+m:mm,:)),mean(mdata(1+m:mm,:)),jd2,mean(xd(1+m:mm,1)));
       fprintf(fid2,'\n');
    
 fprintf(fid2,'\n\n');
 fprintf(fid2,'误差变化：\n');
 ne=length(eee)
for i=1:ne
fprintf(fid2,'%4d %12.6f\n',i,eee(i));
end
  fclose(fid2);