function [x,x_debias,objective,times,debias_start,mses,taus]= ...
GPSR_BB(y,A,tau,varargin)
%
% GPSR_BB version 5.0, December 4, 2007
%
% This function solves the convex problem
% arg min_x = 0.5*|| y - A x ||_2^2 + tau || x ||_1
% using the algorithm GPSR-BB, described in the following paper
%
% "Gradient Projection for Sparse Reconstruction: Application
% to Compressed Sensing and Other Inverse Problems"
% by Mario A. T. Figueiredo, Robert D. Nowak, Stephen J. Wright,
% Journal of Selected Topics on Signal Processing, December 2007
% (to appear).
%
%
% -----------------------------------------------------------------------
% Copyright (2007): Mario Figueiredo, Robert Nowak, Stephen Wright
%
% GPSR is distributed under the terms
% of the GNU General Public License 2.0.
%
% Permission to use, copy, modify, and distribute this software for
% any purpose without fee is hereby granted, provided that this entire
% notice is included in all copies of any software which is or includes
% a copy or modification of this software and in all copies of the
% supporting documentation for such software.
% This software is being provided "as is", without any express or
% implied warranty. In particular, the authors do not make any
% representation or warranty of any kind concerning the merchantability
% of this software or its fitness for any particular purpose."
% ----------------------------------------------------------------------
%
%
% Please check for the latest version of the code and paper at
% www.lx.it.pt/~mtf/GPSR
%
% ===== Required inputs =============
%
% y: 1D vector or 2D array (image) of observations
%
% A: if y and x are both 1D vectors, A can be a
% k*n (where k is the size of y and n the size of x)
% matrix or a handle to a function that computes
% products of the form A*v, for some vector v.
% In any other case (if y and/or x are 2D arrays),
% A has to be passed as a handle to a function which computes
% products of the form A*x; another handle to a function
% AT which computes products of the form A'*x is also required
% in this case. The size of x is determined as the size
% of the result of applying AT.
%
% tau: usually, a non-negative real parameter of the objective
% function (see above). It can also be an array, the same
% size as x, with non-negative entries; in this case,
% the objective function weights differently each element
% of x, that is, it becomes
% 0.5*|| y - A x ||_2^2 + tau^T * abs(x)
%
% ===== Optional inputs =============
%
%
% 'AT' = function handle for the function that implements
% the multiplication by the conjugate of A, when A
% is a function handle. If A is an array, AT is ignored.
%
% 'StopCriterion' = type of stopping criterion to use
% 0 = algorithm stops when the relative
% change in the number of non-zero
% components of the estimate falls
% below 'ToleranceA'
% 1 = stop when the relative
% change in the objective function
% falls below 'ToleranceA'
% 2 = stop when the norm of the difference between
% two consecutive estimates, divided by the norm
% of one of them falls below toleranceA
% 3 = stop when LCP estimate of relative
% distance to solution
% falls below 'ToleranceA'
% 4 = stop when the objective function
% becomes equal or less than toleranceA.
% 5 = stop when the norm of the difference between
% two consecutive estimates, divided by the norm
% of one of them falls below toleranceA
% Default = 3.
%
% 'ToleranceA' = stopping threshold; Default = 0.01
%
% 'Debias' = debiasing option: 1 = yes, 0 = no.
% Default = 0.
%
% 'ToleranceD' = stopping threshold for the debiasing phase:
% Default = 0.0001.
% If no debiasing takes place, this parameter,
% if present, is ignored.
%
% 'MaxiterA' = maximum number of iterations allowed in the
% main phase of the algorithm.
% Default = 10000
%
% 'MiniterA' = minimum number of iterations performed in the
% main phase of the algorithm.
% Default = 5
%
% 'MaxiterD' = maximum number of iterations allowed in the
% debising phase of the algorithm.
% Default = 200
%
% 'MiniterD' = minimum number of iterations to perform in the
% debiasing phase of the algorithm.
% Default = 5
%
% 'Initialization' must be one of {0,1,2,array}
% 0 -> Initialization at zero.
% 1 -> Random initialization.
% 2 -> initialization with A'*y.
% array -> initialization provided by the user.
% Default = 0;
%
% 'Monotone' = enforce monotonic decrease in f, or not?
% any nonzero -> enforce monotonicity
% 0 -> don't enforce monotonicity.
% Default = 1;
%
% 'Continuation' = Continuation or not (1 or 0)
% Specifies the choice for a continuation scheme,
% in which we start with a large value of tau, and
% then decrease tau until the desired value is
% reached. At each value, the solution obtained
% with the previous values is used as initialization.
% Default = 0
%
% 'ContinuationSteps' = Number of steps in the continuation procedure;
% ignored if 'Continuation' equals zero.
% If -1, an adaptive continuation procedure is used.
% Default = -1.
%
% 'FirstTauFactor' = Initial tau value, if using continuation, is
% obtained by multiplying the given tau by
% this factor. This parameter is ignored if
% 'Continuation' equals zero.
% Default = such that the first tau is equal to
% 0.5*max(abs(AT(y))).
%
% 'True_x' = if the true underlying x is passed in
% this argument, MSE evolution is computed
%
% 'AlphaMin' = the alphamin parameter of the BB method.
% (See the paper for details)
% Default = 1e-30;
%
% 'AlphaMax' = the alphamax parameter of the BB method.
% (See the paper for details)
% Default = 1e30;
%
% 'Verbose' = work silently (0) or verbosely (1)
%
% ===================================================
% ============ Outputs ==============================
% x = solution of the main algorithm
%
% x_debias = solution after the debiasing phase;
% if no debiasing phase took place, this
% variable is empty, x_debias = [].
%
% objective = sequence of values of the objective function
%
% times = CPU time after each iteration
%
% debias_start = iteration number at which the debiasing
% phase started. If no debiasing took place,
% this variable is returned as zero.
%
% mses = sequence of MSE values, with respect to True_x,
% if it was given; if it was not given, mses is empty,
% mses = [].
% ========================================================
% test for number of required parametres
if (nargin-length(varargin)) ~= 3
error('Wrong number of required parameters');
end
% flag for initial x (can take any values except 0,1,2)
Initial_X_supplied = 3333;
% Set the defaults for the optional parameters
stopCriterion = 3;
tolA = 0.01;
tolD = 0.0001;
debias = 0;
maxiter = 10000;
maxiter_debias = 500;
miniter = 5;
miniter_debias = 5;
评论21