Carlos Milovic^{1,2}, Berkin Bilgic^{3}, Bo Zhao^{3}, Christian Langkammer^{4}, Julio Acosta-Cabronero^{5}, and Cristian Tejos^{1,2}

QSM requires

**PURPOSE**

**METHODS**

Given the nonlinear functional: $$$argmin_{\chi}\frac{1}{2}\parallel w(e^{id\star m\chi}-e^{i\phi})\parallel^2_2$$$ with d the dipole kernel, χ the susceptibility distribution, w a magnitude dependent weight, m an inverted ROI binary mask and ϕ the acquired phase. In the context of the alternating direction method of multipliers (ADMM)^{8,9}, we first perform variable splitting by introducing an auxiliary variable $$$z=DF\chi$$$ to decouple the equation system into different subproblems and iterate between solutions until convergence. The subproblem for χ becomes: $$$argmin_{\chi} \frac{\mu}{2}\parallel F^HDFm\chi-z+s\parallel^2_2$$$, where s is the associated Lagrangian multiplier and μ is a Lagrangian weight. We introduce the variable $$$ z_2 = m\chi $$$ to further decouple the problem leading to: $$$argmin_{\chi,z_2}\frac{\mu}{2}\parallel F^HDFz_2-z+s\parallel^2_2+\frac{\mu_2}{2}\parallel m\chi-z_2+s_2\parallel^2_2$$$, which has closed-form solutions: $$$\chi=m(z_2-s_2)$$$ and $$$Fz_2=\frac{\mu D^HF(z-s)+\mu_2F(m\chi+s_2)}{\mu D^HD+\mu_2}$$$. Finally, the nonlinear subproblem z can be solved iteratively as in Milovic, et al^{10} with an initial estimate of the unwrapped phase consisting of an ROI-specific least-squares fitting procedure that solves for the step function, $$$\kappa d\star m+\phi_0$$$, with κ and ϕ_{0} constants. Thus, the nonlinear PDF operation can be seen as a perturbation from such preconditioned initial state.

The proposed algorithm—empirically set to converge for 1,500 iterations with μ=0.0018—was compared with original PDF and Laplacian boundary value (LBV) methods—both initialized with Laplacian and iterative-Laplacian unwrapped phases11 in the following experiments:

A) COSMOS^{12 }ground-truth (3D spoiled-GRE on a 3T Siemens Tim Trio with 1.06-mm isotropic voxels, 240×196×120, TE/TR=25/35ms, 12 head orientations), from which a noisy field map was forward simulated with a widespread distribution of external sources (Fig. 1).

B) In vivo data (3D spoiled-GRE on a 3T Philips system, 352×352×170, 0.6×0.6×1.0 mm^{3}, TE/TR=[7.2,13.4,19.6,25.8,32]/44ms, α=17°) processed with an in-house temporal unwrapping algorithm and subsequent magnitude-weighted nonlinear multi-echo combination.

QSM was calculated from all pre-filtered field maps using the FANSI toolbox^{10} with μ_{2} set to 1.0.

**RESULTS**

A) Resulting field maps, RMSE scores and difference maps are shown in Fig. 2. nPDF returned the lowest RMSE values and a less structured residual than the original PDF method (Fig. 3). RMSE scores for QSM (not shown) relative to COSMOS were: LVB=70.3%, PDF=61.5% and nPDF=58.6%, respectively. A stability test on the ADMM μ parameter revealed less than 0.5% RMSE variation for a μ range spanning an order of magnitude.

B) Fig. 4 illustrates the resulting local field estimates from the three background filtering algorithms investigated in the present study. Qualitatively, although all methods incompletely removed the background field, nPDF yielded the least prominent residual from dipole fields originated from cranial cavities.

**DISCUSSION & CONCLUSION**

1. Li L, Leigh JS. Quantifying arbitrary magnetic susceptibility distributions with MR. Magn. Reson. Med. 2004; 51: 1077–1082.

2. Li L, Leigh JS. High-precision mapping of the magnetic field utilizing the harmonic function mean value property. J. Magn. Reson. 2001; 148: 442–448.

3. Zhou D, Liu T, Spincemaille P, Wang Y. Background field removal by solving the Laplacian boundary value problem. NMR in Biomedicine. 2014;27(3):312–319.

4. Li W, Avram AV, Wu B, Xiao X, Liu C. Integrated Laplacian-based phase unwrapping and background phase removal for quantitative susceptibility mapping. NMR Biomed. 2014;27(2):219-227.

5. Schweser F, Deistung A, Lehr BW, Reichenbach JR. Quantitative imaging of intrinsic magnetic tissue properties using MRI signal phase: an approach to in vivo brain iron metabolism? Neuroimage 2011;54:2789–2807.

6. Wu B, Li W, Guidon A, Liu C. Whole brain susceptibility mapping using compressed sensing. Magn Reson Med. 2012;67(1):137–147.

7. 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(9):1129–1136.

8. Bilgic B, Chatnuntawech I, Langkammer C, Setsompop K. Sparse Methods for Quantitative Susceptibility Mapping. Wavelets and Sparsity XVI, SPIE 2015. doi: 10.1117/12.2188535

9. Zhao B, Setsompop K, Ye H, Cauley S. and Wald LL. Maximum likelihood reconstruction for magnetic resonance fingerprinting, IEEE Trans Med Imaging 2016;35:1812-1823.

10. Milovic C, Bilgic B, Zhao B, Acosta-Cabronero J, Tejos C. A Fast Algorithm for Nonlinear QSM Reconstruction with Variational Penalties. Proceedings of the Fourth International Workshop on MRI Phase Contrast & Quantitative Susceptibility Mapping. 2016;1:132.

11. Schofield MA, Zhu Y. Fast phase unwrapping algorithm for interferometric applications. Opt Lett. 2003;28(14):1194–1196.

12. Liu T, Spincemaille P, De Rochefort L, Kressler B, Wang Y. Calculation of susceptibility through multiple orientation sampling (COSMOS): A method for conditioning the inverse problem from measured magnetic field map to susceptibility source image in MRI. Magnetic Resonance in Medicine. 2009;61(1):196–204.

Fig.
1: COSMOS-brain experiment. a) Ground truth susceptibility, b)
magnitude, c) forward simulated phase, d) Noise-corrupted, wrapped
phase (i.e. input data), e) Laplacian unwrapping, f)
Iterative-Laplacian unwrapping.

Fig. 2: Results for the COSMOS-brain experiment. a) Simulated, noise-corrupted local field. b), d), f), h) and j) are the local fields from all algorithms evaluated (in parenthesis the unwrapping method used as pre-processing). c), e), g), i) and k) are the corresponding difference maps. In the top right corners are displayed the RMSE scores for each reconstruction.

Fig. 3 Contrast-enhanced version of Fig. 2e and 2k, respectively.

In
vivo local field results for all algorithms investigated. Input phase
is the result of temporal unwrapping and subsequent nonlinear
multi-echo combination.