admm.bp: Basis Pursuit
In ADMM: Algorithms using Alternating Direction Method of Multipliers

Description Usage Arguments Value Iteration History Examples

View source: R/admm.bp.R

For an underdetermined system, Basis Pursuit aims to find a sparse solution that solves

\textrm{min}_x ~ \|x\|_1 \quad \textrm{s.t} \quad Ax=b

which is a relaxed version of strict non-zero support finding problem. The implementation is borrowed from Stephen Boyd's MATLAB code.

admm.bp(
  A,
  b,
  xinit = NA,
  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
`xinit`	a length-n vector for initial value
`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

x: a length-n solution vector
history: 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,

objval: object (cost) function value
r_norm: norm of primal residual
s_norm: norm of dual residual
eps_pri: feasibility tolerance for primal feasibility condition
eps_dual: 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.

## generate sample data
n = 30
m = 10

A = matrix(rnorm(n*m), nrow=m)    # design matrix
x = c(stats::rnorm(3),rep(0,n-3)) # coefficient
x = base::sample(x)
b = as.vector(A%*%x)              # response

## run example
output  = admm.bp(A, b)
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)