magni.cs.reconstruction.gamp.channel_initialisation module¶
Module providing utility functions for initilisation of in- and output channels in the Generalised Approximate Message Passing (GAMP) algorithm.
Routine listings¶
- get_em_bg_amp_initialisation(problem_params, method=’vila’)
- Get initial parameters for EM Bernoulli-Guassian AMP.
- rho_se(delta, zeta, resolution=1000)
- Return the theoretical noiseless LASSO phase transition.
-
magni.cs.reconstruction.gamp.channel_initialisation.
get_em_bg_amp_initialisation
(problem_params, method='vila')[source]¶ Get initial parameters for EM Bernoulli-Guassian AMP.
If the initialisation method is vila then the scheme from [1] is used. If it is krzakala then the scheme from [2] is used.
Parameters: - problem_params (dict) – The problem parameters used to compute the initialisation.
- method (str) – The initialisation method to use.
Returns: - tau (float) – The initial sparsity level.
- sigma_sq (float) – The initial AWGN noise level.
- theta_tilde (float) – The initial Gaussian prior variance.
See also
magni.cs.reconstruction.gamp._output_channel.wrap_calculate_using_AWGN()
- Related output channel.
magni.cs.reconstruction.gamp._input_channels.wrap_calculate_using_iidsGB()
- Related input channel.
Notes
Independently of the choice of method, the problem_params are:
- y: the measurements
- A: the system matrix
If method is vila, one must also specify:
- SNR: the signal-to-noise ratio
References
[1] J. P. Vila and P. Schniter, “Expectation-Maximization Gaussian-Mixture Approximate Message Passing”, IEEE Transactions on Signal Processing, 2013, vol. 61, no. 19, pp. 4658-4672, Oct. 2013. [2] F. Krzakala, M. Mezard, F. Sausset, Y. Sun, and L. Zdeborova, “Probabilistic reconstruction in compressed sensing: algorithms, phase diagrams, and threshold achieving matrices”, Journal of Statistical Mechanics: Theory and Experiment, vol. P08009, pp. 1-57, Aug. 2012. Examples
For example, get the “vila” initialisation for a SNR of 100
>>> import numpy as np >>> from magni.cs.reconstruction.gamp.channel_initialisation import ... get_em_bg_amp_initialisation as get_init >>> np.random.seed(6012) >>> A = np.random.randn(20, 40) >>> y = np.random.randn(20, 1) >>> problem_params = {'A': A, 'y': y, 'SNR': 100} >>> init_tup = get_init(problem_params, method='vila') >>> [round(float(elem), 3) for elem in init_tup] [0.193, 0.01, 0.123]
or get the corresponding “krzakala” initialisation
>>> del problem_params['SNR'] >>> init_tup = get_init(problem_params, method='krzakala') >>> [round(float(elem), 3) for elem in init_tup] [0.05, 1.0, 0.479]
-
magni.cs.reconstruction.gamp.channel_initialisation.
rho_se
(delta, zeta, resolution=1000)[source]¶ Return the theoretical noiseless LASSO phase transition.
Parameters: - delta (float) – The under sampling ratio.
- zeta ({1, 2}) – The “problem” to get the phase transition for.
- resolution (int) – The resolution used in the brute force optimisation.
Returns: rho_se (float) – The phase transition value.
Notes
The theoretical noiseless LASSO phase transition is computed based on eq. 5 in [3]. A simple brute force optimisation with the specified resolution of that expression is used to find the phase transition. The “problems”, for which the phase transition may be computed, are:
- Sparse nonnegative vectors
- Sparse signed vectors
References
[3] 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, pp. 18914-18919, Nov. 2009. Examples
For example, find a phase transition value for a sparse signed vector:
>>> from magni.cs.reconstruction.gamp.channel_initialisation import rho_se >>> round(float(rho_se(0.19, 2)), 3) 0.238
or find the corresponding value for a sparse nonnegative vector
>>> round(float(rho_se(0.19, 1)), 3) 0.318
and a few more examples
>>> round(float(rho_se(0.0, 1)), 3) 0.0 >>> round(float(rho_se(0.0, 2)), 3) 0.0 >>> round(float(rho_se(0.5, 1)), 3) 0.558 >>> round(float(rho_se(0.5, 2)), 3) 0.386 >>> round(float(rho_se(1.0, 1)), 3) 0.95 >>> round(float(rho_se(1.0, 2)), 3) 0.95