# ISTA Algorithm¶

class fastmat.algorithms.ISTA(fmatA, **kwargs)

Bases: fastmat.algorithms.Algorithm.Algorithm

Iterative Soft Thresholding Algorithm

Definition and Interface: For a given matrix $$A \in \mathbb{C}^{m \times N}$$ with $$m \ll N$$ and a vector $$b \in \mathbb{C}^m$$ we approximately solve

$\min\limits_{ x \in \mathbb{C}^N}\Vert{ A \cdot x - b}\Vert^2_2 + \lambda \cdot \Vert x \Vert_1,$

where $$\lambda > 0$$ is a regularization parameter to steer the trade-off between data fidelity and sparsity of the solution.

>>> # import the packages
>>> import numpy.linalg as npl
>>> import numpy as np
>>> import fastmat as fm
>>> import fastmat.algorithms as fma
>>> # define the dimensions and the sparsity
>>> n, k = 512, 3
>>> # define the sampling positions
>>> t = np.linspace(0, 20 * np.pi, n)
>>> # construct the convolution matrix
>>> c = np.cos(2 * t) * np.exp(-t ** 2)
>>> C = fm.Circulant(c)
>>> # create the ground truth
>>> x = np.zeros(n)
>>> x[np.random.choice(range(n), k, replace=0)] = 1
>>> b = C * x
>>> # reconstruct it
>>> ista = fma.ISTA(C, numLambda=0.005, numMaxSteps=100)
>>> y = ista.process(b)
>>> # test if they are close in the
>>> # domain of C
>>> print(npl.norm(C * y - b))


We solve a sparse deconvolution problem, where the atoms are harmonics windowed by a gaussian envelope. The ground truth $$x$$ is build out of three pulses at arbitrary locations.

Note

The proper choice of $$\lambda$$ is crucial for good perfomance of this algorithm, but this is not an easy task. Unfortunately we are not in the place here to give you a rule of thumb what to do, since it highly depends on the application at hand. Again, consult [1] for any further considerations of this matter.

 [1] Amir Beck, Marc Teboulle, “A Fast Iterative Shrinkage-Thresholding Algorithm for Linear Inverse Problems”, SIAM Journal on Imaging Sciences, 2009, Vol. 2, No. 1 : pp. 183-202
Parameters: fmatA : fm.Matrix the system matrix arrB : np.ndarray the measurement vector numLambda : float, optional the thresholding parameter; default is 0.1 numMaxSteps : int, optional maximum number of steps; default is 100 np.ndarray solution array
__init__(fmatA, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

softThreshold(arrX, numAlpha)

Do a soft thresholding step.