test_l1_constraint.m

%% Basis Pursuit with Douglas Rachford
% Test for DR algorithm for L1 minimization (BP).
% We do here a compressed sensing resolution
% (random matrix).

%%
% Add the toolbox.

addpath('../');
addpath('../toolbox/');

%% 
% Dimensionality of the signal and number of measurements.

n = 200;
p = n/4;

%%
% Sensing matrix.

A = randn(p,n);

%%
% Measurements.

y = randn(p,1);

%%
% We aim at solving 

%%
% |min_{A*x=y} norm(x,1)|

%%
% This can be rewriten as the minimization of |F(x)+G(x)|
% where |F=norm(x,1)| and |G=i_{A*x=y}| is the indicator function.


%%
% The proximity operator of the L1 norm is the soft thresholding.

ProxF = @(x,tau)perform_soft_thresholding(x, tau);

%%
% The proximity operator of the indicator of |A*x=y| is the orthogonal
% projection on A*x=y.

pA = A'*(A*A')^(-1);
ProxG = @(x,tau)x + pA*(y-A*x);

%%
% Create a function to record the values of F and the constraint at each iteration.

F = @(x)norm(x,1);
Constr = @(x)1/2*norm(y-A*x)^2;
options.report = @(x)struct('F', F(x), 'Constr', Constr(x));

%%
% Run the algorithm. 

options.gamma = 5;
options.niter = 5000;
[x,R] = perform_dr(zeros(n,1), ProxF, ProxG, options);

%%
% Display the solution. At convergence, it should be of sparsity |p|.

clf;
plot(x);
axis tight;

%%
% Retrieve the F and constraint function values.

f = s2v(R,'F');
constr = s2v(R,'Constr');

%%
% Display.

clf;
subplot(2,1,1);
plot(f(2:end));
axis tight; title('Objective');
subplot(2,1,2);
plot(constr(2:end));
axis tight; title('Constraint');