Nuclear paramagnetism‐induced MR frequency shift and its implications for MR‐based magnetic susceptibility measurement

To investigate the 1H spin contribution (0.004 parts per million (ppm)) to the water magnetic susceptibility and discuss its implications for high‐precision phase mapping and tissue susceptibility measurement.


INTRODUCTION
Modern high-field MRI allows determination of local proton Larmor frequency shift with a very high precision. For example, in a 3 Tesla (T) gradient-echo (GRE) image with an echo time (TE) of 30 ms and signal-to-noise ratio (SNR) of 50, the image phase can be determined with an uncertainty of 1/50 rad % 1 , which corresponds to the Larmor frequency shift of 1 /30 ms % 0.1 Hz. This is less than 1 part per billion (ppb) of the center frequency at 3T.
Such high precision is partly responsible for the success of quantitative susceptibility mapping (QSM), which solves the magnetic dipolar field inverse problem to reconstruct tissue magnetic susceptibility distribution from multi-echo GRE phase data (1). Recent progress in QSM has allowed spatially resolved, in vivo measurement of tissue magnetic susceptibility with a very high precision (2)(3)(4)(5). For example, tissue susceptibility measurements of the human brain are routinely reported with error bars on the order of 10 ppb (5). A recent ex vivo mouse brain scan at 9.4T (6) exhibited a high-resolution susceptibility map in the range of À50 to 50 ppb, with apparent (visual) susceptibility resolution better than 5 ppb. Finally, susceptibility tensor imaging exploits the orientation-dependent apparent susceptibility of myelinated axons, with a peak-to-peak variation of about 25 ppb (7), to reconstruct the susceptibility tensor map of the white matter.
It appears that recent advances in MR-based susceptibility measurement have reached the point where magnetic susceptibility of hydrogen nuclear paramagnetism, on the order of 4 ppb in water, may start to be significant. So far this well-established but small nuclear contribution to the tissue magnetic susceptibility (8) has not been much discussed in the experimental MRI literature. Despite its small magnitude, the fact that MRI relies on active manipulation of nuclear spin polarization implies that the nuclear contribution to the total tissue susceptibility is dynamically changing and might therefore affect the accuracy of MRI-based susceptibility measurements if not properly accounted for.
In this work, we experimentally demonstrate the effect of nuclear spin polarization on the MR frequency shift. The shift depends on the tip angle and the shape of the (nuclear) susceptibility boundary in a way that is consistent with theoretical predictions. Our observation implies that high-precision B 0 mapping and susceptibility quantification can be affected by the dynamic and spatial variations of the proton spin longitudinal polarization, which in turn are influenced by sequence parameters and tissue longitudinal relaxation time (T 1 ).

Nuclear Paramagnetism and Tissue Susceptibility
Tissue magnetic susceptibility (xÞ comprises the electronic (x elec Þ and the nuclear (x nucl ) components and can be expressed as Here, x elec includes the orbital diamagnetic and spin paramagnetic contributions. For most biological tissue, x elec ¼ À962 parts per million (ppm) (8). On the other hand, x nucl is due to the paramagnetic alignment of nuclear spins with an applied magnetic field. Such alignment gives rise to the equilibrium magnetization M 0 that is used in standard MRI. For 1 H in water at room temperature, x nucl is given by where m 0 ¼ 4p Â 10 À7 H=m, k B ¼ 1:38 Â 10 À23 J=K, T ¼ 293 K, n spin ¼ 110 mol=l, and m p ¼ 1:41 Â 10 À26 Am 2 is the proton magnetic moment. We note that x nucl is inversely proportional to the temperature, such that at the physiological temperature (¼ 310 KÞ it is slightly reduced to 3:87 ppb. In this work, we assume In MRI, the nuclear spin alignment and the spin isocromat's polar angle a with respect to the main magnetic field dynamically change during a pulse sequence execution. Taking this into account, we can rewrite x nucl as the projection of the equilibrium (maximum) susceptibility onto the longitudinal axis, where x eq nucl ¼ 4:09 ppb. The cosa term captures the dependence of tissue susceptibility on the nuclear spin tip angle. Equation [3] represents nuclear susceptibility just after the excitation before magnetization recovery in a sequence with a long repetition time (TR > 5 T 1 ). If TR is short, the longitudinal magnetization reaches a steady state that effectively reduces the equilibrium magnetization as (9) This results in the proportional change in x nucl in a steady state acquisition, Equation [5] reduces to Equation [3] when TR=T 1 ! 1.

Resonance Frequency Shift
A uniformly magnetized object exhibits MR frequency shift proportional to the magnetization, with a spatial distribution that depends on the shape of the object. In this work, we considered doped water phantoms in a thin, flat square-shaped container with its face normal or parallel to B 0 . When normal to B 0 , the phantom approximates a thin oblate spheroid with a demagnetization factor that approaches 1, which means the magnetization M of the object cancels the self-induced part of the H-field such that the total macroscopic magnetic field B mac equals the applied field: In MRI, the resonance frequency is determined by B mac minus 2m 0 M=3, the self field of a small hypothetical sphere (Lorentz sphere) carved out of the magnetized medium at the position of the spin (10-12): Using m 0 M ¼ xm 0 H ¼ xB 1 þ OðxÞ % xB, which is valid up to an error on the order of x $ 10 À5 in biological tissue, the susceptibility-induced field shift in an infinitely thin, flat phantom normal to B 0 is given by On the other hand, an infinitely thin disc with its face parallel to B 0 has a vanishing demagnetization factor. In this case, the susceptibility-induced field shift is larger than the above by m 0 M ¼ xB 0 , namely, Equations [6] and [8] can also be derived from the magnetic field inside a uniformly magnetized infinite slab (13). Demagnetization factors of various geometrical solids are listed in, for example, (14).

Tip Angle-Dependent MR Frequency Shift
From Equations [1], [3], [7] and [8], we find that the nuclear magnetization-dependent B 0 shifts for a thin phantom in the two orientations are given by x eq nucl Á cosa [9] and dB jj 0;nucl nucl Á cosa; [10] in the limit of long TR.

Virtual Phantom: Slice Orientation Dependence
To the extent that we are concerned with the nuclear spin effects, a virtual thin phantom can be created simply by exciting nuclear spins on a slice in a 3D object. A thin slice in the axial plane with tip angle a will exhibit tip angle-dependent frequency shift in the slice according to Equation [9]. Similarly, an excited slice in the coronal (or sagittal) plane will behave like a thin parallel phantom with a shift given by Equation [10].

Effect of Finite Thickness
Equations [7][8][9][10] correspond to limiting cases of an infinitely thin phantom (slice). For a phantom with a finite thickness, the susceptibility-induced field shift can be computed numerically by summing up the dipolar fields from all the voxels in the phantom. The magnitude of the frequency shifts at the center of the phantom is generally smaller than those given by Equations [9] and [10] (Table 1).

METHODS
The tip angle-dependent frequency shift was measured from the free-induction-decay (FID) signals from two different thin square phantoms at 3T and 9.4T. Additionally, nuclear spin-dependent MR frequency shift was examined from 2D B 0 maps obtained at 3T on a spherical phantom with different tip angles and slice orientations.
In each experiment, the 90 tip angle was calibrated by searching for the maximum in the FID amplitude from the whole phantom as the transmit radiofrequency (RF) coil voltage was increased from zero. For both the 3T and 9.4T experiments, the standard body (volume) coil of the scanner was used for transmission. Due to the thin geometry of the square phantoms and relatively small size of the sphere (10-cm diameter, in a whole-body scanner), B 1 field inhomogeneity was expected to be low. This was verified through a modified double-angle (triple angle) method for B 1 mapping, where the ratio of two gradient echo images with nominal flip angles 30 and 90 was converted to the actual flip angle map, using the relationship: sin3a=sina ¼ 2cos2a þ 1. The B 1 inhomogeneity (standard deviation divided by the mean) from this method was less than 10% in the imaged slices.

3T FID Measurement
The first phantom FID measurements were performed on a 3T scanner (Magnetom Prisma, Siemens Healthineers, Erlangen, Germany). The phantom consisted of a squareshaped acrylic container with inner dimensions of 15 Â15 Â2 cm 3 and wall thickness of 2 mm, filled with physiological saline doped with 0.25 mM of Gd-DO3Abutrol (Gadavist, Bayer Healthcare, Whippany, NJ, USA). The bulk magnetic susceptibility of the acrylic material was separately measured to be À10.0 6 0.2 ppm from its effect on the B 0 profile of a nearby spherical phantom. The measured T 1 of the phantom was 700 ms. The phantom was seated inside a 64-channel head-neck RF receiver coil from the scanner vendor with its face parallel or perpendicular to the main magnetic field. A hard pulse of duration 0.5 ms excited the spins with a range of tip angles, from 10 to 90 in 10 steps, after which the FID signal was acquired for 512 ms. The dwell time of data acquisition was 0.488 ms. Successive FID measurements were temporally separated using a sufficient repetition time (TR ¼ 5 s) to ensure full magnetization recovery. The entire experiment was repeated four to eight times to average out random noise. The complex FID signals from individual receiver channels were summed after subtracting the phase of the first time point in order to account for the channel-dependent receiver phase. The combined FID signal was then zeropadded by a factor of 512, apodized by a rectangular window with a 200 ms cutoff time, and Fourier-transformed to determine the peak frequency. The frequency bin size was 1=ð512 Â 0:512 sÞ ¼ 3.8 mHz.

T FID Measurement
The second phantom was constructed similarly but with smaller inner dimensions of 5 Â 5 Â 0.5 cm 3 . The phantom was filled with 1 mM Gd-doped saline with T 1 ¼ 190 ms. The phantom was scanned in a 9.4T scanner (Biospec 94/30, Bruker, Billerica, MA, USA) with a vendorprovided transmit-receive volume coil with its face parallel or perpendicular to B 0 . The FID measurement and analysis were performed in a way similar to the 3T FID experiment, with the following differences: hardpulse duration ¼ 0.2 ms, FID acquisition time ¼ 400 ms, dwell time ¼ 0.1 ms, TR ¼ 3 s, cutoff time ¼ 60 ms, frequency bin size ¼ 1=ð512 Â 0:4 sÞ ¼ 4.9 mHz.
B 0 Mapping at 3 T A spherical phantom with an inner diameter of 10 cm, filled with water containing 5.11 mM of CuSO 4 and 10 g/ l agar, was scanned for B 0 maps at 3T at two different slice orientations (axial, coronal) and tip angles (a ¼ 30 , 90 ). The measured T 1 of the phantom was 330 ms. A 20-channel head-neck receiver coil from the scanner vendor was used for the scan. B 0 maps were calculated from the voxel-by-voxel phase difference between the two echoes in a double-echo gradient echo sequence with the following parameters: TR ¼ 1 s, TE ¼ 10/30 ms, field of view ¼ 140 mm, grid size ¼ 92 Â 92, slice thickness ¼ 5 mm. For each slice orientation, the B 0 map at a ¼ 30 was subtracted from the one at 90 in order to measure the tip angle dependence of the B 0 map. Because the ratio between TR and T 1 was moderate (TR=T 1 ¼ 3:03), in this experiment Equation [5] is the more appropriate expression for the tip angle-dependent nuclear susceptibility. However, for the two particular tip angles used, Equation [3] and Equation [5] do not give significantly

Finite Thickness Simulation
The finite thickness of the phantoms alters the geometric coefficients (À2/3 and 1/3) in Equations [9] and [10]. This effect was simulated by summing up magnetic dipolar fields from all the voxels in the phantom magnetized according to the voxel's equilibrium nuclear susceptibility. The summation was carried out in the Fourier domain without aliasing artifacts following the method of ref. (15), using MatLab (MathWorks, Natick, MA, USA) on a laptop computer with 8 GB of RAM and Intel Core i5 CPU at 1.7 GHz. Table 1 compares the simulated nuclear spin-induced B 0 shifts at the center of the actual and virtual phantoms with those of an ideal, infinitely thin phantom (Eqs. [9] and [10]). At both phantom orientations, the magnitudes of the B 0 shift of the actual phantom are 14% (5 cm phantom) to 18% (15 cm phantom) lower than the ideal case. For the 5 mm-thick slice (virtual phantom), the discrepancy was reduced to 7.8% due to the thickness-towidth ratio being closer to the ideal, infinitely thin case.

Finite-Thickness B 0 Simulation
In what follows, the finite-thickness phantom frequency shifts of Table 1, times cosa, were used for validation of the experimental tip angle-dependent frequency shifts.
Phase Shift in the Time Domain Figure 1 shows examples of the real part of the normalized FID signals from the 15-cm phantom oriented with its face parallel to B 0 at a ¼ 10 and 90 . The signals are shown modulated with a carrier frequency of 21 Hz. Due to large initial transients (signal decay), each FID signal was normalized to unity at t ¼ 25 ms for clarity of display. It is observed that the phase of the a ¼ 90 FID signal gradually lags behind the a ¼ 10 signal, indicating a negative frequency shift. The observed shift in this plot is approximately (À1/20 cycles)/350 ms ¼ À0.14 Hz. This is compared with the theoretical shift due to the nuclear spin polarization, which, using the relevant coefficient in Table 1, equals 0:1383 Â ðcos90 À cos10 Þ ¼ À0:1362 ½Hz: We note that, although the agreement between the theory and the experiment is satisfactory, in this unaveraged illustrative example a high degree of agreement is somewhat fortuitous because T 1 relaxation and magnet drift were not accounted for (see Discussions).

Frequency Shift Versus Tip Angle
The plots of the FID frequency shifts as a function of the tip angle are shown in Figure 2. The data (marker with solid line) represent either 4 (Figs. 2a,b,d) or 8 (Fig. 2c) averages. For both phantom orientations and B 0 strengths, the data closely follow the theoretical prediction (dashed line). The only adjustable parameter for the theoretical curve was the constant offset in the vertical direction, determined to minimize the least-square distance to the data. The data plots clearly show the opposite-sign frequency shifts at different phantom orientations (Eqs. [9] and [10]). The magnitude of the shift was approximately doubled at the perpendicular orientation compared to the parallel orientation, and as expected was roughly tripled at each orientation as B 0 increased from 3 to 9.4T. The scatter of the data was slightly larger at 9.4T compared to 3T. The lower SNR at 9.4T was attributed to poorer filling factor of the phantom in the RF coil and the shorter signal lifetime at 9.4T. , defined as the change of the B 0 map in Hz when the tip angle increased from 30 to 90 on axial and coronal slices in the spherical phantom at 3T. Because in this measurement, B 0 at a ¼ 90 was referenced to B 0 at a ¼ 30 , static effects such as local oscillator frequency offset and B 0 shim errors were effectively cancelled out. No adjustable parameters were used in generating the plots of Figure 3. The sign of the B 0 change, which is positive (negative) for the axial (coronal) plane, agrees with the nuclear susceptibility effect of Equations [9] and [10]. The magnitudes of the shift are somewhat smaller than the theoretical predictions based on the equilibrium shifts for a 5-mm slice ( Table 1, virtual phantom), obtained as À0:3103 Â ðcos90 À cos30 Þ ¼ 0:2687 ½Hz and 0:1551 Âðcos90 À cos30 Þ ¼ À0:1343 ½Hz for an axial and a coronal slice, respectively. Apart from random noise, systematic errors such as imperfect slice profiles, inhomogeneous tip-angle distribution, and magnet drift could be responsible for the observed discrepancy. Repeated measurements indicated reproducibility of df

DISCUSSION
In this paper, we have demonstrated that out of the total (electronic þ nuclear) equilibrium magnetic susceptibility of water, the contribution from the 1 H spins can be separated out through MR frequency shift that depends on the tip angle and slice orientation. Shift in the NMR frequency due to nuclear spin dipolar fields is well documented in hyperpolarized samples such as laser- polarized liquid 129 Xe (16). To our knowledge, a similar effect from 1 H in water in a clinical MRI system has not been published before. The susceptibility and B 0 changes produced by nuclear paramagnetism in water are small, about 1/2,000 compared to those of the electrons. Such differences, however, may become significant as one pushes the quantitative accuracy of tissue susceptibility measurement. At 3T, for example, the 64 ppb swing of the magnetic susceptibility caused by aligned/antialigned proton magnetization in water corresponds to local magnetization changes of m 0 M ¼ 60:012 mT, with corresponding Larmor frequencies of 60:5 Hz. Such frequencies are within the frequency resolution of a typical local frequency shift map obtained by multi-echo gradient echo sequences used in QSM.
As demonstrated in our virtual phantom experiment, nuclear susceptibility-induced resonance shift is affected by transient polarization boundaries created by slice selection. In this respect, a 3D acquisition may be preferred to a 2D acquisition if the nuclear susceptibility effects are to be avoided in, for example, volumetric phase imaging. This is because a 3D acquisition minimizes creating abrupt changes in longitudinal magnetization in the imaged object. For a 2D acquisition, interleaved slice ordering could provide more consistent nuclear spin boundaries across the slices than sequential, gapless ordering. 1 H longitudinal magnetization may affect local frequency shift measurements in two ways: 1. Through proton-density differences among different tissue types and at the tissue-air boundary, in much the same way as the electronic susceptibility affects local B 0 distribution near the tissue boundaries. 2. Through nonuniform tip angles and longitudinal magnetization changes during a sequence execution.
If the proton tip angle remains small, and the magnetization recovery is nearly complete (TR > > T 1 ), the first effect can be made inconsequential by absorbing the equilibrium proton magnetic susceptibility as part of the tissue susceptibility. For example, in water the additional 4 ppb can be viewed as part of the tissue susceptibility that QSM aims to measure. In this case, it is important to note that the low flip-angle condition needs to be applied in a consistent way when comparing results from different studies.
The second effect can be more problematic in experiments involving large tip angles (90 and larger), if susceptibility resolution better than 4 ppb is desired, due to dynamically changing nature of nuclear spin polarization. The sources of nuclear polarization variations (spatial or temporal) include: nonuniform excitation due to B 0 and B 1 inhomogeneity, inversion and saturation pulses, slice profile, and longitudinal relaxation including steady-state formation. In principle, all of these factors can affect phase-based susceptibility measurements. In practice, the magnitude of the susceptibility errors involved is quite small: in most QSM experiments with susceptibility resolution of about 0.01 ppm and utilizing low tip angles, the nuclear paramagnetism effect is likely to be insignificant. Our results should nevertheless be taken as a demonstration of a potential limiting factor for the accuracy of MR-based susceptibility quantification when nuclear effects are ignored.
It is instructive to estimate realistic nuclear susceptibility contrast in brain tissue under representative experimental conditions. In a spoiled gradient echo acquisition, steady state formation results in effective reduction of the longitudinal magnetization, as in Equation [5]. Assuming T 1 ¼ 800 ms for the white matter at 3T and TR¼ 50 ms, the Ernst angle to maximize the white matter signal is a ¼ cos À1 ðe ÀTR=T1 Þ ¼ 20 . Under these conditions, the longitudinal magnetization in the steady state for the white matter (with equilibrium M 0 ¼ 6:5 Â 10 À3 ½A=m) (17) and gray matter (T 1 ¼ 1300 ms, M 0 ¼ 7:3 Â 10 À3 ½A=m) (17) becomes M ss;w ¼ 3:34 Â 10 À3 ½A=m and M ss;g ¼ 2:85 Â 10 À3 ½A=m, respectively. These correspond to the nuclear susceptibility m 0 M ss =B 0 of 1.4 ppb (white) and 1.2 ppb (gray), which are substantially reduced from the thermal equilibrium values and whose difference is likely undetectable.
There are several limitations in our experiments affecting the precision of the measured data. They include magnet drift, temperature effect, T 1 relaxation during the measurement, non-ideal slice profiles, and tip angle nonuniformity. A clinical scanner is known to exhibit measurable magnet B 0 drift over the time scale of minutes. An exemplary B 0 drift time course can be seen in (18), where B 0 drift on the order of 2 Hz/64 MHz/25 min ¼ 1.25 ppb/min can be gleaned, which translates to 0.16 Hz/min at 3T. In our case, during the FID experiments in both the 9.4T and the 3T scanners, we have regularly measured the FID frequencies at a low flip angle (between 10 and 30 ) every 1 to 2 minutes. Between measurements, the system's center frequency was assumed to drift linearly and was determined by interpolation. The estimated drift values were then subtracted from the tip angledependent FID frequency data. In most cases, the observed B 0 drift was much less than 1 ppb/min and was therefore not dominating over the nuclear spin effect studied.
It is conceivable that under certain conditions a tip angle-dependent MR frequency shift can be caused by RF energy-dependent sample temperature rise, coupled with the proton resonance frequency shift of 0.01 ppm/ deg (19). However, such scenario does not explain our experimental result that the tip-angle dependence changed sign with the phantom orientation and thus is unlikely to affect the main conclusion of our work.
In using Equations [3] and [5], we have ignored T 1 relaxation after the excitation. This is justified if the FID measurement time or the TE is much shorter than T 1 . For example, in the 9.4T FID experiment, T 1 ¼ 190 ms and the FID cutoff time for the frequency calculation was 60 ms. For 90 excitation, at 60 ms the longitudinal polarization would have recovered from 0 to 0:27 M 0 , with corresponding change of the tip angle from 90 to 74.3 . Although this is significant, it should be noted that the amplitude of the FID from the whole phantom decayed rapidly in the same time period, so the calculated overall frequency was more heavily weighted by the earlier part of the FID when a remained closer to 90 . We have not attempted systematic analysis of dynamic frequency changes in relation to T 1 relaxation. We note, however, that in Figure 2 the uncertainty in a due to T 1 relaxation was greater for larger tip angles because of the faster rate of recovery when the magnetization is farther away from the equilibrium. This partly motivated us to report the data in Figure 2 referenced to a ¼ 10 , instead of a ¼ 90 where ideally the shift vanishes, and also not to anchor the theory curve to the data at a ¼ 90 .
In the present work, we have not thoroughly investigated the measurement errors due to T 1 relaxation, slice profiles, and tip-angle inhomogeneities. These effects could be best analyzed through a full Bloch simulation including relaxation. Independent and simultaneous probing of the sample temperature and the magnet drift could further isolate systematic errors from the true nuclear spin effects. We note that gradient eddy current is not likely to have affected the tip-angle dependence of dB 0 , provided it was not large enough to affect the tipangle distribution after slice-selective excitation in the virtual phantom experiment. Despite the above limitations, we have observed dB 0 changes as a function of the phantom orientation and tip angle whose magnitude and polarity could be explained reasonably well by the theory based on the nuclear spin paramagnetism.

CONCLUSION
In conclusion, we have demonstrated experimentally that nuclear paramagnetism can produce tip angle-dependent frequency shift in MR experiments. Due to the dynamic change of nuclear longitudinal polarization during a scan acquisition, such shift can potentially affect quantitative accuracy of susceptibility measurements based on MR phase imaging. It is recommended that the 1 H spin longitudinal polarization be noted when high-precision (ppb-level) susceptibility measurements are reported and compared with one another.