Daniel@0: function [h, hdata] = glmhess_weighted(net, x, t, eso_w, hdata)
Daniel@0: %GLMHESS Evaluate the Hessian matrix for a generalised linear model.
Daniel@0: %
Daniel@0: %	Description
Daniel@0: %	H = GLMHESS(NET, X, T) takes a GLM network data structure NET, a
Daniel@0: %	matrix X of input values, and a matrix T of target values and returns
Daniel@0: %	the full Hessian matrix H corresponding to the second derivatives of
Daniel@0: %	the negative log posterior distribution, evaluated for the current
Daniel@0: %	weight and bias values as defined by NET. Note that the target data
Daniel@0: %	is not required in the calculation, but is included to make the
Daniel@0: %	interface uniform with NETHESS.  For linear and logistic outputs, the
Daniel@0: %	computation is very simple and is  done (in effect) in one line in
Daniel@0: %	GLMTRAIN.
Daniel@0: %
Daniel@0: %	See also
Daniel@0: %	GLM, GLMTRAIN, HESSCHEK, NETHESS
Daniel@0: %
Daniel@0: %	Copyright (c) Ian T Nabney (1996-9)
Daniel@0: 
Daniel@0: % Check arguments for consistency
Daniel@0: errstring = consist(net, 'glm', x, t);
Daniel@0: if ~isempty(errstring);
Daniel@0:   error(errstring);
Daniel@0: end
Daniel@0: 
Daniel@0: ndata = size(x, 1);
Daniel@0: nparams = net.nwts;
Daniel@0: nout = net.nout;
Daniel@0: p = glmfwd(net, x);
Daniel@0: inputs = [x ones(ndata, 1)];
Daniel@0: 
Daniel@0: if nargin == 4
Daniel@0:    hdata = zeros(nparams);	% Full Hessian matrix
Daniel@0:    % Calculate data component of Hessian
Daniel@0:    switch net.outfn
Daniel@0:     
Daniel@0:    case 'softmax'
Daniel@0:     bb_start = nparams - nout + 1;	% Start of bias weights block
Daniel@0:     ex_hess = zeros(nparams);	% Contribution to Hessian from single example
Daniel@0:     for m = 1:ndata
Daniel@0:       X = x(m,:)'*x(m,:);
Daniel@0:       a = diag(p(m,:))-((p(m,:)')*p(m,:)); 
Daniel@0:       a=eso_w(m,1)*a;
Daniel@0:       ex_hess(1:nparams-nout,1:nparams-nout) = kron(a, X);
Daniel@0:       ex_hess(bb_start:nparams, bb_start:nparams) = a.*ones(net.nout, net.nout);
Daniel@0:       temp = kron(a, x(m,:));
Daniel@0:       ex_hess(bb_start:nparams, 1:nparams-nout) = temp;
Daniel@0:       ex_hess(1:nparams-nout, bb_start:nparams) = temp';
Daniel@0:       hdata = hdata + ex_hess;
Daniel@0:     end
Daniel@0:     
Daniel@0:     otherwise
Daniel@0:       error(['Unknown activation function ', net.actfn]);
Daniel@0:     end
Daniel@0: end
Daniel@0: 
Daniel@0: [h, hdata] = hbayes(net, hdata);