FISTA Algorithm

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

Bases: Algorithm

Fast Iterative Shrinking-Thresholding Algorithm (FISTA)

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
>>> fista = fma.FISTA(C, numLambda=0.005, numMaxSteps=100)
>>> y = fista.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

fmatAfm.Matrix

the system matrix

arrBnp.ndarray

the measurement vector

numLambdafloat, optional

the thresholding parameter; default is 0.1

numMaxStepsint, optional

maximum number of steps; default is 100

Returns

np.ndarray

solution array

__init__(fmatA, **kwargs)

softThreshold(arrX, numAlpha)

Do a soft thresholding step.