1799

A Bayesian approach for diffusional kurtosis imaging
Eizou Umezawa1, Daichi Ishihara2, Yasunari kono3, Toshiaki Nakai4, and Ryoichi Kato1

1School of Health Sciences, Fujita Health University, Aichi, Japan, 2Department of Radiology, Nagoya City University Hospital, 3Washimi Orthopedic Surgery, 4Department of Radiology, Fujita Health University Hospital

Synopsis

We propose a Bayesian approach for improving the accuracy of diffusional kurtosis imaging in a small number of data acquisitions. Gaussian-approximated prior distributions are made from primary maximum-likelihood estimation (MLE). The approach was tested using a healthy volunteer data in which a part of signals was replaced with simulated glioma signals. Although the approach does not yield further improvement when MLE has a certain degree of accuracy, the approach has effect to reduce large misestimations and did not cause false shrinkage of dispersions that the sample parameters inherently have. The approach reduces the burden of data acquisition.

Purpose

Diffusional kurtosis is a promising quantity to be biomarkers indicating some diseases state, e.g., brain tumors1. An accurate measurement of the kurtosis requires a larger number of diffusion weighted image (DWI) acquisition. We propose a method for improving the accuracy in a small number of data acquisitions.

Methods

Diffusional Kurtosis Imaging (DKI) Algorithms

Normalized signal of q-space imaging is expressed as $$$S(\vec{q})/S(0)=\int_{-\infty}^{\infty}P(\vec{R})\exp{(i\vec{q}\cdot\vec{R})}d^3R=1-q_iq_jM_{ij}/2+q_iq_jq_kq_lM_{ijkl}/4!+\cdots$$$, where $$$M_{ij}=\int_{-\infty}^{\infty}P(\vec{R})R_iR_jd^3R$$$, $$$M_{ijkl}=\int_{-\infty}^{\infty}P(\vec{R})R_iR_jR_kR_ld^3R$$$, and $$$\sum_i a_ib_i$$$ is written as $$$a_ib_i$$$. This is the characteristic function of the probability density function $$$P(\vec{R})$$$ for diffusion displacement $$$\vec{R}$$$. Using $$$b$$$ value,

$$\frac{S(b,\vec{n})}{S(0)}=1-bn_in_jD_{ij}+\frac{1}{6}b^2n_in_jn_kn_lD_{ijkl}+\cdots\hspace{1.7cm}(1),$$

where $$$D_{ij}=M_{ij}/(2t_d)$$$ and $$$D_{ijkl}=M_{ijkl}/(2t_d)^2$$$ with $$$t_d=\Delta-\delta/3$$$, and $$$\delta$$$/$$$\Delta$$$ are the duration/separation of motion probing gradient (MPG). The logarithm of the characteristic function is the cumulant generating function:

$$\log{\left\{\frac{S(b,\vec{n})}{S(0)}\right\}}=-bn_in_jD_{ij}+\frac{1}{6}b^2n_in_jn_kn_lV_{ijkl}+\cdots\hspace{1cm}(2),$$

where

$$V_{ijkl}=D_{ijkl}-D_{ij}D_{kl}-D_{jk}D_{li}-D_{kl}D_{ij}.$$

The kurtosis along a direction $$$\vec{n}$$$ is $$$K(\vec{n})=n_in_jn_kn_lV_{ijkl}/(n_in_jD_{ij})^2$$$.

We used two kinds of DKI algorithms: the characteristic function method (CFM)2 and the cumulant generating function method (CGFM). The CFM uses Eq. (1) for fitting, and the CGFM, Eq. (2)3. We truncated the series of Eqs. (1) and (2) after the second-order terms. Kurtosis maps were made for mean, axial and radial kurtosis (K) defined in Ref. 4. To address systematic errors contained in DKI parameters we implement certain calibration procedures.

Bayesian Estimation

The likelihood function is

$$p(\vec{B}|\vec{X})\propto\exp{\left\{-\frac{1}{2}(A\vec{X}-\vec{B})^T\Sigma^{-1}(A\vec{X}-\vec{B})\right\}},$$

where $$$\vec{X}$$$ consists of $$$D_{ij}$$$'s and $$$V_{ijkl}$$$'s, $$$\vec{B}$$$ is composed of measured signals for given MPG directions ($$$\vec{n_i}$$$) and b-values ($$$b_i$$$), $$$A$$$ is a matrix composed of the $$$n_{ij}$$$'s and $$$b_i$$$'s, and $$$\Sigma^{-1}$$$ is a diagonal matrix. First, we estimate $$$D_{ij}$$$'s and $$$V_{ijkl}$$$'s by the maximum-likelihood estimation (MLE) using $$$p(\vec{B}|\vec{X})$$$:

$$\vec{X}_0=(A^T\Sigma^{-1} A)^{-1}A^T\Sigma^{-1}\vec{B}.$$

In this step, $$$\Sigma^{-1}$$$ was determined by the propagation law of error. We made a distribution for $$$\vec{X_0}$$$ in a slice by conducting the MLE for each voxel, and read off the median, $$$\langle\vec{X_0}\rangle$$$ and normalized interquartile ranges, $$$\Sigma_p^{1/2}$$$ of the distribution. We then made Gaussian prior distribution:

$$p(\vec{X})\propto \exp{\left\{-\frac{1}{2}(\vec{X}-\langle\vec{X_0}\rangle)^T\Sigma_p^{-1}(\vec{X}-\langle\vec{X_0}\rangle)\right\}}. $$

The Bayesian estimation was implemented using the posterior distribution:

$$p(\vec{X}|\vec{B})\propto p(\vec{B}|\vec{X})p(\vec{X}).$$

Estimates of $$$\vec{X}$$$ are expectation values:

$$\vec{X}_{\rm estimate}=(A^T\Sigma^{-1} A+\Sigma_p^{-1})^{-1}(A^T\Sigma^{-1}\vec{B}-\Sigma_p^{-1}\langle\vec{X_0}\rangle).$$

To calculate $$$\vec{X}_{\rm estimate}$$$, overall factor of $$$\Sigma^{-1}$$$ is necessary. We estimated it from the residuals obtained in MLE at the first step.

Data Acquisition and Simulated Glioma Signals

A healthy male volunteer was scanned on 3T Philips Ingenia system with TE/TR=128/12864 ms, NSA=1, isotropic voxel of 2×2×2 mm3. DWIs were acquired along 32 gradient directions with b=0,200,400,1000,1500,2000,2500 s/mm2. The data were thinned out as appropriate to test our approach. We made simulated DWI signals for gliomas using the bi-exponential (exp.) model and replaced a part of acquired signals with them. For the model parameters we referred Refs. 5 and 1, and made 20 signals for high and low grade gliomas by fluctuating the parameters. True values of the DKI parameters for the sample are shown in Fig. 3.

Results

We only show the results obtained by CGFM for b=0, 1000, 2500 s/mm2 and 16 gradient directions in this abstract. Kurtosis map and the numbers of the outlier-voxels were shown in Figs. 1 and 2, respectively. Imaging qualities were clearly improved by the Bayesian estimation. Medians and interquartile ranges of estimated parameters for the simulated glioma samples are shown in Fig. 4, and the normalized values to the contralateral normal-appearing white matter are shown in Fig. 5. All medians in Figs. 4 and 5 were significantly different between low and high grades (P<0.05). Overlaps of interquartile ranges between low and high grade gliomas tend to reduce in the Bayesian estimation, and especially for normalized parameters, error bars become smaller than those of MLE. The results not shown in this abstract indicated similar trends. In particular, the poorer the results of MLE are, the more Bayesian estimations become effective.

Discussion

The Bayesian approach has effect to reduce large misestimation and does not cause false shrinkage of dispersions that the sample parameters inherently have6. The approach is, however, less able to reduce small estimation error, and thus when MLE has a certain degree of accuracy, the approach does not yield further improvement. The phenomena are attributed to the Gaussian assumption for the prior distributions. If the actual distribution has multi-peak, the Gaussian-approximated distribution becomes broad as contains the peaks inside. The fact explains the above mentioned phenomena. Further improvement may be possible by a more rigorous treatment for the prior distribution with trading off the merit of simplification of the calculation due to the Gaussian approximation.

Conclusions

The Bayesian approach based on the Gaussian approximation is effective to avoid large estimation error in DKI analysis and to reduce the burden of data acquisition.

Acknowledgements

We are grateful to Hiroshi Toyama, MD, of Fujita Health University, and Yutaka Kinomura, RT, of Fujita Health University Hospital for use of facilities and constant encouragement. We thank Mitsuru Masumoto, RT, and Takashi Fukuba, RT, of Fujita Health University Hospital supporting the experiments involving MRI equipment. We also wish to acknowledge Kazuhiro Murayama, MD, Masato Abe, MD, and Masayuki Yamada, PhD, of Fujita Health University for fruitful discussions and comments.

References

1. Van Cauter, Veraart J, Sijbers J, et al. Gliomas: diffusion kurtosis MR imaging in grading. Radiology 2012;263(2):492-501.

2. Umezawa E, Yoshikawa M, Yamaguchi K, et al. q-Space imaging using small magnetic field gradient. Magn Reson Med Sci. 2006;5(4):179-189.

3. Jensen JH, Helpern JA, Ramani A, et al. Diffusional kurtosis imaging: the quantification of non-gaussian water diffusion by means of magnetic resonance imaging. Magn Reson Med. 2005;53(6):1432-1440.

4. Tabesh A, Jensen JH, Ardekani BA, et al. Estimation of tensors and tensor-derived measures in diffusional kurtosis imaging. Magn Reson Med. 2011;65(3):823-836.

5. Maier SE, Sun Y, Mulkern RV. Diffusion imaging of brain tumors. NMR Biomed. 2010;23(7):849-864.

6. Orton MR, Collins DJ, Koh DM. Improved intravoxel incoherent motion analysis of diffusion weighted imaging by data driven Bayesian modeling. Magn Reson Med. 2014;71(1):411-420.

Figures

Figure 1: Kurtosis maps of brain in a healthy male volunteer. Simulated glioma signals are of low grade. The maximum-likelihood estimation (MLE) and Bayesian estimation were used with the cumulant generating function method for b=0, 1000, 2500 s/mm2 and 16 gradient directions.

Figure 2: The numbers of the outlier-voxels in kurtosis maps of Fig. 1. The modified Thompson tau test with significance level of 10-4 was used to detect the outliers.

Figure 3: Medians and interquartile ranges of true DKI parameters for simulated glioma sample. The sample number is 20. The bi-exponential model was used to generate the signals of the low and high grade glioma (LGG and HGG). Units of Mean diffusivity (D) are 10-3 mm2/s.

Figure 4: Medians and interquartile ranges of estimated parameters for simulated low and high grade glioma (LGG and HGG) samples embedded in the actual data. The cumulant generating function method was used for b=0, 1000, 2500 s/mm2 and 16 gradient directions. Units of Mean diffusivity (D) are 10-3 mm2/s.

Figure 5: Medians and interquartile ranges of normalized estimated parameters for simulated low and high grade glioma samples embedded in the actual data. The parameters are normalized to the contralateral normal-appearing white matter. The cumlant generating function method was used for b=0, 1000, 2500 s/mm2 and 16 gradient directions.

Proc. Intl. Soc. Mag. Reson. Med. 25 (2017)
1799