wolffd@0: function h = hesschek(net, x, t)
wolffd@0: %HESSCHEK Use central differences to confirm correct evaluation of Hessian matrix.
wolffd@0: %
wolffd@0: %	Description
wolffd@0: %
wolffd@0: %	HESSCHEK(NET, X, T) takes a network data structure NET, together with
wolffd@0: %	input and target data matrices X and T, and compares the evaluation
wolffd@0: %	of the Hessian matrix using the function NETHESS and using central
wolffd@0: %	differences with the function NETERR.
wolffd@0: %
wolffd@0: %	The optional return value H is the Hessian computed using NETHESS.
wolffd@0: %
wolffd@0: %	See also
wolffd@0: %	NETHESS, NETERR
wolffd@0: %
wolffd@0: 
wolffd@0: %	Copyright (c) Ian T Nabney (1996-2001)
wolffd@0: 
wolffd@0: w0 = netpak(net);
wolffd@0: nwts = length(w0);
wolffd@0: h = nethess(w0, net, x, t);
wolffd@0: 
wolffd@0: w = w0;
wolffd@0: hcent = zeros(nwts, nwts);
wolffd@0: h1 =  0.0; h2 =  0.0; h3 =  0.0; h4 = 0.0;
wolffd@0: epsilon = 1.0e-4;
wolffd@0: fprintf(1, 'Checking Hessian ...\n\n');
wolffd@0: for k = 1:nwts;
wolffd@0:   for l = 1:nwts;
wolffd@0:     if(l == k)
wolffd@0:       w(k) = w0(k) + 2.0*epsilon;
wolffd@0:       h1 = neterr(w, net, x, t);
wolffd@0:       w(k) = w0(k) - 2.0*epsilon;
wolffd@0:       h2 = neterr(w, net, x, t);
wolffd@0:       w(k) = w0(k);
wolffd@0:       h3 = neterr(w, net, x, t);
wolffd@0:       hcent(k, k) = (h1 + h2 - 2.0*h3)/(4.0*epsilon^2);
wolffd@0:     else
wolffd@0:       w(k) = w0(k) + epsilon;
wolffd@0:       w(l) = w0(l) + epsilon;
wolffd@0:       h1 = neterr(w, net, x, t);
wolffd@0:       w(k) = w0(k) - epsilon;
wolffd@0:       w(l) = w0(l) - epsilon;
wolffd@0:       h2 = neterr(w, net, x, t);
wolffd@0:       w(k) = w0(k) + epsilon;
wolffd@0:       w(l) = w0(l) - epsilon;
wolffd@0:       h3 = neterr(w, net, x, t);
wolffd@0:       w(k) = w0(k) - epsilon;
wolffd@0:       w(l) = w0(l) + epsilon;
wolffd@0:       h4 = neterr(w, net, x, t);
wolffd@0:       hcent(k, l) = (h1 + h2 - h3 - h4)/(4.0*epsilon^2);
wolffd@0:       w(k) = w0(k);
wolffd@0:       w(l) = w0(l);
wolffd@0:     end
wolffd@0:   end
wolffd@0: end
wolffd@0: 
wolffd@0: fprintf(1, '   analytical    numerical       delta\n\n');
wolffd@0: temp = [h(:), hcent(:), (h(:) - hcent(:))];
wolffd@0: fprintf(1, '%12.6f  %12.6f  %12.6f\n', temp');