wolffd@0: function jac = rbfjacob(net, x)
wolffd@0: %RBFJACOB Evaluate derivatives of RBF network outputs with respect to inputs.
wolffd@0: %
wolffd@0: %	Description
wolffd@0: %	G = RBFJACOB(NET, X) takes a network data structure NET and a matrix
wolffd@0: %	of input vectors X and returns a three-index matrix G whose I, J, K
wolffd@0: %	element contains the derivative of network output K with respect to
wolffd@0: %	input parameter J for input pattern I.
wolffd@0: %
wolffd@0: %	See also
wolffd@0: %	RBF, RBFGRAD, RBFBKP
wolffd@0: %
wolffd@0: 
wolffd@0: %	Copyright (c) Ian T Nabney (1996-2001)
wolffd@0: 
wolffd@0: % Check arguments for consistency
wolffd@0: errstring = consist(net, 'rbf', x);
wolffd@0: if ~isempty(errstring);
wolffd@0:   error(errstring);
wolffd@0: end
wolffd@0: 
wolffd@0: if ~strcmp(net.outfn, 'linear')
wolffd@0:   error('Function only implemented for linear outputs')
wolffd@0: end
wolffd@0: 
wolffd@0: [y, z, n2] = rbffwd(net, x);
wolffd@0: 
wolffd@0: ndata = size(x, 1);
wolffd@0: jac = zeros(ndata, net.nin, net.nout);
wolffd@0: Psi = zeros(net.nin, net.nhidden);
wolffd@0: % Calculate derivative of activations wrt n2
wolffd@0: switch net.actfn
wolffd@0: case 'gaussian'
wolffd@0:   dz = -z./(ones(ndata, 1)*net.wi);
wolffd@0: case 'tps'
wolffd@0:   dz = 2*(1 + log(n2+(n2==0)));
wolffd@0: case 'r4logr'
wolffd@0:   dz = 2*(n2.*(1+2.*log(n2+(n2==0))));
wolffd@0: otherwise
wolffd@0:    error(['Unknown activation function ', net.actfn]);
wolffd@0: end
wolffd@0: 
wolffd@0: % Ignore biases as they cannot affect Jacobian
wolffd@0: for n = 1:ndata
wolffd@0:   Psi = (ones(net.nin, 1)*dz(n, :)).* ...
wolffd@0:     (x(n, :)'*ones(1, net.nhidden) - net.c');
wolffd@0:   % Now compute the Jacobian
wolffd@0:   jac(n, :, :) =  Psi * net.w2;
wolffd@0: end