function xest = uamp(yn, sparsity, A, niter)

    [M, N] = size(A);
    tau = sqrt(2*log10(M/sparsity));

    xest = zeros(N, 1);
    yp = yn;
    z = yp;
    eta = @(x, beta)(x./abs(x).*(abs(x)-beta).*(abs(x)-beta>0));
    for iter = 1:niter 
        sigma = norm(z,2)/sqrt(M);
        xest = eta(xest+A'*z, tau*sigma);
        tmp = sum(abs(xest)>0);
        z = yp - A*xest + (tmp/M)*z; 
    end

end
