Quantitative Susceptibility Mapping (QSM) estimates the underlying tissue magnetic susceptibility $\chi$ from the gradient echo (GRE) phase signal. Conventional QSM algorithms involve sequential solution of multiple ill-posed post-processing steps. Recently, QSM algorithms that directly relate the unprocessed GRE phase to the unknown $\chi$ have been proposed to improve the susceptibility estimation by preventing error propagation through successive operations in conventional QSM^1-5. In this work, we propose a QSM model that combines the benefits from using a single-step processing, total generalized variation regularization (TGV)⁶, and multiple spherical mean value (SMV) kernels^7-9. A fast solver for the proposed single-step model, which enables simple analytical solutions for all of the optimization steps, is also developed. We demonstrate the performance of the proposed method on three in vivo experiments acquired using the rapid 3D-EPI¹⁰ and Wave-CAIPI trajectories¹¹. Herein, we also demonstrate the multi-echo capability of Wave-CAIPI sequence for the first time, and introduce an automated, phase-sensitive coil sensitivity estimation scheme.

The underlying $\chi$ is estimated from the unprocessed GRE phase $\phi$ by solving the following optimization problem

$$\min_{\chi}{\frac{1}{2}\sum_i\parallel{M_{i}(h_{i}*d*\chi)-M_{i}(h_{i}*\Psi(\phi))}\parallel^2_2+TGV^2_\alpha(\chi)}$$

where $\Psi$ is the Laplacian unwrapping operator, $d$ is the spatial dipole kernel, and $h_{i}$ is the SMV kernel for the reliable phase region indicated by the mask $M_{i}$. Large SMV kernels are used to avoid residual error amplification, whereas smaller SMV kernels are used near the object boundary to prevent eliminating cortical phase information^7-9. ${TGV}^2_\alpha$ is the second-order TGV.

The first term combines phase unwrapping, background phase removal using multiple SMV kernels, and dipole inversion. The second term imposes prior information by promoting a piecewise-smooth solution. By using the convolution theorem and discrete form of $TGV^2_\alpha$, we obtain

$$\min_{\chi,v}{\frac{1}{2}\sum_i\parallel{M_{i}F^{-1}H_{i}DF\chi-M_{i}F^{-1}H_{i}F\Psi(\phi)}\parallel^2_2+\alpha_1\parallel G\chi-v\parallel_1+\alpha_0\parallel\mathcal{E}(v)\parallel_1}$$

where $F$ is the discrete Fourier transform, $D=1/3-k_z^2/(k_x^2+k_y^2+k_z^2)$ is the dipole kernel in k-space, $H_{i}$ is $h_{i}$ in k-space, $G$ is the 3D gradient operator, and $\mathcal{E}$ is a symmetrized derivative. By introducing auxiliary variables $z_0$, $z_1$ , and $z_2$ with consensus constraints ($z_0=\mathcal{E}(v)$, $z_1=G\chi-v$, and $z_{2,i}=H_iDF\chi$), applying alternating direction method of multipliers (ADMM)¹², and exploiting the circulant structure of the matrices, analytical solutions are developed for all of the optimization steps.

Experiments

Using in vivo 3D EPI, Multi-Echo Wave-CAIPI at 3T, and high-resolution Wave-CAIPI data at 7T, we compared the performances of four QSM methods:

(i) SMV background filtering, followed by the dipole inversion with $l_2$-quadratic smoothing regularization defined as $\parallel G\chi\parallel_2^2$, (V-SHARP L2),

(ii) SMV background filtering, followed by TGV-regularized dipole inversion (V-SHARP TGV),

(iii) Single-step QSM with $l_2$-quadratic smoothing regularization (Single-Step L2), and

(iv) Single-step QSM with TGV regularization (Proposed: Single-Step TGV).

This work also introduces the multi-echo extension of Wave-CAIPI that fully utilizes the repetition time to sample multiple echoes, which are combined to improve both magnitude and phase signal-to-noise ratios (Figure 1). In addition, we develop a novel fully-automated coil sensitivity estimation based on ESPIRiT¹³ and virtual body coil concept¹⁴ (described in Figure 2).

For all the reconstruction methods, the GRE phase was unwrapped using Laplacian-based phase unwrapping^15-16. For V-SHARP L2 and Single-Step L2, the MATLAB pcg function was used as a solver with the residual tolerance of 0.1%, and 14 SMV kernels were used with their radii ranging from 1 to 14 voxels. For V-SHARP TGV and Single-Step TGV, ADMM was used with the stopping criterion that iterations were terminated when the change in the solution between consecutive iterations was less than 1%, and 5 SMV kernels were used with their radii ranging from 1 to 5 voxels.Figure 3 through 5 display the reconstructed susceptibility maps from the total field map using the different QSM algorithms. The proposed method had lowest level of dipole artifacts and background phase contamination while preserving the edges.The proposed method better mitigated the dipole artifacts and improved homogeneity of image contrast compared to other QSM methods evaluated thanks to its three main components: (i) the single-step processing that prevents error propagation through successive operations, (ii) multiple SMV kernels that permit high fidelity background removal while retaining the cortical phase information, and (iii) TGV regularization that promotes a piecewise-smooth solution without staircasing artifacts. The reconstruction time was minimized by using ADMM with variable splitting in both the data and regularization terms¹⁷ and exploiting the structure of the matrices. We also introduced an automated sensitivity estimation scheme that relies on a rapid calibration acquisition for phase-sensitive coil combination. This permitted high quality parallel imaging with Wave-CAIPI, which was extended to acquire multiple echoes for improved SNR. Combination of Multi-Echo Wave-CAIPI and single-step QSM with TGV regularization thus enabled high-quality susceptibility mapping with efficient acquisition and reconstruction.