{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Reconstruction of the diffusion signal with the WMTI model\n\nDKI can also be used to derive concrete biophysical parameters by applying\nmicrostructural models to DT and KT estimated from DKI. For instance,\nFieremans et al. [Fierem2011]_ showed that DKI can be used to\nestimate the contribution of hindered and restricted diffusion for well-aligned\nfibers - a model that was later referred to as the white matter tract integrity\nWMTI technique [Fierem2013]_. The two tensors of WMTI can be also\ninterpreted as the influences of intra- and extra-cellular compartments and can\nbe used to estimate the axonal volume fraction and diffusion extra-cellular\ntortuosity. According to previous studies [Fierem2012]_ [Fierem2013]_,\nthese latter measures can be used to distinguish processes of axonal loss from\nprocesses of myelin degeneration.\n\nIn this example, we show how to process a dMRI dataset using\nthe WMTI model.\n\nFirst, we import all relevant modules:\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import numpy as np\nimport matplotlib.pyplot as plt\nimport dipy.reconst.dki as dki\nimport dipy.reconst.dki_micro as dki_micro\nfrom dipy.core.gradients import gradient_table\nfrom dipy.data import get_fnames\nfrom dipy.io.gradients import read_bvals_bvecs\nfrom dipy.io.image import load_nifti\nfrom dipy.segment.mask import median_otsu\nfrom scipy.ndimage import gaussian_filter" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As the standard DKI, WMTI requires multi-shell data, i.e. data acquired from\nmore than one non-zero b-value. Here, we use a fetcher to download a\nmulti-shell dataset which was kindly provided by Hansen and Jespersen\n(more details about the data are provided in their paper [Hansen2016]_).\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fraw, fbval, fbvec, t1_fname = get_fnames('cfin_multib')\n\ndata, affine = load_nifti(fraw)\nbvals, bvecs = read_bvals_bvecs(fbval, fbvec)\ngtab = gradient_table(bvals, bvecs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For comparison, this dataset is pre-processed using the same steps used in the\nexample for reconstructing DKI (see `example_reconst_dki`).\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# data masking\nmaskdata, mask = median_otsu(data, vol_idx=[0, 1], median_radius=4, numpass=2,\n autocrop=False, dilate=1)\n\n# Smoothing\nfwhm = 1.25\ngauss_std = fwhm / np.sqrt(8 * np.log(2))\ndata_smooth = np.zeros(data.shape)\nfor v in range(data.shape[-1]):\n data_smooth[..., v] = gaussian_filter(data[..., v], sigma=gauss_std)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The WMTI model can be defined in DIPY by instantiating the\n'KurtosisMicrostructureModel' object in the following way:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "dki_micro_model = dki_micro.KurtosisMicrostructureModel(gtab)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Before fitting this microstructural model, it is useful to indicate the\nregions in which this model provides meaningful information (i.e. voxels of\nwell-aligned fibers). Following Fieremans et al. [Fieremans2011]_, a simple way\nto select this region is to generate a well-aligned fiber mask based on the\nvalues of diffusion sphericity, planarity and linearity. Here we will follow\nthese selection criteria for a better comparison of our figures with the\noriginal article published by Fieremans et al. [Fieremans2011]_. Nevertheless,\nit is important to note that voxels with well-aligned fibers can be selected\nbased on other approaches such as using predefined regions of interest.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Diffusion Tensor is computed based on the standard DKI model\ndkimodel = dki.DiffusionKurtosisModel(gtab)\ndkifit = dkimodel.fit(data_smooth, mask=mask)\n\n# Initialize well aligned mask with ones\nwell_aligned_mask = np.ones(data.shape[:-1], dtype='bool')\n\n# Diffusion coefficient of linearity (cl) has to be larger than 0.4, thus\n# we exclude voxels with cl < 0.4.\ncl = dkifit.linearity.copy()\nwell_aligned_mask[cl < 0.4] = False\n\n# Diffusion coefficient of planarity (cp) has to be lower than 0.2, thus\n# we exclude voxels with cp > 0.2.\ncp = dkifit.planarity.copy()\nwell_aligned_mask[cp > 0.2] = False\n\n# Diffusion coefficient of sphericity (cs) has to be lower than 0.35, thus\n# we exclude voxels with cs > 0.35.\ncs = dkifit.sphericity.copy()\nwell_aligned_mask[cs > 0.35] = False\n\n# Removing nan associated with background voxels\nwell_aligned_mask[np.isnan(cl)] = False\nwell_aligned_mask[np.isnan(cp)] = False\nwell_aligned_mask[np.isnan(cs)] = False" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Analogous to DKI, the data fit can be done by calling the ``fit`` function of\nthe model's object as follows:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "dki_micro_fit = dki_micro_model.fit(data_smooth, mask=well_aligned_mask)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The KurtosisMicrostructureFit object created by this ``fit`` function can then\nbe used to extract model parameters such as the axonal water fraction and\ndiffusion hindered tortuosity:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "AWF = dki_micro_fit.awf\nTORT = dki_micro_fit.tortuosity" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "These parameters are plotted below on top of the mean kurtosis maps:\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "MK = dkifit.mk(0, 3)\n\naxial_slice = 9\n\nfig1, ax = plt.subplots(1, 2, figsize=(9, 4),\n subplot_kw={'xticks': [], 'yticks': []})\n\nAWF[AWF == 0] = np.nan\nTORT[TORT == 0] = np.nan\n\nax[0].imshow(MK[:, :, axial_slice].T, cmap=plt.cm.gray,\n interpolation='nearest', origin='lower')\nim0 = ax[0].imshow(AWF[:, :, axial_slice].T, cmap=plt.cm.Reds, alpha=0.9,\n vmin=0.3, vmax=0.7, interpolation='nearest', origin='lower')\nfig1.colorbar(im0, ax=ax.flat[0])\n\nax[1].imshow(MK[:, :, axial_slice].T, cmap=plt.cm.gray,\n interpolation='nearest', origin='lower')\nim1 = ax[1].imshow(TORT[:, :, axial_slice].T, cmap=plt.cm.Blues, alpha=0.9,\n vmin=2, vmax=6, interpolation='nearest', origin='lower')\nfig1.colorbar(im1, ax=ax.flat[1])\n\nfig1.savefig('Kurtosis_Microstructural_measures.png')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ ".. figure:: Kurtosis_Microstructural_measures.png\n :align: center\n\n Axonal water fraction (left panel) and tortuosity (right panel) values\n of well-aligned fiber regions overlaid on a top of a mean kurtosis all-brain\n image.\n\n\n## References\n\n.. [Fierem2011] Fieremans E, Jensen JH, Helpern JA (2011). White matter\n characterization with diffusion kurtosis imaging. NeuroImage\n 58: 177-188\n.. [Fierem2012] Fieremans E, Jensen JH, Helpern JA, Kim S, Grossman RI,\n Inglese M, Novikov DS. (2012). Diffusion distinguishes between\n axonal loss and demyelination in brain white matter.\n Proceedings of the 20th Annual Meeting of the International\n Society for Magnetic Resonance Medicine; Melbourne, Australia.\n May 5-11.\n.. [Fierem2013] Fieremans, E., Benitez, A., Jensen, J.H., Falangola, M.F.,\n Tabesh, A., Deardorff, R.L., Spampinato, M.V., Babb, J.S.,\n Novikov, D.S., Ferris, S.H., Helpern, J.A., 2013. Novel\n white matter tract integrity metrics sensitive to Alzheimer\n disease progression. AJNR Am. J. Neuroradiol. 34(11),\n 2105-2112. doi: 10.3174/ajnr.A3553\n.. [Hansen2016] Hansen, B, Jespersen, SN (2016). Data for evaluation of fast\n kurtosis strategies, b-value optimization and exploration of\n diffusion MRI contrast. Scientific Data 3: 160072\n doi:10.1038/sdata.2016.72\n\n.. include:: ../links_names.inc\n\n" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.13" } }, "nbformat": 4, "nbformat_minor": 0 }