function w = logisticRegression(X,y)
% w = logistic_regression(X,y)
% X is your data with the constant term added, i.e.
% X(i,:) is [1 xi1 xi2...]
% y is your labels (column vector)
% Returns w_hat

[N,d] = size(X);

w = zeros(d,1);

% initial values for control variables
cont = 1;
err=inf;

t=1; % iteration counter

while (cont),
  % yhat denotes P(y=1|x,w) with logistic model
  yhat = g(X*w); 
  
  % you will need to write a function gradient 
  grad = gradient(X,y,yhat);
  
  % compute the Hessian
  Z = repmat(yhat.*(1-yhat),1,d).*X;
  hess = -(eye(d) + Z'*X); % Hessian (second derivatives)
  
  wold = w;
  w = wold - inv(hess)*grad; % Newton-Raphson step
  
  % error function
  errold=err;
  err=-sum(ourlog(y,yhat)+ourlog(1-y,1-yhat));
  cont = norm(err-errold)>1e-2; % stopping criterion
  
  fprintf(2,'%d: change in w = %.6f  err=%.6f\n',t,norm(w-wold),err);
  t=t+1;
end;


function l=ourlog(a,b)
% we use the convention that 0*log(0) is 0.
f0=find(b==0);
f1=find(b>0);
l=zeros(size(b));
l(f0)=0;
l(f1)=a(f1).*log(b(f1));
 
function grad = gradient(X,y,yhat)
grad  = ((y-yhat)'*X)';  % gradient

function p=g(x)
% the logistic function. You may want to make it a separate file
p=1./(1+exp(-x));