compsense.algorithms.TwIST_raw

compsense.algorithms.TwIST_raw(y, A, tau, psi_function=None, phi_function=None, lam1=0.0001, alpha=0, beta=0, AT=None, stop_criterion=1, tolA=0.01, debias=0, tolD=0.001, maxiter=1000, miniter=5, maxiter_debias=200, miniter_debias=5, init=0, enforce_monotone=True, sparse=True, true_x=None, verbose=True)

Two-step Iterative Shrinkage/Thresholding Algorithm for Linear Inverse Problems. Solves the regularization problem

arg min_x = 0.5*|| y - A x ||_2^2 + \tau \phi( x ),

where A is a generic matrix and phi() is a regularizarion function such that the solution of the denoising problem

\Psi_{\tau}(y) = arg min_x = 0.5*|| y - x ||_2^2 + \tau \phi( x ),

is known.

For further details about the TwIST algorithm, see the paper:

J. Bioucas-Dias and M. Figueiredo, "A New TwIST: Two-Step
Iterative Shrinkage/Thresholding Algorithms for Image 
Restoration",  IEEE Transactions on Image processing, 2007.

and

J. Bioucas-Dias and M. Figueiredo, "A Monotonic Two-Step 
Algorithm for Compressive Sensing and Other Ill-Posed 
Inverse Problems", submitted, 2007.
Parameters :

y : array,

1D vector or 2D array (image) of observations.

A : {array, function handle},

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 Av, 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 Ax; another handle to a function AT which computes products of the form A^Tx is also required in this case. The size of x is determined as the size of the result of applying AT.

tau : float,

regularization parameter, usually a non-negative real parameter of the objective function (see above).

psi_function : function handle, optional

handle to denoising function (the default is soft threshold)

phi_function : function handle, optional

handle to regularizer needed to compute the objective function. (the default = ||x||_1)

lam1 : float, optional (default=0.04)

parameter of the TwIST algorithm: Optimal choice: lam1 = min eigenvalue of :A^T*A. If min eigenvalue of A^T*A equals 0, or unknwon, set lam1 to a value much smaller than 1. Rule of Thumb:

  • lam1=1e-4 for severyly ill-conditioned problems
  • lam1=1e-2 for mildly ill-conditioned problems
  • lam1=1 for A unitary direct operators

Note

If max eigenvalue of :A^T*A > 1, the algorithm may diverge. This is to be avoided by taking one of the follwoing measures:

  1. Set enforce_monotone=True (default)
  2. Solve the equivalenve minimization problem

min_x = 0.5*|| (y/c) - (A/c) x ||_2^2 + (tau/c^2) \phi( x ),

where c > 0 ensures that max eigenvalue of :(A^TA/c^2) \leq 1.

alpha : float, optional (default=calculated as function of lam1)

parameter alpha of TwIST (see ex. (22) of the paper)

beta : float, optional (default=calculated as function of lam1)

parameter beta of twist (see ex. (23) of the paper)

AT : function handle, optional

function that implements the multiplication by the conjugate of A, when A is a function handle. If A is an array, AT is ignored.

stop_criterion : {0, 1, 2, 3}, optional (default=0)

type of stopping criterion to use

  • stop_criterion=0 algorithm stops when the relative change in the number of non-zero components of the estimate falls below tolA
  • stop_criterion=1 stop when the relative change in the objective function falls below tolA
  • stop_criterion=2 stop when the relative norm of the difference between two consecutive estimates falls below tolA
  • stop_criterion=3 stop when the objective function becomes equal or less than tolA.

tolA : float, optional (default=0.01)

stopping threshold.

debias : bool, optional (default=False)

debiasing option

Note

Debiasing is an operation aimed at the computing the solution of the LS problem

arg min_x = 0.5*|| y - A^T x ||_2^2

where :A^T is the submatrix of A obatained by deleting the columns of A ``corresponding of components of ``x set to zero by the TwIST algorithm

tolD : float, optional (default=0.0001)

stopping threshold for the debiasing phase. If no debiasing takes place, this parameter, is ignored.

maxiter : int, optional (default=1000)

maximum number of iterations allowed in the main phase of the algorithm.

miniter : int, optional (default=5)

minimum number of iterations performed in the main phase of the algorithm.

maxiter_debias : int, optional (default=5)

maximum number of iterations allowed in the debising phase of the algorithm.

miniter_debias : int, optional (default=5)

minimum number of iterations to perform in the debiasing phase of the algorithm.

init : {0, 1, 2, array}, optional (default=0)

must be one of

  • init=0 Initialization at zero.
  • init=1 Random initialization.
  • init=2 initialization with :A^Ty.
  • init=array initialization provided by the user.

enforce_monotone : bool, optional (default=True)

enforce monotonic decrease in f.

sparse : bool, optional (default=True)

Accelarates the convergence rate when the regularizer \Phi(x) is sparse inducing, such as :||x||_1.

true_x : array, optional (default=None)

if the true underlying x is passed in this argument, MSE evolution is computed

Verbose : bool, optional (default=False)

work silently (False) or verbosely (True)

Returns :

x : array,

solution of the main algorithm

x_debias : array,

solution after the debiasing phase; if no debiasing phase took place, this variable is [].

objective : array,

sequence of values of the objective function

times : arraym

CPU time after each iteration

debias_start : int,

iteration number at which the debiasing phase started. If no debiasing took place, this variable is returned as zero.

mses : array,

sequence of MSE values, with respect to true_x, if it was given; if it was not given, mses is [].

max_svd : float,

inverse of the scaling factor, determined by TwIST, applied to the direct operator (A/max_svd) such that every IST step is increasing.

Previous topic

compsense.algorithms.TwIST

Next topic

License for pycompsense

This Page