Diffusion Kurtosis Imaging (DKI)[1] aims to approximate the diffusion weighted signal in a more precise manner by accounting for the leading non-Gaussian diffusion effects. Although the kurtosis expression originates from Taylor expansion of the log diffusion signal around $$$b=0$$$[2], normally b-values up to $$$b ∼ 2-3 \mathrm{ms /μ m}^2$$$ are used. Practically, it is necessary to ensure sufficient signal attenuation to quantify the second order term, while still remaining within the appropriate regime of the cumulant expansion without hitting the noise floor. In DKI imaging studies, a wide range of gradient strengths is used, which affects the estimated diffusivity and kurtosis parameters[2,3].

DKI is being increasingly used both as a biomarker, and as a starting point for microstructural modelling[4,5,6]. Hence, there is a need to assess the range of validity of the kurtosis expansion and how the accuracy of the estimated parameters depends on b-values and fitting strategy. A high level of discrepancy between the analytically evaluated kurtosis and the fit results was previously established[7], but here we evaluate the variability and error of the kurtosis parameter for different fitting techniques using diverse data.

The standard DKI approximation reads$$\log\left(\frac{S\left(b,\mathbf{\hat{n}}\right)}{S_{0}}\right)=-b\mathbf{\hat{n}}^{\top}\mathbf{D}\mathbf{\hat{n}}+\frac{b^{2}\bar{D}^{2}}{6}\mathbf{\hat{n}}\mathbf{\hat{n}}^{\top}\mathbf{W}\mathbf{\hat{n}}\mathbf{\hat{n}}^{\top}$$and the mean of the kurtosis tensor may be calculated through[8,9]$$\bar{W}=\frac{1}{5}\left(W_{xxxx}+W_{yyyy}+W_{zzzz}+2W_{xxyy}+2W_{xxzz}+2W_{yyzz}\right)=\frac{1}{5}\mathrm{Tr}\left(W\right)$$To assess the accuracy of mean kurtosis estimation, a microstructural model for which kurtosis can be calculated is necessary. Here the biexponential and a biophysical model[10] were used to provide the ground truths. The latter model is described by $$S=\left(1-\nu\right)S_{h}+\nu\cdot S_{c}$$ where $$$S_{h}\left(b,\mathbf{\hat{n}}\right)=\exp\left(-b\mathbf{\hat{n}}^{\top}\mathbf{D}\mathbf{\hat{n}}\right)$$$ and $$$S_{c}$$$ is the expansion of the signal from a population of cylindrically symmetric diffusion tensors in spherical harmonics. Neurite densities derived from this model were found to be highly correlated with histology findings[11]. In this study we implemented fitting routines that have their data requirements met by shorter scanning protocols[3,12]: weighted linear least squares fit (LLSQ), non-linear least squares fit (NLSQ), linear and non-linear constrained objective function minimisation (CLLSQ,CNLSQ), directionwise fit (signal along each direction fit linearly, and then diffusion and kurtosis tensors are constructed) where higher terms of the cumulant expansion are absent (DW) or present ($$$b^{3}$$$:DW+, $$$b^{3}$$$,$$$b^{4}$$$:DW++).

Human data was acquired in normal volunteers on a Siemens Trio 3T with 32 channel head coil and a double spin echo DW EPI sequence. Data consisted of $$$b=0$$$ image and 14 shells with 33 directions at $$$b=0.2−3.0\mathrm{ms /μm}^2$$$, 2.5mm isotropic, TR/TE=4300/103ms.

Fixed rat brain data acquired as in [7], was fit with the biophysical[10] and biexponential models using NLLS to provide ground truth $$$\bar{W}_\mathrm{th}$$$ computed from the estimated model parameters. These parameters were then used to generate synthetic signals at 12 directions and 14 equispaced b-values so that $$$0<b_{i}\leq b_{\mbox{max}},\,i=1\ldots14$$$($$$b_{\mathrm{max}},\mbox{SNR}$$$ are varied) and were subsequently fit with Eq1. The fit results were then compared to $$$\bar{W}_\mathrm{th}$$$ as a function of $$$b_{\mathrm{max}}$$$, reported in terms of the normalised absolute difference between the experimental $$$\bar{W}_\mathrm{exp}$$$ and $$$\bar{W}_\mathrm{th}$$$, henceforth referred as an error.

The fit to clinical data revealed that kurtosis varies substantially as the maximum b-value changes (Figures 1-2).

We observe that the degree of the variation depends significantly on the fitting techniques. Our findings are confirmed with analysis of ex-vivo spinal cord (data will be presented). These results are in agreement with [7] and show that the results of kurtosis evaluation depend on $$$b_{\mbox{max}}$$$ and algorithm choice.

We also examined a number of fitting techniques using simulated data (Figure 3). The fit with one higher term of cumulant expansion performs best for the grey matter tissue for 'clinical' $$$b_{\mbox{max}}\sim2-2.5\mbox{ms}/\mbox{μm}^{2}$$$. However we should take into consideration that the standard deviation of the error is considerably higher for this type of fit. A high error in evaluation of kurtosis was also observed with the biexponential model as a ground truth, for grey and especially white matter. In this case there is no obvious advantage of any specific fitting technique for $$$b_{\mbox{max}}\sim2-2.5\,\mbox{ms}/\mbox{μm}^{2}$$$.We demonstrated that estimation of kurtosis parameters is very sensitive to the choice of fitting technique and maximum b-value. Even the most successful strategy could not enhance the estimation accuracy above the 20% for the biophysical model in grey matter or biexponential in white matter. The measurement bias for each one of the techniques evidently depends on the microstructural model generating the signal and tissue types, and in particular influences the choice of the maximum b-value. Our results suggest caution when using kurtosis to draw conclusions about tissue structure or estimating microstructural parameters[4,5] and even comparing studies that employ different methods. This does not subtract from the value of the kurtosis model in providing sensitive biomarkers.