Carlos Milovic^{1}, Kristian Bredies^{2}, Christian Langkammer^{3}, and Karin Shmueli^{4}

^{1}University College London, London, United Kingdom, ^{2}Institute of Mathematics and Scientific Computing, University of Graz, Graz, Austria, ^{3}Department of Neurology, Medical University of Graz, Graz, Austria, ^{4}Department of Medical Physics and Biomedical Engineering, University College London, London, United Kingdom

Single-step quantitative susceptibility mapping (QSM) algorithms simplify the processing pipeline and promise to be more robust against background fields than traditional two-step methods but they often underestimate tissue susceptibilities. Here, we propose a highly efficient gradient descent Tikhonov-regularized proximal solver and a highly accurate ADMM TV-regularized proximal solver to improve the accuracy of two Laplacian-based single-step methods. Our solvers outperformed current single-step methods and showed in-vivo performance very similar to traditional two-step methods. This will simplify QSM processing pipelines, allowing further automation in future, although more research is needed to improve robustness against noise and boundary-conditions-related artifacts.

$$Laplacian-model:argmin_χ\frac{1}{2}\parallel m(\triangledown^2(d*χ)-\triangledown^2Φ)\parallel^2_2+Ω(χ)$$

where $$$m$$$ is the ROI binary mask, $$$\triangledown^2$$$ is the Laplacian operator, $$$d$$$ is the dipole kernel convolving the susceptibility distribution χ, Φ is the total field map and Ω(χ) is a suitable regularizer. Alternatively, it was proposed

$$Graz-model:argmin_{χ,Ψ}\frac{1}{2}\parallelΨ\parallel^2_2 +\frac{σ}{2}\parallel m(\triangledown^2(d*χ)-\triangledown^2Φ-\triangledown^2Ψ)\parallel^2_2+Ω(χ)$$

Ψ is theorized to contain phase discrepancies, such as the background fields. Here σ is a penalty weight. In these two models, the Laplacian and dipole kernel convolution are merged into a single wave-like operator: $$$\frac{1}{3}\partial_x^2+\frac{1}{3}\partial_y^2-\frac{2}{3}\partial_z^2$$$.

We propose two optimization solvers and corresponding regularizers for these two models. First, similarly to NDI

To minimize boundary-condition-related artifacts (as shown in Figure 3 for the proximal Graz-model solvers) in the proximal GD-Tik solvers we used a dynamic Tikhonov weight: $$$λ_n=2^{-n/5}$$$ where $$$n$$$ is the current iteration number. A 3-voxel ROI erosion was performed in all cases.

To assess the accuracy of these algorithms, we used a synthetic realistic total field simulation based on the 2019 QSM Reconstruction Challenge (RC2) phantom

1. Shmueli K. Chapter 31 - Quantitative Susceptibility Mapping. Advances in Magnetic Resonance Technology and Applications, Academic Press, 2020(1):819-838. doi:10.1016/B978-0-12-817057-1.00033-0

2. Schweser, F., Robinson, S. D., de Rochefort, L., Li, W., and Bredies, K. (2017) An illustrated comparison of processing methods for phase MRI and QSM: removal of background field contributions from sources outside the region of interest. NMR Biomed., 30: e3604. doi: 10.1002/nbm.3604.

3. Langkammer C, Bredies K, Poser B, Barth M, Reishofer G, Fan AP, Bilgic B, Fazekas F, Mainero C, Ropele S. Fast quantitative susceptibility mapping using 3D EPI and total generalized variation. Neuroimage 2015;111:622-630.

4. Chatnuntawech I, McDaniel P, Cauley SF, Gagoski BA, Langkammer C, Martin A, Grant PE, Wald LL, Setsompop K, Adalsteinsson E, Bilgic B. Single-step quantitative susceptibility mapping with variational penalties. NMR Biomed. 2017;30:e3570. doi:10.1002/nbm.3570

5. Liu T, Zhou D, Spincemaille P, Yi W. Differential approach to quantitative susceptibility mapping without background field removal. In: Proceedings of the 22nd Annual Meeting of ISMRM, Milan, Italy. ; 2014:p597.

6. Liu Z, Kee Y, Zhou D, Wang Y, Spincemaille P. Preconditioned total field inversion (TFI) method for quantitative susceptibility mapping. Magn Reson Med. 2017;78:303-315. doi:10.1002/mrm.26331

7. Sun H, Ma Y, MacDonald ME, Pike GB. Whole head quantitative susceptibility mapping using a least-norm direct dipole inversion method. Neuroimage. 2018;179:166-175. doi: 10.1016/j.neuroimage.2018.06.036.

8. Milovic C, Bilgic B, Zhao B, Langkammer C, Tejos C, Acosta-Cabronero J. Weak-harmonic regularization for quantitative susceptibility mapping. Magn Reson Med. 2019;81:1399-1411.

9. Polak, D, Chatnuntawech, I, Yoon, J, et al. Nonlinear dipole inversion (NDI) enables robust quantitative susceptibility mapping (QSM). NMR Biomed. 2020;e4271. doi:10.1002/nbm.4271

10. Milovic C and Shmueli K. Automatic, Non-Regularized Nonlinear Dipole Inversion for Fast and Robust Quantitative Susceptibility Mapping. 29th International Conference of the ISMRM, 2021:p3982.

11. Milovic C, Bilgic B, Zhao B, Acosta-Cabronero J, Tejos C. Fast nonlinear susceptibility inversion with variational regularization. Magn Reson Med. 2018;80:814-821. doi:10.1002/mrm.27073

12. Marques JP, Meineke J, Milovic C, Bilgic B, Chan K-S, Hedouin R, van der Zwaag W, Langkammer C, and Schweser F. QSM Reconstruction Challenge 2.0: A Realistic in silico Head Phantom for MRI data simulation and evaluation of susceptibility mapping procedures. Magn Reson Med. 2021;86: 526– 542 doi:10.1002/mrm.28716

13. Liu T, Wisnieff C, Lou M, Chen W, Spincemaille P, Wang Y. Nonlinear formulation of the magnetic field to source relationship for robust quantitative susceptibility mapping. Magn Reson Med. 2013;69:467-476. doi:10.1002/mrm.24272

14. Milovic C, Tejos C, and Irarrazaval P. Structural Similarity Index Metric setup for QSM applications (XSIM). 5th International Workshop on MRI Phase Contrast & Quantitative Susceptibility Mapping, Seoul, Korea, 2019.

15. Karsa A, Punwani S, and Shmueli K. The effect of low resolution and coverage on the accuracy of susceptibility mapping. Magn Reson Med. 2019; 81: 1833– 1848

16. Liu T, Khalidov I, de Rochefort L, Spincemaille P, Liu J, Tsiouris AJ, Wang Y. A novel background field removal method for MRI using projection onto dipole fields (PDF). NMR Biomed. 2011;24:1129-36.

Figure
1.
Laplacian-based reconstructions lead to severe susceptibility
underestimation for the GD-Tik solver (top). The spatial frequency
analysis (right column) shows strong attenuation of low and medium
frequencies near the dipole-kernel-related-cone at the magic angle,
and mild attenuation elsewhere. The Graz model promotes results with
further suppression of high frequencies, which prevents noise
amplification. ADMM-TV solvers are more accurate, but also show
attenuated low and medium frequencies near the cone at the magic
angle.

Figure
2.
The proximal GD-Tik (Graz model shown) is effective in recovering low
and medium frequencies. Using a Tikhonov regularization (left column)
prevents the recovery of low and medium frequencies, leading to
underestimation. Non-regularized, coefficients are quickly recovered
after a few iterations far from the magic cone, and slowly recovered
in its proximity (right column). This acts as an implicit
regularization, and iterations should be stopped before noise is
amplified.

Figure
3.
Proximal solvers are more prone to amplify artifacts originating at
the boundaries (e.g. see arrows), due to noise or other ill-posed
boundary conditions. ROI mask erosion can be used to mitigate this
effect, but this should be further investigated in future research.

Figure
4.
Optimal reconstructions (i.e. those optimizing the RMSE) of the total
field map based on the RC2 dataset with additional background fields.
GrazTGV solution with (α1,
α0)=(0.0015, 0.0005) and 10000 iterations is shown as comparison.

Figure
5.
In-vivo reconstructions, using the regularization parameters that
gave optimal reconstructions of the RC2 total field map dataset. For
comparison, NDI (non-regularized and automatically stopped^{10})
and FANSI results are also shown (top), using PDF for background
removal.

DOI: https://doi.org/10.58530/2022/2467