function [alpha,w,b]=reestimate_X(X,Y,alpha0, alpha_prune, gamma_a,gamma_b)
%% The "alpha" from this method has length k*(C-1)
%% this is because each term in "w" has a corresponding alpha
    
  % check dimension of inputs
  [n,k]=size(X);
  C=max(Y);
  assert(length(Y)==n);
  assert(length(alpha0)==k*(C-1));
  
  % proceed as usual
%  alphaprune=1E3;
  alpha=alpha0;
  w=zeros(size(alpha));
  
  numiter=0;
  maxnumiter=20;
  terminate=0;
  eps_tol=1E-5;
  
  while ~terminate
    tick=cputime;
    numiter=numiter+1;
    alpha_old=alpha;

    
    inind=(abs(alpha)<alpha_prune);   %% ALPHA may become negative!!
    fprintf('iter=%i\tnum pruned features=%i\n',numiter, sum(~inind));
    [w,b, invH]=logist_X(X,Y,alpha,inind);
%    [w,b, invH]=logist_X(X,Y,alpha,inind,w);
    dlogdetH=diag(invH);
    
    gamma=1-alpha(inind).*dlogdetH;
    alpha(inind)=(gamma+2*gamma_a)./(w(inind).^2+2*gamma_b);
%    alpha=abs(alpha);  %% CLAMP ALPHA TO BE POSTIIVE!
    alpha(alpha<0)=inf;
    
    if (norm(1./(alpha+eps)-1./(alpha_old+eps))<eps_tol | numiter >= maxnumiter)
      terminate=1;
    end
  end

  % last run to ensure wb agrees with alpha
  inind=(abs(alpha)<alpha_prune);   %% ALPHA may become negative!!
  [w,b, invH]=logist_X(X,Y,alpha,inind,w);

  w(~inind)=0;