X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?p=CreaPhase.git;a=blobdiff_plain;f=octave_packages%2Fnnet-0.1.13%2F__calcjacobian.m;fp=octave_packages%2Fnnet-0.1.13%2F__calcjacobian.m;h=25470a2eff328d68e2ae62819d6cc9e756df275f;hp=0000000000000000000000000000000000000000;hb=c880e8788dfc484bf23ce13fa2787f2c6bca4863;hpb=1705066eceaaea976f010f669ce8e972f3734b05 diff --git a/octave_packages/nnet-0.1.13/__calcjacobian.m b/octave_packages/nnet-0.1.13/__calcjacobian.m new file mode 100644 index 0000000..25470a2 --- /dev/null +++ b/octave_packages/nnet-0.1.13/__calcjacobian.m @@ -0,0 +1,283 @@ +## Copyright (C) 2006 Michel D. Schmid +## +## +## This program is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 2, or (at your option) +## any later version. +## +## This program is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with this program; see the file COPYING. If not, see +## . + +## -*- texinfo -*- +## @deftypefn {Function File} {}@var{Jj} = __calcjacobian (@var{net},@var{Im},@var{Nn},@var{Aa},@var{vE}) +## This function calculates the jacobian matrix. It's used inside the +## Levenberg-Marquardt algorithm of the neural network toolbox. +## PLEASE DO NOT USE IT ELSEWEHRE, it proparly will not work! +## @end deftypefn + +## Author: Michel D. Schmid + + +function [Jj] = __calcjacobian(net,Im,Nn,Aa,vE) + + ## comment: + ## - return value Jj is jacobi matrix + ## for this calculation, see "Neural Network Design; Hagan, Demuth & Beale page 12-45" + + + ## check range of input arguments + error(nargchk(5,5,nargin)) + + ## get signals from inside the network + bias = net.b; + + ## calculate some help matrices + mInputWeight = net.IW{1} * Im; + nLayers = net.numLayers; + for i=2:nLayers + mLayerWeight{i,1} = net.LW{i,i-1} * Aa{i-1,1}; + endfor + + ## calculate number of columns and rows in jacobi matrix + ## firstly, number of columns + a = ones(nLayers+1,1); # +1 is for the input + a(1) = net.inputs{1}.size; + for iLayers = 1:nLayers + a(iLayers+1) = net.layers{iLayers}.size; + endfor + nColumnsJacobiMatrix = 0; + for iLayers = 1:nLayers + nColumnsJacobiMatrix = (a(iLayers)+1)*a(iLayers+1) + nColumnsJacobiMatrix; + endfor + ## secondly, number of rows + ve = vE{nLayers,1}; + nRowsJacobiMatrix = length(ve(:)); + + + ## FIRST STEP ----------------------------------------------------- + ## calculate the neuron outputs without the transfer function + ## - n1_1 = W^1*a_1^0+b^1: the ^x factor defined the xth train data set + ## the _x factor defines the layer + ## ********** this datas should be hold in Nn + ## ********** should be calculated in "__calcperf" + ## ********** so Nn{1} means hidden layer + ## ********** so Nn{2} means second hidden layer or output layer + ## ********** and so on ... + ## END FIRST STEP ------------------------------------------------- + + ## now we can rerange the signals ... this will be done only for + ## matrix calculation ... + [nRowsError nColumnsError] = size(ve); + errorSize = size(ve(:),1); # this will calculate, if only one row + # of errors exist... in other words... two rows will be reranged to + # one row with the same number of elements. + rerangeIndex = floor([0:(errorSize-1)]/nRowsError)+1; + nLayers = net.numLayers; + + for i = 1:nLayers + Nn{i,1} = Nn{i,1}(:,rerangeIndex); + Aa{i,1} = Aa{i,1}(:,rerangeIndex); + [nRows nColumns] = size(Nn{i,1}); + bTemp = bias{i,1}; + bias{i,1} = repmat(bTemp,1,nColumns); + bias{i,1} = bias{i,1}(:,rerangeIndex); + endfor + mInputWeight = mInputWeight(:,rerangeIndex); + for i=2:nLayers + mLayerWeight{i,1} = mLayerWeight{i,1}(:,rerangeIndex); + endfor + Im = Im(:,rerangeIndex); + + ## define how the errors are connected + ## ATTENTION! this happens in row order... + numTargets = net.numTargets; + mIdentity = -eye(numTargets); + cols = size(mIdentity,2); + mIdentity = mIdentity(:,rem(0:(cols*nColumnsError-1),cols)+1); + errorConnect = cell(net.numLayers,1); + startPos = 0; + for i=net.numLayers + targSize = net.layers{i}.size; + errorConnect{i} = mIdentity(startPos+[1:targSize],:); + startPos = startPos + targSize; + endfor + + ## SECOND STEP ---------------------------------------------- + ## define and calculate the derivative matrix dF + ## - this is "done" by the two first derivative functions + ## of the transfer functions + ## e.g. __dpureline, __dtansig, __dlogsig and so on ... + + ## calculate the sensitivity matrix tildeS + ## start at the end layer, this means of course the output layer, + ## the transfer function is selectable + + ## for calculating the last layer + ## this should happen like following: + ## tildeSx = -dFx(n_x^x) + ## use mIdentity to calculate the number of targets correctly + ## for all other layers, use instead: + ## tildeSx(-1) = dF1(n_x^(x-1))(W^x)' * tildeSx; + + for iLayers = nLayers:-1:1 # this will count from the last + # layer to the first layer ... + n = Nn{iLayers}; # nLayers holds the value of the last layer... + ## which transfer function should be used? + if (iLayers==nLayers) + switch(net.layers{iLayers}.transferFcn) + case "radbas" + tildeSxTemp = __dradbas(n); + case "purelin" + tildeSxTemp = __dpurelin(n); + case "tansig" + n = tansig(n); + tildeSxTemp = __dtansig(n); + case "logsig" + n = logsig(n); + tildeSxTemp = __dlogsig(n); + otherwise + error(["transfer function argument: " net.layers{iLayers}.transferFcn " is not valid!"]) + endswitch + tildeSx{iLayers,1} = tildeSxTemp .* mIdentity; + n = bias{nLayers,1}; + switch(net.layers{iLayers}.transferFcn) + case "radbas" + tildeSbxTemp = __dradbas(n); + case "purelin" + tildeSbxTemp = __dpurelin(n); + case "tansig" + n = tansig(n); + tildeSbxTemp = __dtansig(n); + case "logsig" + n = logsig(n); + tildeSbxTemp = __dlogsig(n); + otherwise + error(["transfer function argument: " net.layers{iLayers}.transferFcn " is not valid!"]) + endswitch + tildeSbx{iLayers,1} = tildeSbxTemp .* mIdentity; + endif + + if (iLayers double of rows in the + ## jacobi matrix ... 3 targets --> three times the number of rows like + ## with one target and so on. + + ## now calculate jacobi matrix + ## to do this, define first the transposed of it + ## this makes it easier to calculate on the "batch" way, means all inputs + ## at the same time... + ## and it makes it easier to use the matrix calculation way.. + + JjTrans = zeros(nRowsJacobiMatrix,nColumnsJacobiMatrix)'; # transposed jacobi matrix + + ## Weight Gradients + for i=1:net.numLayers + if i==1 + newInputs = Im; + newTemps = tildeSx{i,1}; + gIW{i,1} = copyRows(newTemps,net.inputs{i}.size) .* copyRowsInt(newInputs,net.layers{i}.size); + endif + if i>1 + Ad = cell2mat(Aa(i-1,1)'); + newInputs = Ad; + newTemps = tildeSx{i,1}; + gLW{i,1} = copyRows(newTemps,net.layers{i-1}.size) .* copyRowsInt(newInputs,net.layers{i}.size); + endif + endfor + + for i=1:net.numLayers + [nRows, nColumns] = size(Im); + if (i==1) + nWeightElements = a(i)*a(i+1); # n inputs * n hidden neurons + JjTrans(1:nWeightElements,:) = gIW{i}(1:nWeightElements,:); + nWeightBias = a(i+1); + start = nWeightElements; + JjTrans(start+1:start+nWeightBias,:) = tildeSbx{i,1}; + start = start+nWeightBias; + endif + if (i>1) + nLayerElements = a(i)*a(i+1); # n hidden neurons * n output neurons + JjTrans(start+1:start+nLayerElements,:)=gLW{i}(1:nLayerElements,:); + start = start + nLayerElements; + nLayerBias = a(i+1); + JjTrans(start+1:start+nLayerBias,:) = tildeSbx{i,1}; + start = start + nLayerBias; + endif + endfor + Jj = JjTrans'; + ## END THIRD STEP ------------------------------------------------- + + +#======================================================= +# +# additional functions +# +#======================================================= + + function k = copyRows(k,m) + # make copies of the ROWS of Aa matrix + + mRows = size(k,1); + k = k(rem(0:(mRows*m-1),mRows)+1,:); + endfunction + +# ------------------------------------------------------- + + function k = copyRowsInt(k,m) + # make copies of the ROWS of matrix with elements INTERLEAVED + + mRows = size(k,1); + k = k(floor([0:(mRows*m-1)]/m)+1,:); + endfunction + +# ===================================================================== +# +# END additional functions +# +# ===================================================================== + +endfunction