Multi-Echo Quantitative Susceptibility Mapping for Strategically Acquired Gradient Echo (STAGE) Imaging

Purpose: To develop a method to reconstruct quantitative susceptibility mapping (QSM) from multi-echo, multi-flip angle data collected using strategically acquired gradient echo (STAGE) imaging.

Methods: The proposed QSM reconstruction algorithm, referred to as “structurally constrained Susceptibility Weighted Imaging and Mapping” scSWIM, performs an ℓ₁ and ℓ₂ regularization-based reconstruction in a single step. The unique contrast of the T1 weighted enhanced (T1WE) image derived from STAGE imaging was used to extract reliable geometry constraints to protect the basal ganglia from over-smoothing. The multi-echo multi-flip angle data were used for improving the contrast-to-noise ratio in QSM through a weighted averaging scheme. The measured susceptibility values from scSWIM for both simulated and in vivo data were compared to the: original susceptibility model (for simulated data only), the multi orientation COSMOS (for in vivo data only), truncated k-space division (TKD), iterative susceptibility weighted imaging and mapping (iSWIM), and morphology enabled dipole inversion (MEDI) algorithms. Goodness of fit was quantified by measuring the root mean squared error (RMSE) and structural similarity index (SSIM). Additionally, scSWIM was assessed in ten healthy subjects.

Results: The unique contrast and tissue boundaries from T1WE and iSWIM enable the accurate definition of edges of high susceptibility regions. For the simulated brain model without the addition of microbleeds and calcium, the RMSE was best at 5.21ppb for scSWIM and 8.74ppb for MEDI thanks to the reduced streaking artifacts. However, by adding the microbleeds and calcium, MEDI’s performance dropped to 47.53ppb while scSWIM performance remained the same. The SSIM was highest for scSWIM (0.90) and then MEDI (0.80). The deviation from the expected susceptibility in deep gray matter structures for simulated data relative to the model (and for the in vivo data relative to COSMOS) as measured by the slope was lowest for scSWIM + 1%(−1%); MEDI + 2%(−11%) and then iSWIM −5%(−10%). Finally, scSWIM measurements in the basal ganglia of healthy subjects were in agreement with literature.

Conclusion: This study shows that using a data fidelity term and structural constraints results in reduced noise and streaking artifacts while preserving structural details. Furthermore, the use of STAGE imaging with multi-echo and multi-flip data helps to improve the signal-to-noise ratio in QSM data and yields less artifacts.

Introduction

Magnetic resonance imaging (MRI) offers many different contrast mechanisms. Today, it is possible to obtain magnetic susceptibility maps, $\mathrm{\chi }\left(\overrightarrow{r}\right)$, of the human brain (and other parts of the body) that show the underlying tissue susceptibility distribution. Quantitative susceptibility mapping (QSM) data are reconstructed from phase information, which represents the magnetic field variations caused by the magnetization of an object in the presence of an external magnetic field (Haacke et al., 2015). The resulting susceptibility maps can be used to assess bleeding (Bilgic et al., 2012), calcium deposits (Deistung et al., 2013; Chen et al., 2014) and oxygen saturation (Haacke et al., 2010). The knowledge of the susceptibility source and the quantity of either iron or calcium can help improve the diagnosis of neurodegenerative diseases such as multiple sclerosis, Parkinson’s disease, stroke, Sturge-Weber syndrome and traumatic brain injury (Haacke et al., 2015) to name a few.

Extracting the susceptibility, χ, from Gradient Recalled Echo (GRE) phase data is an ill-posed problem because the dipole kernel has zeroes along a conical surface and, therefore, under-samples k-space (Haacke et al., 2015). Many studies have attempted to solve this problem. A fast and direct method to reconstruct χ is the Thresholded K-space Division (TKD) approach (Wharton et al., 2010) that uses a threshold to ignore the smaller values near the zeroes in the inversion process. However, the TKD reconstructed susceptibility map suffers from streaking artifacts and underestimates χ. An alternative approach referred to as iterative Susceptibility Weighted Imaging and mapping (iSWIM) has been used to fill in the missing parts of k-space to overcome these artifacts (Tang et al., 2013). This was accomplished by constraining the susceptibility values in regions with high susceptibility. However, the final images are still noisy in regions of uniform susceptibility. A better approach, in theory, but one that requires multiple scans, is the Calculation Of Susceptibility through Multiple Orientation Sampling (COSMOS) (Liu et al., 2009). This method utilizes the phase images from multiple orientations to stabilize the inversion process and remove the singularities by weighted linear least squares. This method is usually used as a gold standard in the evaluation of any QSM reconstruction method.

A number of other approaches use regularization techniques with different a priori information to reconstruct the susceptibility. Although, these methods are computationally more expensive than TKD approaches, the reconstruction times are still reasonable, and they are designed to smooth over regions that have homogeneous susceptibilities. For example, morphology enabled dipole inversion (MEDI) exploits the structural consistency between χ and the magnitude image in the form of an ℓ₁-norm (Liu et al., 2012). However, this constraint can cause errors in regions where there are inconsistencies between the magnitude images and the susceptibility maps. Homogeneity Enabled Incremental Dipole Inversion (HEIDI) (Schweser et al., 2012) is another method that uses structural information from both magnitude and phase images to correct this issue. An alternative approach, structural feature based collaborative reconstruction (SFCR) (Bao et al., 2016), argues that the edge information from either magnitude or phase images does not reflect all the structural features in χ and the reconstructed image suffers from over-smoothed edges. The key steps in this approach are to include a structural feature-based ℓ₁-norm constraint and a voxel fidelity-based ℓ₂-norm constraint. This allows both edges and small objects to be recovered while still minimizing artifacts. Furthermore, most of these methods find the total field through a linear fitting of multi-echo phase data. However, the inclusion of long echo times can lead to blooming artifact, an increase in signal loss at the edges of the object and, potentially, an underestimation of χ.

Strategically acquired gradient echo imaging (STAGE) is a rapid multi-contrast multi-parametric imaging approach that employs two fully flow compensated double-echo GRE scans using low and high flip angles (FAs) relative to the Ernst angle of white matter. It provides not only a variety of qualitative images such as the T1weighted enhanced (T1WE) image, but also provides multiple quantitative information such as ${\text{R}}_2^{*}$, T1, and susceptibility maps (Chen et al., 2018b; Wang et al., 2018; Haacke et al., 2020). The T1WE image is generated from the combination of two GRE scans with low and high FAs (Chen et al., 2018b) where the radiofrequency (RF) transmit field variation is corrected (Wang et al., 2018). When compared with conventional T1W or T2^∗W images, the T1WE images derived from STAGE have improved contrast between cortical gray matter and white matter, and between deep gray matter and white matter (Chen et al., 2018b). The improved contrast in the T1WE image benefits structural segmentation. STAGE has also become more broadly tested for a number of neurodegenerative diseases (Haacke et al., 2020). Therefore, in this study, we propose a “structurally constrained Susceptibility Weighted Imaging and Mapping” (scSWIM) method that reconstructs the susceptibility using multiple echo, multiple flip angle STAGE data. Similar to SFCR, scSWIM utilizes the structural information from both magnitude data and the susceptibility maps but in a single step. The scSWIM approach specifically uses the enhanced contrast available in STAGE imaging to define prior information about the edges of the white matter and gray matter. In this paper, we introduce scSWIM, evaluate it on simulated data and test it on in vivo brain data.

Materials and Methods

Calculating the Susceptibility From an L1 and L2 Norm Cost Function

Based on Maxwell’s equations, the relationship between the phase image φ(r) (obtained from a 3D GRE imaging approach) and susceptibility χ(r) in ppm (parts per million) can be written as (Haacke and Reichenbach, 2011):

$\mathrm{\phi }\left(r\right)=\mathrm{\gamma }{B}_0\mathit{TE}{F}^{-1}\left\{D\left(k\right)F\left\{\mathrm{\chi }\left(r\right)\right\}\right\},(1)$

where r, B_o and TE are the voxel position vector in image domain, the main magnetic field strength (in T) and the echo time, respectively; γ = 2.675×10⁸rad/s/T is the gyromagnetic ratio; F and F⁻¹ denote the Fourier and inverse Fourier transform operators, respectively; and D(k) is the Fourier transform of the unit dipole kernel at the position k = [k_x, k_y, k_z] in k-space and is defined as:

$D\left(k\right)=\frac{k_x^2+k_y^2-2k_z^2}{3\left(k_x^2+k_y^2+k_z^2\right)}=\frac{1}{3}-\frac{k_z^2}{{\left|k\right|}^2}.(2)$

The objective function of scSWIM is similar to the S-step of SFCR (Bao et al., 2016) with changes in constraint definitions and is given as:

$\begin{array}{c}f\left(\mathrm{\chi }\left(r\right)\right)=\frac{1}{2}{||W\left(r\right)\left(F^{-1}D\left(k\right)F\mathrm{\chi }\left(r\right)-\mathrm{\delta }B\left(r\right)\right)||}_2^2\hfill \\ +{\mathrm{\lambda }}_1{||P(r)G\mathrm{\chi }(r)||}_1+\frac{{\mathrm{\lambda }}_2}{2}{||R(r)\mathrm{\chi }(r)||}_2^2,\hfill \end{array}(3)$

and the final solution for the susceptibility is given by:

${\mathrm{\chi }}_{\mathit{scSWIM}}\left(r\right)=\underset{\mathrm{\chi }(r)}{argmin}f\left(\mathrm{\chi }\left(r\right)\right),(4)$

where δB(r) = φ(r)/(γB₀TE) and W in the data fidelity term is a weighting matrix proportional to the image magnitude that defines the reliability of the magnetic field shift for each voxel and G denotes the gradient operator.

In the S-step of the SFCR method, the edge matrix, P, is a binary mask that is derived from the initial susceptibility, $\widehat{\mathrm{\chi }}$ (where for convenience we have now dropped the dependence on r). This initial $\widehat{\mathrm{\chi }}$ (which is reconstructed from the first regularized minimization step of the SFCR, called the M-step) is based on an objective function that is similar to Eq. [3] but its constraints are based on the magnitude image. Also, R in the S-step of the SFCR method is a fidelity mask where voxels with high signal-to-noise ratio (SNR) are mapped to zero, low SNR to one and voxels corresponding to susceptibility artifact to two. However, the choice of R, P and the starting input are different for scSWIM as described below.

In scSWIM, we replaced the SFCR first regularized minimization (M-step) with iSWIM (Tang et al., 2013) since it has no smoothing, provides an initial susceptibility map with sharp vessels and the reconstruction times are short. Then, in the l₁ regularization term of Eq. [3], we used the edge matrix, P, which is the binary mask that is derived from the product of the thresholded gradients of the T1WE image, P_T1WE, and the initial susceptibility map,$P_{\widehat{\mathrm{\chi }}}$:

$P_{T1\mathit{WE},i}=\{\begin{array}{c}\hfill 0,\left|G_i\mathrm{\rho }\right|\ge {\mathrm{\mu }}_1\hfill \\ \hfill 1,\left|G_i\mathrm{\rho }\right|<{\mathrm{\mu }}_1\hfill \end{array}\mathit{\hspace{1em}}\text{and}\mathit{\hspace{1em}}{P}_{\widehat{\mathrm{\chi }},i}=\{\begin{array}{c}\hfill 0,\left|G_i\widehat{\mathrm{\chi }}\right|\ge {\mathrm{\mu }}_2\hfill \\ \hfill 1,\left|G_i\widehat{\mathrm{\chi }}\right|<{\mathrm{\mu }}_2\hfill \end{array},(5)$

where G_i denotes the gradient operator, i is an indicator to the x, y or z directions, and ρ denotes the T1WE image. Both μ₁ and μ₂ are threshold values chosen to be 2.5 times the noise level of the derivatives of ρ and $\widehat{\mathrm{\chi }}$, respectively, in order to maintain the edges of the gray/white matter, veins and other structures in the brain. Essentially, P_T1WE excludes the edges of the white matter and gray matter and $P_{\widehat{\mathrm{\chi }}}$ excludes the edges of the vessels and basal ganglia including the globus pallidus (GP), caudate nucleus (CN), putamen (PT), thalamus (THA), substantia nigra (SN), and red nucleus (RN) and P = P_T1WE × $P_{\widehat{\mathrm{\chi }}}$.

In the l₂ regularization term, we have used a structural matrix R to protect voxels in the regions of high susceptibility, such as veins and basal ganglia structures, from being over-smoothed while still smoothing other regions. The matrix R is generated from the normalized T1WE image excluding the regions detected in the R_DGM (where DGM stands for “deep gray matter”) and $R_{\widehat{\mathrm{\chi }}}$ masks defined next. The R_DGM mask is calculated using an atlas-based segmentation method developed in-house (Wang et al., 2019). This method segments the deep gray matter structures from the high flip angle magnitude image (T1W), initial susceptibility map $\widehat{\mathrm{\chi }}$ and STAGE T1 weighted data and T1 maps. The $R_{\widehat{\mathrm{\chi }}}$ mask is generated from the method used in Tang et al. (2013) by applying a threshold to the homodyne filtered $\widehat{\mathrm{\chi }}$ map. Finally, the constants λ₁ and λ₂ are found using the L-curve approach (Hansen and O’Leary, 1992).

The single-echo scSWIM approach just described was then adopted to handle the multiple echo, multiple flip angle STAGE data. For this purpose, iSWIM was used as the initial input into scSWIM for the low flip angle, short echo data (FA_LTE₁). Then for the other three echoes (FA_HTE₁, FA_LTE₂, and FA_HTE₂), the reconstructed scSWIM from the previous echo was used as the initial guess for processing the scSWIM of the next echo. Finally, an averaged scSWIM was generated by an $R_2^{*}$-based weighted average of the individual echo scSWIMs (χ_i):

$\mathrm{\chi }=\sum _i{w}_i^2{\mathrm{\chi }}_i/\sum _i{w}_i^2,(6)$

where $w_i={\mathit{TE}}_i{e}^{-{\mathit{TE}}_i{R}_2^{*}}$ and $R_2^{*}$ is from the STAGE data and is created from averaging the $R_2^{*}$ maps from each of the flip angle images (Chen et al., 2018b; Wang et al., 2018):

$R_2^{*}=\frac{1}{{\mathit{TE}}_1-{\mathit{TE}}_2}\ln \left({\mathrm{\rho }}_2/{\mathrm{\rho }}_1\right),(7)$

where ρ₁ and ρ₂ are the magnitudes of the first (TE₁) and second (TE₂) echoes, respectively.

This multi-echo approach has three advantages: first, each echo can be reviewed; second, the weighted scSWIM will have a better SNR; and third, loss of tissues associated with the use of a phase quality control map (especially at longer echoes) will be, to a large degree, replaced with the shorter echo scSWIM value. This weighting automatically ensures that wherever there is a measured susceptibility from one echo it will contribute to the final QSM result (while echoes with zeroes will not make a contribution). Figure 1 shows the block diagram of the proposed multi-echo, multi-flip angle scSWIM processing steps for STAGE.

FIGURE 1

Figure 1. Block diagram of multi-echo, multi-flip angle scSWIM for STAGE imaging. Here, φ, $\widehat{\mathrm{\chi }}$ denote the phase and initial estimate of the susceptibility map from the multi-echo R2^∗ weighted iSWIM, respectively. FA_L and FA_H denote the double-echo low and high flip angles scans of STAGE imaging, respectively.

Simulated Data

The 3D isotropic susceptibility model developed in this laboratory (Buch et al., 2012) was used to test the algorithm. This model includes the general structures of the human brain such as gray matter (GM), white matter (WM), basal ganglia and midbrain structures [PT, GP, CN, THA, RN, SN, and crus cerebri (CC)], cerebrospinal fluid (CSF) and the major veins. The susceptibility values for these structures are summarized in the first row of Table 1.

TABLE 1

Table 1. Susceptibility, T₁ relaxation time, and relative proton density (ρ₀) values for different structures in the simulated brain model.

To test the performance of the reconstruction in the presence of cerebral microbleeds (CMB) or calcium deposits (CaD), two spherical objects with susceptibility values (radius) of 1000 ppb (5 mm) and 3000 ppb (3 mm), respectively, were added to the frontal white matter and two spherical objects with susceptibility values of −1000 ppb (5 mm) and −3000 ppb (3 mm) were added to the posterior white matter. Additionally, one spherical object with a radius of 3 mm with a susceptibility of −3000 ppb was added to the model to mimic the pineal gland (PG). The values for CMBs were taken from our experience in the field of traumatic brain injury and stroke where we usually see CMBs with susceptibilities as large as 1000 ppb but on occasion higher values up to 2000 ppb and 3000 ppb have been seen so we used both 1000 ppb and 3000 ppb to test the metal of the method. For the CaD, the values are around −3000 ppb but can range lower and slightly higher than this as the calcium is highly diamagnetic (Buch et al., 2015).

This final susceptibility model, χ_ideal, was used to generate the magnitude and phase images using the STAGE imaging parameters: FA = 6^o/24^o, TE₁ = 7.5/8.75ms, TE₂ = 17.5/18.75ms, and TR = 25ms. The phase images were simulated from the forward model in Equations [1] and [2] at B₀ = 3T. To create the magnitude images, first an R2^∗ map was generated from χ_ideal using the relationship R2^∗ = 20/s + 0.125χ (Ghassaban et al., 2019b) assuming R2^∗ = 40/s for CMB, PG, and CaD objects. Then, the magnitude image was calculated using the Ernst equation (Brown et al., 2014). The proton density and T₁ relaxation times for different brain structures are summarized in Table 1 while they are assumed to be zero for CMB, PG, and CaD objects. These values were adopted from the literature (Lee et al., 2006; Brown et al., 2014) or manually measured from the in vivo STAGE PD-map and T1-map.

Finally, Gaussian noise was added to the complex signal to produce an SNR of 10:1. The reconstructed susceptibility map using the proposed method was compared with the TKD, iSWIM, and MEDI methods. The original simulated susceptibility model (χ_ideal) was used as the gold standard to measure the performance of each method using RMSE and SSIM as measures of goodness of fit (Wang et al., 2004) where SSIM = 1 corresponds to the perfect structural similarity while SSIM = 0 indicates no similarity between the two images.

In vivo Data

The proposed scSWIM method was tested on two sets of in vivo datasets that are discussed below. All subjects involved in this study signed a consent form to be part of this IRB approved research.

Single Test Case Including COSMOS

The in vivo MRI data for this single test case was acquired from a 29-year old male volunteer on a 3T Siemens scanner (Siemens Healthcare, Erlangen, Germany) at Wayne State University. The imaging parameters were: 6° and 24° for the low and high flip angle scans with TR = 25 ms, TE₁ = 6.5/7.5 ms, TE₂ = 17.5/18.5 ms, bandwidth: 277 Hz/pixel, and GRAPPA = 2. The matrix size, voxel resolution, and FOV were 384 × 288 × 104, 0.67 × 0.67 × 1.33 mm³, and 256 × 192 × 139mm³, respectively. The total scan time for the high-resolution STAGE is about 10 min. For the purpose of generating COSMOS, two additional orientations with the same imaging parameters were collected for this subject. The reconstructed susceptibility map using the proposed scSWIM method was compared with those from the TKD, iSWIM and MEDI methods and compared to COSMOS as the reference image.

Evaluation on a Set of Healthy Human Subjects

Additionally, we tested scSWIM for ten healthy subjects acquired using a Siemens 3T Prisma scanner with lower resolution compared to the above-mentioned in vivo case. The imaging parameters were the same for the sample used above in the simulated data except the matrix size, voxel resolution, and FOV were 384 × 144 × 64, 0.67 × 1.33 × 2 mm³ (interpolated to 0.67 × 0.67 × 2 mm³) and 256 × 192 × 128 mm³, respectively, TE₁ = 7.5/8.5 ms, and a bandwidth of 240 Hz/pixel. The total scan time for this resolution is about 5 min.

Data Pre-processing

The entire processing pipeline was implemented in MATLAB (The Mathworks, Inc., Natick, MA, United States) on a workstation with Windows 10, Intel CPU i7-3770 with 4 cores and 16GB RAM. The phase image was first unwrapped using the bootstrapping (Chen et al., 2018a) and quality guided 3D phase unwrapping (Abdul-Rahman et al., 2007) methods in the simulated and in vivo data, respectively.

For the in vivo data, the induced background field from the air/tissue interfaces was removed from the unwrapped phase using the Sophisticated Harmonic Artifact Reduction for Phase data (SHARP) algorithm (Schweser et al., 2011) with a kernel size of 6 pixels. Next, the brain mask for the in vivo data was extracted from the magnitude images using BET (Brain Extraction Tool) (Smith, 2002). Finally, the resulting phase was zero-padded symmetrically in the spatial domain to a matrix size of 256 × 256 × 256 or 512 × 512 × 512 for simulated and in vivo datasets, respectively.

Susceptibility Map Reconstruction

In Eq. [3], the parameters λ₁ and λ₂ were determined by plotting the measured residual errors of the data fidelity and the two regularization terms for each of the individual STAGE scans using the L-curve method (Hansen and O’Leary, 1992). In theory, λ₁ controls the spatial smoothness and λ₂ helps to preserve the high susceptibility regions and small objects such as vessels from being over-smoothed. As mentioned in section “Calculating the Susceptibility From an L1 and L2 Norm Cost Function,” an atlas-based segmentation method developed in-house (Wang et al., 2019) was used to generate the R_DGM mask. This method provided the labeled mask segmenting the right and left subcortical deep gray matter structures from the T1W, STAGE T1WE, T1 map, and $\widehat{\mathrm{\chi }}$. This labeled mask was carefully reviewed and if needed fine-tuned manually (this was done on 3 cases for the GP and SN structures which sometimes were smaller than what would have been drawn manually). If these regions had not been corrected, the algorithm would have smoothed that part of the GP not protected. Finally, the R_DGM mask is generated from binarizing the labeled mask.

Several algorithms were chosen to compare with scSWIM, including TKD, iSWIM, and MEDI. In generating the MEDI results, a regularization parameter of 250 (350) was used for the simulated (in vivo) data. For TKD processing, a threshold of 0.1 was used and iSWIM was performed with 4 iterations. All of these parameters were adjusted to give the lowest RMSE. Additionally, COSMOS was used as the gold standard for the in vivo data. Multi orientation images for the COSMOS data were co-registered using ANTs (Avants et al., 2009, 2012). In the TKD, iSWIM, and scSWIM methods, the final multiple echo, multiple flip angle QSM data were generated using a multi-echo R2^∗-based weighted averaging of the individual QSM images from each echo and each flip angle data. In MEDI, the final QSM was generated by averaging the reconstructed QSM images from the fitted phases in each of the multi-echo low and high flip angle scans.

Quantitative Analysis for Susceptibility Map

For the quantitative analysis of the data, the susceptibility mean and standard deviation were found from the entire 3D structure of interest. In the simulated model, all the structures of interest were measured automatically (since we know the location of each structure). For the in vivo data, the susceptibility of the midbrain structures were also automatically measured since they have been determined in creating the R_DGM masks for the boundaries of these structures as described earlier. On the other hand, the susceptibility of the CSF, WM, and major veins [SSV and internal cerebral vein (ICV)] were measured manually by tracing the ROIs on the QSM data using SPIN (SpinTech, Inc., Bingham Farms, MI. United States). The manual tracing was performed in the axial view for CSF and WM, but veins were traced in the sagittal view for easier localization. A linear regression model was used to compare the measured susceptibility values from each reconstruction method with those from the susceptibility model and COSMOS to assess the accuracy of midbrain structures in the simulated and in vivo data, respectively.

Results

Simulated Data

By comparing the P and R masks for the simulated data (discussed in Section “Simulated Data”) and also the first and second regularization terms, and for the purpose of bringing the two terms to the same order, we set λ₁ = 0.005λ₂. This is further reviewed in the Discussion section. Based on this assumption and simulations in the human brain model, λ₂ = {6.81,1.47,3.16,1.00}×10⁻³ provided the best results in terms of residual errors for the four different scans (FA_LTE₁, FA_HTE₁, FA_LTE₂, and FA_HTE₂), respectively (see Figures 2A–C for FA_HTE₁). A comparison of scSWIM with TKD, iSWIM, and MEDI along with their absolute errors and structural similarity maps relative to the simulation model are shown in Figure 3. In the simulated data (Figures 3A–C), we have used the exact known edge and structural matrices from χ_ideal to create P_ideal (Figure 3D) and R_ideal (Figure 3E). The TKD results (Figures 3F–J) show severe streaking artifacts while the iSWIM results have much less streaking (Figures 3K–O). MEDI does an excellent job (Figures 3P–T) as does scSWIM (Figures 3U–Y) in reproducing the model with minimal artifacts and noise. In both these last two reconstructions, the streaking artifact is highly reduced compared to both TKD and iSWIM and the images look much better in terms of SNR. However, MEDI does not resolve the streaking artifact from the CMB, pineal gland, or calcified objects with higher susceptibility values.

FIGURE 2

Figure 2. Determination of the scSWIM regularization parameter λ₂ in the simulated (A–C) and in vivo (D–F) data for the higher flip angle, short echo (FA_HTE₁) scan using L-curve method. The curves in the first column show the log-log L-curve. The curvature and RMSE/residual error plot vs λ₂ values are displayed in the third column. The optimal values (shown by the red circle) for the scSWIM at FA_HTE₁ scan were determined to be λ₂ = 1.47×10⁻³ andλ₂ = 1.47×10⁻⁴ for the simulated and in vivo data, respectively, where λ₁ was set equal to 0.005λ₂. This process is repeated for the other scans (FA_LTE₁, FA_LTE₂, and FA_HTE₂) and the optimal parameters were selected.

FIGURE 3

Figure 3. Depiction of multi-echo, multi-flip angle QSM images using different methods for the simulated data. This figure shows the orthogonal views of the susceptibility model (A–C), and reconstructed QSM images from TKD (F–H), iSWIM (K–M), MEDI (P–R) and scSWIM (U–W) along with the scSWIM constraints P_ideal (D) and R_ideal (E). The cerebral microbleeds (CMB), pineal gland (PG) and calcium deposits (CaD) are labeled on the model (A). Streaking artifacts are indicated by the arrows. The last two columns show the corresponding susceptibility absolute error map (I,N,S,X) and structural similarity map (J,O,T,Y) for the different methods. In this simulated data, scSWIM provides better reconstruction with less artifacts, less error and higher similarity relative to the numerical model. Please note the complements of the P and R mask are shown in this figure (D,E) for better visualization.

In the simulated data with (or without) CMBs, PG and CaDs, the RMSE for TKD, iSWIM, MEDI, and scSWIM were 32.91 (22.09), 24.61 (18.21), 47.53 (8.74) and 5.01 (5.21) ppb, respectively. Also, the SSIM index was measured as 0.52 (0.59), 0.62 (0.63), 0.80 (0.86) and 0.90 (0.91) for TKD, iSWIM, MEDI, and scSWIM, respectively, for these two conditions. Based on these results, scSWIM has the lowest error and the highest similarity to the model compared to the other methods. The measured susceptibility values in different structures are summarized in Table 2 showing that the measured susceptibilities in the midbrain structures for both MEDI and scSWIM are closer to the expected susceptibilities in the model while scSWIM has smaller standard deviations. The measured susceptibilities of the straight sinus vein, calcium deposition and CMB objects show that scSWIM provides the most accurate results in these structures as well.

TABLE 2

Table 2. Measured susceptibility values (mean ± standard deviation) in ppb for different structures in the reconstructed QSM images using TKD, iSWIM, MEDI, and scSWIM methods for the simulated human dataset along with the reference values.

In vivo Data

Based on the L-curve analysis using the single high resolution human in vivo data (discussed in Section “Single Test Case Including COSMOS”) and by assuming λ₁ = 0.005λ₂ for the purpose of bringing two regularization terms to the same scale, λ₂ = {1, 1.47, 1.00, 1.00}×10⁻⁴ provided the best results in terms of residual errors for FA_LTE₁, FA_HTE₁, FA_LTE₂, and FA_HTE₂, respectively (see Figures 2D–F). The structural terms used in the scSWIM cost function are illustrated in Figure 4. Specifically, Figures 4A–D show the edge and structural matrices P (includes P_x, P_y, and P_z) and R. The binary matrix P, excludes the extracted edges from the enhanced T1-weighted and initial susceptibility while the binary mask R excludes the deep gray matter structures, vessels and other high susceptibility regions (the complement of P and R masks are shown in the figure for better visualization). Figures 4E–H show the conventional T1-weighted (Figure 4E) and T1WE (Figure 4H) from STAGE and their corresponding extracted edges (final P representation of extracted edges in three directions). It can be seen visually that the contrast between gray matter and white matter of the T1WE is higher than the conventional T1W image and its corresponding edge matrix, P_T1WE, provides more information about the edge.

FIGURE 4

Figure 4. Illustration of scSWIM constraints and comparison of constraints extracted from conventional T1W and STAGE T1WE for the single high-resolution in vivo data. The first row shows the scSWIM structural constraints for the single high-resolution in vivo data: edge matrix, P, in the x, y, and z directions (A–C), and the structural matrix, R (D). The second row shows the advantage of extracting the constraints from STAGE versus conventional GRE data: conventional T1W (G), STAGE T1WE (H), and the extracted edges (product of three directions) from conventional T1W (F) and STAGE T1WE (G). As seen, (G) provides more information about the white and gray matter edges (white arrow) and is less noisy than (F). Please note the complement of the P and R mask is shown in this figure for better visualization.

Figure 5 shows three orthogonal views of the reconstructed multi-echo, multi-flip angle susceptibility images for this high-resolution human data set using the TKD (Figures 5A–C), iSWIM (Figures 5D–F), MEDI (Figures 5G–I), scSWIM (Figures 5J–L), and COSMOS (Figures 5M–O) methods. It can be seen in these images that scSWIM has less noise while the sharpness of the vessels and other brain structures are well-preserved. MEDI also provides a smooth reconstruction but in the regions that are close to the veins there are still some remaining artifacts. The measured susceptibility values in different structures are summarized in Table 3.

FIGURE 5

Figure 5. Depiction of multi-echo, multi-flip angle QSM images using different methods for the single high-resolution in vivo data. This figure shows three orthogonal views of the reconstructed multi-echo, multi-flip angle susceptibility maps from TKD (A–C), iSWIM (D–F), MEDI (G–I), scSWIM (J–L), and COSMOS (M–O) for the single high-resolution in vivo data. All of the images are displayed with the same window/level settings. White arrows show streaking artifacts. The SNR and image quality are best in the scSWIM images while the sharpness of the vessels and other brain structures are preserved.

TABLE 3

Table 3. Measured susceptibility values (mean ± standard deviation) in ppb for different structures in the reconstructed QSM images using TKD, iSWIM, MEDI, scSWIM, and COSMOS methods for the multi-echo, multi-flip angle in the single high-resolution in vivo data.

The structural terms used in the scSWIM cost function for two selected healthy subjects with lower resolution data (discussed in Section “Evaluation on a Set of Healthy Human Subjects”) are illustrated in Figure 6. Here it can be seen that the edges are still well preserved with this in vivo STAGE approach. Figure 7 shows the reconstructed multi-echo, multi-flip angle susceptibility images using TKD (Figures 7A,E), iSWIM (Figures 7B,F), MEDI (Figures 7C,G) and scSWIM (Figures 7D,H) methods for two examples of this data. There are artifacts around the basal ganglia structures and larger veins in the TKD, iSWIM and MEDI (shown with white arrows). Furthermore, in the second slice (Figures 7E–H), the PG looks dilated in MEDI compared to the other methods (marked by a red arrow). Table 4 summarizes the averaged measured susceptibility values (mean ± standard deviation) in the reconstructed QSM images from the four different methods for the ten healthy subjects.

FIGURE 6

Figure 6. Illustration of scSWIM structural constraints for the two healthy subjects from the low-resolution dataset. The scSWIM structural constraints, edge matrix, P, in the x, y, and z, and structural matrix, R, are shown for the low- resolution in vivo data from a 62-year old healthy subject (A–D) and a 54-year old healthy subject (E–H). Please note the complement of the P and R mask is shown in this figure for better visualization.

FIGURE 7

Figure 7. Depiction of multi-echo, multi-flip angle QSM images using different methods for the two healthy subjects from the low-resolution dataset. Multi-echo, multi-flip angle susceptibility maps from TKD (A,E), iSWIM (B,F), MEDI (λ = 350) (C,G), and scSWIM (D,H) are shown for the two healthy subjects from Figure 6. The artifacts around the basal ganglia and larger veins in the TKD, iSWIM, and MEDI are shown by the white arrows. In the second row (E–H), the pineal gland looks dilated in MEDI compared to other methods (red arrow).

TABLE 4

Table 4. Averaged susceptibility values (mean ± standard deviation) in ppb for midbrain structures in the reconstructed QSM images using TKD, iSWIM, MEDI, and scSWIM for ten healthy subjects from the low-resolution in vivo dataset from a Siemens 3T PRISMA scanner.

Figure 8 shows the correlation between the zero-referenced estimated susceptibility for deep gray matter structures from different reconstruction methods with the actual susceptibility from the numerical model for the simulated data and reconstructed COSMOS for the in vivo data. The measured CSF susceptibility for each method is used to zero-reference the measurements. Among these methods, scSWIM (in blue color) has the closest values to the reference image in both datasets. The slope of scSWIM is 1.01(0.99) while TKD, iSWIM and MEDI are 0.89(0.78), 0.95(0.90), and 1.02(0.89) for simulated (and in vivo) data, respectively. The correlation coefficients in all methods are close to one and p-values to zero.

FIGURE 8

Figure 8. This figure shows the correlation of the susceptibilities of different basal ganglia structures (bilateral, that is the average of left and right) in the reference image with the ones in the reconstructed images using different methods in the simulated data (A) and in vivo data (B) [TKD (black), iSWIM (green), MEDI (red), and scSWIM (blue)]. All methods correlated well with iron content but scSWIM provided the best result relative to the correct absolute susceptibility. The dashed pink line corresponds to the line of identity between the individual reconstruction method and the reference susceptibility model and COSMOS for simulated and in vivo data, respectively.

The current implementation of scSWIM for a single echo, converges in less than 3 and 5 iterations for the simulated and in vivo data, respectively. Each iteration consists of a minimization process that uses a preconditioned conjugate gradient solver. For the zero-padded in vivo data with a matrix size of 512×512×128, the total processing time for each single-echo scSWIM is currently 2∼5 minutes depending on the number of iterations using a Windows 10, Intel CPU i7-3770 with 4 cores and 16GB RAM.

Discussion

The quantitative and qualitative analysis on both simulated data and human in vivo data showed that the reconstructed TKD suffers from streaking artifacts and underestimates the susceptibility values of deep gray matter and veins. The streaking artifact is reduced in iSWIM by using constraints from high susceptibility structures, but the final image is still noisy in the homogeneous regions. Thanks to the use of an ℓ₁-norm regularization MEDI creates high SNR results. However, some streaking artifacts remain in regions where magnitude data is inconsistent with the susceptibility map. On the other hand, scSWIM uses both ℓ₁ and ℓ₂regularization terms to protect edges and structures while also allowing smoothing to increase SNR in regions without structure and it successfully reduces streaking artifacts leading to less noise and faithful estimates of the susceptibility. Furthermore, scSWIM outperforms other methods in reconstructing the susceptibility map in the presence of CMBs and CaDs with high susceptibilities. In simulated data, both microbleeds with susceptibilities of 1000 ppb and 3000 ppb and calcium objects with susceptibilities of −1000 ppb and −3000 ppb were reconstructed accurately using scSWIM compared to other methods. Also, in scSWIM, the standard deviation of the measured susceptibilities (Table 2) in all structures even in the CMB or CaD with the highest susceptibility values are much lower than other methods showing the strengths of this multi-echo approach. Although, MEDI provides a smooth QSM image under normal circumstances it appears to have trouble in reconstructing the data in the presence of high susceptibilities such as seen with the CMB and CaD in the simulated model and for the pineal gland in the in vivo data (which appeared dilated compared to that in scSWIM). This could be due to the fact that MEDI uses phase fitting across multiple echoes and high susceptibilities can cause both signal loss at the edge of the object and severe aliasing at longer echoes. Furthermore, in the in vivo data, one could observe slight streaking with MEDI around the large veins that could be due to the inconsistency between the magnitude and susceptibility data.

The in vivo results for scSWIM showed average susceptibilities for the ten healthy subjects very close to the reported values in the literature (Ghassaban et al., 2019a). Also, the measured susceptibilities in the reconstructed COSMOS (Table 3) were not as close to scSWIM and MEDI as one would have hoped because it can contained errors due to registration of the different orientation data and noise in the data. The registration error is higher and more noticeable in the regions near the surface of the brain. Luckily most of the regions of interest (the deep gray matter) in this paper are near the core of the brain where the registration error is smaller therefore this central region can still be used as a baseline to compare different methods.

Structural Constraints in scSWIM

The cost function of scSWIM includes two regularization terms. The ℓ₁-norm regularization term is based on a P mask to penalize the noisy non-edge pixels and the ℓ₂-norm regularization term is based on the R mask that prevents smoothing in the excluded high susceptibility regions. If the pre-processing fails to extract the edges of a true structure, then the P mask will penalize and smooth them. On the other hand, if R fails to exclude a high susceptibility structure, the streaking artifacts from this structure will remain and its mean susceptibility will be reduced due to smoothing. This is because the R mask protects the structures of high susceptibility from being over smoothed by the ℓ₁-norm regularization term. The overall performance of the cost function works well when the edges and structures are best defined.

Optimal Parameter Selection for scSWIM

In the regularization-based approaches, there is always a trade-off between obtaining accurate susceptibility values, reducing streaking artifacts, and increasing SNR. In scSWIM, the λ₁ parameter controls the spatial smoothness by applying the sparsity constraint on the gradient of the susceptibility map. The larger the λ₁, the smoother the non-edge regions will be for both the background and basal ganglia (basically increasing the SNR). On the other hand, λ₂ also controls smoothing the background but protects the objects defined by the R mask. Smaller λ₂ reduces the effect of the regularization term and increases the effect of the data fidelity term and the streaking artifact will not be handled as well. On the other hand, larger λ₂ will increase the effect of the regularization term and reduce the effect of the data fidelity term and will result in an over-smoothed image where the background such as WM and GM and smaller objects would be washed out.

Therefore, the challenging part of scSWIM is to find the optimal parameters to keep sharp edges, smooth where appropriate, and satisfy the data fidelity condition. However, finding optimal values for more than one parameter in regularization problems is still a difficult problem. With the admission of sub-optimality, we assumed that the ratio of λ₁ and λ₂ is fixed. For this purpose, we compared the P and R masks and also the first and second regularization terms and observed that λ₁ = 0.005λ₂ will bring the two terms to the same order. The final step was to determine the optimal value for λ₂. This was accomplished using the L-curve approach that plots the residual data fidelity versus the regularization for different regularization parameters and selecting the value that results in the maximum curvature. For multi-echo, multi-flip angle scSWIM, the L-curves were analyzed for each individual scan separately and the optimal λ₂ values selected accordingly.

Multi-Echo, Multi-Flip Angle scSWIM

As mentioned before, STAGE imaging uses double-flip angle, double-echo GRE scans. The multi-echo, multi-flip angle scSWIM, or STAGE scSWIM is generated by an $R_2^{*}$-based weighted averaging of the individual echo scSWIM data sets. Besides having higher SNR in the STAGE scSWIM results, each individual scSWIM dataset can be reviewed separately if desired. It would be of interest to compare the QSM results with those from the R2^∗ maps or even the T1maps given that iron can affect the T1 of tissue. Recently, there has been more interest in multi-contrast quantitative mapping in diseases such as Parkinson’s disease and dementia where a more systemic quantitative approach is being taken with 3D data. Iron has played a key role in these studies not just in the basal ganglia but also in the hippocampus, motor cortex and cortical gray matter in general.

More importantly, the final STAGE scSWIM will keep regions that have been removed by the phase quality control map at longer echo times. An alternate approach would be reconstructing QSM from the linear fit to the phase as done in MEDI. However, regions of high susceptibility phase aliasing can be severe and phase fitting may not be successful. Furthermore, severe loss of signal in and around the object (blooming artifacts) will occur for high susceptibilities that will result in a significantly under-estimated susceptibility. The use of shorter echo times and the weighting factors can favor the short echo data replacing the long echo data when the susceptibilities are very high as in the case of the CMBs and CaD as shown in the results section.

STAGE uses the conventional SWI with two flip angles and is effectively available at any site that can run 3D GRE imaging. It is a 5 min scan (2.5 min for each flip angle) that provides eight qualitative and seven quantitative clinically useful images such as T1maps, spin density maps, QSM, R2^∗, B1 field corrections and etc. Although, the high resolution STAGE scan time may take longer (∼10 min), using a compressed sense factor of 3 to 4 the scan times can be brought back to a time frame of 5 to 7 min. The proposed scSWIM method achieved the best results when processing double-echo, double-flip angle STAGE data by using the derived T1WE images to extract reliable geometry constraints, but it can also be performed on a single-echo T1W SWI dataset.

Conclusion

In this paper, we have proposed a constraint based QSM reconstruction algorithm scSWIM which uses STAGE and iSWIM inputs to reconstruct the susceptibility map from multiple flip angle, multiple echo data. The results show for both simulated and in vivo human brain data that streaking artifacts are suppressed, and SNR is increased. Further, the measured susceptibilities are accurate relative to the brain model used and scSWIM works well even for regions with high susceptibility such as microbleeds and calcifications.

Data Availability Statement

The data analyzed in this study is subject to the following licenses/restrictions: The data is currently not available. The data is from a collaboration with a clinical site. Requests to access these datasets should be directed to Nmrimaging@aol.com.

Ethics Statement

The studies involving human participants were reviewed and approved by Ethics Committee of Wayne State University and Ethics Committee of The First Affiliated Hospital of Zhengzhou University. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

SG and SL designed, executed, and implemented the research project under the supervision of EMH. SB, YC, YW, and MJ helped in designing the experiments. CZ and BW helped with data acquisition. JC was the senior supporter of the data acquisition. SG wrote the first draft of the manuscript. EMH, SL, and SB critically discussed the results and reviewed the manuscript. EMH, TW, YC, and NK provided feedback and helped shape the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This research was supported in part by the NIH/NHLBI grant R44HL145826. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.

Conflict of Interest

SG, MJ, and EMH are employees of the Magnetic Resonance Innovations Inc., Bingham Farms, MI, United States. YW and EMH are employees of The MRI Institute for Biomedical Research, Bingham Farms, MI, United States.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnins.2020.581474/full#supplementary-material

References

Abdul-Rahman, H. S., Gdeisat, M. A., Burton, D. R., Lalor, M. J., Lilley, F., and Moore, C. J. (2007). Fast and robust three-dimensional best path phase unwrapping algorithm. Appl. Opt. 46, 6623–6635. doi: 10.1364/AO.46.006623

PubMed Abstract | CrossRef Full Text | Google Scholar

Avants, B. B., Tustison, N., and Song, G. (2009). Advanced normalization tools (ANTS). Insight J. 2, 1–35.

Google Scholar

Avants, B. B., Tustison, N. J., Song, G., Cook, P. A., Klein, A., and Gee, C. (2012). A reproducible evaluation of ANTs similarity metric performance in brain image registration. 54, 2033–2044. doi: 10.1016/j.neuroimage.2010.09.025.A

CrossRef Full Text | Google Scholar

Bao, L., Li, X., Cai, C., Chen, Z., and Van Zijl, P. C. M. (2016). Quantitative susceptibility mapping using structural feature based collaborative reconstruction (SFCR) in the human brain. IEEE Trans. Med. Imaging 35, 2040–2050. doi: 10.1109/TMI.2016.2544958

PubMed Abstract | CrossRef Full Text | Google Scholar

Bilgic, B., Pfefferbaum, A., Rohl, T., Sullivan, E. V., and Adalsteinsson, E. (2012). NeuroImage MRI estimates of brain iron concentration in normal aging using quantitative susceptibility mapping. Neuroimage 59, 2625–2635. doi: 10.1016/j.neuroimage.2011.08.077

PubMed Abstract | CrossRef Full Text | Google Scholar

Brown, R. W., Cheng, Y. C. N., Haacke, E. M., Thompson, M. R., and Venkatesan, R. (2014). Magnetic Resonance Imaging: Physical Principles and Sequence Design. Hoboken, NY: Wiley-Blackwell.

Google Scholar

Buch, S., Liu, S., Haacke, E. M., and Neelavalli, J. (2012). “Simulated 3D brain model to study the phase behavior of brain structures,” in Proceedings of the 20th Annual Meeting of ISMRM Melbourne, 2332.

Google Scholar

Buch, S., Liu, S., Ye, Y., Cheng, Y. C. N., Neelavalli, J., and Haacke, E. M. (2015). Susceptibility mapping of air, bone, and calcium in the head. Magn. Reson. Med. 73, 2185–2194. doi: 10.1002/mrm.25350

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, W., Zhu, W., Kovanlikaya, I., Kovanlikaya, A., Liu, T., Wang, S., et al. (2014). Intracranial calcifications and hemorrhages: characterization with quantitative susceptibility mapping. Radiology 270, 496–505. doi: 10.1148/radiol.13122640

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., Liu, S., Kang, Y., and Haacke, E. M. (2018a). “A rapid, robust multi-echo phase unwrapping method for quantitative susceptibility mapping (QSM) using strategically acquired gradient echo (STAGE) data acquisition,” in Proceedings of the SPIE 10573, Medical Imaging 2018: Physics of Medical Imaging (Washington, DC: SPIE), doi: 10.1117/12.2292951

CrossRef Full Text | Google Scholar

Chen, Y., Liu, S., Wang, Y., Kang, Y., and Haacke, E. M. (2018b). STrategically acquired gradient echo (STAGE) imaging, part I: creating enhanced T1 contrast and standardized susceptibility weighted imaging and quantitative susceptibility mapping. Magn. Reson. Imaging 46, 130–139. doi: 10.1016/j.mri.2017.10.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Deistung, A., Schweser, F., Wiestler, B., Abello, M., Roethke, M., Sahm, F., et al. (2013). Quantitative susceptibility mapping differentiates between blood depositions and calcifications in patients with glioblastoma. PLoS One 8:e57924. doi: 10.1371/journal.pone.0057924

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghassaban, K., He, N., Sethi, S. K., Huang, P., and Chen, S. (2019a). Regional high iron in the substantia nigra differentiates parkinson’s disease patients from healthy controls. Front. Aging Neurosci. 11:106. doi: 10.3389/fnagi.2019.00106

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghassaban, K., Liu, S., Jiang, C., and Haacke, E. M. (2019b). Quantifying iron content in magnetic resonance imaging. Neuroimage 187, 77–92. doi: 10.1016/j.neuroimage.2018.04.047

PubMed Abstract | CrossRef Full Text | Google Scholar

Haacke, E. M., Chen, Y., Utriainen, D., Wu, B., Wang, Y., Xia, S., et al. (2020). STrategically acquired gradient echo (STAGE) imaging, part III: technical advances and clinical applications of a rapid multi-contrast multi-parametric brain imaging method. Magn. Reson. Imaging 65, 15–26. doi: 10.1016/j.mri.2019.09.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Haacke, E. M., Liu, S., Buch, S., Zheng, W., Wu, D., and Ye, Y. (2015). Quantitative susceptibility mapping: current status and future directions. Magn. Reson. Imaging 33, 1–25. doi: 10.1016/j.mri.2014.09.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Haacke, E. M., and Reichenbach, J. R. (2011). Susceptibility Weighted Imaging in MRI: Basic Concepts and Clinical Applications. Hoboken, NY: Wiley Blackwell.

Google Scholar

Haacke, E. M., Tang, J., Neelavalli, J., and Cheng, Y. C. N. (2010). Susceptibility mapping as a means to visualize veins and quantify oxygen saturation. J. Magn. Reson. Imaging 32, 663–676. doi: 10.1002/jmri.22276.Susceptibility

CrossRef Full Text | Google Scholar

Hansen, P. C., and O’Leary, D. P. (1992). The use of the L-curve in the regularization of discrete Ill-posed problems. SIAM J. Sci. Comput. 14, 1487–1503. doi: 10.1137/0914086

CrossRef Full Text | Google Scholar

Lee, C., Baker, E., and Thomasson, D. (2006). Normal regional T1 and T2 relaxation times of the brain at 3T. Proc. Intl. Soc. Magn. Reson. Med. 14:2006.

Google Scholar

Liu, J., Liu, T., De Rochefort, L., Ledoux, J., Khalidov, I., Chen, W., et al. (2012). Morphology enabled dipole inversion for quantitative susceptibility mapping using structural consistency between the magnitude image and the susceptibility map. Neuroimage 59, 2560–2568. doi: 10.1016/j.neuroimage.2011.08.082

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, T., Spincemaille, P., De Rochefort, L., Kressler, B., and Wang, Y. (2009). 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. Magn. Reson. Med. 61, 196–204. doi: 10.1002/mrm.21828

PubMed Abstract | CrossRef Full Text | Google Scholar

Schweser, F., Deistung, A., Lehr, B. W., and Reichenbach, J. R. (2011). Quantitative imaging of intrinsic magnetic tissue properties using MRI signal phase: an approach to in vivo brain iron metabolism? Neuroimage. 54, 2789–2807. doi: 10.1016/j.neuroimage.2010.10.070

PubMed Abstract | CrossRef Full Text | Google Scholar

Schweser, F., Sommer, K., Deistung, A., and Reichenbach, J. R. (2012). Quantitative susceptibility mapping for investigating subtle susceptibility variations in the human brain. Neuroimage 62, 2083–2100. doi: 10.1016/j.neuroimage.2012.05.067

PubMed Abstract | CrossRef Full Text | Google Scholar

Smith, S. M. (2002). Fast robust automated brain extraction. Hum. Brain Mapp. 17, 143–155. doi: 10.1002/hbm.10062

PubMed Abstract | CrossRef Full Text | Google Scholar

Tang, J., Liu, S., Neelavalli, J., Cheng, Y. C. N., Buch, S., and Haacke, E. M. (2013). Improving susceptibility mapping using a threshold-based k-space/image domain iterative reconstruction approach. Magn. Reson. Med. 69, 1396–1407. doi: 10.1038/jid.2014.371

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Y., Chen, Y., Utriainen, D. T., and Haacke, E. M. (2019). “Automatic segmentation of deep grey matter structures for iron quantification,” in Proceedings of the 27th Annual Meeting of ISMRM Montréal, QC.

Google Scholar

Wang, Y., Chen, Y., Wu, D., Wang, Y., Sethi, S. K., Yang, G., et al. (2018). STrategically Acquired gradient echo (STAGE) imaging, part II: correcting for RF inhomogeneities in estimating T1 and proton density. Magn. Reson. Imaging 46, 140–150. doi: 10.1016/j.mri.2017.10.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., Bovik, A. C., Sheikh, H. R., and Simoncelli, E. P. (2004). Image quality assessment: From error visibility to structural similarity. IEEE Trans. Image Process. 13, 600–612. doi: 10.1109/TIP.2003.819861

PubMed Abstract | CrossRef Full Text | Google Scholar

Wharton, S., Schäfer, A., and Bowtell, R. (2010). Susceptibility mapping in the human brain using threshold-based k-space division. Magn. Reson. Med. 63, 1292–1304. doi: 10.1002/mrm.22334

PubMed Abstract | CrossRef Full Text | Google Scholar