Synopsis

We describe a new algorithm for model-based MRI reconstruction with separate magnitude and phase regularization. The algorithm, named NAPALM, combines the existing proximal alternating linearized minimization (PALM) algorithm for nonsmooth and nonconvex optimization with Nesterov's acceleration and adaptive gradient (AdaGrad) acceleration methods. Results demonstrate the advantages of NAPALM over existing state-of-the-art algorithms.

Introduction

Many model-based MRI reconstruction methods assume that data acquisition is modeled as $$\mathbf{b}=\mathbf{A}\mathbf{f}+\mathbf{n},$$ where $$$\mathbf{b}$$$ is the measured data, $$$\mathbf{A}$$$ is the encoding matrix, $$$\mathbf{f}$$$ is the vector of image voxel values, and $$$\mathbf{n}$$$ is noise. Based on this forward model, image reconstruction is frequently formulated as $$\hat{\mathbf{f}}=\arg\min_{\mathbf{f}\in\mathbb{C}^N}\frac{1}{2}\|\mathbf{A}\mathbf{f}-\mathbf{b}\|_2^2+R(\mathbf{f}),$$ where $$$R(\cdot)$$$ is a regularization penalty.

In MRI (and in several other imaging modalities), the image $$$\mathbf{f}$$$ is complex-valued. The magnitude image and phase image often represent manifestations of different physical phenomenon, and have distinct spatial characteristics that are not always easily captured by a single unified regularization penalty. This has caused some authors to consider separate regularization of the magnitude and phase of the image [1-8]. These works often consider optimization problems of the form:$$\{\hat{\mathbf{m}},\hat{\mathbf{p}}\}=\arg\min_{\mathbf{m},\mathbf{p}\in\mathbb{R}^N}\frac{1}{2}\|\mathbf{A}(\mathbf{m}\odot e^{i\mathbf{p}})-\mathbf{b}\|_2^2+R_1(\mathbf{m})+R_2(e^{i\mathbf{p}}) \triangleq F(\mathbf{m},e^{i\mathbf{p}}),$$ where $$$\mathbf{m}$$$ and $$$\mathbf{p}$$$ represent the magnitude and phase of the image $$$\mathbf{f}$$$ with $$$\mathbf{f}=\mathbf{m}\odot e^{i\mathbf{p}}$$$, and $$$R_1(\cdot)$$$ and $$$R_2(\cdot)$$$ are corresponding regularization penalties. This formulation has the advantage that regularization can be customized independently for the magnitude and phase. However, this optimization problem is never convex and difficult to optimize.

Most existing methods use Alternating Minimization (AM), which iteratively alternates between estimating $$$\mathbf{m}$$$ and estimating $$$\mathbf{p}$$$ [2,4-8]. Estimating $$$\mathbf{m}$$$ (while holding $$$\mathbf{p}$$$ fixed) results in a standard regularized least-squares problem, and there are many efficient algorithms for this problem. However, estimating $$$\mathbf{p}$$$ (while holding $$$\mathbf{m}$$$ fixed) is nontrivial. Existing methods have frequently relied on the nonlinear conjugate gradient (NCG) algorithm [2,4-6,8] or phase cycling (PhsCyc) [7] to solve the phase subproblem.

In this work, we propose a new algorithm (named NAPALM), which is based on combining Nesterov's acceleration method [9], adaptive gradient (AdaGrad) methods [10,11], and the proximal alternating linearized minimization (PALM) algorithm [12].

Methods

The recent PALM algorithm [12] is the foundation of NAPALM. PALM is a simple and computationally efficient algorithm for solving nonconvex and nonsmooth optimization problems, and has global convergence properties. To make PALM even faster, we combine it with Nesterov's acceleration approach (a popular momentum-based approach for improving convergence speed) [9] and the AdaGrad approach [10,11] (which is also becoming increasingly popular for improving convergence speed, and rescales different components of the gradient based on the observed geometry of the cost function). Pseudo-code for NAPALM is shown in Fig. 1. NAPALM has three tuning parameters: the weight $$$\beta$$$ for the weighted average step in AdaGrad and two global step sizes $$$\gamma_{1},\gamma_{2}$$$ [10-12]. We use the recommended setting for $$$\beta=0.99$$$ [11] and manually tuned $$$\gamma_{1}$$$ and $$$\gamma_{2}$$$ around the recommended value $$$0.01$$$ [11] to make sure NAPALM has a good trade-off between accuracy and speed.

NAPALM was compared against AM+NCG [2,4-6,8] and AM+PhsCyc [7] in two different scenarios. The first scenario involves denoising an MRI brain image [8], with a total variation (TV) penalty applied to the magnitude and a quadratic smoothness penalty applied to the phase. The gold standard and noisy images are shown in Fig. 2. The second scenario involves compressed sensing and parallel imaging reconstruction of an MRI brain image [5-7], with $$$\ell_1$$$ Daubechies-4 wavelets penalties applied independently to the magnitude and phase. Reconstruction used SENSE parallel imaging [13] with 8 channels, an acceleration factor of 7.7, and Poisson disc undersampling. The gold standard images and the sampling mask are shown in Fig. 3.

Results

Convergence curves for both scenarios are shown in Fig. 4. As can be seen, NAPALM converges faster than AM+PhsCyc and AM+NCG in both cases. NAPALM and AM+NCG were observed to both be stable and converge to similar local minima. AM+PhsCyc was stable for compressed sensing (the scenario it was designed for [7]), but was unstable for denoising.

NAPALM is built from PALM [12] and two acceleration methods (Nesterov's method [9] and AdaGrad [10,11]). Additional experiments were performed to analyze the relative contribution of the acceleration methods. Results are shown in Fig. 5, and suggest that Nesterov's acceleration may have a bigger contribution than AdaGrad, although it's reasonable to use both together.

Conclusion

Acknowledgements

References

[1] M. Cetin and W. C. Karl. "Feature-enhanced synthetic aperture radar image formation based on nonquadratic regularization." IEEE Trans Image Process, vol. 10, pp. 623-631, 2001.

[2] J. A. Fessler and D. C. Noll. "Iterative image reconstruction in MRI with separate magnitude and phase regularization." Proc. IEEE Int Symp Biomed Imaging, 2004, pp. 209-212.

[3] A. Tuysuzoglu, J. M. Kracht, R. O. Cleveland, M. Cetin, and W. C. Karl. "Sparsity driven ultrasound imaging." J Acoust Soc Am, vol. 131, pp. 1271-1281, 2012.

[4] J. P. Haldar, Z. Wang, G. Popescu, and Z.-P. Liang. "Deconvolved spatial light interference microscopy for live cell imaging." IEEE Trans Biomed Eng, vol. 58, pp. 2489-2497, 2011.

[5] F. Zhao, D. C. Noll, J. F. Nielsen, and J. A. Fessler. "Separate magnitude and phase regularization via compressed sensing." IEEE Trans Med Imaging, vol. 31, pp. 1713-1723, 2012.

[6] M V. W. Zibetti and A. R. De Pierro. "Improving compressive sensing in MRI with separate magnitude and phase priors." Multidim Syst Signal Process, vol. 28, pp. 1109-1131, 2017.

[7] F. Ong, J. Y. Cheng, and M. Lustig. "General phase regularized reconstruction using phase cycling." Magn Reson Med, vol. 80, pp. 112-125, 2018.

[8] J. P. Haldar, Q. Fan, and K. Setsompop. "Whole-brain quantitative diffusion MRI at 660 $$$\mu$$$m resolution in 25 minutes using gSlider-SMS and SNR-enhancing joint reconstruction." Proc. ISMRM, 2016, p. 102.

[9] Y. Nesterov. "Introductory lectures on convex optimization." Springer, 2004.

[10] J. Duchi, E. Hazan, and Y. Singer. "Adaptive subgradient methods for online learning and stochastic optimization." J Mach Learn Res, vol. 12, pp. 2121-2159, 2011.

[11] D. P. Kingma and J. L. Ba. Adam: A method for stochastic optimization. ArXiv e-prints, Dec. 2014. arXiv: 1412.6980 [cs.LG].

[12] J. Bolte, S. Sabach, and M. Teboulle. "Proximal alternating linearized minimization for nonconvex and nonsmooth problems." Math Program, vol. 146, pp.459-494, 2014.

[13] K. P. Pruessmann, M. Weiger, M. B. Scheidegger, and P. Boesiger. "SENSE: Sensitivity encoding for fast MRI." Magn Reson Med, vol. 42, pp. 952-962, 1999.