Correction of Gradient Nonlinearity Bias in Quantitative Diffusion Parameters of Renal Tissue with Intra Voxel Incoherent Motion.

Spatially non-uniform diffusion weighting bias due to gradient nonlinearity (GNL) causes substantial errors in apparent diffusion coefficient (ADC) maps for anatomical regions imaged distant from magnet isocenter. Our previously-described approach allowed effective removal of spatial ADC bias from three orthogonal DWI measurements for mono-exponential media of arbitrary anisotropy. The present work evaluates correction feasibility and performance for quantitative diffusion parameters of the two-component IVIM model for well-perfused and nearly isotropic renal tissue. Sagittal kidney DWI scans of a volunteer were performed on a clinical 3T MRI scanner near isocenter and offset superiorly. Spatially non-uniform diffusion weighting due to GNL resulted both in shift and broadening of perfusion-suppressed ADC histograms for off-center DWI relative to unbiased measurements close to isocenter. Direction-average DW-bias correctors were computed based on the known gradient design provided by vendor. The computed bias maps were empirically confirmed by coronal DWI measurements for an isotropic gel-flood phantom. Both phantom and renal tissue ADC bias for off-center measurements was effectively removed by applying pre-computed 3D correction maps. Comparable ADC accuracy was achieved for corrections of both b-maps and DWI intensities in presence of IVIM perfusion. No significant bias impact was observed for IVIM perfusion fraction.

Introduction

Recent multicenter oncology trials have evaluated quantitative diffusion-weighted imaging (DWI) as a radiological marker of tumor malignancy and response to therapy (1–3). The underlying physical principle for this technology is that oncogenic processes and therapeutic interventions induce regional changes in cellularity of the imaged tissue that can be detected and quantified by mean (isotropic) diffusivity (4, 5). In clinical applications outside of the brain, tissues with low fractional anisotropy are typically assessed by combining 3 orthogonal DWI acquisitions as a function of diffusion gradient weighting, quantified by a b-value to provide a mean diffusivity measure of the tissue. The optimal number of acquired b-values depends on the diffusion model utilized to appropriately characterize tissue diffusivity (6–8).

The default measure of mean diffusivity in current clinical trials is the apparent diffusion coefficient (ADC), which assumes monoexponential signal decay with increasing b-values (4, 9–11). Advanced body oncology trials are designed to allow for multiexponential DWI signal decay, either because of true multicomponent diffusion or perfusion effects, such as intravoxel incoherent motion (IVIM) (7, 8, 12, 13). For IVIM, the typically derived metrics include perfusion-suppressed ADC values and perfusion fraction. Characterization and minimization of technical errors in diffusion metrics is imperative for standardizing DWI measurements so that meaningful and consistent clinical trial results can be obtained to further establish the diagnostic and clinical response value of DWI-derived biomarkers (14, 15).

Recent multisite DWI phantom studies (16–18) have revealed the major sources of technical errors that confounded ADC metrics originating from diffusion weighting (DW) bias caused by spatially dependent deviations from the nominal b-value for off-center anatomic locations. In contrast, excellent reproducibility was demonstrated for ADC measurements acquired at the magnet isocenter (16–18) using a temperature-controlled (ice-water) phantom (variability <3%). Multi-institutional phantom studies conducted across the National Institute of Health-founded Quantitative Imaging Network (19) confirmed that gradient nonlinearity (GNL) is a main contributor to spatial DW bias and variability in off-center ADC measurements across clinical magnetic resonance imaging (MRI) platforms (17). This platform-dependent bias was shown to stem from nonuniform DW that resulted from GNL (20–22) and ranged from −55% to +25% depending on the anatomic location and gradient system design (17). In fact, detected GNL bias accounted for ∼95% of the observed absolute ADC error on a single MRI platform and resulted in an average 20% variation across MRI scanners.

Our previous work has shown that the bulk of the ADC error resulting from spatial GNL bias could be effectively removed for monoexponential diffusion medium of arbitrary anisotropy using 3 orthogonal DWI measurements (21, 23). The proposed ADC correction framework was based on the rotation of the system nonlinearity tensor into the acquired DWI frame, where system GNL tensor characteristics were obtained empirically. In this work, the DW bias correction was tested for IVIM diffusion in (nearly isotropic) renal tissue on a clinical scanner with GNL characteristics provided by the vendor (24, 25). The theoretically predicted DW bias contribution resulting from GNL was validated by the ADC measurements on an isotropic flood phantom. The effect of GNL bias correction via the elimination of error from either DWI intensities or b-values was compared for the perfusion-suppressed ADC and perfusion fraction.

Methodology

The experimental design for this study was tailored to illustrate the feasibility of GNL correction in the presence of IVIM. The renal tissue was chosen as a model IVIM medium (24, 25) for its known high (∼20%) perfusion fraction, relatively low anisotropy (fractional anisotropy <0.3), and because we could select substantial tissue regions of interest (ROIs) with reasonably uniform parametric maps. A large isotropic gel phantom was prepared to empirically confirm spatial GNL characteristics of the scanner within the imaged (torso-sized) volume. The DWI acquisition was optimized to improve the signal-to-noise ratio such that the random measurement errors (DWI SD) in the studied b-value range were lower than the predicted systematic GNL bias at a chosen spatial location. All acquired data were stored in Digital Image Communication in Medicine (DICOM) format (26), and data analysis was automated using routines developed in MATLAB 7 (MathWorks, Natick, MA).

DWI of Gel-Flood Phantom

An isotropic flood-DWI phantom was prepared in a large 300 × 380 × 150-mm3 container using 1.8% weight gelatin (Gelita USA, Sioux City, IA) with 12.5 L of tap water. Coronal DWI scans of the phantom were acquired on a 3-T Philips (Best, the Netherlands) Ingenia MRI scanner with large field of view (FOV = 480 × 480 mm2) using three b values (0, 500, 1000), with DWI directions along primary magnet axis (LAB) and 8 excitations per b value. Other relevant scan parameters were as follows: retention time/echo time (TR/TE) = 4.0/0.066 s; 21 slices; slice thickness/gap = 4/1 mm; in-plane resolution = 5 × 5 mm; and pixel bandwidth = 2686 Hz. The acquired DWI had a signal-to-noise ratio >20 for the highest b (1000).

DWI Acquisition and IVIM Analysis for Renal Tissue

Sagittal DWI scans of an IVIM renal tissue (volunteer consented according to local institutional review board guidelines) were performed on a 3T Philips Ingenia MRI scanner near the isocenter and offset superiorly by 120 mm using a 32-channel torso phase-array coil. Five b values (0, 100, 200, 500, and 800) were acquired for 2 sets of orthogonal DWI directions: ; and . Two distinct gradient direction scenarios were used to empirically test for GNL bias dependence on DWI orientation both for individual directions and the trace. The offsets and angles for sagittal FOV = 375 × 375 mm2 were fixed to 0, whereas the table with the volunteer was physically moved from SI∼0 to SI∼120 mm, keeping the initial landmark and ignoring table position. Other acquisition parameters were as follows: TR = 4.0 s; TE (LAB/OVP) = 0.0937/0.08 s; 11 slices; slice thickness = 5.5 mm; in-plane resolution = 1.67 mm; and pixel bandwidth = 2583 Hz. Eight free-breathing single-shot echo-planar imaging (SS-EPI) dynamics were acquired and stored individually and then coregistered for each slice using home-built 2D full-affine transformation (allowing in-plane scale, shear, rotation, and translation) before averaging for each DWI direction and b-value. Coregistration efficiency (for removing the breathing artifact) was visually evaluated from difference images with respect to b = 0 with and without coregistration.

The perfusion-suppressed ADC component in the presence of IVIM was obtained as a slope of monoexponential fit for b > 100 values (24, 25) of each pixel for log-trace-DWI (direction-average image) intensity ratios relative to the b = 0 image. The perfusion fraction was derived as an intercept of the linear fit. The original spatial ADC bias error off-center (SI∼120 mm) was measured as the deviation from the “true” reference ADC at SI∼0 mm for the same anatomy. The anatomic slice with the most uniform parametric map close to RL∼0 mm was selected as a reference. The slice ROI was defined by manually tracing the kidney border on the T2-weighted (b = 0) image. The ROI edges were defined to exclude edge artifacts that resulted from susceptibility gradients near phantom-container walls or residual misregistration for kidney anatomy. The ADC histograms were binned with the step of 0.01 between 0.5 and 3.5 (×10−3 mm2/s), while for perfusion fraction histograms, a bin size of 0.005 between 0.05 and 0.7 was used. All histograms were smoothed with a 3-point moving-average. The ROI histogram statistics were characterized by mean and SD.

Systematic Bias Prediction and Correction

System nonlinearity tensor (20) was constructed using gradient design (spherical harmonics) coefficients provided by the vendor. The Frobenius norm of the biased matrix normalized to the nominal b value at the isocenter was used to generate (static) bias corrector maps for each gradient direction (21) on a Cartesian grid sampled every 5 mm within a 360-mm FOV. As defined, the correction factors for each pixel were dimensionless and positive, with an allowed range between 0 and 1 for negative GNL and between 1 and 2 for positive GNL ( at the isocenter, where GNL was absent). For experimental data, assuming a nearly isotropic medium, a single direction-average corrector map was constructed for each orthogonal DWI -schema (LAB and OVP) and interpolated according to DICOM header information on imaged volume and resolution. Because of the cylindrical symmetry of the GNL model for the horizontal-bore system, the predicted corrector maps were symmetric around the SI along the AP vs the RL direction (coronal vs sagittal slices).

The corrector was then applied pixel-by-pixel to yield corrected DWI intensities or b maps to derive an unbiased ADC (21, 23). A corrected ADC map was derived from pixel-by-pixel correction (21, 23) of trace-DWI image intensities (S-correction) or of b-maps (b-correction) . (Here, image intensities acquired without DW were denoted as S, whereas S referred to biased DWI intensities.) The effect of correction on log-intensity dependence on the b-value (utilized to derive ADC and perfusion fraction) was different for S versus b correction. Numerically, for each spatial location, S-correction scaled biased log-intensities by inverse correction factor (with unaltered b-values) versus b-correction resulting in a direct multiplication of nominal b-values by correction factors (with preserved intensities). Note also that the b-correction for the isotropic medium was equivalent to the direct correction of the “measured” ADC map by . The correction efficiency was assessed by comparing histogram statistics (mean and SD) before and after correction for the reference ROI. The effect of both correction scenarios on the slope and intercept of linear regression fit was directly visualized for the mean ROI intensities as a function of b > 100.

Results

Figure 1 illustrates how the nonuniformity in DW (b-value) is directly reflected in the measured ADC map for a coronal slice through the isotropic gel phantom (Figure 1A). The apparent b value is symmetrically lower SI (negative GNL) and higher RL (positive GNL) than nominal (isocenter), leading to correspondingly under- or overestimated ADC values (Figure 1A). The color-bar scale in Figure 1A (right) reports on the observed b-value bias range between 0.8 and 1.2 with respect to the nominal value at the isocenter. Within a relatively large ROI (220 × 240 mm2), such nonuniformity resulted in artificial broadening of the ADC histogram that was accompanied by a shift of the mean ADC value (Figure 1B). Knowledge of specific gradient design information allowed for the deterministic prediction of GNL bias and effective removal of nonuniformity in the ADC map (Figure 1C). Effective bias removal was demonstrated by narrowing the ADC histogram down to a measurement uncertainty of ±2.5% and shifting its mean to the isocenter reference value (Figure 1B). The observed ADC bias (Figure 1A), normalized to the isocenter reference value, agreed with the predicted by GNL model for the scanner (see Methodology), with ROI pixel-by-pixel difference falling within 3%. As expected for the cylindrically symmetric GNL model, the bias measured along the SI for coronal phantom orientation (Figure 1A) also agreed with that predicted for the sagittal orientation (Figure 2). Furthermore, the predicted average corrector maps (assuming isotropic medium) were identical for OVP versus LAB DWI orientations.

Figure 1.

(A) Measured ADC nonuniformity (left) for a coronal slice (AP offset of 70 mm) through a flood-phantom is corrected (C) using the vendor-provided GNL model. (B) The wide (2SD∼20%) histogram bias is reduced by correction down to the measurement error (2SD∼5%) in (C). The DW-bias color bar in (A) and top-axis scale in (B) reflect the relative deviation from the nominal b-value at the isocenter (r0). In the absence of GNL bias, this ratio is a unity (b(r)=b(r)).

Figure 2.

(A) Predicted direction-average (LAB or OVP) DW-bias maps across kidney ROIs at scanned SI offsets (z = 0 and 120 mm). The same bias is predicted for LAB and OVP DWI orientations. The gray scale bar shows the range of predicted DW bias with respect to the nominal b-value. Scale value of 1 corresponds to absent GNL bias (uniform DW) predicted for the z = 0 mm reference. (B) Predicted width of DW-bias histogram, color-coded as their corresponding ROIs in (A), changes depending on scan position consistent with the phantom measurements in Figure 1.

Figure 2A illustrates DW bias expected across kidney ROIs at 2 locations measured in this work. The corresponding bias histograms in Figure 2B show how steep GNL along the SI near z = 120 mm results in the broad and shifted ROI histogram compared to the reference at z = 0 mm. Figure 3, A and D, illustrate that the observed bias for the perfusion-suppressed ADC maps at z = 120 mm was consistent with the one predicted from the system GNL model (Figure 2) and virtually independent of DWI direction schema as expected for the nearly isotropic medium. The strong nonuniform bias gradient along the SI was evident across the kidney parametric map at the z = 120-mm superior offset location (Figure 3, A and D). Similar to the phantom data in Figure 1C, the uniformity of the kidney ADC map was nominally restored after GNL b-correction (Figure 3, C and F), closely reproducing the ADC of the reference parametric maps acquired near the isocenter (RL∼0 mm) (Figure 3, B and E).

Figure 3.

Observation (A and D) and b-correction (C and F) of GNL bias in perfusion-suppressed ADC for the parametric maps of the uniform sagittal slice (RL∼0 mm) through the kidney at the 120-mm superior offset (A, D, C, and F) versus isocenter references (B and E) for (A–C) and (D–F) DWI. The common scale of the quantitative parametric ADC maps is given by the color bar.

The efficiency of bias correction is further quantified by the changes observed for ADC histograms of kidney ROIs in Figure 4A. ADC histograms for all ROI pixels of a uniform renal tissue slice were narrower in the vicinity of the isocenter (green) compared with the superior offset for both LAB and OVP DWI orientations (solid and dotted magenta). The steep DW nonuniformity bias across kidney ROIs observed in Figure 3, A and D, resulted in additional (nonbiological) broadening of the corresponding ADC histograms (Figure 4A, solid and dotted magenta). The mean ADC value for the reference histogram (Figure 4A, green) was ∼20% higher than mean ADCs at z = 120 mm either for OVP (dotted magenta) or LAB (solid magenta) DWI schema. With a similar initial bias resulting from GNL and identical corrector maps (Figure 2A), the effect of correction was similar independent of DWI orientation. The example of corrected ROI histograms for LAB DWI is shown in Figure 4A. The original mean bias of 20% (Figure 4, solid magenta) for ADC (z = 120 mm) was reduced to <2.5% after correcting the GNL bias either in DWI intensities (S-correction; Figure 4A, orange) or in b-values (b-correction; Figure 4A, blue), nearly matching the unbiased reference histogram (green) for ADC (z = 0 mm). The bias correction by DWI intensity route apparently slightly overcorrected the ADC histogram, shifting it to somewhat higher values (orange trace) relative to the reference (green trace).

Figure 4.

(A) Perfusion-suppressed ADC histograms are shown for kidney ROIs corresponding to Figure 3A (solid magenta, LAB DWI) and D (dotted magenta, OVP DWI) before correction, and to Figure 3B (green histogram) and C (blue histogram) for the reference (isocenter) map and b-value-corrected map of LAB DWI, respectively. The orange trace in (A) corresponds to a histogram of the corrected ADC map (not shown) achieved via S-correction for LAB DWI intensity. (B) Mean ROI log-intensity signal (symbols) and fit (lines) are plotted as a function of the b-value for LAB DWI at the isocenter (green), before (magenta) and after bias correction via b values (blue pluses), and DWI intensities (orange crosses). Note the horizontal shift of data points with respect to measured (biased, magenta) signal after b-correction vs vertical shift after S-correction. The error bar of the mean reference signal at the highest b-value (green) reflects the 2SD of the corresponding log-intensity within the kidney slice ROI. The figure labels are color-coded to mark correspondence between the histogram ROIs in (A) and mean data values shown in (B).

The slope error of linear fit for ROI-mean log-intensity dependence on the b-value shown in Figure 4B (magenta) is effectively corrected either for b-values (blue pluses) or DWI signal S-intensities (orange crosses). As expected from the corresponding correction formalism, b-correction scaled biased data points (magenta circles) horizontally along the b-axis, whereas S-correction scaled them vertically. Because bias is a multiplicative factor for b-correction, the observed difference between the biased and corrected data location along the b-axis was larger for higher b-values. Both correction methods brought corrected data closer to the reference fit line (Figure 4B, green). The correction efficiency was similar by either method within measurement and fit uncertainty, as is evident from the proximity of the experimental and corrected data points to the fit lines. Although the original GNL bias and bias correction have a noticeable effect on the slope (ADC = 1.9 vs 1.5 × 10−3 mm2/s) of the fit lines in Figure 4B, the effect on their intercept (∼0.23, perfusion fraction) was barely detectible. An insignificant change (<2%) resulting from the bias for the IVIM perfusion fraction (absolute intercept value ∼0.2) was likewise confirmed from observing the corresponding parametric maps and their histograms (without correction) in Figure 5. Both offset (Figure 5A) and reference (Figure 5C) parametric perfusion maps produce nominally overlapping histograms (Figure 5D, magenta and green) without GNL correction. For completeness, we confirmed that, like b-correction (applied directly to ADC maps), S-correction (Figure 5B) did not introduce significant sporadic bias into the perfusion maps and corresponding ROI histograms (Figure 5D).

Figure 5.

IVIM perfusion fraction maps for sagittal slice (RL∼0 mm) acquired with OVP DWI at a superior offset of 120 mm (A) and obtained after S-correction (B) show negligible bias compared to the isocenter reference (C). Similar results are observed for their corresponding ROI histograms in (D). The figure labels are color-coded to mark correspondence between the histogram ROIs in (A–C).

Discussion

The purpose of this study was to illustrate that GNL bias correction is feasible for quantitative diffusion metrics of the IVIM model. Although the current study design was tailored to emphasize the effect of GNL bias, similar effects may be observed routinely for large FOVs typical of body DWI applications. Consistent with previous findings (16–18), GNL causes nonuniform DW that follows a spatially static pattern for a given system independent of the nominal b-value and is readily predictable from the deterministic gradient system design and applied DWI directions (20–22). For the isotropic medium, direction-average bias is independent of the DWI schema. This further simplifies deriving spatial DW correction maps. In the multisite trial setting, the static GNL corrector maps would need to be calculated once for a specific system model and could then be applied to any acquired DWI scan according to its DICOM geometry.

In this study, DW bias observed across a sagittal slice through the kidney followed the trends predicted from the system gradient design for both LAB and OVP DWI orientations. The spatially dependent bias resulted in significant (∼20%) nonuniformity error for the perfusion-suppressed ADC map off-center (z = 120 mm) compared to the same anatomy close to the isocenter. DW nonuniformity as a result of GNL shifted and broadened kidney ADC histograms. For multisite clinical trials, such system and location specific errors would likely lead to significant technical variability that would confound a population-wide analysis of predictive power for obtained ADC metrics if left uncorrected (14, 15). The observed bias was independent of the nominal b-value and the DWI direction schema and mainly affected the ADC derived from the slope of the linear fit for the perfusion-suppressed component of the IVIM diffusion. As expected, the IVIM perfusion fraction derived from the fit intercept was immune to the b-value bias. DW bias correctors reduced the ADC nonuniformity and mean error to <3% for diffusion in the presence of perfusion. Similar to anisotropic monoexponential medium (23), the same correction efficiency (within measurement and fit uncertainty) was achieved for ADC extracted using either corrected trace DWI intensities or corrected b-values. These corrections did not introduce sporadic bias into the perfusion fraction maps.

Although the described correction removes the errors related to GNL-induced bias in DW, it does not address other sources of non-GNL bias error in quantitative diffusion metrics or geometric distortions. Although these bias sources are presumed to have a minor effect on b-values for most horizontal-bore systems (17, 18, 22), including the scanner in this study, they may need to be treated differently depending on severity when observed for a specific scanner model. To reduce the variability of quantitative IVIM diffusion metrics derived in a multisite clinical trial setting, the system-specific correction of GNL bias is best performed for perfusion-suppressed ADC maps before a combined population-wide analysis.

In conclusion, significant DW nonuniformity bias at off-center locations results both in shifting and broadening of perfusion-suppressed ADC histograms for renal tissue. For this well-perfused, nearly isotropic tissue, ADC bias for off-center measurements could be effectively removed by applying direction-average DW bias correctors based on known gradient design specifications. Comparable performance was achieved using corrected DWIs, b-values, or ADC maps independent of DWI orientation. No significant bias impact was observed for IVIM perfusion fraction with or without correction. The demonstrated system-specific correction of GNL bias in DW for off-center anatomy is feasible for clinical trials that utilize quantitative parametric maps based on the IVIM diffusion model.