comparison toolboxes/FullBNT-1.0.7/netlab3.3/gpgrad.m @ 0:e9a9cd732c1e tip

first hg version after svn

+function g = gpgrad(net, x, t)
+%GPGRAD	Evaluate error gradient for Gaussian Process.
+%
+%	Description
+%	G = GPGRAD(NET, X, T) takes a Gaussian Process data structure NET
+%	together  with a matrix X of input vectors and a matrix T of target
+%	vectors, and evaluates the error gradient G. Each row of X
+%	corresponds to one input vector and each row of T corresponds to one
+%	target vector.
+%
+%	See also
+%	GP, GPCOVAR, GPFWD, GPERR
+%
+
+%	Copyright (c) Ian T Nabney (1996-2001)
+
+errstring = consist(net, 'gp', x, t);
+if ~isempty(errstring);
+  error(errstring);
+end
+
+% Evaluate derivatives with respect to each hyperparameter in turn.
+ndata = size(x, 1);
+[cov, covf] = gpcovar(net, x);
+cninv = inv(cov);
+trcninv = trace(cninv);
+cninvt = cninv*t;
+
+% Function parameters
+switch net.covar_fn
+
+  case 'sqexp'		% Squared exponential
+    gfpar = trace(cninv*covf) - cninvt'*covf*cninvt;
+
+  case 'ratquad' 	% Rational quadratic
+    beta = diag(exp(net.inweights));
+    gfpar(1) = trace(cninv*covf) - cninvt'*covf*cninvt;
+    D2 = (x.*x)*beta*ones(net.nin, ndata) - 2*x*beta*x' ... 
+      + ones(ndata, net.nin)*beta*(x.*x)';
+    E = ones(size(D2));
+    L = - exp(net.fpar(2)) * covf .* log(E + D2); % d(cn)/d(nu)
+    gfpar(2) = trace(cninv*L) - cninvt'*L*cninvt;
+
+  otherwise
+    error(['Unknown covariance function ', net.covar_fn]);
+end
+
+% Bias derivative
+ndata = size(x, 1);
+fac = exp(net.bias)*ones(ndata);
+gbias = trace(cninv*fac) - cninvt'*fac*cninvt;
+
+% Noise derivative
+gnoise = exp(net.noise)*(trcninv - cninvt'*cninvt);
+
+% Input weight derivatives
+if strcmp(net.covar_fn, 'ratquad')
+  F = (exp(net.fpar(2))*E)./(E + D2);
+end
+
+nparams = length(net.inweights);
+for l = 1 : nparams
+  vect = x(:, l);
+  matx = (vect.*vect)*ones(1, ndata) ... 
+	- 2.0*vect*vect' ... 
+	+ ones(ndata, 1)*(vect.*vect)';
+  switch net.covar_fn
+    case 'sqexp'	% Squared exponential
+      dmat = -0.5*exp(net.inweights(l))*covf.*matx;
+      
+    case 'ratquad'	% Rational quadratic
+      dmat = - exp(net.inweights(l))*covf.*matx.*F;
+    otherwise
+      error(['Unknown covariance function ', net.covar_fn]);
+  end
+
+  gw1(l) = trace(cninv*dmat) - cninvt'*dmat*cninvt;
+end
+
+g1 = [gbias, gnoise, gw1, gfpar];
+g1 = 0.5*g1;
+
+% Evaluate the prior contribution to the gradient.
+if isfield(net, 'pr_mean')
+  w = gppak(net);
+  m = repmat(net.pr_mean, size(w));
+  if size(net.pr_mean) == [1 1]
+    gprior = w - m;
+    g2 = gprior/net.pr_var;
+  else
+    ngroups = size(net.pr_mean, 1);
+    gprior = net.index'.*(ones(ngroups, 1)*w - m);
+    g2 = (1./net.pr_var)'*gprior;
+  end
+else
+  gprior = 0;
+  g2 = 0;
+end
+
+g = g1 + g2;