{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Crossing invariant fiber response function with FORECAST model\n\nWe show how to obtain a voxel specific response function in the form of\naxially symmetric tensor and the fODF using the FORECAST model from\n[Anderson2005]_ , [Kaden2016]_ and [Zucchelli2017]_.\n\nFirst import the necessary modules:\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import numpy as np\nimport matplotlib.pyplot as plt\nfrom dipy.reconst.forecast import ForecastModel\nfrom dipy.viz import actor, window\nfrom dipy.data import fetch_hbn, get_sphere\nimport nibabel as nib\nimport os.path as op\nfrom dipy.core.gradients import gradient_table" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Download and read the data for this tutorial. Our implementation of FORECAST\nrequires multi-shell `data.fetch_hbn()` provides data that was acquired using\nb-values of 1000 and 2000 as part of the Healthy Brain Network study\n[Alexander2017]_ and was preprocessed and quality controlled in the HBN-POD2\ndataset [RichieHalford2022]_.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "data_path = fetch_hbn([\"NDARAA948VFH\"])[1]\ndwi_path = op.join(\n data_path, \"derivatives\", \"qsiprep\", \"sub-NDARAA948VFH\",\n \"ses-HBNsiteRU\", \"dwi\")\n\nimg = nib.load(op.join(\n dwi_path,\n \"sub-NDARAA948VFH_ses-HBNsiteRU_acq-64dir_space-T1w_desc-preproc_dwi.nii.gz\"))\n\ngtab = gradient_table(\n op.join(dwi_path,\n\"sub-NDARAA948VFH_ses-HBNsiteRU_acq-64dir_space-T1w_desc-preproc_dwi.bval\"),\n op.join(dwi_path,\n\"sub-NDARAA948VFH_ses-HBNsiteRU_acq-64dir_space-T1w_desc-preproc_dwi.bvec\"))\n\ndata = np.asarray(img.dataobj)\n\nmask_img = nib.load(\n op.join(dwi_path,\n\"sub-NDARAA948VFH_ses-HBNsiteRU_acq-64dir_space-T1w_desc-brain_mask.nii.gz\"))\n\nbrain_mask = mask_img.get_fdata()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let us consider only a single slice for the FORECAST fitting\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "data_small = data[:, :, 50:51]\nmask_small = brain_mask[:, :, 50:51]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Instantiate the FORECAST Model.\n\n\"sh_order\" is the spherical harmonics order used for the fODF.\n\ndec_alg is the spherical deconvolution algorithm used for the FORECAST basis\nfitting, in this case we used the Constrained Spherical Deconvolution (CSD)\nalgorithm.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fm = ForecastModel(gtab, sh_order=6, dec_alg='CSD')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Fit the FORECAST to the data\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "f_fit = fm.fit(data_small, mask_small)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Calculate the crossing invariant tensor indices [Kaden2016]_ : the parallel\ndiffusivity, the perpendicular diffusivity, the fractional anisotropy and the\nmean diffusivity.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "d_par = f_fit.dpar\nd_perp = f_fit.dperp\nfa = f_fit.fractional_anisotropy()\nmd = f_fit.mean_diffusivity()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Show the indices and save them in FORECAST_indices.png.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig = plt.figure(figsize=(6, 6))\nax1 = fig.add_subplot(2, 2, 1, title='parallel diffusivity')\nax1.set_axis_off()\nind = ax1.imshow(d_par[:, :, 0].T, interpolation='nearest',\n origin='lower', cmap=plt.cm.gray)\nplt.colorbar(ind, shrink=0.6)\nax2 = fig.add_subplot(2, 2, 2, title='perpendicular diffusivity')\nax2.set_axis_off()\nind = ax2.imshow(d_perp[:, :, 0].T, interpolation='nearest',\n origin='lower', cmap=plt.cm.gray)\nplt.colorbar(ind, shrink=0.6)\nax3 = fig.add_subplot(2, 2, 3, title='fractional anisotropy')\nax3.set_axis_off()\nind = ax3.imshow(fa[:, :, 0].T, interpolation='nearest',\n origin='lower', cmap=plt.cm.gray)\nplt.colorbar(ind, shrink=0.6)\nax4 = fig.add_subplot(2, 2, 4, title='mean diffusivity')\nax4.set_axis_off()\nind = ax4.imshow(md[:, :, 0].T, interpolation='nearest',\n origin='lower', cmap=plt.cm.gray)\nplt.colorbar(ind, shrink=0.6)\nplt.savefig('FORECAST_indices.png', dpi=300, bbox_inches='tight')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ ".. figure:: FORECAST_indices.png\n :align: center\n\n **FORECAST scalar indices**.\n\n\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Load an ODF reconstruction sphere\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "sphere = get_sphere('repulsion724')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Compute the fODFs.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "odf = f_fit.odf(sphere)\nprint('fODF.shape (%d, %d, %d, %d)' % odf.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Display a part of the fODFs\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "odf_actor = actor.odf_slicer(odf[30:60, 30:60, :], sphere=sphere,\n colormap='plasma', scale=0.6)\nscene = window.Scene()\nscene.add(odf_actor)\nwindow.record(scene, out_path='fODFs.png', size=(600, 600), magnification=4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ ".. figure:: fODFs.png\n :align: center\n\n **Fiber Orientation Distribution Functions, in a small ROI of the brain**.\n\n## References\n\n.. [Anderson2005] Anderson A. W., \"Measurement of Fiber Orientation\n Distributions Using High Angular Resolution Diffusion Imaging\", Magnetic\n Resonance in Medicine, 2005.\n\n.. [Kaden2016] Kaden E. et al., \"Quantitative Mapping of the Per-Axon Diffusion\n Coefficients in Brain White Matter\", Magnetic Resonance in Medicine,\n 2016.\n\n.. [Zucchelli2017] Zucchelli E. et al., \"A generalized SMT-based framework for\n Diffusion MRI microstructural model estimation\", MICCAI Workshop on\n Computational DIFFUSION MRI (CDMRI), 2017.\n\n.. [Alexander2017] Alexander LM, Escalera J, Ai L, et al. An open resource for\n transdiagnostic research in pediatric mental health and learning\n disorders. Sci Data. 2017;4:170181.\n\n.. [RichieHalford2022] Richie-Halford A, Cieslak M, Ai L, et al. An\n analysis-ready and quality controlled resource for pediatric brain\n white-matter research. Scientific Data. 2022;9(1):1-27.\n\n\n.. include:: ../links_names.inc\n\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 }