function [f,D,E,N,lambda] = rpca_admm(X,Lambda,tau,rho,homo)
% minimize lambda*norm_nuc(D)+norm(E,{2,1})
% subject to
% X = D+E+N
% norm(N,'fro')<= tau
% lambda is chosen from the given candidate vector Lambda, stop if D is
% rank-deficient by 1
% f is regressor taken from the null space of D
% X is noisy data
% E is outlier bias part
% N is noise bias part
% mu update: mu{k+1} = mu{k}*rho

[d,n] = size(X);
for lambda = Lambda
    D = X;
    E = X*0;
    N = E;
    Y = E;
    mu = d*n/norm(D(:),1);
    
    for iter = 1:3000
        Z = Y/mu+X-E-N;
        [U,S,V] = svd(Z);
        s = diag(S);
        s = max(s-lambda/mu,0);
        D = U(:,1:length(s))*diag(s)*V(:,1:length(s))';
        Z = Y/mu+X-D-N;
        if homo==1
            Z(end,:) = 0;
        end
        E = bsxfun(@times,max(sqrt(sum(Z.^2,1))-0.15/mu,0),normc(Z));
        Z = Y/mu+X-D-E;
        if homo==1
            Z(end,:) = 0;
        end
        if norm(Z,'fro')>tau
            Z = Z/norm(Z,'fro')*tau;
        end
        N = Z;
        Y = Y+mu*(X-D-E-N);
        if norm(X-D-E-N,'fro')/norm(X,'fro')<=1e-8
            break;
        end
        mu = min(rho*mu,1e3);
    end
    
    [~,S,V] = svd(D');
    s = diag(S);
    if s(end)<=1e-4
        break;
    end
end
f = V(:,d);