aprod2ann.m

书籍“Regularization tools for training large feed-forward neural networks using Automatic Differentiat

function [Af] =aprod2ANN(mode,m,n,e,iw,rw,w1,b1,Output_Layer2,Deriv_Layer2,...
	df1,w2,b2,Output_Layer3,Deriv_Layer3,df2,P,T)
%
%Call: [Af] =aprod2ANN(mode,m,n,e,iw,rw,w1,b1,Output_Layer2,Deriv_Layer2,...
	%df1,w2,b2,Output_Layer3,Deriv_Layer3,df2,P,T)
%This routine computes A*f or A'*f, where A is the Jacobian corresponding
%to a FFANN with one hidden layer.
%It should normally be used with the routine LSQR written by Paige
%and Saunders
%Note that if mode == 1 the input parameter e should be a column vector.
%If mode == 2 it is assumed that e is a Thrd*q matrix corresponding to the
%output layer.

f=e(:);
% If mode == 1 aprod2ANN returns y=A*f
% If mode == 2 aprod2ANN returns y=A'*f
%
[Fst,q]=size(P);	% # nodes in Input layer (first layer)
[Snd,r]=size(w1);	% # nodes in 1'st hidden layer (second layer)
if Fst ~= r		% q is # samples presented to the inputs
   error('In aprod2ANN: Fst .ne. r');
end
[Thrd,r]=size(w2);	% # nodes in output layer (third layer)
if Snd ~= r
   error('In aprod2ANN: Snd .ne. r')
end
%[Fth,r]=size(w3);	% # nodes in output layer (forth layer)

Om=Fst*Snd+Snd;		% # weights in first hidden layer
Nm=Om+Snd*Thrd+Thrd;	% # total weights
%Mm=Nm+Thrd*Fth+Fth;	% # total weights
lenf=max(size(f));
if mode == 1
  if Nm ~= lenf
     error('In aprod2ANN: Nm .ne. lenf (mode == 1)');
  else
     Af=zeros(q*Thrd,1);
  end
elseif (q*Thrd) ~= lenf
      error('In aprod2ANN: q*Thrd .ne. lenf (mode == 2)');
    else
      Af=zeros(Nm,1);
end
%
if mode == 1
TEMP=ones(1,q);
%Compute A*f by using forward mode
%time1=cputime;
%unpack the vector elements such that same kind of computation as
%in the function evaluation can be done.

%First hidden layer
%for k=1:Fst,
   %r1=(k-1)*Snd;
   %w1(:,k)=f(r1+1:r1+Snd);
%end
w1=reshape(f(1:Fst*Snd),Snd,Fst);
b1=f(Fst*Snd+1:Fst*Snd+Snd);
O_L2=b1*TEMP+w1*P;
S2=Deriv_Layer2 .* O_L2;

%Output layer
%for k=1:Snd,
   %r1=Om+(k-1)*Thrd;
   %w2temp(:,k)=f(r1+1:r1+Thrd);
%end
w2temp=reshape(f(Om+1:Om+Snd*Thrd),Thrd,Snd);
b2=f(Om+Snd*Thrd+1:Om+Snd*Thrd+Thrd);
O_L3=b2*TEMP+w2*S2+w2temp*Output_Layer2;
S3=Deriv_Layer3 .* O_L3;
Af=S3(:);

else % Now mode should be 2
   %Compute A'*f
   d2=feval(df2,Output_Layer3,e);
   d1=feval(df1,Output_Layer2,d2,w2);
   [dw1,db1]=learnbp(P,d1,1);
   [dw2,db2]=learnbp(Output_Layer2,d2,1);
   Af=[dw1(:), db1(:); dw2(:); db2(:)];
end %of if mode == 1
end