Quantitative Susceptibility Mapping of Magnetic Quadrupole Moments

Cho, Junghun; Zhou, Dong; Kee, Youngwook; Spincemaille, Pascal; Wang, Yi

doi:https://doi.org/10.1155/2019/7174937

Concepts in Magnetic Resonance Part A, Bridging Education and Research

On this page

Abstract Introduction Materials and Methods Results Discussion Conclusion Appendix Data Availability Disclosure Conflicts of Interest Acknowledgments Supplementary Materials References Copyright Related Articles

Research Article | Open Access

Volume 2019 | Article ID 7174937 | https://doi.org/10.1155/2019/7174937

Quantitative Susceptibility Mapping of Magnetic Quadrupole Moments

Junghun Cho,¹Dong Zhou,²Youngwook Kee,²Pascal Spincemaille,²and Yi Wang^1,2

Academic Editor: Rongchun Zhang

Received30 Apr 2019

Revised03 Jul 2019

Accepted14 Jul 2019

Published06 Aug 2019

Abstract

We modeled the magnetic field up to the quadrupole term to investigate not only the average susceptibility (dipole), but also the susceptibility distribution (quadrupole) contribution. Expanding the magnetic field up to the order provides the quadrupole (: monopole, : dipole). Numerical simulations were performed to investigate the quadrupole contribution with subvoxel nonuniformity. Conventional dipole and our dipole + quadrupole models were compared in the simulation, the phantom and human brain. Furthermore, the quadrupole field was compared with the anisotropic susceptibility field in the dipole tensor model. In a nonuniformity case, numerical simulations showed a nonnegligible quadrupole field contribution. Our study showed a difference between the two methods in the susceptibility map at the edges; both the phantom and human studies showed sharper structural edges with the dipole + quadrupole model. Quadrupole moments showed contrast mainly at the structural boundaries. The quadrupole moment field contribution was smaller but nonnegligible compared to the anisotropic susceptibility contribution. Nonuniform and uniform source distributions can be separately considered by quadrupole expansion, which were mixed together in the dipole model. In the presence of nonuniformity, the susceptibility maps may be different between the two models. For a comprehensive field model, the quadrupole might need to be considered along with susceptibility anisotropy and microstructure effects.

1. Introduction

Quantitative susceptibility mapping (QSM) is a magnetic resonance imaging (MRI) technique that deconvolves the nonlocal field to reveal local tissue magnetic properties [1–7]. A wide range of clinical applications of QSM have been demonstrated, including calcifications and hemorrhages [8, 9], iron overload [10–14], cerebral microbleeds [15], Parkinson disease [16], and oxygen consumption [17–20]. The dipole field model, usually used in QSM, has the following implicit assumption: each voxel is considered uniform, and the corresponding magnetic field induced by each voxel is approximated as the field due to a single dipole. Under this assumption, the ill-posed field-source inverse problem has been well solved using multiple sample orientation [21, 22], or tissue structure guided regularization with the clinical single orientation [1, 2, 5, 23–25].

It is well known that dipole-approximation is accurate only when the observation point is far from the sources [5, 21]. When close to the sources, a higher order multipole field needs to be considered [26, 27], which can be substantial when the susceptibility distribution inside the voxel is not uniform. P-space imaging has been proposed to study neural architecture caused by quadrupoles [28]. However the definition of the quadrupoles in the p-imaging is the field inhomogeneity not the inhomogeneity of the susceptibility within a voxel. We introduce here the classical image-space quadrupole formulation [26, 27] that is explicitly defined by the susceptibility inhomogeneity at the voxel scale. This model expansion might improve the QSM accuracy and provide the information regarding subvoxel source distribution. We validate this quadrupole inclusion in QSM using numerical simulation and imaging of gadolinium balloon phantoms and human brains.

2. Materials and Methods

2.1. Theory

For an isotropic tissue with a volume magnetic susceptibility distribution under the polarization of the main MRI field , the tissue gains a magnetic dipole density pointing along the or z direction. We define the zeroth moment (dipole, ) and first moment (quadrupole ) of susceptibility in a voxel occupying space as [27] (appendix):

For a uniform susceptibility distribution, the voxel quadrupole moment = 0. This definition of quadrupole by the first moment of susceptibility distribution clearly reflects the nonuniform susceptibility, especially asymmetric source distribution, on a voxel scale, typically on the millimeter scale in human imaging. Microscopic variations at the molecular (~nm) and cellular (~ μm) scales are averaged out with a null contribution from the quadrupole moments.

Then, with accounting for the Lorentz sphere correction, the z-component of the tissue magnetic field ) is sensitized in the MR signal phase at voxel , and can be expressed in the following discretized form:where represents Cartesian coordinates. is the dipole kernel (, and is the quadruple kernel component [29].

The corresponding k-space expression iswhere the dipole and quadrupole kernels in k-space are and the dipole and quadrupole kernels are illustrated in Figure 1.

Figure 1

The dipole and quadrupole kernels in image- (top) and k-space (bottom). The kernels in image-space are based on (3) and (4), and those in k-space are based on (6) and (7). In image-space, all the kernels are real. The dipole kernel (d) decays as , and the quadrupole kernel () decays as . The dipole kernel, therefore, affects a wider range than the quadrupole kernel. In k-space, the dipole kernel (D) is real. The quadrupole kernel () is purely imaginary and proportional to . The magnitude of the quadrupole kernels, hence, increases as the k-coordinate increases.

The formulation above assumes that the main field is along z. If the main field is along , the dipole kernel (), quadrupole kernel (), and quadrupole moment () are changed correspondingly as follows:where R is the rotation matrix,

To solve the four unknowns in (5), multiple orientation measurements are needed.where is the local field in the jth orientation in k-space, and denotes the jth rotation. The first term on the right of (10) represents the dipole field dependence on rotation, and the second term represents the quadrupole field dependence on rotation.

2.2. Numerical Simulations

To investigate the quadrupole contribution to the field, numerical simulation was performed (Figure 2). In a high resolution (matrix size = [240 × 240 × 240], voxel size = [1 × 1 × 1 mm]), a dipole field was generated by a spherical susceptibility input source with continuous boundary (Input in Figure 2 and (11)). This is considered as ground truth field (black curve in Figure 2). Based on this high-resolution susceptibility input, the susceptibility and quadrupole moments (, , , and ) at lower resolutions (larger voxels) were calculated via (1) in image-space. Then, the dipole and dipole + quadrupole field were calculated with the low-resolution susceptibility and quadrupole moments, respectively (red and blue curve in Figure 2) via (5). The two fields were compared to the ground truth field. Note that the downsampling means making the voxel larger, not downsampling the high-resolution field. By increasing the downsampling factors (making the voxel size larger), we expect that the effect of the nonuniform source distribution within a voxel will increase. We used 2, 4, and 6 as downsampling factors (making voxel 2³, 4³, and 6³ times larger). where is the distance from the origin, is the inner sphere radius (12mm), and is the outer sphere radius (20mm).

Figure 2

Numerical simulation of the quadrupole contribution to the magnetic field. The field profile of the gray line in the Input χ (top left) is plotted as blue (dipole + quadrupole), red (dipole), and black line (ground truth) in the field graphs (right). Increasing the downsampling factor, effectively equivalent to increasing voxel size, the dipole + quadrupole field (blue) remains closer to the ground truth (black). This is especially notable for downsampling factors 4 and 6, where the dipole + quadrupole field (blue) captures the field peak at the boundary of the input sphere, whereas the dipole field (red) shows a noticeable drop.

To investigate the effect of quadrupole inclusion in model fitting, a second simulation was performed using a numerical brain phantom from the Cornell MRI research group (http://weill.cornell.edu/mri/pages/qsm.html) (Figure 3). The multi-echo gradient echo signal was simulated at a high resolution by the dipole model, and then the signal was downsampled to lower resolutions to mimic actual experiments with larger voxels. At the lower resolutions, the susceptibility was reconstructed with the dipole and dipole + quadrupole model, respectively, to compare which model provides more accurate susceptibility map. The detailed steps are as follows. First, the in the dataset was upsampled to a higher resolution 480 × 480 × 480 by zero-padding k-space. The resulting susceptibility map was then considered as the ground truth. To simulate a high-resolution gradient echo signals, the resolution was set to 1 × 1 × 1 mm³ and twenty directions (taken from the 12 acquired in the human subject 1 and 8 from the phantom acquisition as detailed below) were used. For each direction (k), was rotated to set the direction from to to mimic an actual experiment. High-resolution fields were then generated with the fixed direction 0 0 1 by convolving each rotated susceptibility map with the dipole kernel. Subsequently, gradient echo signals were set to M for each TE (2.6:2.6:28.6 ms) and each orientation, with . Then, downsampling was performed on these signals by taking the average over 2 × 2 × 2 (downsampling factor 2) or 4 × 4 × 4 (downsampling factor 4) voxel blocks. The downsampled signal was then used to fit the lower resolution field, which was then rotated back ( direction:). Finally, the lower resolution field was fitted with the dipole and dipole + quadrupole models ((12), with the identity field SNR weighing, = 1), respectively.

2.3. Phantom Experiment

A 1% agarose phantom was constructed with five balloons: one filled with tap water and four filled with gadolinium solutions with susceptibility values of 3.26, 1.63, 0.82, and 0.41 ppm. For data acquisition, all experiments were conducted on a 3.0T clinical Siemens scanner with multi-echo gradient echo sequence. The resolution was set to 1 × 1 × 1 mm³, TR to 27 ms, bandwidth to 260 Hz/pixel, and flip angle to 15°. Four TEs were used: 6.2, 11.5, 16.9, and 22.2 ms. The phantom was rotated with respect to the main magnetic field. Eight different rotations were acquired at the COSMOS optimal angles (0, 60, 120, and 270°) in the xy and yz planes (B₀ in z-axis) [21]

2.4. Human Subject Study

Human subject study was performed with approval from our institutional review board. Two volunteers, a 25 year old female (subject 1) and a 32 year old male (subject 2), were recruited to be scanned on a 3T clinical MR scanner (General Electric Excite HD; GE Healthcare, Waukesha, WI, USA) with a transmit/receive head coil to allow for rotation. For subject 1, a total of 12 orientations were acquired using a combination of neutral, left-, right-, forward-, and backward-leaning orientations in the supine and prone positions (shown in Figure 8). For each orientation, a gradient echo sequence (GRE) with the following imaging parameters was performed: 11 echoes with TE spacing of 2.64 ms, TR of 46.94 ms, and 0.9 × 0.9 × 1.5 mm³ resolution. A 2D echo planar imaging (EPI) dual spin echo DTI sequence was used to acquire diffusion tensor data in the neutral position with a resolution of 2×2×2.4 mm³, with the following parameters: 33 directions, b value of 1000 s/mm², 22cm FOV, acquisition matrix of 110×110, 2.4 mm slice thickness, 17s TR, 85.3 ms TE, and bandwidth of 1953.12 Hz/pixel. For subject 2, a total of 8 orientations were acquired. For each orientation, a GRE with the following imaging parameters was performed: 11 echoes with TE spacing of 2.97 ms, TR of 36.1 ms, and 1.2 × 1.2 × 1.2 mm³ resolution. A 2D echo planar imaging (EPI) dual spin echo DTI sequence was used to acquire diffusion tensor data in the neutral position with a resolution of 0.86×0.86×2.6 mm³, with the following parameters: 33 directions, b value of 1000 s/mm², 22 cm FOV, acquisition matrix of 256×256, 2.6 mm slice thickness, 17s TR, 107.4 ms TE, and bandwidth of 1953.12 Hz/pixel. DTI results are registered to QSM via the FSL FLIRT algorithm [30, 31].

2.5. Data Processing

For the phantom and human data postprocessing, phase unwrapping and background field removal were used to obtain the local field for each orientation [32]. The local fields were coregistered into the first orientation, and the rotation matrices () were obtained using the FSL FLIRT algorithm [30, 31] with 6 degrees of freedom (3 rotation + 3 translation). Then, the susceptibility and quadrupole moment were fitted by the dipole (D) or dipole + quadrupole (D+Q) model as follows. where is the field SNR weighting [5], indicates Fourier transform, and represents registered local fields. N is the total number of directions.

For the dipole (D) model,

For the dipole + quadrupole (D+Q) model,

We applied a right preconditioning technique [33] to improve convergence speed and reduce error propagation in the solution for the quadrupole model. This is similar to recent work in quantitative susceptibility mapping that deals with a large dynamic range [34], and QSM-based cerebral metabolic rate of oxygen (CMRO2) mapping that obtains two parameters in different scales, oxygen extract rate and non-blood susceptibility, simultaneously [18]. The system of equations was refined as , where and . is the preconditioner defined as , where are scaling factors chosen such that the elements of have a similar order of magnitude. , where (x, y, z direction) and and are voxel size and gradient of susceptibility in direction, respectively. The dipole model of susceptibility was used to calculate the preconditioning factors. We set this preconditioning constant based on the assumption that the susceptibility profile is linear as in , which provides the susceptibility of and the quadrupole moment of with (1) and (2). We used the conjugate gradient (CG) algorithm to solve (12). The maximum number of iterations was set to 40, similar to COSMOS [21].

In the second numerical simulation (the numerical brain phantom), we measured the sharpness at the white matter boundary. The white matter mask was constructed by thresholding the QSM value at -0.1 ppm. The boundary mask was constructed by subtracting the dilated white matter mask by the eroded mask (~ 2 voxel thickness) (Supporting Figure S1). For lower resolutions, the boundary mask was downsampled by cropping k-space, followed by thresholding at 0.5, resulting in boundary mask . To measure the sharpness at the structure boundary in the susceptibility maps, we calculated the average norm of the susceptibility gradient in the boundary region, , where is the voxel index and is the total number of voxels in the boundary region. A higher value indicates sharper edges.

In the phantom study, we performed ROI analysis and estimated the boundary sharpness in susceptibility maps. For the ROI analysis, the balloons were segmented by thresholding the norm of every echo magnitude image, and eroded 2 voxels. The susceptibility of each gadolinium balloon was measured relative to the average susceptibility of a water balloon. To obtain the boundary mask, we dilated and eroded the balloon masks, and the dilated mask was subtracted by the eroded mask, ~ 2 voxel thickness (Supporting Figure S1).

In the human study, we first used the same method described above in the numerical simulation to measure the sharpness at the white matter boundary on susceptibility maps. The white matter mask was constructed by thresholding the DTI fractional anisotropy (FA) value at 0.3. The boundary mask was again constructed by subtracting the dilated white matter mask by the eroded mask (~ 2 voxel thickness). Second, we measured the contribution of the quadrupole moments to the magnetic field using a forward model calculation, starting from the results in subject 1 with in the z direction, and compared it with χ₂₃, χ₁₃, and χ₃₃ obtained following a conventional susceptibility tensor imaging (STI) reconstruction [35–37]: where and . A conjugate gradient solver with a maximum of 40 iterations was used to solve this system.

Third, we investigated the possibility of susceptibility anisotropy and white matter microstructure contributions to the quadrupole moment. A scatter plot between quadrupole moment and white matter fiber tract orientation by the independent DTI measurement was obtained, and linear regression was subsequently performed. Fourth, we investigated the effect of the number of orientations on χ, , , and maps in subject 1. The optimization was performed with 4, 6, 8, and 12 orientations. 4 was set as the minimum number of B0 orientations required for the model inversion of D+Q, which has four unknowns. Additionally, 6 and 8 B0 orientations were used to investigate how much the accuracy of result (χ, , , ) increases as orientations are added. For a given number (equal to 4, 6 or 8), we considered all the possible groups of selections, e.g., groups for , and, for each group, computed the average inner product of the directions within that group. A large average inner product for a group indicates that the directions in that group tend to be parallel and are likely to provide a high condition number in inverse problem. A small average inner product group indicates that the directions in that group tend to be orthogonal to each other and present a small condition number. We then selected the 3 groups of orientations with the minimum, maximum, and mean average inner products. The relative result difference between N = 4, 6, or 8 orientations and 12 orientations, , was calculated in the three groups. Then, the mean and standard deviation were obtained. The same stopping criteria of a maximum of 40 iterations were used for all cases.

3. Results

3.1. Numerical Simulations

Figure 2 shows the ground truth (black), dipole (red), and dipole + quadrupole fields (blue) at different downsampling factors, i.e., different degrees of nonuniformity. These fields are plotted along the z-axis from the center of the object (gray line in the input χ in Figure 2). The same ground truth (black) is plotted at all downsampling factors. There is no significant difference in the fields when the downsampling factor is 2, but as it is increased to 4 and 6, the dipole + quadrupole field remains closer to the ground truth than the dipole field. Notably, the field peak at the object boundary is captured by the dipole + quadrupole field (blue), whereas a noticeable drop is observed in the dipole field (red).

Figure 3 shows the reconstructed susceptibility maps in lower resolutions (downsampling factors 2 and 4) with the dipole (D) and dipole + quadrupole (D+Q) models. D+Q provides a sharper susceptibility map than the D model. For a downsampling factor of 2, the average norm of the susceptibility gradient in the white matter boundary region was 6.78×10⁻³, 1.32×10⁻², and 1.32×10⁻² ppm/mm for D, D+Q, and the ground truth, respectively. For a downsampling factor of 4, it was 3.60×10⁻³, 6.51×10⁻³, and 7.33×10⁻³ ppm/mm for D, D+Q, and the ground truth, respectively (Table 1). The quadrupole moments () become greater in magnitude for the downsampling factor of 4 as compared to a factor of 2.

3.2. Phantom Experiment

Figure 4 shows the susceptibility () and quadrupole moments () of the phantom for the D and D+Q models. There are no significant differences in maps between two models except for the balloon edges: For balloons with susceptibility values of 3.26, 1.63, 0.82, and 0.41 ppm, the susceptibility values from the D model are 3.05±0.03, 1.43±0.02, 0.75±0.01, and 0.36±0.01 ppm, respectively, and the values from the D+Q model are 3.04±0.04, 1.43±0.02, 0.76±0.01, and 0.35±0.01 ppm, respectively. The susceptibility map shows slightly sharper balloon edges with the D+Q model (Figure 4 and Table 1); the average norm of the susceptibility gradient in the highest susceptibility balloon boundary is 1.36 for D+Q and 1.31 ppm/mm for the D model (Table 1). The quadrupole moments are apparent only at the balloon boundaries.

3.3. Human Study

Figure 5 shows the susceptibility () and quadrupole moments () of human brains for the D+Q model. In both subjects, the D+Q model (χ_q) shows sharper boundaries between gray and white matter than the D model (χ_d). Consequently, detailed structures are shown more clearly (indicated by yellow arrows in susceptibility maps). The average norm of the susceptibility gradient in the white matter boundary region is 1.60×10⁻² for D+Q and 1.13×10⁻² ppm/mm for the D model in subject 1; and 1.58×10⁻² for D+Q and 1.19×10⁻² ppm/mm for D model in subject 2 (Table 1).

Figure 5

Comparison between dipole and dipole + quadrupole model in the human brain. The susceptibility map from the dipole + quadrupole model shows sharper edges between white and gray matter than the map from the dipole model for both subjects 1 and 2. Consequently, the detailed structures are more clearly shown by the dipole + quadrupole model, as indicated by the yellow arrows (left). The quadrupole moment generally looks obvious at the brain structure boundaries (yellow arrow in axial Q₁), but lacks contrast at some boundaries (red arrow in axial Q₁). Also, in Q₂, the boundary between the Putamen and Globus pallidus (yellow arrows in axial Q₂) is bright on the left brain side, whereas right brain side appears dark.

Quadrupole moments are generally nontrivial at the brain structure boundaries (yellow arrow in axial ), but lack contrast at some boundaries (red arrow in axial ). They also show nonzero contrast at nonboundary regions. In addition, in axial , the same brain structures show opposite contrast between the left and right brain side, as seen in the boundary between the Putamen and Globus pallidus (yellow arrows in axial ). By comparing the quadrupole moment maps between subjects 1 and 2, the contrast strength appears different even though the signs are matched (Figure 5). For instance, compared to subject 2, subject 1 shows a smaller contrast in and , but greater contrast in .

Figure 6 shows the field generated by each of the quadrupole moments in the D+Q model and the field generated by the anisotropic susceptibility tensor components (left column). Additionally, a histogram of each field contribution across the brain is shown in the right column. The quadrupole contributions (red) to the field were substantially smaller than the isotropic susceptibility contribution (black), and smaller but nonnegligible compared to the off-diagonal anisotropy tensor contributions (blue) in magnitude. The standard deviation of the fields generated by , , , , , and were 0.0020, 0.0020, 0.0039, 0.0055, 0.0062, and 0.0119 ppm, respectively.

Figure 7 shows the scatter plot between the quadrupole moment and white matter fiber tract orientation as measured by DTI. The quadrupole moment did not show a significant linear dependence on orientation (R² = 6.4×10⁻⁴, 8.7×10⁻⁵, 4.0×10⁻⁴ for , , and, respectively).

Figure 8

The effect of the number of orientations on χ, Q₁, Q₂, and Q₃ maps in subject 1. Optimization was performed with 4, 6, 8, and 12 different orientations. As the number of orientation increased, both susceptibility and quadrupole moment maps appeared more similar to those obtained using 12 orientations. The χ, Q₁, Q₂, and Q₃ maps were obtained from the group of directions with the minimum average inner product (, , and for = 8, 6, and 4 orientations, respectively). The average and standard deviation of the relative result difference () over three different groups of orientations (the maximum, mean, and minimum average inner product within that group) were shown below each map.

Figure 8 shows the effect of the number of included orientations on the , , , and maps in subject 1. As the number of orientation increased, both susceptibility and quadrupole moment maps appeared more similar to those obtained using 12 orientations. With , and 8 orientations, the difference with respect to the 12 orientation result was 56.3 ± 14.2%, 24.9 ± 11.9%, and 15.0 ± 5.6%, respectively. The quadrupole moment difference was 71.8 ± 15.4%, 50.6 ± 15.1%, and 32.1 ± 11.3% for , and 8 orientations, respectively.

4. Discussion

The quadrupole formulated in the image-space can be added to the standard dipole model of QSM to account for nonuniform susceptibility distribution in a voxel. Our data confirm that the quadrupole terms are substantial in regions of nonuniform susceptibility, such as the gray-white matter interface. Also, the field contribution by the quadrupole moment in our dipole + quadrupole (D+Q) model was not negligible compared to the anisotropic susceptibility field in the dipole tensor model (Figure 6). This suggests that the quadrupole expansion might need to be considered along with susceptibility tensor imaging and microstructure effect to construct a comprehensive field model.

The standard quadrupole definition from classic electrodynamics is used in this work [26, 27]. In the previous p-space imaging, multipole expansion was also used to capture spatial variations of the phase or magnetic field in a voxel [28]. The major difference is that the quadrupole here is defined on the spatial inhomogeneity of the susceptibility in a voxel, while p-imaging is defined on the field inhomogeneity. P-imaging data allows reconstruction of the water proton magnetization distribution within a voxel, the phase of which allows estimation of the field distribution within the voxel. If QSM is applied on this phase to generate susceptibility distribution in the voxel, we can then compute the classic quadrupole moment as defined by (1). Because a strong susceptibility source in a neighboring voxel can also generate field variation, p-imaging without applying QSM to field distribution within a voxel can be confounded by sources in neighboring voxels or blooming artifacts.

The susceptibility values estimated by the dipole (D) and dipole + quadrupole (D+Q) models were not different for Gd doped water balloons (Figure 4), an expected result due to the uniform distribution of paramagnetic Gd inside the balloons. The edges in the D+Q map were sharper than those in the D map (Figure 4 and Table 1). In the human study, the D+Q susceptibility map showed sharper edges between white and gray matter as compared to the D susceptibility map (Figure 5 and Table 1). The numerical brain simulation showed that for several image resolutions the inclusion of the quadrupole could provide a sharper susceptibility map in D+Q as compared to the D model (Figure 3), similar to the sharpness of the ground truth. This indicates that partial volume artifact in QSM may be alleviated by inclusion of quadrupole moment.

The quadrupole moments appeared to be substantial at structural boundaries in both the phantom and human brains, as expected. However, nonzero contrast was also observed at nonboundary regions (Figure 5), which might reflect actual subvoxel susceptibility variations. In white matter, variations in myelin concentration from anterior to posterior corpus callosum have been observed in recent myelin water fraction [18, 38] and histology [39, 40] studies, thus implying that the quadrupole moment within white matter may be nonzero. Further studies are needed to investigate myelin nonuniformity and quadrupole moments in depth. Another possible cause is the effects of microstructures and susceptibility anisotropy in white matter. To investigate their effects on the estimated quadrupole moments in white matter, we attempted to find a relationship between these moments and the local fiber orientation independently obtained using DTI analysis, but found no significant correlation (Figure 7). This suggests that the quadrupole moment does not result from fitting to the fields by the microstructure or susceptibility anisotropy in a significant manner.

The definition of the quadrupole moment (1) explains the positive/negative signs of the quadrupole moments as the susceptibility value increases/decreases along an axis (Figures 4 and 5). For example, the at the boundary of the Putamen and Globus pallidus showed a positive contrast on the left brain side, but a negative contrast on the right brain side (yellow arrows in axial , Figure 5), indicating that the Putamen has a lower susceptibility than Globus pallidus. Within the boundary voxel on the left brain side, the subvoxel susceptibility increases as y increases, thus providing a positive quadrupole moment, based on (1). On the contrary, the subvoxel susceptibility decreases as y increases on the right brain side, which provides a negative quadrupole moment. In terms of the quadrupole moment comparison between subjects 1 and 2, subject 1 showed smaller contrast in and but greater contrast in . This is likely due to voxel size differences between the two subjects because the quadrupole moment is a function of voxel size. For instance, if the susceptibility profile within a voxel is linear, the quadrupole moment is linearly proportional to the square of voxel linear dimension (1). The data from subject 1 has smaller in-plane voxel size (0.9 vs 1.2 mm) but a larger voxel size in the z direction (1.5 vs 1.2 mm), which is constant with smaller contrast in and , and greater contrast in . Additionally, 12 and 8 directions were used for subjects 1 and 2, respectively, which might result in more accurate maps for subject 1.

Despite the advantage of quadrupole inclusion, the multiple orientations in data acquisition needed to determine quadruple moments are difficult to perform on research volunteers and impractical to perform in clinical practice. Fewer orientations, e.g., 4 or 6, could generate maps, but they would be relatively inaccurate (Figure 8). Tissue structural information may be used to regularize the inverse problem for (5), conceptually similar to the morphology enabled dipole inversion [1, 41] and the DTI guided susceptibility tensor imaging [37, 42]. Also, only isotropic (scalar) susceptibility distributions are considered here. A tensor model that could cope with existing field models such as the generalized Lorentz tensor approach might be needed [35, 37, 43, 44]. In addition, larger datasets, e.g., more human subjects and animal disease models, are needed to confirm the observations in this study.

5. Conclusion

In summary, tissue magnetic susceptibility sources can be expressed in a multipole expansion, and the dipole model in current QSM may be improved by including the quadrupole moments, defined as the first moment of susceptibility distribution within a voxel. In regions with nonuniform susceptibility distribution such as gray-white matter interfaces, the quadrupole moments can be substantial and their field contribution is nonnegligible compared to susceptibility anisotropy field contributions. The quadrupole moment might need to be considered for a comprehensive field model.

Appendix

For derivation of (2), we start with the “H” field that has a scalar potential solution to the magnetostatics problem [26].

This has singularity at the origin where the measurement point () and source location () are identical. To facilitate calculus on , the origin has to be excluded [29].

where H is the Heaviside function defined by the following.

The term can be expanded up to the order with Taylor expansion in powers of with the localized source condition [26]. This leads to (A.6) with the along .

where the dipole interaction, , and the susceptibility, , are as follows.

The quadrupole interaction, , and the quadrupole moments, , are as follows.

The field z-component can be obtained as follows, based on (A.6).

TheB field in (2) can be obtained with . The limit to 0 was also taken to the and terms in (A.9) [29].

Data Availability

Phantom and human subject data come from Wang lab [41] and are available on request. Reconstruction MATLAB codes regarding local field estimation are available at http://pre.weill.cornell.edu/mri/pages/qsm.html. D and D+Q model reconstruction codes are available on request.

Disclosure

An earlier version of this study has been presented in a meeting at ISMRM Annual Meeting & Exhibition, Singapore, 2016.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this article.

Acknowledgments

This study was supported by the National Institutes of Health under award numbers R01NS072370, R01NS090464, R01NS095562, and S10OD021782. We thank Kelly McCabe Gillen, Ph.D., for her assistance in manuscript editing.

Supplementary Materials

Supporting Figure S1. The boundary mask for measuring sharpness. For the numerical brain phantom, the white matter mask was constructed by thresholding the QSM value at -0.1 ppm. For the phantom, the balloons were segmented by thresholding the norm of every echo magnitude image. For the human subjects, the white matter mask was constructed by thresholding the DTI FA value at 0.3. The boundary mask was constructed by subtracting the dilated mask by the eroded mask (~2 voxel thickness). (Supplementary Materials)

References

L. de Rochefort, T. Liu, B. Kressler et al., “Quantitative susceptibility map reconstruction from MR phase data using bayesian regularization: validation and application to brain imaging,” Magnetic Resonance in Medicine, vol. 63, no. 1, pp. 194–206, 2010.
View at: Publisher Site | Google Scholar
B. Kressler, L. de Rochefort, T. Liu, P. Spincemaille, Q. Jiang, and Y. Wang, “Nonlinear regularization for per voxel estimation of magnetic susceptibility distributions from MRI field maps,” IEEE Transactions on Medical Imaging, vol. 29, no. 2, pp. 273–281, 2010.
View at: Publisher Site | Google Scholar
K. Shmueli, J. A. de Zwart, P. van Gelderen, T.-Q. Li, S. J. Dodd, and J. H. Duyn, “Magnetic susceptibility mapping of brain tissue in vivo using MRI phase data,” Magnetic Resonance in Medicine, vol. 62, no. 6, pp. 1510–1522, 2009.
View at: Publisher Site | Google Scholar
S. Wharton, A. Schäfer, and R. Bowtell, “Susceptibility mapping in the human brain using threshold-based k-space division,” Magnetic Resonance in Medicine, vol. 63, no. 5, pp. 1292–1304, 2010.
View at: Publisher Site | Google Scholar
J. Liu, T. Liu, L. de Rochefort et al., “Morphology enabled dipole inversion for quantitative susceptibility mapping using structural consistency between the magnitude image and the susceptibility map,” NeuroImage, vol. 59, no. 3, pp. 2560–2568, 2012.
View at: Publisher Site | Google Scholar
F. Schweser, A. Deistung, B. W. Lehr, and J. R. Reichenbach, “Quantitative imaging of intrinsic magnetic tissue properties using MRI signal phase: An approach to in vivo brain iron metabolism?” NeuroImage, vol. 54, no. 4, pp. 2789–2807, 2011.
View at: Publisher Site | Google Scholar
J. Li, S. Chang, T. Liu et al., “Reducing the object orientation dependence of susceptibility effects in gradient echo MRI through quantitative susceptibility mapping,” Magnetic Resonance in Medicine, vol. 68, no. 5, pp. 1563–1569, 2012.
View at: Publisher Site | Google Scholar
W. Chen, W. Zhu, I. Kovanlikaya et al., “Intracranial calcifications and hemorrhages: characterization with quantitative susceptibility mapping,” Radiology, vol. 270, no. 2, pp. 496–505, 2014.
View at: Publisher Site | Google Scholar
A. Deistung, F. Schweser, B. Wiestler et al., “Quantitative susceptibility mapping differentiates between blood depositions and calcifications in patients with glioblastoma,” PLoS ONE, vol. 8, no. 3, Article ID e57924, 2013.
View at: Publisher Site | Google Scholar
W. Chen, S. A. Gauthier, A. Gupta et al., “Quantitative susceptibility mapping of multiple sclerosis lesions at various ages,” Radiology, vol. 271, no. 1, pp. 183–192, 2014.
View at: Publisher Site | Google Scholar
C. Langkammer, T. Liu, M. Khalil et al., “Quantitative susceptibility mapping in multiple sclerosis,” Radiology, vol. 267, no. 2, pp. 551–559, 2013.
View at: Publisher Site | Google Scholar
A. I. Blazejewska, A. M. Al-Radaideh, S. Wharton et al., “Increase in the iron content of the substantia nigra and red nucleus in multiple sclerosis and clinically isolated syndrome: A 7 Tesla MRI study,” Journal of Magnetic Resonance Imaging, vol. 41, no. 4, pp. 1065–1070, 2015.
View at: Publisher Site | Google Scholar
X. Li, D. M. Harrison, H. Liu et al., “Magnetic susceptibility contrast variations in multiple sclerosis lesions,” Journal of Magnetic Resonance Imaging, vol. 43, no. 2, pp. 463–473, 2016.
View at: Publisher Site | Google Scholar
T. Christen, H. Schmiedeskamp, M. Straka, R. Bammer, and G. Zaharchuk, “Measuring brain oxygenation in humans using a multiparametric quantitative blood oxygenation level dependent MRI approach,” Magnetic Resonance in Medicine, vol. 68, no. 3, pp. 905–911, 2012.
View at: Publisher Site | Google Scholar
T. Liu, K. Surapaneni, M. Lou, L. Cheng, P. Spincemaille, and Y. Wang, “Cerebral microbleeds: burden assessment by using quantitative susceptibility mapping,” Radiology, vol. 262, no. 1, pp. 269–278, 2012.
View at: Publisher Site | Google Scholar
D. Bulte, M. Kelly, M. Germuska et al., “Quantitative measurement of cerebral physiology using respiratory-calibrated MRI,” NeuroImage, vol. 60, no. 1, pp. 582–591, 2012.
View at: Publisher Site | Google Scholar
J. Zhang, T. Liu, A. Gupta, P. Spincemaille, T. D. Nguyen, and Y. Wang, “Quantitative mapping of cerebral metabolic rate of oxygen (CMRO₂) using quantitative susceptibility mapping (QSM),” Magnetic Resonance in Medicine, vol. 74, no. 4, pp. 945–952, 2015.
View at: Publisher Site | Google Scholar
J. Zhang, D. Zhou, T. D. Nguyen, P. Spincemaille, A. Gupta, and Y. Wang, “Cerebral metabolic rate of oxygen (CMRO 2 ) mapping with hyperventilation challenge using quantitative susceptibility mapping (QSM),” Magnetic Resonance in Medicine, vol. 77, no. 5, pp. 1762–1773, 2017.
View at: Publisher Site | Google Scholar
A. P. Fan, B. Bilgic, L. Gagnon et al., “Quantitative oxygenation venography from MRI phase,” Magnetic Resonance in Medicine, vol. 72, no. 1, pp. 149–159, 2014.
View at: Publisher Site | Google Scholar
J. D. Dickson, T. W. Ash, G. B. Williams et al., “Quantitative BOLD: The effect of diffusion,” Journal of Magnetic Resonance Imaging, vol. 32, no. 4, pp. 953–961, 2010.
View at: Publisher Site | Google Scholar
T. Liu, P. Spincemaille, L. de Rochefort, B. Kressler, and Y. Wang, “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, vol. 61, no. 1, pp. 196–204, 2009.
View at: Publisher Site | Google Scholar
R. G. Wise, A. D. Harris, A. J. Stone, and K. Murphy, “Measurement of OEF and absolute CMRO2: MRI-based methods using interleaved and combined hypercapnia and hyperoxia,” NeuroImage, vol. 83, pp. 135–147, 2013.
View at: Publisher Site | Google Scholar
F. Schweser, A. Deistung, K. Sommer, and J. R. Reichenbach, “Toward online reconstruction of quantitative susceptibility maps: Superfast dipole inversion,” Magnetic Resonance in Medicine, vol. 69, no. 6, pp. 1581–1593, 2013.
View at: Publisher Site | Google Scholar
J. Tang, S. Liu, J. Neelavalli, Y. C. N. Cheng, S. Buch, and E. M. Haacke, “Improving susceptibility mapping using a threshold-based K-space/image domain iterative reconstruction approach,” Magnetic Resonance in Medicine, vol. 69, no. 5, pp. 1396–1407, 2013.
View at: Publisher Site | Google Scholar
B. Wu, W. Li, A. Guidon, and C. Liu, “Whole brain susceptibility mapping using compressed sensing,” Magnetic Resonance in Medicine, vol. 67, no. 1, pp. 137–147, 2012.
View at: Publisher Site | Google Scholar
J. D. Jackson, Classical Electrodynamics, Wiley, New York, NY, USA, 1999.
View at: MathSciNet
J. M. Blatt and V. F. Weisskopf, Theoretical Nuclear Physics, Dover Publications, New York, NY, USA, 1991.
C. Liu and W. Li, “Imaging neural architecture of the brain based on its multipole magnetic response,” NeuroImage, vol. 67, pp. 193–202, 2013.
View at: Publisher Site | Google Scholar
R. P. Kanwal, Generalized Functions Theory and Technique: Theory and Technique, Springer Science & Business Media, 2012.
M. Jenkinson and S. Smith, “A global optimisation method for robust affine registration of brain images,” Medical Image Analysis, vol. 5, no. 2, pp. 143–156, 2001.
View at: Publisher Site | Google Scholar
M. Jenkinson, P. Bannister, M. Brady, and S. Smith, “Improved optimization for the robust and accurate linear registration and motion correction of brain images,” NeuroImage, vol. 17, no. 2, pp. 825–841, 2002.
View at: Publisher Site | Google Scholar
T. Liu, I. Khalidov, L. de Rochefort et al., “A novel background field removal method for MRI using projection onto dipole fields (PDF),” NMR in Biomedicine, vol. 24, no. 9, pp. 1129–1136, 2011.
View at: Publisher Site | Google Scholar
M. Benzi, “Preconditioning techniques for large linear systems: a survey,” Journal of Computational Physics, vol. 182, no. 2, pp. 418–477, 2002.
View at: Publisher Site | Google Scholar | MathSciNet
Z. Liu, Y. Kee, D. Zhou, Y. Wang, and P. Spincemaille, “Preconditioned total field inversion (TFI) method for quantitative susceptibility mapping,” Magnetic Resonance in Medicine, vol. 78, no. 1, pp. 303–315, 2017.
View at: Publisher Site | Google Scholar
C. Liu, “Susceptibility tensor imaging,” Magnetic Resonance in Medicine, vol. 63, no. 6, pp. 1471–1477, 2010.
View at: Publisher Site | Google Scholar
C. Liu, W. Li, B. Wu, Y. Jiang, and G. A. Johnson, “3D fiber tractography with susceptibility tensor imaging,” NeuroImage, vol. 59, no. 2, pp. 1290–1298, 2012.
View at: Publisher Site | Google Scholar
C. Wisnieff, T. Liu, P. Spincemaille, S. Wang, D. Zhou, and Y. Wang, “Magnetic susceptibility anisotropy: Cylindrical symmetry from macroscopically ordered anisotropic molecules and accuracy of MRI measurements using few orientations,” NeuroImage, vol. 70, pp. 363–376, 2013.
View at: Publisher Site | Google Scholar
D. Hwang, D.-H. Kim, and Y. P. Du, “In vivo multi-slice mapping of myelin water content using T2^* decay,” NeuroImage, vol. 52, no. 1, pp. 198–204, 2010.
View at: Publisher Site | Google Scholar
N. Stikov, J. S. Campbell, T. Stroh et al., “In vivo histology of the myelin g-ratio with magnetic resonance imaging,” NeuroImage, vol. 118, pp. 397–405, 2015.
View at: Publisher Site | Google Scholar
F. Aboitiz, A. B. Scheibel, R. S. Fisher, and E. Zaidel, “Fiber composition of the human corpus callosum,” Brain Research, vol. 598, no. 1-2, pp. 143–153, 1992.
View at: Publisher Site | Google Scholar
Y. Wang and T. Liu, “Quantitative susceptibility mapping (QSM): decoding MRI data for a tissue magnetic biomarker,” Magnetic Resonance in Medicine, vol. 73, no. 1, pp. 82–101, 2015.
View at: Publisher Site | Google Scholar
X. Li, D. S. Vikram, I. A. Lim, C. K. Jones, J. A. Farrell, and P. C. van Zijl, “Mapping magnetic susceptibility anisotropies of white matter in vivo in the human brain at 7T,” NeuroImage, vol. 62, no. 1, pp. 314–330, 2012.
View at: Publisher Site | Google Scholar
D. A. Yablonskiy and A. L. Sukstanskii, “Generalized Lorentzian Tensor Approach (GLTA) as a biophysical background for quantitative susceptibility mapping,” Magnetic Resonance in Medicine, vol. 73, no. 2, pp. 757–764, 2015.
View at: Publisher Site | Google Scholar
D. A. Yablonskiy and A. L. Sukstanskii, “Effects of biological tissue structural anisotropy and anisotropy of magnetic susceptibility on the gradient echo MRI signal phase: theoretical background,” NMR in Biomedicine, vol. 30, no. 4, p. e3655, 2017.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2019 Junghun Cho et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

1406

Downloads

1598

Citations