function z = rbfpreimg2(varargin)% RBFPREIMG2 RBF pre-image problem by Gradient optimization.%% Synopsis:%  z = rbfpreimg2(model)%  z = rbfpreimg2(model,options)%%  Description:%   z = rbfpreimg2(model) it uses gradient method to solve %     the pre-image problem for the Radial Basis Function (RBF) %     kernel. The function 'fminunc' of the Matlab Optimization %     toolbox is exploited for 1D search along the gradient %     direction.%%   z = rbfpreimg2(model,options) use to specify the control%     parameters of the gradient optimization.% % Input:%  model [struct] Kernel expansion:%   .Alpha [num_data x 1] Weight vector.%   .sv.X [dim x num_data] Vectors determining the kernel expansion.%   .options.arg [1x1] Argument of the RBF kernel (see 'help kernel').%%  options [struct] Control parameters:%   .min_improvement [1x1] Minimal allowed improvement of the objective %     function in two consecutive steps (default 1e-3).%  options.start_t [1x1] Starting value of the 1D search procedure %   (default 1e-3).%% Output:%  z [dim x 1] Found preimage.%% See also %  RBFPREIMG, RBFPREIMG3, RSRBF, KPCAREC.%% About: Statistical Pattern Recognition Toolbox% (C) 1999-2003, Written by Vojtech Franc and Vaclav Hlavac% <a href="http://www.cvut.cz">Czech Technical University Prague</a>% <a href="http://www.feld.cvut.cz">Faculty of Electrical Engineering</a>% <a href="http://cmp.felk.cvut.cz">Center for Machine Perception</a>% Modifications:% 10-jun-2004, VF% 03-dec-2003, VF% process input arguments%----------------------------------if nargin > 2,  z=foo(varargin{1},varargin{2},varargin{3},varargin{4},varargin{5});  return;else  model = varargin{1};  if nargin < 2, options=[]; else options = c2s( varargin{2}); end
  if ~isfield(options,'min_improvement'), options.min_improvement = 1e-3; end  if ~isfield(options,'start_t'), options.start_t = 1e-3; end  if ~isfield(options,'attempts'), options.attempts = 10; endend[dim,num_sv]=size(model.sv.X);ker = 'rbf';arg = model.options.arg;iXi = sum( model.sv.X.^2)';s2 = arg^2;% Selection of the starting point out of the model.sv.X.% The point in which is the objective function minimal is taken.% Minimum over 50 randomly drawn points is used.%--------------------------------------------------------------rand_inx = randperm( num_sv );rand_inx = rand_inx(1:min([num_sv,50]));Z = model.sv.X(:,rand_inx);fval = kernel(Z,model.sv.X,ker,arg)*model.Alpha(:);fval = -fval.^2;[dummy, inx ] = min( fval );z = Z(:,inx );% Gradient descent optimization  %--------------------------------------change=inf;opt=optimset('display','off','Diagnostics','off','LargeScale','off');warning off MATLAB:divideByZero;while change > options.min_improvement,      % compute gradient   dotp = kernel( model.sv.X,z,ker,arg ).*model.Alpha(:);   dz = z*sum(dotp) - model.sv.X*dotp;   % auxiciliary variables   zXi = model.sv.X'*z;   dzXi = model.sv.X'*dz;      Ai = -(1/(2*s2)) * (iXi - 2*zXi + z'*z);   Bi = -(1/s2) * (z'*dz - dzXi);   C = -(1/(2*s2)) * (dz'*dz);      % 1D-search to determine the size of the step in the gradient direction   [t,fval] = fminunc('rbfpreimg2',options.start_t,opt,model.Alpha,Ai,Bi,C);   
   old_z = z;   z = z + dz*t;         change = sum((z-old_z).^2);endreturn;%---------------------------------function f=foo(t,Alpha,Ai,Bi,C)% f = Alpha(:)'*exp(Ai + Bi*t + C*t^2);f = -f^2;return;% EOF