Extrapolation-Based References Improve Motion and Eddy-Current Correction of High B-Value DWI Data: Application in Parkinson’s Disease Dementia

Purpose Conventional motion and eddy-current correction, where each diffusion-weighted volume is registered to a non diffusion-weighted reference, suffers from poor accuracy for high b-value data. An alternative approach is to extrapolate reference volumes from low b-value data. We aim to compare the performance of conventional and extrapolation-based correction of diffusional kurtosis imaging (DKI) data, and to demonstrate the impact of the correction approach on group comparison studies. Methods DKI was performed in patients with Parkinson’s disease dementia (PDD), and healthy age-matched controls, using b-values of up to 2750 s/mm2. The accuracy of conventional and extrapolation-based correction methods was investigated. Parameters from DTI and DKI were compared between patients and controls in the cingulum and the anterior thalamic projection tract. Results Conventional correction resulted in systematic registration errors for high b-value data. The extrapolation-based methods did not exhibit such errors, yielding more accurate tractography and up to 50% lower standard deviation in DKI metrics. Statistically significant differences were found between patients and controls when using the extrapolation-based motion correction that were not detected when using the conventional method. Conclusion We recommend that conventional motion and eddy-current correction should be abandoned for high b-value data in favour of more accurate methods using extrapolation-based references.


Introduction
Diffusion tensor imaging (DTI) and diffusion tensor tractography are widely used neuroimaging techniques that yield information on the tissue microstructure [1,2]. These techniques are based on the analysis of data acquired with relatively low diffusion weighting using a b-value of approximately 1000 s/mm 2 . Acquiring data using higher b-values, i.e., with b above 2000 s/ mm 2 enables the use of diffusional kurtosis imaging (DKI) [3], which is a technique that complements DTI and is sensitive to the within-voxel variability of diffusion tensors [4]. High angular resolution diffusion imaging (HARDI) for improved tractography in regions of crossing fibres also benefits from the use of high b-values [5,6]. Moreover, high b-value acquisitions serves as the basis for microstructure imaging, where the data is analysed using biophysical models in order to estimate physiologically important parameters such as the axon diameter or the cell-membrane permeability [7]. From a clinical perspective, high b-value data may enable, for example, improved tumour grading and better characterization of stroke lesions [8,9].
Prior to analysis of high b-value data, it must be corrected for subject motion and eddy-current distortions induced by the diffusion encoding gradients [10,11]. A simple approach to this problem is to register each diffusion-weighted volume to a reference volume acquired without diffusion encoding (b = 0 s/mm 2 ). To account for both motion and eddy-current induced distortions, a full affine transform can be used with 12 degrees of freedom representing translation, rotation, shearing, and scaling along or around the x, y, and z axes [12]. Approaches other than the full affine transform can be used to reduce the complexity of the transform [13], or to fully model the distortion [11,14]. Regardless of the technique in use, the transform parameters are optimised using an objective function, which is commonly based on mutual information [15]. The choice of mutual information over a simpler metric, such as the sum-of-squares of the image difference, is motivated by the markedly different contrast in images acquired with and without diffusion encoding.
We will refer to the method where a full affine transform is used to register diffusionweighted volumes to a non diffusion-weighted reference by optimizing the mutual information as the 'conventional method'. Although mutual information metrics are designed to allow registration of images with different contrasts, the conventional method suffers from poor accuracy [16], mainly due to the large difference in contrast between low and high b-value images. This difference is most pronounced where a rim of CSF surrounds the brain, which is common in elderly subjects and patients with cerebral atrophy. This rim of CSF is clearly visible in non diffusion-weighted images, but is completely attenuated in high b-value images, which results in an inward shift of the apparent brain outline (Fig 1). This shift may induce errors in conventional motion and eddy-current correction of high b-value data, since image registration algorithms tend to match borders. Extrapolation-based correction, where reference volumes with appropriate high b-value contrast are extrapolated from low b-value data, was suggested by Ben-Amitay et al as a potential solution [16]. The benefit of such a method is that low b-value data can be corrected using the conventional method with a precision sufficient for accurate extrapolation of undistorted reference volumes. Ben-Amitay et al used standard DTI analysis on the low b-value data, and adapted the results for the CHARMED model in order to extrapolate to high b-value shells. However, standard DTI performs poorly in regions where tissue is mixed with free water such as CSF [17], which may have a negative impact on the extrapolation.
In this study, we compared conventional and extrapolation-based motion and eddy-current correction methods in a group of patients with Parkinson's disease dementia (PDD) and healthy controls. We employed DTI as the starting point for extrapolations, and compared the CHARMED-based extrapolation with a novel method that employs partial-volume correction for CSF. Comparisons comprised the magnitude of registration errors between low and high b- The challenge: Motion and eddy-current correction of diffusion MRI data in the elderly brain. Top row: Non diffusion-weighted MR-image (left), and diffusion-weighted images, encoded using a b-value of 1000 s/mm 2 (middle) and 2750 s/mm 2 (right) and averaged across multiple directions. Note that the rim of CSF surrounding the anterior part of the brain visualized in the non-diffusion weighted image to the left is absent in the diffusion weighted images. Bottom: Normalized plot of the logarithm of the MR signal as a function of position along the lines indicated in the images at the top. The anterior position where the MR signal turns into background is shifted posteriorly when comparing the zero b-value profile (blue) to the high b-value profile (orange). If the correction employs a non diffusion-encoded reference, the difference in contrast may cause the high b-value data to be erroneously scaled in the antero-posterior direction to fit the low b-value signal outline. doi:10.1371/journal.pone.0141825.g001 Improved Motion and EC Correction by Extrapolation-Based References value data, and we demonstrated the impact of the choice of motion and eddy-current correction method on the accuracy of tractography. Finally, we tested whether the choice of correction method can have an impact on the conclusion from a group comparison by comparing DTI and DKI parameters between PDD patients and healthy controls.

Image acquisition
Multi-shell diffusion MRI data were acquired in a group of patients with PDD (mean age ± SD; 74±7 years, n = 11) and age-matched controls (age 70±4 years, n = 27). Imaging was performed on a Siemens Skyra 3T scanner equipped with a 32 channel head coil. The diffusion MRI protocol comprised 99 DWI volumes, where the choice of b-values and encoding directions was inspired by Poot et al [18]. Three volumes were acquired with b = 0 s/mm 2 , and the remaining volumes were acquired using b-values of 250, 500, 1000 and 2750 s/mm 2 distributed over 6, 6, 20, and 64 directions, respectively. A single-shot spin-echo with EPI read-out was used, with the following settings: TR = 8100 ms, TE = 103 ms, voxel size = 2.3×2.3×2.3 mm 3 , FOV = 294×294 mm 2 , iPAT = 2, and partial Fourier factor = 6/8. The imaging volume comprised 52 contiguous axial slices adjusted to include the whole cerebrum. Total acquisition time was 14 minutes.

Ethics statement
The study, including the consent procedure, was approved by the Regional Ethical Review Board ("Regionala Etikprövningsnämnden i Lund"), Lund, Sweden (number 2011-277), and was performed in accordance with the Helsinki Declaration. Patients and healthy controls were informed of the study content in both oral and written form. Informed consent was taken in written form. All participants were judged to have the capacity to consent, gave written informed consent prior to participation, and underwent cognitive testing and neurological examination by a medical doctor. Patients diagnosed with PDD met the Clinical Diagnostic Criteria for Dementia Associated with PD according to Emre et al [19]. Controls exhibited no parkinsonian symptoms or cognitive deficits.

Conventional motion and eddy-current correction
Conventional motion and eddy-current correction was performed by registering each volume to the first volume acquired with b = 0 s/mm 2 . The registration utilized an affine transformation with 12 degrees of freedom [12], and was performed using ElastiX [20]. As the objective function, we employed mutual information, and the final image was interpolated with a bspline of order 3. The Jacobian determinant of the affine transformation matrix was applied to rescale the image intensity [12].
For completeness, we also performed conventional correction using the 'eddy_correct' tool in FSL, which registers the DWI volumes to the first volume acquired with b = 0 s/mm 2 , using an affine transform implemented in FLIRT (www.fmrib.ox.ac.uk/fsl).

Extrapolation-based correction
Extrapolation-based correction was performed by registering each volume to a reference extrapolated from a model fitted to low b-value data that had first been corrected using conventional correction. The registration employed the same settings as for the conventional correction with ElastiX. Two different extrapolation methods were investigated, as described below. Each set of reference volumes was extrapolated from low b-value data (b 1000 s/mm 2 ). As a starting point for the extrapolation, we used DTI defined by where B is the b-matrix, '< Á, Á >' denotes the tensor inner product, and D init is the initial tensor resulting from the model fit. In previous work, the diffusion tensor model itself has been used to interpolate references for low b-value shells [21]. This approach is, however, not suitable for extrapolation of high b-value references since DTI is not accurate for b-values above approximately 1000 s/mm 2 . Here we investigate the use of the CHARMED model for extrapolation, as proposed by Ben Amitay [16], and a novel model where we correct for partial volume effects from CSF, as CSF is known to corrupt diffusion tensor estimates [22]. Other procedures that have been explored for creating references include, for example, averaging of the diffusion-weighted images in each shell [23] the use of CSF suppression by fluid attenuated inversion recovery (FLAIR) [24], were not included in this investigation. CHARMED-based extrapolation. Extrapolation of references for registration has previously been performed using the CHARMED model [16]. This model is based on two tensors, modelling hindered and restricted diffusion, with signal fractions given by f r and f h (see Appendix 1 for the full model). The value of f r can be approximated from the fractional anisotropy (FA) of D init according to f r = FA/0.78 and f h from f h = 1 -f r [16]. Calculation of the tensor for the restricted component requires an estimate of the axon radius, but how to calculate this parameter from D init is not fully specified in [16]. We therefore assumed the effective radius to be zero, which is a reasonable assumption for most of the white matter and most experimental protocols [7,25]. The restricted tensor D r was thus calculated according to D r = AD r uu T , where AD r = 1.7 μm 2 /ms is the axial diffusivity of the restricted component, and u is a column vector of unity length that was obtained from the principal diffusion direction of D init . The hindered tensor was approximated by the initial tensor, D h = D init , according to [16].
Using this model for extrapolation, we found that the edge between extrapolated signal and background depended on the diffusion encoding direction (compare Fig 3B and 3C, two lower rows), in particular where the CSF-filled subarachnoid space is wide. The dependency was caused by an assumption in the CHARMED model that enforced highly anisotropic diffusion of the slowly diffusing component. This assumption is, however, not valid in grey matter where the voxel-scale diffusion appears nearly isotropic. Another reason for the discrepancy is that in the cortical GM region, the signal from tissue is mixed with that from CSF. This partial-volume effect is not accounted for in the initial DTI fit, and resulted in an extrapolation with lower image intensity than expected.
CSF-corrected extrapolation. To improve the signal prediction, we developed an extrapolation approach that corrected for partial volume effects from CSF and that did not assume high anisotropy by default. This resulted in references more similar to the acquired data (compare Fig 3B and 3D).
Partial volume correction for CSF was performed by separating the initial diffusion tensor D init into two components representing tissue and CSF. Assuming an isotropic tensor for CSF Comparison of image contrasts. A: Acquired without diffusion encoding. B: Acquired using a high b-value. C: Extrapolated using the CHARMED model. D: Extrapolated using the proposed method with CSF correction. The two bottom rows show magnifications of anterior segments of the brain for two diffusion encoding directions. The yellow lines show the outline of the brain at a fixed position across all of the images. Note how the anterior rim of CSF is completely attenuated in the high b-value image (compare A and B). Both extrapolation methods yielded images with gross contrast similar to the acquired images. However, the CHARMED model introduced a shift in the outline of the brain for some encoding directions (compare outlines between bottom two rows in column C). For the proposed method, the outline of the brain did not vary substantially with encoding direction (column D). with a diffusivity of D CSF , we calculated the tissue tensor D tissue from the relation where D init is obtained from fitting the tensor to the low b-value shell, I is the unit tensor, and f CSF is the signal fraction of CSF. Three assumptions were required to solve this equation for D tissue and f CSF . First, we assumed D CSF = 2.1 μm 2 /ms. This value is lower than the expected value of 3 μm 2 /ms, but we found that a lower value made the procedure more robust. Second, we assumed all tissue to have equal mean diffusivity, i.e., that Tr(D tissue )/3 = D tissue = 0.8 μm 2 / ms, where Tr(Á) denotes the trace. This assumption is reasonable in both white and grey matter where the mean diffusivity is in the range 0.8-0.9 μm 2 /ms during adulthood [26,27]. Third, we assumed that f CSF can be approximated from Values of f CSF were then constrained to the interval from zero to unity, after which D tissue was calculated from Eq 2. Note that the values on f CSF obtained through this procedure are not appropriate for any quantitative comparisons [22], and are thus only used in the extrapolation procedure. A second correction of D tissue was introduced to increase the accuracy of the extrapolation; voxels with a mean diffusivity below D min = 0.3 μm 2 /ms were modified according to : To avoid introduction of artificial diffusional anisotropy in grey matter, we based the extrapolation solely on D 0 tissue . To model slow diffusion, we employed the stretched-exponential model for the extrapolation [28]. This model is based on the assumption that there is a distribution of diffusion coefficients present in the voxel, where the width of the distribution is modulated by the parameter α. We believe this model to be useful for extrapolation, since there is pronounced orientation dispersion in most voxels of the brain, and orientation dispersion induce a rather wide distribution of diffusion coefficients [29]. We fixed α to 0.8, which is representative for cortex, but slightly high for white matter [28]. The final model was given by Comparing correction results Both the CHARMED model and the CSF-corrected extrapolation model were used to extrapolate reference volumes for all DWI volumes acquired, which were then used for motion and eddy-current correction (see flowchart in Fig 2).
To determine the accuracy of the motion and distortion correction procedures, we calculated FA maps from both the low and the high b-value part of the data and registered the latter to the former using an affine transform. Our hypothesis was that a bias-free correction should result in zero translation, shearing, rotation, and scaling between the two FA maps. We assumed that the overall contrast of the FA maps, important to the registration, were robust to slight misregistrations within the low and the high b-value subset.
We also analysed the time series of the transform parameters from the corrections, averaged across all subjects. For each volume acquired, the affine registration yields parameters associated with translation, rotation, scaling, and skewing along and around the x, y, and z axes (a i , i = 1-12). Some of the registration parameters depend on the eddy-currents, which in turn depend on the diffusion-encoding gradient, g = (g x , g y , g z ) [30]. In order to assess the degree to which eddy-currents may explain variation in these parameters, they were each modelled as a linear function of g for volume j according to [30] where β 0 , β 1 , β 2 , and β 3 can be obtained by linear regression. We performed the regression and determined the coefficient of determination (R2) that represents the amount of variation between volumes explained by the regression. High values of R2 indicate that variation in transform parameters are explained by eddy-currents. Low values of R2 are expected if eddycurrents do not affect the transform parameter, which should be the case e.g. for rotations around and translations along the x and z-directions [30].

DTI and DKI analysis
Diffusion data were used to calculate FA-maps by fitting the DTI model to the low b-value data (b 1000 s/mm 2 ). Maps of the mean diffusivity (MD) and the mean of the kurtosis tensor (MW) were obtained by fitting the following model to the geometrical average of the MR signal from the different shells [31,32] SðbÞ where S 0 , MD, MW and f CSF were free variables and D CSF was fixed to 3.1 μm 2 /ms. Note that the model in Eq 7 is an extension of the conventional DKI model that also accounts for partial volume effects from CSF. In order to assess the effect of the motion and eddy-current correction method on DKI parameters, we obtained histograms of MD, FA, and MW in regions of interest defined to encompass the white matter, i.e. all voxels with FA above 0.4 and where at least 14 out of 26 neighbouring voxels also fulfilled the FA condition.

Tractography and group comparisons
One of the main objectives of this study was to evaluate the effect of the motion and eddy-current correction method on the outcome of a group comparison. For this purpose, we performed HARDI tractography to extract the cingulum and the anterior thalamic projection tracts, and compared DTI and DKI parameters between the PDD patients and the healthy controls. Our choice of tracts was inspired by a study where DKI of the cingulum bundle was suggested as a means to improve the ability to diagnose Parkinson's disease [33].
After motion and eddy-current correction, the high b-value shell with 64 encoding directions was processed using MRtrix [34]. The analysis employed constrained spherical deconvolution (CSD) to model multiple fibre orientations in each voxel, with a maximum harmonic order of 8 [35]. Probabilistic fibre-tracking was performed using a step-size of 0.5 mm, and terminated in regions where FA was below 0.1.
Three segments of the cingulum were extracted according to the scheme proposed by Jones et al [36], i.e., the subgenual, retrosplenial and parahippocampal segments. The tract segments were cropped so that only well-specified parts were retained. The subgenual streamlines were cropped at the level of the inferior genu and at the middle of the corpus callosum. Retrosplenial tracts were cropped at the middle of the corpus callosum and at the inferior border of the splenium. Parahippocampal segments were cropped at the inferior border of the splenium and at a level corresponding to the middle of the pons. Tracts projecting anteriorly from the thalamus were also extracted, followed by a cropping procedure that removed segments of the tracts entering the thalamus or proceeding anteriorly to the genu of the corpus callosum. In all cases, regions of interest (ROIs) were placed in MNI-space and projected back to native space using warp fields calculated by FNIRT (www.fmrib.ox.ac.uk/fsl) [37], applied to the FA volumes. Seeding was performed in a region fully encompassing the region in which tracks were expected to be found in order to address the path-length dependency of tractography [38].
Quantification of DKI parameters was adapted to minimize the influence of partial volume effects from CSF. We calculated the weighted parameter average for each streamline in the tract, according to where y is the parameter of interest and w i = 1 -(f CSF ) i with f CSF obtained from Eq 7, and i index positions along the tract. This procedure was most important in the thalamic projection tracts, which pass close to the ventricles, and occasionally included streamlines passing through voxels with a high CSF content where the parameter values were less well determined than in pure WM. Values of FA, MD and MW were then averaged across all streamlines, and tested for differences between controls and patients with PDD while correcting for age using linear regression, implemented in MATLAB (MathWorks, Natick, MA). The significance threshold was set to 0.05.

Accuracy of the correction
The accuracy of conventional and extrapolation-based correction was evaluated by comparing FA maps calculated from the low and high b-value data. While the overall contrast in the FA maps did not differ much between these maps, substantial shifts of up to two voxels were observed between the FA map from the low and high b-value data set when using conventional correction. This discrepancy was quantified by registering the FA volumes from the high b-value data set to those from the low b-value data set ( Table 1). The conventional approach resulted in substantial translation and scaling regardless of whether the registration was performed with ElastiX or FSL. For data corrected using extrapolated references, the average transform parameters were close to the expected values of zero for both methods, except for a scaling of approximately 1% along the z-direction when using the CHARMED model for extrapolation. The superior performance of the extrapolation-based methods is also supported by inspection of the translation, rotation, scaling and skew parameters, averaged across all subjects ( Fig  5 and S1 Fig). Nonzero values can be due to eddy-currents, subject motion, or errors caused by the registration. Parameters correcting only for subject motion, i.e., rotations and translations in the x-and z-directions, respectively, were expected to average to zero in the population. By contrast, most parameters from the conventional methods showed substantial shifts from zero in the high b-value volumes (Fig 5 and S1 Fig, red dots). Parameters from the extrapolationbased corrections (Fig 5 and S1 Fig, blue dots), however, followed expectations and were scattered around zero for all except four parameters: rotation around the x-axis, translation along the y-axis and skew along the y and z-axis. The variation in these parameters were generally explained well by regressing the gradient amplitudes onto the data, where R2 was 0.62, 0.44,  . Notably, MW was clearly higher in the data corrected using the extrapolation-based method compared to the conventional method. We believe this discrepancy can be attributed to the erroneous scaling of high b-value volumes by the conventional method (based on Elas-tiX). This effect is confirmed in Fig 7, which shows histograms of MD, MW and FA in the white matter for the conventional method and the two extrapolation-based methods. Again, MW is clearly higher for the two extrapolation-based approaches, and slightly higher for the CSF-based approach than for the CHARMED-based approach. For MD, the CHARMEDbased approach yield lower values than the CSF-based approach.

Effects on tractography
The effect of misregistration of high b-value volumes was examined by visualising a coronal cross section of the retrosplenial part of the cingulate bundle obtained from high b-value Table 1. Accuracy of motion and eddy-current correction methods. The table shows mean (standard deviation) of rotation, translation, scale and skew parameters obtained by registering FA-volumes calculated using high b-values to those calculated using low b-values. Top two sub tables show results from conventional motion and eddy-current correction (C-MEC) using FSL and ElastiX whereas the bottom two subpanels show results from CSF-corrected and CHARMED-based extrapolation-based motion and eddy-current correction (EB-MEC).

C-MEC (FSL)
Axis Improved Motion and EC Correction by Extrapolation-Based References tractography ( Fig 8A). The conventional correction (using ElastiX) produced tracks that exhibited a gross misalignment between the tracks and the expected position of the cingulate bundle ( Fig 8B). The tracts were closer to the expected position for data corrected with the CHARMED-based extrapolation method (Fig 8C), but only the CSF-corrected extrapolationbased correction produced tracks that were reliably confined within the intended structure ( Fig 8D).

Effects on group comparisons
The main objective of this study, to assess the impact of the choice of motion and eddy-current correction method on group comparisons, such as comparing a disease group to controls, was pursued by comparing DKI parameters between PDD patients and healthy controls in the tracts shown in Fig 9. When data were corrected using the conventional method (based on ElastiX), significant differences were detected only for FA within the parahippocampal part of the cingulum (Table 2). However, when data was corrected using the CSF-corrected extrapolation-based method, FA differed significantly between the diagnostic groups in the retrosplenial part of the cingulum, and MW differed significantly in the thalamic projection tracts and in the retrosplenial part of the cingulum ( Table 2, Fig 10). In this aspect, the CHARMED-based method and the CSF-corrected method yielded similar results. The increase in sensitivity can be explained by increased statistical power mediated by a reduction in parameter variability, in particular for FA and MW when using extrapolation-based registration; the standard deviation of MW was reduced by up to 50% in the subgenual part of the cingulum in the PDD patients. The FA and MW were also generally higher in the data corrected using the extrapolation-based method.  doi:10.1371/journal.pone.0141825.g008

Discussion
Conventional motion and eddy-current correction, where each diffusion-weighted volume is registered to a non diffusion-weighted reference, suffers from poor accuracy for high b-value data. We show that motion and eddy-current correction of high b-value data can be improved significantly by extrapolating reference images from low b-value data. In group comparisons, improved correction resulted in reduces variability in the DKI metrics and may thus lead to increased statistical power and thereby a higher sensitivity for detection of group-wise differences.
We also compared two methods for extrapolation, where one was similar to that of Ben-Amitay et al. [16] and based on the CHARMED framework [39], and the other was a novel approach where we corrected for partial volume effects with CSF (Eqs 2-4). The overall results were similar, but correction for CSF yielded registrations that were visible improved ( Table 1, Fig 8 and S2 Fig). The CSF corrections assume prior knowledge of the mean diffusivity of tissue and CSF, which can be readily estimated from data or found in the literature. Although these priors may be different in, for example, tumours, ischemia, or neurodegeneration, we do not expect the method to be overly sensitive to such deviations due to the use of mutual information in the image registration, which makes it robust to local deviations in image contrast.
The main reason why extrapolation-based methods can be expected to be superior to conventional methods is the substantially different contrast in high b-value versus non diffusionweighted volumes (Figs 1 and 3), as also reported by other groups [16,40]. For low b-values, all methods apart from the FSL-based conventional correction performed similarly well (e.g., compare the parameters for low b-values in Fig 5 and S1 Fig), probably because the CSF signal is not fully attenuated and thereby represented in the mutual information. However, benefits of using model-based references also for low b-values have been reported [21]. At higher b-values, CSF is fully attenuated; the MR-signal in CSF at b = 0, 1000 and 3000 s/mm 2 is 100%, 5%, and 0.01%. Corresponding values in WM, assuming encoding perpendicular to the structure, are 100%, 90%, and 70%. Even at b = 28000 s/mm 2 , approximately 20% of the WM signal can remain [41]. This indicates that regardless of whether the image comparison metric takes contrast differences into account or not, adequate alignment of the rim of CSF around the brain is virtually impossible in high b-value volumes where the CSF signal is completely attenuated (Fig 1).
We would like to emphasize that optimal motion and eddy-current correction of DWI data is a question open to discussion and further research. For example, the use of a global affine transformation neglects slice-to-slice variations of the eddy-current distortions [11]. However, affine registrations have been reported to be sufficient for whole-brain correction [30]. Another aspect of correction that we, and others, neglect is the interaction between susceptibility distortions and motion. This may result in additional smoothing because the effective direction in which the susceptibility effects distort the image will depend on rotations of the head. Optimal correction of susceptibility distortion requires that data is acquired with reversed polarity of the phase encoding gradients in the EPI sequence [42]. For such data, and with gradient directions spread out across the whole sphere, a Gaussian-process approach has shown promising results [43]. Unfortunately, standard acquisition protocols yield data unsuitable for this approach since the EPI is acquired with unipolar phase encoding gradients and gradient directions are spread out across a half sphere only (e.g. see instructions for 'eddy' from FSL). Addressing the motion and eddy-current correction of high b-value data acquired with standard protocols is important because high b-value diffusion data are gaining popularity, warranting more research into this topic. Here, we show that the step from conventional to Table 2. Effect of motion and eddy-current correction method on group comparison. Results are presented as mean (standard deviation), and were for the left pair of columns obtained using conventional motion and eddy-current correction (C-MEC) whereas the right pair of columns show values obtained using extrapolation-based motion and eddy-current correction. Significant differences were observed more frequently between the healthy controls and PDD patients when using the extrapolation-based correction. Results were similar but not identical when each hemisphere was compared separately.  Improved Motion and EC Correction by Extrapolation-Based References extrapolation-based correction is of utmost importance for the accuracy of the registration (Table 1), while the benefit of improving the extrapolation method (CSF-corrected vs CHARMED) is moderate but visible (Fig 8 and S2 Fig), and possibly important e.g. for cortical analyses.

C-MEC (ElastiX
Determining the quality of a registration is a non-trivial task, as discussed in the context of non-rigid registrations by Crum et al. [44]. We found that simple visual inspection of the DWI series as a cine loop was not sufficient to detect the subtle shifts reported in Table 1. Inspection of the average registration parameters may raise suspicion of suboptimal registration. For example, the y and z scaling for high b-value image volumes in Fig 5A appear unreasonable, but such an observation cannot singlehandedly be considered as proof of erroneous registration. A quantitative test of the accuracy of the registration can, however, be performed by registering FA volumes calculated from high b-value data to those calculated from low b-values  (Table 1) because of the high contrast similarity in the two images (Fig 4). This comparison revealed a substantial scaling of the high b-value volumes, for data corrected using the conventional method, of up to 2.8%, while the CSF-corrected extrapolation-based correction reduced the erroneous y-scaling to 0.7%. In the z-direction, the CHARMED-based method resulted in a z-scaling of 1.0%, compared to 0.4% for the CSF-corrected method, likely due to the accumulation of CSF in superior part of the subarachnoid space and the smaller field of view in the zdirection compared to the y-direction (compare Fig 2 and S2 Fig). While both extrapolationbased methods generally yielded superior results compared to the conventional methods, some residual errors were observed. For example, high b-value volumes were unexpectedly rotated by approximately 1 degree around the x-axis (Fig 5B, top left plot). A potential explanation is the difficulty to distinguish a weak skew from a small rotation.
Erroneous scaling and rotation may cause secondary errors. True eddy-currents cause scaling of the image volume, which affects the image intensity. This can be corrected for by rescaling the image intensity by the Jacobian determinant of the affine transformation matrix [12]. True rotations are also represented in the transformation matrix, which can be used to correct the b-matrix for improved tractography [12]. Both of these corrections may deteriorate rather than improve the quality of the data, if the eddy-current and motion correction procedure is inaccurate. Erroneous scaling and rotation was most prominent when using the conventional method, but are to a lesser degree also present when using extrapolation-based methods. Scaling errors were also larger for data corrected using the CHARMED-based method compared to the CSF-corrected approach (Table 1). Erroneous image rescaling when using the conventional method may explain the large effect on MW observed in the parameter maps (cf Fig 6D and  6E), since this metric is sensitive to the signal amplitude of high b-value volumes. Table 2 shows an increase in the average MW of approximately 30% and a reduction of its standard deviation of approximately 20% when using the CSF-corrected extrapolation-based method instead of the conventional method. Similar numbers apply for the FA, whereas the effect on MD is negligible. The result for FA and MD is probably due to the misalignment of tracts for the conventional method (Fig 8). Similarly, large effects of tract misalignment may be expected in maps with high contrast within the white matter (FA), but not in maps with low contrast (MD). We expect that the reduced standard deviation in these metrics will yield higher statistical power in group comparisons based on, for example, DKI [45] and filter exchange imaging [46].
The choice of motion and eddy-current correction method can impact the outcome of group comparisons, as we demonstrated by comparing data from patients with PDD to data from healthy elderly controls (Table 2, Fig 10). Elderly subjects tend to have a pronounced rim of CSF around the brain, and the PDD group is expected to exhibit substantial cerebral atrophy [47]. In PDD, alterations in the tissue microstructure have been observed by DTI in the substantia nigra and the putamen [48][49][50][51], and in the genu of the corpus callosum, the superior longitudinal fasciculus and the cingulum [49]. The DTI findings have not been substantially lateralised [49,52], which is why we chose to average across the hemispheres. A large fraction of patients with PDD show, in addition to α-synuclein pathology, amyloid pathology similar to that of Alzheimer's disease (AD), which is associated to structures relevant for episodic memory function such as the hippocampus and the anterior nucleus of the thalamus that connect to the hippocampus via the cingulum [53,54]. Effects of PDD on the anterior thalamic tracts and the retrosplenial part of the cingulum bundle is thus not surprising, but the increase in FA of the PDD patients in the parahippocampal cingulum is unexpected. However, FA is a sensitive but not very specific parameter. Elevated FA has, for example, been reported in AD patients in regions of crossing fibres, due to selective degeneration of one fibre population [55]. In the left parahippocampal cingulum, ageing in an elderly group unexpectedly resulted in elevated FA, which was hypothesized to be the result of reduced sprouting of hippocampal axons [56]. These unexpected findings are the result of inherent limitations of the DTI model, which can be addressed by new imaging protocols that enable imaging of microscopic anisotropy [4,29,57]. Finally, we note that diffusion MRI data should always be interpreted with care, particularly since residual effects of subject motion may yield spurious group differences [58].

Conclusion
Conventional motion and eddy-current correction, where diffusion-weighted volumes are registered to a non diffusion-weighted volume, is inadequate for high b-value data, due to the fundamentally different contrast between the volumes being registered. In particular, the complete attenuation of CSF signal in volumes acquired with b-values above approximately 1500 s/mm 2 renders a challenge for the conventional approach. Registration to extrapolated references that have a contrast similar to the acquired data results in substantial improvements in registration accuracy and improved geometrical correspondence between parameter maps based on low and high b-value data. This is of importance for DKI as well as for high b-value tractography. In data processed with the extrapolation-based method, mean kurtosis was approximately 30% higher than in the data processed with the conventional method, probably due to the erroneous scaling induced by the conventional method. The CSF-corrected extrapolation approach yielded more accurate registration results than the CHARMED-based approach. Motion and eddy-current correction has a significant influence on DKI parameters, and thus quantitative comparisons between different DKI studies should take the choice of correction method into account. When applied to a comparison between patients with PDD and healthy controls, the extrapolation-based method yielded a higher statistical power for detecting the diffusion changes, compared to the conventional method. In summary, extrapolation-based methods are strongly recommended for motion and eddy-current correction of high b-value data.
where n is the direction of the fibre, AD r is the intra-axonal axial diffusivity represented by D || in Eq 7 of [59]. The time-dependent radial diffusivity RD r is given by where R is the radius of the fibre, τ is half the echo time, and D ? is the intrinsic diffusivity in the direction perpendicular to the fibres.
Supporting Information