FRANABP - Frame Analysis Basis Pursuit

Usage

c = franabp(F,f)
c = franabp(F,f,lambda,C,tol,maxit)
[c,relres,iter,frec,cd] = franabp(...)

Input parameters

Output parameters

Description

c = franabp(F,f) solves the basis pursuit problem

for a general frame F using SALSA (Split Augmented Lagrangian Srinkage algorithm) which is an appication of ADMM (Alternating Direction Method of Multipliers) to the basis pursuit problem.

The algorithm given F and f and parameters C \(>0\), lambda \(>0\) (see below) acts as follows:

Initialize c,d
repeat
  v <- soft(c+d,lambda/C) - d
  d <- F*(FF*)^(-1)(f - Fv)
  c <- d + v
end

When compared to other algorithms, Fc = f holds exactly (up to a num. prec) in each iteration.

For a quick execution, the function requires analysis operator of the canonical dual frame F*(FF*)^(-1). By default, the function attempts to call framedual to create the canonical dual frame explicitly. If it is not available, the conjugate gradient method algorithm is used for inverting the frame operator in each iteration of the algorithm. Optionally, the canonical dual frame object or an anonymous function acting as the analysis operator of the canonical dual frame can be passed as a key-value pair 'Fd',Fd see below.

Optional positional parameters (lambda,C,tol,maxit)

lambda

A parameter for weighting coefficients in the objective function. For lambda~=1 the basis pursuit problem changes to

\begin{align*} \arg\,\min_c \ ||\lambda c||_1 \\\\ \text{subject to } Fc = f \end{align*}

lambda can either be a scalar or a vector of the same length as c (in such case the product is carried out elementwise). One can obtain length of c from length of f by frameclength. framecoef2native and framenative2coef will help with defining weights specific to some regions of coefficients (e.g. channel-specific weighting can be achieved this way). The default value of lambda is 1.

C

A step parameter of the SALSA algorithm. The default value of C is the upper frame bound of F. Depending on the structure of the frame, this can be an expensive operation.

tol

Defines tolerance of relres which is a norm or a relative difference of coefficients obtained in two consecutive iterations of the algorithm. The default value 1e-2.

maxit

Maximum number of iterations to do. The default value is 100.

Other optional parameters

Key-value pairs:

'Fd',Fd

A canonical dual frame object or an anonymous function acting as the analysis operator of the canonical dual frame.

'printstep',printstep

Print current status every printstep iteration.

Flag groups (first one listed is the default):

'print','quiet'

Enables/disables printing of notifications.

'zeros','frana'

Starting point of the algorithm. With 'zeros' enabled, the algorithm starts from coefficients set to zero, with 'frana' the algorithm starts from c=frana(F,f).

Returned arguments:

[c,relres,iter] = franabp(...) returns the residuals relres in a vector and the number of iteration steps done iter.

[c,relres,iter,frec,cd] = franabp(...) returns the reconstructed signal from the coefficients, frec (this requires additional computations) and a coefficients cd minimising the ||c||_2 norm (this is a byproduct of the algorithm).

The relationship between the output coefficients frec and c is given by

frec = frsyn(F,c);

And cd and f by

cd = frana(framedual(F),f);

Examples:

The following example shows how franabp produces a sparse representation of a test signal greasy still maintaining a perfect reconstruction:

f = greasy;
% Gabor frame with redundancy 8
F = frame('dgtreal','gauss',64,512);
% Solve the basis pursuit problem
[c,~,~,frec,cd] = franabp(F,f);
% Plot sparse coefficients
figure(1);
plotframe(F,c,'dynrange',50);

% Plot coefficients obtained by applying an analysis operator of a
% dual Gabor system to *f*
figure(2);
plotframe(F,cd,'dynrange',50);

% Check the reconstruction error (should be close do zero).
% frec is obtained by applying the synthesis operator of frame *F*
% to sparse coefficients *c*.
norm(f-frec)

% Compare decay of coefficients sorted by absolute values
% (compressibility of coefficients)
figure(3);
semilogx([sort(abs(c),'descend')/max(abs(c)),...
sort(abs(cd),'descend')/max(abs(cd))]);
legend({'sparsified coefficients','dual system coefficients'});