magni.cs.reconstruction.amp._algorithm module

Module providing the core Approximate Message Passing (AMP) algorithm.

Routine listings

run(y, A)

Run the AMP reconstruction algorithm.

See also

magni.cs.reconstruction.amp._config

Configuration options.

magni.cs.reconstruction.amp.threshold_operator

Threshold operators.

magni.cs.reconstruction.amp.stop_criterion

Stop criteria.

magni.cs.reconstruction.amp.util

Utilities.

Notes

The default configuration of the AMP algorithm provides the Donoho, Maleki, Montanari (DMM) AMP from [1] using the soft threshold with the “residual” thresholding strategy from [2].

References

[1]	D. L. Donoho, A. Maleki, and A. Montanari, “Message-passing algorithms for compressed sensing”, Proceedings of the National Academy of Sciences of the United States of America, vol. 106, no. 45, p. 18914-18919, Nov. 2009.

[2]	A. Montanari, “Graphical models concepts in compressed sensing” in Compressed Sensing: Theory and Applications, Y. C. Eldar and G. Kutyniok (Ed.), Cambridge University Press, ch. 9, pp. 394-438, 2012.

magni.cs.reconstruction.amp._algorithm.run(y, A)

Run the AMP reconstruction algorithm.

Parameters:

• y (ndarray) – The m x 1 measurement vector.
• A (ndarray or magni.utils.matrices.{Matrix, MatrixCollection}) – The m x n matrix which is the product of the measurement matrix and the dictionary matrix.

Returns:

• alpha (ndarray) – The n x 1 reconstructed coefficient vector.
• history (dict, optional) – The dictionary of various measures tracked in the GAMP iterations.

See also

magni.cs.reconstruction.amp._config()

Configuration options.

magni.cs.reconstruction.amp.threshold_operator()

Threshold operators.

magni.cs.reconstruction.amp.stop_criterion()

Stop criteria.

magni.cs.reconstruction.amp.util()

Utilities.

Notes

Optionally, the algorithm may be configured to save and return the iteration history along with the reconstruction result. The returned history contains the following:

• alpha_bar : Coefficient estimates (the reconstruction coefficients).
• MSE : solution mean squared error (if the true solution is known).
• threshold_parameters : The state of the threshold parameters.
• stop_criterion : The currently used stop criterion.
• stop_criterion_value : The value of the stop criterion.
• stop_iteration : The iteration at which the algorithm stopped.
• stop_reason : The reason for termination of the algorithm.

Examples

For example, recovering a vector from noiseless measurements using AMP with soft thresholding

>>> import numpy as np
>>> from magni.cs.reconstruction.amp import run, config
>>> from magni.cs.reconstruction.amp.util import theta_mm
>>> np.random.seed(seed=6028)
>>> k, m, n = 10, 200, 400
>>> A = 1 / np.sqrt(m) * np.random.randn(m, n)
>>> alpha = np.zeros((n, 1))
>>> alpha[:k] = np.random.normal(scale=2, size=(k, 1))
>>> np_printoptions = np.get_printoptions()
>>> np.set_printoptions(suppress=True, threshold=k+2)
>>> alpha[:k + 2]
array([[ 1.92709461],
       [ 0.74378508],
       [-3.2418159 ],
       [-1.32277347],
       [ 0.90118   ],
       [-0.19157262],
       [ 0.82855712],
       [ 0.24817994],
       [-1.43034777],
       [-0.21232344],
       [ 0.        ],
       [ 0.        ]])
>>> y = A.dot(alpha)
>>> threshold_params = {'threshold_level_update_method': 'residual',
... 'theta': theta_mm(float(m) / n), 'tau_hat_sq': 1.0}
>>> config['threshold_parameters'] = threshold_params
>>> alpha_hat = run(y, A)
>>> alpha_hat[:k + 2]
array([[ 1.92612707],
       [ 0.74185284],
       [-3.24179313],
       [-1.32200929],
       [ 0.90089208],
       [-0.19097308],
       [ 0.82658038],
       [ 0.24515825],
       [-1.42980997],
       [-0.2111469 ],
       [ 0.00002815],
       [ 0.00047293]])
>>> np.sum(np.abs(alpha - alpha_hat) > 1e-2)
0

or recover the same vector returning a history comparing the pr. iteration solution to the true vector

>>> config['report_history'] = True
>>> config['true_solution'] = alpha
>>> alpha_hat, history = run(y, A)
>>> alpha_hat[:k + 2]
array([[ 1.92612707],
       [ 0.74185284],
       [-3.24179313],
       [-1.32200929],
       [ 0.90089208],
       [-0.19097308],
       [ 0.82658038],
       [ 0.24515825],
       [-1.42980997],
       [-0.2111469 ],
       [ 0.00002815],
       [ 0.00047293]])
>>> np.array(history['MSE']).reshape(-1, 1)[1:11]
array([[ 0.01273323],
       [ 0.0053925 ],
       [ 0.00276334],
       [ 0.00093472],
       [ 0.0004473 ],
       [ 0.00021142],
       [ 0.00009618],
       [ 0.00004371],
       [ 0.00002275],
       [ 0.00001053]])

or recover the same vector using the median based update method

>>> config['report_history'] = False
>>> threshold_params['threshold_level_update_method'] = 'median'
>>> config['threshold_parameters'] = threshold_params
>>> alpha_hat = run(y, A)
>>> alpha_hat[:k + 2]
array([[ 1.92577079],
       [ 0.7413161 ],
       [-3.24131125],
       [-1.32044146],
       [ 0.89998407],
       [-0.18981235],
       [ 0.82584016],
       [ 0.24445332],
       [-1.42934708],
       [-0.21043221],
       [ 0.        ],
       [ 0.        ]])
>>> np.sum(np.abs(alpha - alpha_hat) > 1e-2)
0

>>> np.set_printoptions(**np_printoptions)

magni.cs.reconstruction.amp._algorithm._get_amp_initialisation(y, A)

Return an initialisation of the AMP algorithm.

Parameters:

• y (ndarray) – The m x 1 measurement vector.
• A (ndarray or magni.utils.matrices.{Matrix, MatrixCollection}) – The m x n matrix which is the product of the measurement matrix and the dictionary matrix.

Returns:

init (dict) – The initialisation of the the AMP algorithm.