Least Absolute Shrinkage and Selection Operator
LASSO, or L1-regularized regression, is an optimization problem to solve
\textrm{min}_x ~ \frac{1}{2}\|Ax-b\|_2^2 + λ \|x\|_1
for sparsifying the coefficient vector x. The implementation is borrowed from Stephen Boyd's MATLAB code.
admm.lasso( A, b, lambda = 1, rho = 1, alpha = 1, abstol = 1e-04, reltol = 0.01, maxiter = 1000 )
A |
an (m\times n) regressor matrix |
b |
a length-m response vector |
lambda |
a regularization parameter |
rho |
an augmented Lagrangian parameter |
alpha |
an overrelaxation parameter in [1,2] |
abstol |
absolute tolerance stopping criterion |
reltol |
relative tolerance stopping criterion |
maxiter |
maximum number of iterations |
a named list containing
a length-n solution vector
dataframe recording iteration numerics. See the section for more details.
When you run the algorithm, output returns not only the solution, but also the iteration history recording following fields over iterates,
object (cost) function value
norm of primal residual
norm of dual residual
feasibility tolerance for primal feasibility condition
feasibility tolerance for dual feasibility condition
In accordance with the paper, iteration stops when both r_norm and s_norm values
become smaller than eps_pri and eps_dual, respectively.
Tibshirani R (1996). “Regression Shrinkage and Selection via the Lasso.” Journal of the Royal Statistical Society. Series B (Methodological), 58(1), 267–288. ISSN 00359246.
## generate sample data
m = 50
n = 100
p = 0.1 # percentange of non-zero elements
x0 = matrix(Matrix::rsparsematrix(n,1,p))
A = matrix(rnorm(m*n),nrow=m)
for (i in 1:ncol(A)){
A[,i] = A[,i]/sqrt(sum(A[,i]*A[,i]))
}
b = A%*%x0 + sqrt(0.001)*matrix(rnorm(m))
## set regularization lambda value
lambda = 0.1*base::norm(t(A)%*%b, "F")
## run example
output = admm.lasso(A, b, lambda)
niter = length(output$history$s_norm)
history = output$history
## report convergence plot
opar <- par(no.readonly=TRUE)
par(mfrow=c(1,3))
plot(1:niter, history$objval, "b", main="cost function")
plot(1:niter, history$r_norm, "b", main="primal residual")
plot(1:niter, history$s_norm, "b", main="dual residual")
par(opar)Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.