Characterization and Correction of Time-Varying Eddy Currents for Diffusion MRI
Characterization and Correction of Time-Varying Eddy Currents for Diffusion MRI
Jake J. Valsamis , Paul I. Dubovan , Corey A. Baron Department of Medical Biophysics, Schulich School of Medicine & Dentistry, Western University Centre for Functional and Metabolic Mapping (CFMM), Robarts Research Institute, London, Ontario, Canada * JV and PD are co-first authors
Corresponding Author
Jake J. Valsamis Email: [email protected]
Keywords diffusion MRI, OGSE, eddy currents, validation, field monitoring, brain
List of Abbreviations dMRI, diffusion MRI; PGSE, pulsed gradient spin echo; Δ eff , effective diffusion time; OGSE, oscillating gradient spin echo; B , static magnetic field; ECC, eddy current compensation; TCEDDY, time-constant eddy current correction algorithm; TVEDDY, time-varying eddy current correction algorithm; MSE, mean squared error; NUFFT, nonuniform FFT; ΔD, diffusion dispersion; FM, field monitoring. Funding Information
Natural Sciences and Engineering Research Council of Canada (NSERC), Grant Number: RGPIN-2018-05448; Canada Research Chairs, Number: 950-231993; Canada First Research Excellence Fund to BrainsCAN; NSERC CGS M program; Ontario Graduate Scholarship program
Word Count
Submitted to NMR in Biomedicine
Abstract Summary
Diffusion MRI (dMRI) generally suffers from eddy currents induced by strong diffusion gradients, which introduce artefacts that can impair subsequent diffusion metric analysis. Existing correction techniques assume that eddy currents do not decay during data acquisition, which is generally effective for traditional Pulsed Gradient Spin Echo (PGSE) diffusion encoding. However, these methods do not necessarily apply to advanced forms of dMRI that require substantial gradient slewing, such as Oscillating Gradient Spin Echo (OGSE). In this work, dynamic field monitoring was used to characterize diffusion gradient induced eddy currents for both PGSE and OGSE to assess the applicability of these assumptions that are typically used. Building from this, the performance of an in-house algorithm (TVEDDY) that models eddy current decay was evaluated by correcting in-vivo PGSE and OGSE brain images and comparing correction quality with conventional correction methods using mean-squared error (MSE) between diffusion weighted images acquired with opposite polarity diffusion gradients. As a ground truth comparison, the spatially varying field dynamics up to third order in space were used in a model based iterative reconstruction to eliminate the eddy current distortions. Time-varying eddy currents were observed for OGSE, which introduced blurring along the phase-encode direction that was not reduced using the traditional approach, but was diminished considerably with TVEDDY and model-based reconstruction. No MSE difference was observed between the conventional approach and TVEDDY for PGSE, but for OGSE TVEDDY resulted in significantly lower MSE than the conventional approach. The field-monitoring-informed model based reconstruction had the lowest MSE for both PGSE and OGSE. Accordingly, this work establishes that it is possible to estimate time-varying eddy currents from the diffusion data itself, which provides substantial image quality improvements for gradient-intensive dMRI acquisitions like OGSE.
Introduction
The diagnostic capabilities of diffusion MRI (dMRI) are largely dependent on the gradients used to encode the behavior of water molecules within biological tissue. Advancements in gradient hardware have the potential to broaden the clinical role of dMRI, however eddy currents induced by these strong gradients corrupt data integrity and degrade diagnostic efficacy. Traditionally, dMRI data is acquired using pulsed gradient spin echo (PGSE), where two gradient lobes with long duration are applied successively. The time between the two gradient lobes and their duration determines the effective diffusion time (Δ eff ), which characterizes how long water molecules are allowed to probe their local environment. In PGSE, the Δ eff is inherently long, which limits sensitivity to varying microstructural length scales. In contrast, oscillating gradient spin echo (OGSE) has been shown to enable short diffusion times through the use of successive short diffusion weighting periods, with diffusion time scaling inversely with OGSE waveform frequency. Variation of the apparent diffusion coefficient (ADC) with OGSE frequency has been observed in the healthy human brain and OGSE has been used in conjunction with PGSE to provide additional insight into acute ischemic stroke. However, the strong oscillating gradients required for successful OGSE acquisitions necessitate the use of rapidly changing, high amplitude gradients. These gradients ultimately induce strong eddy currents that degrade the quality of OGSE data to a greater extent than PGSE.
Accordingly, a robust technique is needed to correct the induced eddy currents that may be unique to this diffusion weighting method. There are a variety of precautionary measures taken to compensate for the effects of eddy currents including shielding of magnetic field coils and pre-emphasis of magnetic field gradients. Despite these precautions, eddy current distortions remain a significant burden to dMRI prompting the need for post-processing techniques that can sufficiently correct eddy current corrupted data. A prominent method for correction of eddy currents induced by the diffusion gradients is the “eddy” algorithm in FSL. FSL eddy and other similar packages correct for eddy currents that are assumed to remain constant with time, which leads to simple affine distortions that are easy to detect and correct. While generally effective for PGSE, the applicability of this assumption has not been verified for advanced forms of dMRI, such as OGSE, which have multiple gradient ramps that run the full span from the negative maximum to positive maximum gradient level. In this work, we utilize dynamic field monitoring (FM) to characterize eddy currents induced by in vivo PGSE and OGSE acquisitions. We then use this characterization to validate a new automatic method that uses the diffusion-weighted data itself to estimate the parameters for an eddy current correction model that, contrary to existing automatic approaches, includes a finite time-constant of eddy current decay. Improved image quality and diffusion dispersion mapping is shown when using this technique compared to FSL eddy. Methods
MR Acquisition Two healthy male subjects were scanned on a Siemens MAGNETOM 7T Plus head-only system (80 mT/m gradient strength and 350 T/m/s slew rate). This study was approved by the Institutional Review Board at Western University and informed consent was obtained before scanning. Diffusion MRI data was acquired in a single scan using standard PGSE (0Hz) and cosine modulated trapezoidal OGSE (40Hz). The remaining parameters were b = 750 s/mm , 4 direction tetrahedral encoding and b = 0 acquisitions with 6 averages each, full Fourier encoding, TE = 124 ms, TR = 11 seconds, field of view = 224 × 224 mm , 2 mm isotropic in-plane resolution, 40 slices (2 mm), and scan time 9.5 minutes. The sequence was prepared similarly to the diffusion dispersion sequence optimization parameters that have been previously determined, except here half of the averages were acquired with negative diffusion gradient polarity. Accordingly, for each diffusion-weighted gradient that was applied, another gradient with the same magnitude was applied in the opposite direction to acquire pairs of images that exhibited opposite polarity eddy current distortions. Eddy Current Characterization and Model-Based Reconstruction An FM system consisting of 16 F NMR field probes connected to an acquisition system (Skope, Zurich Switzerland) was used to monitor the evolution of field dynamics during the EPI readout period for the same diffusion dispersion protocol performed on the healthy patients. The probes were spherically distributed on a scaffold in verified, fixed positions, and the scaffold was positioned near isocenter in the magnet’s frame of reference. Excitation and signal reading of the F probes was handled by inserting a transistor-transistor logic pulse in the MR sequence before the initialization of the readout gradients. The dynamic spatially varying phase data, which is described with 3rd order spherical harmonics, were computed from the FM raw data using vendor provided MATLAB software. Ignoring the effects of static magnetic field (B ) inhomogeneity, the signal measured under the influence of the gradients and diffusion gradient eddy currents is given by: 𝑆 𝑡 = 𝜌 𝑟 exp −𝑖 𝑘 -,/012 𝑡 + 𝑘 -,4556 𝑡 𝑏 - (𝑟) - (1) where S(t) is the signal acquired at time t , 𝜌(𝑟) is the transverse magnetization at voxel position r, 𝑙 counts the 16 spherical harmonic terms (0th to 3rd order), 𝑘 -,/012 is the nominal trajectory applied combined with the 1st order eddy currents induced by the EPI gradients, 𝑘 -,4556 is the phase accrued from diffusion gradient induced eddy currents, and 𝑏 - are the spherical harmonic basis functions. The FM system measures 𝑘 -,/012 (𝑡) + 𝑘 -,4556 (𝑡) , and the 1st order terms along the readout and phase encode directions in a non-diffusion weighted scan were subtracted from the corresponding terms in the diffusion weighted scans to estimate 𝑘 -,4556 . The spherical harmonic data was interpolated to match the number of acquired raw data points along the phase-encode direction, after which each data point was replicated to fill its respective readout line. Accordingly, eddy current variation along the readout direction was ignored, which results in substantial time savings during image reconstructions. The measured zeroth order eddy current term was corrected to account for the eddy current compensation (ECC) that is applied to the raw data by the scanner to counteract predicted phase accumulation using off-line simulation of the pulse sequence. An overview of the reconstruction pipeline is portrayed in Fig. 1a, where the eddy current correction step that occurs after N/2 nyquist ghost correction and ramp sampling regridding along the readout direction is either the model-based method described here or TVEDDY (below). For the model-based approach, the k-space data was corrected for diffusion gradient eddy currents using an iterative conjugate gradient solution for 𝜌(𝑟) using Eq. 1 with 𝑘 -,/012 set to 0 and 𝑘 -,4556 determined as described above. Note that 𝑘 -,/012 is set to 0 because performing Nyquist ghost correction and regridding of ramp sampling along the frequency encode direction implicitly applied 𝑘 -,/012 already at this point of the reconstruction pipeline. All receiver channels were considered separately, and were combined using a SENSE1 reconstruction after eddy current correction. Receiver combination was followed by PCA denoising applied to the complex data, removal of the phase, export to the NIFTI data standard, and application of FSL eddy to correct motion between different diffusion weighted acquisitions, which is not corrected by Eq. 1. Data Driven Correction: TVEDDY Eddy currents lead to the k-space trajectory deviating from the intended trajectory, which primarily manifests as distortions and blurring along the phase encode direction due to its slow k-space traversal. Eddy currents are generally spatially and temporally varying in nature. An approximation of the eddy current fields is given by a Taylor expansion of the spatially and temporally varying terms, and a first order approximation typically accounts for the majority of eddy-current induced phase. As such, eddy current fields produced by diffusion gradients can be accurately defined by a spatially invariant term (B ) and linear field gradient terms (G x , G y , and G z ), all of which are generally time-varying. Each term is often modelled in time as a sum of exponentially decaying functions with discrete time constants. FSL eddy and similar data-driven approaches assume that there is a single, infinite time constant. Our approach, which will be referred to as Time-Varying Eddy (TVEDDY) because it considers the time-varying nature of eddy currents, assumes the eddy current fields approximated by the B and linear gradient eddy currents are each described by: 𝑘 -,4556 = 𝐴 =,- 𝑒𝑥𝑝 −𝑡/𝜏 - + 𝐴 CDE,- (2) where 𝐴 =,- and 𝐴 CDE,- are eddy current amplitudes, 𝜏 - is a finite time-constant, and 𝑙 = {0, 1, 2, 3} corresponds to 𝑏 - (r) = {1, x, y, z}, respectively. It will be assumed henceforth that x, y, and z correspond to the readout, phase encode, and slice-select axes, respectively. For TVEDDY, spherical harmonic terms of 2nd order and higher are ignored. Erroneous bulk distortions (primarily from 𝐴 CDE,- ) and signal smearing (primarily from 𝐴 =,- ) that are the result of the eddy currents induced by diffusion gradients will be equal in magnitude but opposite in direction if gradient polarity is reversed. Accordingly, to improve the ability of our approach to use the data to determine a time constant and two separate eddy current amplitudes, TVEDDY currently requires that each diffusion direction is acquired with both positive and negative net diffusion gradient polarity. The optimization procedure aims to determine the eddy current model parameters that minimize the mean squared error (MSE) between images acquired with positive and negative diffusion gradient polarity via multistart gradient descent (MATLAB). Every iteration, the parameters 𝐴 =,- , 𝐴 CDE,- , and 𝜏 - are used to estimate 𝜌(𝑟) from the raw data non-iteratively using a combination of complex phase shifts and nonuniform FFT (NUFFT) (Fig. 1b). Notably, a non-iterative estimation is possible because only the zeroth and first order terms are considered for TVEDDY, and because the k-space position changes from the diffusion gradient eddy currents are not typically large enough to introduce violations of Nyquist criteria. The eddy current parameters are applied to all slices simultaneously (i.e., all slices was assumed to have the same 𝐴 =,- , 𝐴 CDE,- , and 𝜏 - ), and 𝜏 - was constrained to be between 1.6 ms and 25 ms, since smaller time constants decay to zero within only a few phase-encode lines and larger time constants become close to approximating an infinite time constant (which is already modelled by 𝐴 CDE,- ). During a scan, the subject may move between the acquisition of volumes with opposite polarity diffusion gradients. To overcome this limitation, rigid body motion parameters were jointly determined within the optimization to improve the accuracy of the eddy current parameter estimations. PE shifts were omitted from the motion model to avoid redundancy with the 𝐴 CDE,G term, which also creates bulk PE shifts. Finally, our implementation of TVEDDY allows for manually setting 𝐴 =,- to zero, which results in a constant-time correction similar to FSL eddy; for simplicity, this will be referred to as Time-Constant Eddy (TCEDDY). The primary purpose of TVEDDY is to regrid trajectory errors that occur due to time varying eddy currents, and it does not compensate for subject motion. To correct for motion, FSL eddy was applied after our custom image reconstruction that contains the TVEDDY algorithm. Also, to mitigate the number of interpolations applied to the data, the final iteration of TVEDDY does not apply the motion parameters. Accordingly, our custom reconstruction uses TVEDDY to correct the in-plane eddy current effects without any object-domain interpolation (since NUFFT is used on the complex data), and all object-domain interpolation to correct for motion is applied as a postprocessing step via FSL eddy. FIG. 1. a) Reconstruction Pipeline. TVEDDY, TCEDDY, or FM Correction can be selected to perform eddy current correction. b) TVEDDY Algorithm. Eddy current induced artefacts are corrected starting with k-space data that has been acquired with opposite polarity diffusion gradients. The eddy current parameters 𝐴 =,- , 𝐴 CDE,- , and 𝜏 - , which are applied to the data using phase ramps and a NUFFT (net operator “F”), and translations (T) and rotations (R) are iteratively adjusted until the MSE has converged, where 𝑙 = {0,1,2,3} corresponds to the zeroth and first order eddy currents (Eq. 2). The data sets shown here were synthetically distorted with large eddy current and motion parameters to illustrate the correction process. TVEDDY Validation Validation 1: Algorithm Convergence. The optimization solved by TVEDDY is generally non-convex, and to assess its ability to find the global minimum, 80 pairs of images were synthetically generated from one of the PGSE diffusion weighted 2D slices acquired in vivo . PGSE was chosen for this purpose because time-varying eddy current effects were small for PGSE. Each pair of images was generated assuming eddy currents with opposite polarity and no motion using Eq 2, where the eddy current parameters 𝐴 =,- , 𝐴 CDE,- , and 𝜏 - were randomly generated. After these distortions with known parameters were applied using the TVEDDY forward model (Eq 2; Fig. 1b), noise with a signal to noise ratio of 10 was added to each image. Eddy current correction was performed on a single pair of oppositely distorted slices at a time, and the recovered parameters were correlated against the known parameters to determine if TVEDDY could accurately correct distorted data after the addition of noise. Validation 2: Model Accuracy. Here, synthetic data was generated from a volume of in vivo PGSE slices. Instead of using Eq. 2, eddy current-distorted raw data was generated on a slice-by-slice basis via Eq. 1 using the zeroth and first order (i.e., 𝑙 ={0,1,2,3}) 𝑘 -,4556 data acquired from an OGSE acquisition using FM, which does not necessarily fit the exponential assumption used in Eq. 2. The volume of simulated images was then processed using TCEDDY and TVEDDY, and eddy model parameters were extracted for comparison of the zeroth and first order eddy current deviations predicted by the two models to the ground truth. Accordingly, this simulation tests the accuracy of the exponential assumption used by TVEDDY. Prior to comparison of TVEDDY to the ground truth eddy current data, global background phase shifts common in both diffusion polarities were removed from the FM data. Given that TCEDDY and TVEDDY only capture inverse eddy current effects, any time varying field terms that are consistent between the two acquisitions will not be accounted for. Accordingly, to separate the common background phase effects from eddy-current induced phase profiles in opposite polarity scans, the profiles were decomposed as follows: 𝑘 -,4556 HIJ = 𝑘 - KL 𝑡 + 𝑘 - MN 𝑡 3𝑎 𝑘 -,4556 D4Q = 𝑘 - KL 𝑡 − 𝑘 - MN 𝑡 3𝑏 where 𝑘 -,4556 HIJ and 𝑘 -,4556 D4Q represent the total FM k-coefficient deviations for positive and negative diffusion gradient directions, k [BG] and k [EC] are the predicted distortions resulting from background and eddy current effects respectively. A zeroth order decomposition example is shown in Fig. 2. FIG. 2. FM PGSE zeroth order eddy current phase variations decomposed into a common “background” term that does not change when the diffusion gradient polarity is reversed and an “eddy current” term that reverse polarity along with the applied gradients. Technique Comparisons Reconstructions that used no correction (neither the custom eddy current correction step in Fig. 1a nor FSL eddy), FSL eddy only, TCEDDY, TVEDDY, and model-based reconstruction were compared by computing the MSE in pixel intensity between images acquired with opposing diffusion gradient polarities. For the TCEDDY, TVEDDY, and model-based reconstruction cases, FSL eddy was still applied afterward to correct for subject motion; accordingly, any comparisons investigate the incremental gain of combining these eddy current correction approaches with FSL eddy. For each approach and subject, the MSE was measured for all volume pairs with opposing polarity diffusion gradients, normalized by the mean MSE of each subject’s PGSE volumes without correction, and averaged over all subjects, diffusion directions, and averages (12 total pairs per subject, each for PGSE and OGSE). A mask made using FSL BET was used to exclude voxels outside the brain. Paired t-tests were used to determine statistical significance between the MSE measurements of different techniques (N=24). The Diffusion Dispersion (ΔD) is a dMRI parameter that represents the difference in the apparent diffusion coefficient between PGSE and OGSE according to Eq 4: ∆𝐷 = 1𝑁 U − ln 𝑆 U,C 𝑏 U,C + 1𝑁 UG ln 𝑆 UG,C 𝑏 UG,CX YZ C[\ X Y C[\ (4) where N ω and N ω0 are the number of OGSE and PGSE acquisitions, respectively, 𝑆 U,C and 𝑆 UG,C are the direction-dependent diffusion weighted signals at the frequencies ω > 0 and ω = 0, respectively, and 𝑏 U,C and 𝑏 UG,C are the direction-dependent b-values at frequencies ω > 0 and ω = 0, respectively. ΔD maps were computed after acquired OGSE and PGSE data were corrected with FSL eddy, TCEDDY, TVEDDY and the model-based reconstruction with FM eddy current data.
Results
Generally, PGSE eddy current induced phase increased linearly over the course of the readout, while OGSE exhibited non-linear time dependence (Fig. 3). FIG. 3. FM zeroth and first order 𝑘 -,4556 field dynamics along the 4 diffusion gradient directions for PGSE and OGSE acquisitions. Results shown are after the subtraction of FM b = 0 s/mm dynamics, but prior to the removal of background effects as per Eq. 3. Zeroth order curves (in black) pertain to the left axis and first order curves (in colour) pertain to the right axis, as indicated by the direction of the arrows. The scale is not consistent between panels. TVEDDY Validation 1: Convergence. Synthetically applied eddy current parameters were accurately recovered using the algorithm’s parameter determination capabilities, as indicated in Figure 4 by the good agreement with the expected linear slope of 1 between ground truth and predicted parameters. The slope of the regression line for t , 𝐴 ^ , and 𝐴 CDE is 0.689, 0.939, and 0.915, respectively, which indicates a tendency to underestimate the larger values of t . Bland-Altman analysis revealed little overall bias. FIG. 4. a) TVEDDY’s ability to recover known eddy current parameters. The linear regression line with slope m and standard error SE is shown in comparison to the expected true line with a slope of 1. b) A Bland-Altman plot where the difference between a recovered and applied parameter pair is plotted against the mean of that pair. The solid lines in each plot show the average difference between recovered and applied parameters while the dotted lines represent the interval in which 95% of differences are included. In both panels, t is in units of milliseconds and amplitudes are approximately in units of voxels, where 𝐴 CDE = 1 indicates a distortion that causes a shift of one voxel and 𝐴 ^ = 1 indicates blurring on the order of one voxel. TVEDDY Validation 2: Accuracy. For uniform and relatively linear phase accrual from eddy currents, as typically exhibited by PGSE, the time-constant correction manages to apply phase shifts that reflect the ground-truth phase error well (Fig. 5). However, the temporally varying OGSE eddy currents are not well characterized by the time-constant model TCEDDY. The corrections applied using TVEDDY match well with the linear applied corrections for PGSE, while also providing a closer approximation to the more complex eddy currents produced by OGSE compared to TCEDDY. FIG. 5. Ground truth FM zeroth and first order eddy current terms (i.e., 𝑘 - MN in Eq. 3) compared to the output of TCEDDY and TVEDDY for one of the four diffusion directions acquired. Both approaches perform a good correction for PGSE, but only TVEDDY can at least partially correct for the time varying eddy currents observed for OGSE. Similar results were observed for the other three directions. Figure 6 shows example reconstructed images without correction, after correction with FSL eddy, after correction with TCEDDY and TVEDDY, and after correction using a model-based reconstruction with the full set of FM results. In the original OGSE image, blurring is most noticeable near the edges of the frontal lobe, which are still apparent after correction with FSL eddy and TCEDDY. This artefact is greatly reduced using the time-varying eddy current model, which is qualitatively comparable to the FM correction. FIG. 6. Qualitative comparison between eddy current corrupted OGSE slice without correction and after correction with FSL eddy, TCEDDY, TVEDDY, and FM eddy current modelling. Yellow arrow: blurring artefact corrected by modelling eddy current decay; notably, the blurring artefact is present throughout the slices, but is most conspicuous at the frontal edge of the brain. Decreasing MSE values are observed with increasing model complexity (Fig. 7). Statistically significant changes in MSE (p < 0.01) were observed between each neighboring correction technique in Figure 7, except between TCEDDY and TVEDDY with PGSE. FIG. 7. Quantitative assessment of eddy current correction quality between FSL eddy only (FSL), our time-constant model (TCEDDY), our time-varying model (TVEDDY), and field monitored eddy current correction (FM). For both PGSE and OGSE acquisitions, the MSE was measured between all volumes with inverse distortions of both subjects, normalized by the mean PGSE MSE of each subject without correction, and averaged over all directions. Statistically significant changes in MSE between adjacent correction techniques is indicated by *(p < 0.01). The ΔD maps made using TVEDDY and FM correction show recovered signal voids and improved signal homogeneity of problematic regions in maps made using FSL eddy only, highlighted by yellow circles in both slices of each subject (Fig. 8). FIG. 8. Qualitative comparison between ΔD maps of two slices per subject made with data corrected using FSL eddy only (column 1), data corrected using TVEDDY and FSL eddy (column 2), and data corrected using FM eddy currents and FSL eddy (column 3). Yellow circles highlight problematic regions. The distortions visible in the frontal lobe are due to EPI distortions that stem from B inhomogeneity. Discussion
In this work we have revealed that advanced diffusion weighting approaches may introduce eddy currents that are not compensated well by standard vendor precompensation and introduced a new eddy current correction technique which models eddy current decay. Current methods do not account for this decaying behavior, which primarily leads to image blurring along the phase-encode direction. The high R values and low standard error of the regression line slopes in Fig. 4a), and the average difference between applied and recovered parameters being less than or equal to 0.02 in Fig. 4b), indicate a consistent relationship between recovered and applied parameters. The 𝐴 ^ and 𝐴 CDE regression line slopes close to unity in Fig. 4a) suggest that TVEDDY can accurately recover a wide range of eddy current amplitudes. The analysis of recovered parameters reveals an underestimation between recovered and applied t that tends to increase with increasing t . This finding is likely because at large t the blurring becomes very subtle and difficult to detect. Our current prototype version of TVEDDY correction of a pair of volumes acquired with opposite polarity diffusion gradients requires approximately 1 hour. Thus, a 2D version of TVEDDY that ignores k z was used to assess convergence which enabled a wide range of eddy current parameters to be tested efficiently. Assessing convergence in 3D may result in slightly different performance, however, TVEDDY’s ability to capture eddy current behavior in a volume of data (Fig. 5) suggests it can still converge when k z is included to correct multiple slices. Looking at the FM zeroth order and first order k-coefficient deviations produced by eddy currents, the PGSE sequences generally exhibited linear behaviour that could be modelled using a single exponential term with an infinite decay constant. Conversely, the OGSE scans displayed more complex eddy current evolutions that are not accurately characterized by linear phase accrual. Fig. 5 shows how TCEDDY can accurately model the most basic linear phase accruals but performs poorer with any offsets that deviate from linearity. Since TCEDDY is very similar in approach and performance to FSL eddy, it can be inferred that FSL eddy also struggles with correcting distortions from eddy currents with decay constants on the order of the readout duration. Improvement in distortion correction can be seen with TVEDDY, where including an additional eddy current term with a finite time constant provides additional degrees of freedom to correct for time-varying eddy currents. Overall, this results in correction profiles that more accurately capture the behaviour of applied distortions, which resulted in less blurring in reconstructed OGSE images (Fig. 6) and ΔD maps (Fig. 8), as well as improved similarity between directions acquired with opposite diffusion gradient polarity (Fig. 7). As described in the methods, the displayed FM eddy current-induced profiles in Fig. 5 are not results directly monitored by the system, as the original data was decomposed into two time-varying functions, 𝑘 - KL 𝑡 and 𝑘 - MN 𝑡 , with 𝑘 - MN 𝑡 being the component that reversed polarity with opposite diffusion gradient polarity (Fig. 2). We suspect the background effects arise from temperature increases from the high duty-cycle OGSE diffusion gradients. The temporal zeroth and first order phase accrual originating from the background effects, 𝑘 - KL 𝑡 , were generally linear (Fig. 2), only affecting alignment of images depending on the diffusion gradient direction, while not inducing blurring. Since FSL eddy focuses on correcting linear phase variations, this residual misalignment between non-parallel diffusion directions was likely corrected by FSL in the last step of Fig. 1a. An unexpected result was the oscillatory eddy current phase accumulation for OGSE (Fig. 3). Upon closer investigation, the periodicity of the oscillations is comparable to the period of the OGSE waveforms (~25 ms), suggesting that the gradient oscillations excite oscillatory eddy current modes or that there may be an interaction between residual mechanical vibration of the gradient hardware and eddy current fields. Although TVEDDY provides improved correction by accounting for eddy current decay, it does not capture this oscillating behaviour. Accurately capturing this behaviour would improve correction quality, as seen by the reduced MSE for the field-monitored recon compared to TVEDDY, but the required degrees of freedom for modelling would be high and a robust fitting procedure would likely not be possible. While FM can fully account for these oscillations, it requires dedicated field-monitoring hardware. Alternatively, gradient-impulse response approaches could be implemented, but are also burdened by lengthy calibration scans and it is unclear whether they could predict the oscillatory behaviour observed here. Notably, this type of eddy current behaviour may be vendor or hardware specific. TVEDDY was developed to correct time-varying eddy current distortions that lead to blurring. This requires pairs of images to be acquired with opposing diffusion gradient polarity which will generally come at no cost to scan time as typical OGSE scans already require many signal averages. FSL eddy was used to correct remaining misalignment between diffusion directions that occur from subject motion and the background eddy current terms (Fig. 2). However, it is likely possible to combine the behaviours of TVEDDY and FSL eddy into a single, unified framework. The efficacy of TVEDDY was assessed by determining both quantitatively and qualitatively whether TVEDDY in conjunction with FSL eddy provides any additional benefit to using FSL eddy alone. Given the opposing polarity of distortions in data acquired with opposite polarity diffusion gradients, volume pairs that are not corrupted by eddy current distortions should have very low MSE (Fig. 7). For both PGSE and OGSE volumes, FM correction resulted in the most improved correction while FSL eddy alone performed the worst out of the tested techniques. In OGSE volumes, TVEDDY outperformed TCEDDY for all volumes acquired. In PGSE volumes, there was no significant difference between correction with TVEDDY or TCEDDY, which is consistent with the conventional assumption that modelling eddy current decay is not necessary for PGSE, as the linear phase shifts are adequately handled using a static eddy current model. The reduced MSE in OGSE volumes corrected with TVEDDY and FSL eddy versus FSL eddy alone demonstrates the benefit of accounting for eddy current decay when correcting OGSE data. It is also notable that the MSE of PGSE volumes is typically higher than that of OGSE volumes. OGSE gradients are intrinsically velocity compensated which ensures there is minimal change in cerebrospinal fluid between OGSE volumes. The substantial movement of cerebrospinal fluid between acquisition of PGSE volumes leads to higher MSE than a corresponding OGSE volume pair once eddy currents have been accounted for. The ΔD maps illustrate the impact of eddy current correction on a parametric map that can be computed using combinations of OGSE and PGSE scans. Notably, ΔD is vulnerable to eddy current distortions and blurring that present differently between the PGSE and OGSE scans, which is problematic without a correction like TVEDDY given that only the OGSE scans have significant blurring due to time-varying eddy currents. The images displayed in this work exhibit distortions from B inhomogeneity. Employing parallel imaging would reduce these distortions, but as described by Arbabi et al, the low contrast between PGSE and OGSE volumes necessitates maximizing SNR and also makes ΔD maps particularly sensitive to residual aliasing artefacts that can manifest through parallel imaging. It is likely possible to combine the acquired data with the FM trajectory in a B map informed model-based reconstruction to drastically reduce these distortions, but TVEDDY is designed for use in situations where such hardware is not available and it is thus more relevant to validate it in these realistic scanning conditions. Conclusion
This work presented and validated a new computational method for correcting diffusion gradient eddy currents that outperforms FSL eddy for OGSE acquisitions. The capacity to correct distortions induced by advanced dMRI techniques without the substantial cost of an FM system will promote the development and clinical application of dMRI.
Acknowledgements
The authors would like to thank the Natural Sciences and Engineering Research Council of Canada (NSERC) [grant number RGPIN-2018-05448], Canada Research Chairs [number 950-231993], Canada First Research Excellence Fund to BrainsCAN, the NSERC CGS M program, and the Ontario Graduate Scholarship program for financial support.
References
1. Setsompop K, Kimmlingen R, Eberlein E, et al. Pushing the limits of in vivo diffusion MRI for the Human Connectome Project.
Neuroimage . 2013;80:220-233. 2. Tanner JE, Stejskal EO. Restricted Self - Diffusion of Protons in Colloidal Systems by the Pulsed - Gradient, Spin - Echo Method.
J Chem Phys . 1968;49(4):1768-1777. 3. Baron CA, Beaulieu C. Oscillating gradient spin-echo (OGSE) diffusion tensor imaging of the human brain.
Magn Reson Med . 2014;72(3):726-736. 4. Schachter M, Does MD, Anderson AW, Gore JC. Measurements of restricted diffusion using an oscillating gradient spin-echo sequence.
J Magn Reson . 2000;147(2):232-237. 5. Van AT, Holdsworth SJ, Bammer R. In vivo investigation of restricted diffusion in the human brain with optimized oscillating diffusion gradient encoding.
Magn Reson Med . 2014;71(1):83-94. 6. Arbabi A, Kai J, Khan AR, Baron CA. Diffusion dispersion imaging: Mapping oscillating gradient spin-echo frequency dependence in the human brain.
Magn Reson Med . 2019;(July):1-12.
7. Baron CA, Kate M, Gioia L, et al. Reduction of Diffusion-Weighted Imaging Contrast of Acute Ischemic Stroke at Short Diffusion Times.
Stroke . 2015;46(8):2136-2141. 8. Baron CA, Lebel RM, Wilman AH, Beaulieu C. The effect of concomitant gradient fields on diffusion tensor imaging.
Magn Reson Med . 2012;68(4):1190-1201. 9. Gao F, Shen X, Zhang H, et al. Feasibility of oscillating and pulsed gradient diffusion MRI to assess neonatal hypoxia-ischemia on clinical systems.
J Cereb Blood Flow Metab . Published online August 18, 2020:271678X20944353. 10. Andersson JLR, Sotiropoulos SN. An integrated approach to correction for off-resonance effects and subject movement in diffusion MR imaging.
Neuroimage . 2016;125:1063-1078. 11. Barmet C, De Zanche N, Pruessmann KP. Spatiotemporal magnetic field monitoring for MR.
Magn Reson Med . 2008;60(1):187-197. 12. Wilm BJ, Barmet C, Gross S, et al. Single-shot spiral imaging enabled by an expanded encoding model: Demonstration in diffusion MRI.
Magn Reson Med . 2017;77(1):83-91. 13. Bruder H, Fischer H, Reinfelder HE, Schmitt F. Image reconstruction for echo planar imaging with nonequidistant k-space sampling.
Magn Reson Med . 1992;23(2):311-323. 14. Fessler JA. MODEL-BASED IMAGE RECONSTRUCTION FOR MRI Jeffrey A . Fessler EECS Dept ., University of Michigan.
IEEE Signal Processing . Published online 2010. 15. Sotiropoulos SN, Moeller S, Jbabdi S, et al. Effects of image reconstruction on fiber orientation mapping from multichannel diffusion MRI: reducing the noise floor using SENSE.
Magn Reson Med . 2013;70(6):1682-1689. 16. Schmitt F, Stehling MK, Turner R.
Echo-Planar Imaging: Theory, Technique and Application . Springer Science & Business Media; 2012. 17. Jones DK.
Diffusion MRI . Oxford University Press; 2010. 18. Bodammer N, Kaufmann J, Kanowski M, Tempelmann C. Eddy Current Correction in Diffusion-Weighted Imaging Using Pairs of Images Acquired with Opposite Diffusion Gradient Polarity.
Magn Reson Med . 2004;51(1):188-193. 19. Vannesjo SJ, Haeberlin M, Kasper L, et al. Gradient system characterization by impulse response measurements with a dynamic field camera.
Magn Reson Med . 2013;69(2):583-593.. 2013;69(2):583-593.