mybp.m

自己编写的小波神经网络的代码。已经得到成功的应用

p=[37.07 41	41.25 38.34	38.73 39.13	38.73 40.46	38.12 39.49	38.79 39.72	42.46 42.93	41.93 42.71	42.83 45.55	44.54 44.7 44.04 42.04 45.02 45.67 45.64 46.53 47.18 45.59]';
x1=p(1:2);y1=p(3);
x2=p(2:3);y2=p(4);
x3=p(3:4);y3=p(5);
x4=p(4:5);y4=p(6);
x5=p(5:6);y5=p(7);
x6=p(6:7);y6=p(8);
x7=p(7:8);y7=p(9);
% x=p(8:9);y8=p(10);
% x=p(9:10);y9=p(11);
% x=p(10:11);y10=p(12);
% x=p(11:12);y11=p(13);
% x=p(12:13);y12=p(14);
% x=p(13:14);y=p(15);
% x=p(14:15);y=p(16);
% x=p(15:16);y=p(17);
% x=p(17:18);y=p(19);
% x=p(18:19);y=p(20);
% x=p(19:20);y=p(21);
x=[x1 x2 x3 x4 x5 x6 x7];
y=[y1 y2 y3 y4 y5 y6 y7];
%%%%%%进行输入数据的归一化
[xn,meanx,stdx,yn,meany,stdy]=prestd(x,y);
% [x1n,meanx1,stdx1]=prestd(x1);%列向量
% [x2n,meanx2,stdx2]=prestd(x2);
% [x3n,meanx3,stdx3]=prestd(x3);
% [x4n,meanx4,stdx4]=prestd(x4);
% [x5n,meanx5,stdx5]=prestd(x5);
% [x6n,meanx6,stdx6]=prestd(x6);
% [x7n,meanx7,stdx7]=prestd(x7);
%[x8n,meanx8,stdx8]=prestd(x8);
%[x9n,meanx9,stdx9]=prestd(x9);
%[x10n,meanx10,stdx10]=prestd(x10);
%[x11n,meanx11,stdx11]=prestd(x11);
%xn=[x1n x2n x3n x4n x5n x6n x7n];
I=size(xn,1);P=size(xn,2);
%%%%%%%%%%%%使目标向量归一化
% yn=y./29;
%%%%%%%%%%%%给出学习率、动量因子等参数
eta=0.02;%学习率取的适当可以加快误差的减小
aerfa=0.965;%动量因子？？？？？？？？？？
max_perf_inc=1.04;%????????
lr_dec=0.7;
lr_inc=1.05;c=0.5;%%%%%%%%%伸缩引子
%假设小波函数节点数为：H个；样本数为P;
%输出节点数为:J个；输入节点数为：I个；
H=2*I+1;
J=1;
%%%%%%%%%%%%%%%%%%%%初始化权系数,阈值，小波参数；
whi=rand(I,H);%?????
wjh=rand(H,J);
b1=rand(H,1);
b=rand(H,1);
a=rand(H,1);
b2=rand(J);
%%%%%%%%%%%保存的误差；
for p=1:P
     Err_NetOut=[];
     flag(p)=1;count(p)=0;
     while flag(p)>0  %进入小波神经网络的运行过程
           flag(p)=0;
           count(p)=count(p)+1;
           xhp1(p)=0;
           %%%%%%%%%%%求隐层的输入值
           for h=1:H
               for i=1:I
                   xhp1(p)=xhp1(p)+whi(i,h)*xn(i,p);%求权值和输入值的乘积和
               end
               ixhp(h,1)=xhp1(p)+b1(h,1);%求隐层的输入值
               xhp1(p)=0;
           end
           %%%%%%%%%%%求隐层的输出值
           for h=1:H
               oxhp(h,1)=fai((ixhp(h,1)-b(h,1))/a(h,1));%将小波函数(mallet)作为隐层的传递函数求隐层的输出值
           end
           %%%%%%%%%%%求输出层的输入值
           ixjp1=0;
           for h=1:H
               ixjp1=ixjp1+wjh(h,1)*oxhp(h,1);
           end
           ixjp=ixjp1+b2(1);%求输出层的输入值
           %%%%%%%%%%求输出层的输出值
           ixjp2=c.*ixjp;
           oxjp(p)=fnn(ixjp2);   %求输出层的输出（运用的是sigmoid函数(fnn)加了一个伸缩因子c）??????????????????
           %%%%%%%%%%%求误差并保存
           wuchayy=1/2*sumsqr(oxjp(p)-yn(p));%求输出值与目标值差值的平方和
           Err_NetOut=[Err_NetOut wuchayy];%保存每次的误差;
           %调节学习率eta==================================================
           if count(p)==1
              eta=eta;
           else 
                   if Err_NetOut(count(p))/Err_NetOut(count(p)-1)>max_perf_inc
                         eta=eta*lr_dec;
                   else eta=eta*lr_inc;
                   end
          end    
          %%%%%%%%%%%%%求网络参数的梯度值
          detaj=-c*(oxjp(p)-yn(p))*oxjp(p)*(1-oxjp(p));%在误差函数中对输出层的输入值求偏导并取负值????????
          for h=1:H
              detawjh(h,1)=eta*detaj*oxhp(h,1);%误差函数对隐层和输出层之间的权值的偏导
          end
          detab2=eta*detaj;%误差函数对输出层的阈值的偏导
          for h=1:H
               detah(h,1)=detaj*wjh(h,1)*diffai((ixhp(h,1)-b(h,1))/a(h,1))/a(h,1);%difmallet()函数是mallet（）函数的导数，该行是求误差函数对隐层输入值的偏导
          end
          for h=1:H
              for i=1:I
                  detawhi(i,h)=eta*detah(h,1)*xn(i,p);%求误差函数对输入层和隐层之间的权值的偏导
              end
          end
          detab1=eta*detah;
          detab=-eta*detah;%求误差函数对小波函数中的平移因子的偏导
          for h=1:H
               detaa(h,1)=-eta*detah(h,1)*((ixhp(h,1)-b(h,1))/a(h,1));%求误差函数对隐层的小波函数中的收缩因子的偏导
          end
          detac=ixjp*detaj/c
          %%%%%%%%%%%%%%%%%%引入动量因子aerfa，修正各个系数==========================================
          wjh=wjh+(1+aerfa)*detawjh;
          whi=whi+(1+aerfa)*detawhi;
          a=a+(1+aerfa)*detaa;
          b=b+(1+aerfa)*detab;
          b1=b1+(1+aerfa)*detab1;
          b2=b2+(1+aerfa)*detab2;
          c=c+(1+aerfa)*detac;
          %%%%%%%%%%%%%%%%%%%将输出恢复到原始范围=======================================
          %oxjpn(p)=oxjp(p)*29;
          %%%%%%%%%%%%%%%%%%网络循环计算的判断条件
          if Err_NetOut(count(p))>1e-8
               flag(p)=flag(p)+1;
          else
%                 figure;
%                 plot(Err_NetOut);
%                 title('误差曲线x(p)');
%                 disp('目标达到');
%                 disp(oxjpn(p)); 
%                 disp(count(p));
                break;%跳出while循环
          end
          if count(p)>2000
                figure;
                plot(Err_NetOut);
                title('误差曲线');
                disp('目标未达到');
             %   disp(oxjpn(p)); 
               break;
          end
         %break;  %跳出避免死循环
        % disp(oxjpn(p));
     end
end