% Calculation of gradient and objective for Logistic
% Multi-Class Classifcation.
%
% function [obj,grad] = mcclogistic(v,Y,V,lambda,l,varargin)
% v - vector of parameters [n*p*l,1]
% Y - rating matrix (labels) [n,m]
% V - the feature matrix [m,p]
% lambda - regularization parameter [scalar]
% l - # of labels (1..l)
% obj - value of objective at v [scalar]
% grad - gradient at v [n*p*l,1]
% 
% Written by Jason Rennie, April 2005
% Last modified: Tue Jul 25 15:08:38 2006

function [obj,grad] = mcclogistic(v,Y,V,lambda,l,varargin)
  fn = mfilename;
  if nargin < 5
    error('insufficient parameters')
  end
  % Parameters that can be set via varargin
  verbose = 1;
  % Process varargin
  paramgt;
  
  t0 = clock;
  [n,m] = size(Y);
  p = length(v)./n./l;
  if p ~= floor(p) | p < 1
    error('dimensions of v and Y don''t match l');
  end
  U = reshape(v,n,p,l);
  Z = zeros(n,m,l);
  for i=1:l
    Z(:,:,i) = U(:,:,i)*V';
  end
  obj = lambda.*sum(sum(sum(U.^2)))./2;
  dU = zeros(n,p,l);
  YY = full(Y==0) + Y;
  YI = sub2ind(size(Z),(1:n)'*ones(1,m),ones(n,1)*(1:m),YY);
  ZY = Z(YI);
  for i=1:l
    obj = obj + sum(sum(h(ZY-Z(:,:,i)).*(Y~=i).*(Y>0)));
  end
  ZHP = zeros(n,m);
  for i=1:l
    ZHP = ZHP + hprime(ZY-Z(:,:,i)).*(Y~=i).*(Y>0);
  end
  for i=1:l
    dU(:,:,i) = ((Y==i).*ZHP - (Y~=i).*(Y>0).*hprime(ZY-Z(:,:,i)))*V + lambda.*U(:,:,i);
  end
  grad = dU(:);
  if verbose
    fprintf(1,'lambda=%.2e obj=%.4e grad''*grad=%.4e time=%.1f\n',lambda,obj,grad'*grad,etime(clock,t0));
  end

function [ret] = h(z)
  ret = log(1+exp(-z));

function [ret] = hprime(z)
  ret = -(exp(-z)./(1+exp(-z)));

% ChangeLog
% 7/25/06 - Added varargin, verbose
% 3/23/05 - made calcultions take better advantage of sparseness
% 3/18/05 - fixed bug in objective (wasn't squaring fro norms)
% 3/1/05 - added objective calculation
% 2/23/05 - fixed bug in hprime()
