An Approximate Message Passing Algorithm for Rapid Parameter-Free Compressed Sensing MRI

by   Charles Millard, et al.

For certain sensing matrices, the Approximate Message Passing (AMP) algorithm and more recent Vector Approximate Message Passing (VAMP) algorithm efficiently reconstruct undersampled signals. However, in Magnetic Resonance Imaging (MRI), where Fourier coefficients of a natural image are sampled with variable density, AMP and VAMP encounter convergence problems. In response we present a new approximate message passing algorithm constructed specifically for variable density partial Fourier sensing matrices with a sparse model on wavelet coefficients. For the first time in this setting a state evolution has been observed. A practical advantage of state evolution is that Stein's Unbiased Risk Estimate (SURE) can be effectively implemented, yielding an algorithm with no free parameters. We empirically evaluate the effectiveness of the parameter-free algorithm on simulated data and find that it converges over 5x faster and to a lower mean-squared error solution than Fast Iterative Shrinkage-Thresholding (FISTA).



There are no comments yet.


page 3

page 4


Approximate Message Passing with a Colored Aliasing Model for Variable Density Fourier Sampled Images

The Approximate Message Passing (AMP) algorithm efficiently reconstructs...

Performance Analysis of Approximate Message Passing for Distributed Compressed Sensing

Bayesian approximate message passing (BAMP) is an efficient method in co...

MRI Image Recovery using Damped Denoising Vector AMP

Motivated by image recovery in magnetic resonance imaging (MRI), we prop...

Performance Limits with Additive Error Metrics in Noisy Multi-Measurement Vector Problem

Real-world applications such as magnetic resonance imaging with multiple...

Compressive Imaging using Approximate Message Passing and a Markov-Tree Prior

We propose a novel algorithm for compressive imaging that exploits both ...

Denoising-based Turbo Message Passing for Compressed Video Background Subtraction

In this paper, we consider the compressed video background subtraction p...

Sparse Bayesian Learning Using Approximate Message Passing with Unitary Transformation

Sparse Bayesian learning (SBL) can be implemented with low complexity ba...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

I Introduction

We consider a complex linear regression problem, where complex data vector

is formed of noisy linear measurements of a signal of interest :


where and , where is the identity matrix. Here,

denotes the complex normal distribution with mean

, covariance and white phase. A well-studied approach is to seek a solution of


where is a model-based penalty function. Compressed sensing [8, 7] concerns the reconstruction of from underdetermined measurements . Commonly sparsity in is promoted by solving (2) with for sparse weighting and sparsifying transform .

The Approximate Message Passing (AMP) algorithm [9] is an iterative method that estimates in linear regression problems such as (1). At iteration AMP implements a denoiser on estimate with mean-squared error estimate . For instance, for problems of the form of (2), is the proximal operator associated with :


which is equal to soft thresholding in the case of and orthogonal . For certain sensing matrices and given mild conditions on , AMP’s state evolution guarantees that in the large system limit , , vector behaves as the original signal corrupted by zero-mean white Gaussian noise:



is an iteration-dependant standard deviation. The state evolution of AMP has been proven for real i.i.d. Gaussian measurements in

[4] and i.i.d. sub-Gaussian measurements in [3]. It has also been empirically shown that state evolution holds for uniformly undersampled Fourier measurements of a random artificial signal [9]. When state evolution holds, AMP is known to exhibit very fast convergence. However, for generic , the behavior of AMP is not well understood and it has been noted by a number of authors [21, 6, 15, 12, 22] that it can encounter convergence problems. The recently proposed Vector Approximate Message Passing (VAMP) [20] algorithm broadened the class of measurement matrices for which (4) holds, namely to those matrices that are ‘right-orthogonally invariant’, and was also found to perform very well on certain reconstruction tasks.

I-a Approximate message passing for compressed sensing MRI

In compressed sensing MRI [14], measurements are formed of undersampled Fourier coefficients, so that , where

is a 2D or 3D discrete Fourier transform and

is a undersampling mask that selects the th row of if for sampling set . The signal of interest is a natural image, so typically has a highly non-uniform spectral density that is concentrated at low frequencies. Accordingly, the sampling set

is usually generated such that there is a higher probability of sampling lower frequencies. This work considers an

with elements drawn independently from a Bernoulli distribution with non-uniform probability, such that

. In this variable density setting there are no theoretical guarantees for AMP or VAMP and in practice the behavior of (4) is not observed and the algorithms typically perform poorly.

This letter presents a new method for undersampled signal reconstruction which we term the Variable Density Approximate Message Passing (VDAMP) algorithm. For Fourier coefficients of a realistic image randomly sampled with variable density and orthogonal wavelet we have empirical evidence that a state evolution occurs. Unlike the white effective noise of (4), the

of VDAMP behaves as the ground truth corrupted by zero-mean Gaussian noise with a separate variance for each wavelet subband, such that


where and the effective noise covariance is diagonal so that for a with decomposition scales


where and refer to the variance and dimension of the th subband respectively. We refer to (5) as the colored state evolution of VDAMP.

Selecting appropriate regularisation parameters such as is a notable challenge in real-world compressed sensing MRI applications. We present an approach to parameter-free compressed sensing reconstruction using Stein’s Unbiased Risk Estimate (SURE) [23] in conjunction with VDAMP, building on the work on AMP introduced in [17]. A strength of automatic parameter tuning via SURE is that the it is possible to have a richer regularizer than would usually be feasible for a hand-tuned . We implement a variation on the SureShrink denoiser [10], using a iteration-dependant regularizer that has a separate sparse weighting per subband:


where is the -dimensional column vector of ones. SURE has previously been employed for parameter-free compressed sensing MRI in [13], where the Fast Iterative Shrinking-Thresholding Algorithm (FISTA) algorithm [5] was used with SureShrink in place of the usual shrinkage step. This algorithm is herein referred to as SURE-IT. The effective noise of SURE-IT is highly non-Gaussian, so deviates from a proper theoretical basis for using SURE for threshold selection. To our knowledge, VDAMP is the first algorithm for variable density Fourier sampling of natural images where a state evolution has been observed.

Ii Description of algorithm

For AMP and VAMP, (4) states that the effective noise is white, so can be fully characterised by a scalar . This is appropriate for the kind of uniform measurement matrices and separable, identical sparse signal models that are often encountered in abstract compressed sensing problems. However, when Fourier coefficients of a natural image are sampled the effective noise is colored, so is poorly represented by a scalar [24]. VDAMP models the richer structure of effective noise in the variable density setting with a vector that has one real number per wavelet subband.

The VDAMP algorithm is shown in Algorithm 1. Here, is the diagonal matrix formed of sampling probabilities for . The function refers to some denoiser with a colored effective noise model. The notation in line 7 refers to the (sub)-gradient of the denoiser averaged over subbands, so that for decomposition scales has unique entries, having the same structure as the of (7). In line 8, the notation is used for entry-wise multiplication and for entry-wise division. refers to the entry-wise absolute magnitude of a vector or matrix.

Require: Sensing matrix , orthogonal wavelet transform , probability matrix , measurements , denoiser and number of iterations .

1:  Set and compute
2:  for  do
9:  end for
10:  return  
Algorithm 1 VDAMP

The original AMP paper considers and an i.i.d. Gaussian that is normalized such that . This ensures that

is an unbiased estimate of

. For variable density sampling, the correct normalisation can be achieved by rescaling (1) by . In VDAMP this is manifest in the gradient step of lines 3-4, which feature a crucial weighting by . This provides the correct normalisation in expectation over : , which implies that


for any . Such a rescaling is referred to as ‘density compensation’ in the MRI literature [18, 19], and was used in the original compressed sensing MRI paper with zero-filling to generate a unregularized, non-iterative baseline [14] . Line 5 of Algorithm 1 computes an estimate of the colored effective noise variance . It can be shown that is an unbiased estimate of the expected entry-wise squared error:


for any . We assume both estimators in (8) and (10) concentrate around their expectation, and leave the study of the constraints this imposes on for future works. Note that has unique columns, so for fixed the complexity of VDAMP is governed by and , whose fast implementations have complexity .

Lines 6-8 are the model-based regularization phase from VAMP, but with a colored effective noise model, as in [11]. This phase includes the message passing Onsager correction term, which we have observed leads to the Gaussian effective noise of (5). Line 10 implements an unweighted gradient step that enforces exact consistency of the image estimate with the measurements.

Iii Numerical experiments

Fig. 1: NMSE of VDAMP, FISTA and SURE-IT versus compute time for undersampling factors 4, 6 and 8.

In the experiments presented in this section the denoiser was the complex soft thresholding operator with an automatically tuned subband-dependant threshold. In other words, we used a complex adaption of SURE to approximately solve


where is the entry-wise square root of and is of the form of (7). (11) was solved using a procedure similar to the usual SureShrink but with a version of SURE adapted for effective noise that is complex and colored. Consider a ground truth subband corrupted by complex Gaussian noise with white phase: . For complex soft thresholding with threshold , it can be shown that


is an unbiased estimate of the risk. For each subband the optimal threshold was estimated via


by evaluating (12) for trial thresholds . For large dimension

one would expect by the law of large numbers that cSURE is close to the true risk, and for the threshold to be almost optimal. Since a larger number of decomposition scales

give subbands with lower dimensionality, there is a trade-off between the size of and the quality of threshold selection with cSURE.

We considered the reconstruction of a Shepp-Logan artificially corrupted with complex Gaussian noise with white phase so that . We assumed that was known a priori; in practice it can be well estimated with an empty prescan. All sampling probabilities were generated using the polynomial variable density sampling function from the Sparse MRI package available at We considered a Haar wavelet at decomposition scales, which are referred to as scales 1-4, where scale 1 is the finest and scale 4 is the coarsest. All experiments were conducted in MATLAB 9.4 on a 1.70 GHz Intel i5-8350U processor with 8GB of RAM, and can be reproduced with code available at

Fig. 2: Reconstruction of a 8x undersampled Shepp-Logan after 2 seconds with (b) FISTA (NMSE = -19.3dB) (c) SURE-IT (NMSE = -16.6dB) and (d) VDAMP (NMSE -34.9dB), where absolute values are shown.

Iii-a Comparison with FISTA and SURE-IT

To establish comparative reconstruction quality and convergence speed, VDAMP was compared with FISTA and SURE-IT. For FISTA we used a sparse weighting tuned with an exhaustive search so that the mean-squared error was minimised after 10 seconds. For SURE-IT the mean-squared error estimate was updated using the ground truth: , and (11) with was implemented. All three algorithms were initialised with a vector of zeros.

Three sampling sets were generated at undersampling factors of approximately 4, 6 and 8, and VDAMP, FISTA and SURE-IT were run for 10 seconds. Fig. 1 shows the NMSE as a function of time for each algorithm. at every iteration was calculated so that exact data consistency was ensured, as in line 10 of VDAMP. The mean per-iteration compute time was 0.065s for FISTA, 0.077s for SURE-IT, and 0.091s for VDAMP. Fig. 2 shows the ground truth image, the sampling set, and FISTA and VDAMP reconstruction at undersampling factor 8 after 2s. The entry-wise error is also shown for both algorithms.

Iii-B State evolution and error tracking

Fig. 3: of VDAMP for , and .
Fig. 4:

Normalized quantile-quantile plots against a Gaussian of three subbands of

for , and

in blue, and points along a straight line in red. The real part is plotted in the top and bottom rows and the imaginary is plotted in the middle row. Linearity of the blue points indicates that that the data comes from a Gaussian distribution. Finite dimensional effects causing small deviations from an exact Gaussian are more apparent at coarse scales, where the dimension is smaller.

Fig. 3: of VDAMP for , and .

This section empirically analyses VDAMP in greater depth, continuing to use the illustrative case of undersampling factor 8. In Fig. 4, the absolute value of the residual of is shown for three representative iterations: , and . For the same iterations, Fig. 4 shows quantile-quantile plots of the real parts of three illustrative subbands of : the diagonal detail at scale 1, the horizontal detail at scale 2 and the vertical detail scale 4. These figures provide empirical evidence that the effective noise of VDAMP evolves as (5).

The efficacy of automatic threshold selection with cSURE depends on how accurately the diagonal of from (6) is modelled by . For Fig. 5 shows the ground truth NMSE of the wavelet coefficients at all four scales and the prediction of VDAMP, where the NMSE is per subband: .

Fig. 5: NMSE versus iteration number of for each subband. Lines show the actual NMSE and crosses show the predictions from .

Iv Conclusions

Like the original AMP paper [8], the claim of a state evolution has been substantiated with empirical evidence only. Theoretical work is required to establish the generality of the state evolution observed in these experiments.

In order for the algorithm to be compatible with realistic MRI acquisition protocols, practical developments are required to allow sampling with readout lines and multiple coils. The experiments conducted in this letter indicate that given such developments, VDAMP would have a number of significant advantages over algorithms currently used in the compressed sensing MRI literature. VDAMP’s state evolution provides an informed and efficient way to tune model parameters via SURE, implying that a single algorithm can be used for any image type and variable density scheme without the need for manual adjustment. More degrees of freedom are allowed in the model, enabling higher order prior information such as anisotropy, variability across scales and structured sparsity, without the need to estimate the structure a priori such as in model-based compressed sensing

[2]. The substantial reduction in required compute time observed in these experiments could be clinically valuable as real MRI reconstruction tasks typically have very high dimensionality, so it is often only viable to run a reconstruction algorithm for a limited number of iterations.

It is known that the state evolution of AMP holds for a wide range of denoisers [16]. In a recent survey [1] a number of standard compressed sensing algorithms that leverage image denoisers designed for Gaussian noise were shown to perform well on MRI reconstruction tasks, despite the mismatch between the effective noise and its model. A sophisticated denoiser equipped to deal with wavelet coefficients corrupted with known colored Gaussian noise would be expected to perform well in conjunction with VDAMP.


  • [1] R. Ahmad, C. A. Bouman, G. T. Buzzard, S. Chan, E. T. Reehorst, and P. Schniter (2019-03) Plug and play methods for magnetic resonance imaging. External Links: Link Cited by: §IV.
  • [2] R. G. Baraniuk, V. Cevher, M. F. Duarte, and C. Hegde (2010-04) Model-Based Compressive Sensing. IEEE Transactions on Information Theory 56 (4), pp. 1982–2001. External Links: Link, Document, ISSN 0018-9448 Cited by: §IV.
  • [3] M. Bayati, M. Lelarge, and A. Montanari (2015-04)

    Universality in polytope phase transitions and message passing algorithms

    The Annals of Applied Probability 25 (2), pp. 753–822. External Links: Link, Document, ISSN 1050-5164 Cited by: §I.
  • [4] M. Bayati and A. Montanari (2011-02) The Dynamics of Message Passing on Dense Graphs, with Applications to Compressed Sensing. IEEE Transactions on Information Theory 57 (2), pp. 764–785. External Links: Link, Document, ISSN 0018-9448 Cited by: §I.
  • [5] A. Beck and M. Teboulle (2009-01) A Fast Iterative Shrinkage-Thresholding Algorithm for Linear Inverse Problems. SIAM Journal on Imaging Sciences 2 (1), pp. 183–202. External Links: Link, Document, ISSN 1936-4954 Cited by: §I-A.
  • [6] F. Caltagirone, L. Zdeborova, and F. Krzakala (2014-06) On convergence of approximate message passing. In 2014 IEEE International Symposium on Information Theory, pp. 1812–1816. External Links: Link, ISBN 978-1-4799-5186-4, Document Cited by: §I.
  • [7] E.J. Candes, J. Romberg, and T. Tao (2006-02) Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information. IEEE Transactions on Information Theory 52 (2), pp. 489–509. External Links: Link, Document, ISSN 0018-9448 Cited by: §I.
  • [8] D.L. Donoho (2006-04) Compressed sensing. IEEE Transactions on Information Theory 52 (4), pp. 1289–1306. External Links: Link, Document, ISSN 0018-9448 Cited by: §I, §IV.
  • [9] D. L. Donoho, A. Maleki, and A. Montanari (2009-11) Message-passing algorithms for compressed sensing.. Proceedings of the National Academy of Sciences of the United States of America 106 (45), pp. 18914–9. External Links: Link, Document, ISSN 1091-6490 Cited by: §I.
  • [10] D. L. Donoho and I. M. Johnstone (1995-12) Adapting to Unknown Smoothness via Wavelet Shrinkage. Journal of the American Statistical Association 90 (432), pp. 1200–1224. External Links: Link, Document, ISSN 0162-1459 Cited by: §I-A.
  • [11] A. Fletcher, M. Sahraee-Ardakan, S. Rangan, and P. Schniter (2016-02) Expectation consistent approximate inference: Generalizations and convergence. In IEEE International Symposium on Information Theory - Proceedings, Vol. 2016-Augus, pp. 190–194. External Links: Link, ISBN 9781509018062, Document, ISSN 21578095 Cited by: §II.
  • [12] Q. Guo and J. Xi (2015-04) Approximate Message Passing with Unitary Transformation. External Links: Link Cited by: §I.
  • [13] K. Khare, C. J. Hardy, K. F. King, P. A. Turski, and L. Marinelli (2012-11) Accelerated MR imaging using compressive sensing with no free parameters. Magnetic Resonance in Medicine 68 (5), pp. 1450–1457. External Links: Link, Document, ISSN 07403194 Cited by: §I-A.
  • [14] M. Lustig, D. Donoho, and J. M. Pauly (2007-12) Sparse MRI: The application of compressed sensing for rapid MR imaging. Magnetic Resonance in Medicine 58 (6), pp. 1182–1195. External Links: Link, Document, ISSN 07403194 Cited by: §I-A, §II.
  • [15] J. Ma and L. Ping (2017) Orthogonal AMP. IEEE Access 5, pp. 2020–2033. External Links: Link, Document, ISSN 2169-3536 Cited by: §I.
  • [16] C. A. Metzler, A. Maleki, and R. G. Baraniuk (2016-09) From Denoising to Compressed Sensing. IEEE Transactions on Information Theory 62 (9), pp. 5117–5144. External Links: Link, Document, ISSN 0018-9448 Cited by: §IV.
  • [17] A. Mousavi, A. Maleki, and R. G. Baraniuk (2013-10) Parameterless Optimal Approximate Message Passing. External Links: Link Cited by: §I-A.
  • [18] J. G. Pipe and P. Menon (1999-01) Sampling density compensation in MRI: Rationale and an iterative numerical solution. Magnetic Resonance in Medicine 41 (1), pp. 179–186. External Links: Link, Document, ISSN 0740-3194 Cited by: §II.
  • [19] K. P. Pruessmann, M. Weiger, P. Börnert, and P. Boesiger (2001-10) Advances in sensitivity encoding with arbitrary k-space trajectories. Magnetic Resonance in Medicine 46 (4), pp. 638–651. External Links: Link, Document Cited by: §II.
  • [20] S. Rangan, P. Schniter, and A. K. Fletcher (2019) Vector Approximate Message Passing. IEEE Transactions on Information Theory, pp. 1–1. External Links: Link, Document, ISSN 0018-9448 Cited by: §I.
  • [21] S. Rangan, P. Schniter, and A. Fletcher (2014-06) On the convergence of approximate message passing with arbitrary matrices. In 2014 IEEE International Symposium on Information Theory, pp. 236–240. External Links: Link, ISBN 978-1-4799-5186-4, Document Cited by: §I.
  • [22] S. Rangan, P. Schniter, E. Riegler, A. K. Fletcher, and V. Cevher (2016-12) Fixed Points of Generalized Approximate Message Passing With Arbitrary Matrices. IEEE Transactions on Information Theory 62 (12), pp. 7464–7474. External Links: Link, Document, ISSN 0018-9448 Cited by: §I.
  • [23] C. M. Stein (1981-11) Estimation of the Mean of a Multivariate Normal Distribution. The Annals of Statistics 9 (6), pp. 1135–1151. External Links: Link, Document, ISSN 0090-5364 Cited by: §I-A.
  • [24] P. Virtue and M. Lustig (2017-12) The Empirical Effect of Gaussian Noise in Undersampled MRI Reconstruction. Tomography 3 (4), pp. 211–221. External Links: Link, Document, ISSN 2379-1381 Cited by: §II.