{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Calculate Path Length Map\n\nWe show how to calculate a Path Length Map for Anisotropic Radiation Therapy\nContours given a set of streamlines and a region of interest (ROI).\nThe Path Length Map is a volume in which each voxel's value is the shortest\ndistance along a streamline to a given region of interest (ROI). This map can\nbe used to anisotropically modify radiation therapy treatment contours based\non a tractography model of the local white matter anatomy, as described in\n[Jordan_2018_plm]_, by executing this tutorial with the gross tumor volume\n(GTV) as the ROI.\n\nNOTE: The background value is set to -1 by default\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from dipy.core.gradients import gradient_table\nfrom dipy.data import get_fnames, default_sphere\nfrom dipy.io.gradients import read_bvals_bvecs\nfrom dipy.io.image import load_nifti_data, load_nifti, save_nifti\nfrom dipy.reconst.shm import CsaOdfModel\nfrom dipy.direction import peaks_from_model\nfrom dipy.tracking.stopping_criterion import ThresholdStoppingCriterion\nfrom dipy.tracking import utils\nfrom dipy.tracking.local_tracking import LocalTracking\nfrom dipy.tracking.streamline import Streamlines\nfrom dipy.viz import actor, window, colormap as cmap\nfrom dipy.tracking.utils import path_length\nimport numpy as np\nimport matplotlib as mpl\nfrom mpl_toolkits.axes_grid1 import AxesGrid" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First, we need to generate some streamlines and visualize. For a more complete\ndescription of these steps, please refer to the\n`example_probabilistic_fiber_tracking` and the Visualization of ROI\nSurface Rendered with Streamlines Tutorials.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "hardi_fname, hardi_bval_fname, hardi_bvec_fname = get_fnames('stanford_hardi')\nlabel_fname = get_fnames('stanford_labels')\n\ndata, affine, hardi_img = load_nifti(hardi_fname, return_img=True)\nlabels = load_nifti_data(label_fname)\n\nbvals, bvecs = read_bvals_bvecs(hardi_bval_fname, hardi_bvec_fname)\ngtab = gradient_table(bvals, bvecs)\n\nwhite_matter = (labels == 1) | (labels == 2)\n\ncsa_model = CsaOdfModel(gtab, sh_order=6)\ncsa_peaks = peaks_from_model(csa_model, data, default_sphere,\n relative_peak_threshold=.8,\n min_separation_angle=45,\n mask=white_matter)\n\nstopping_criterion = ThresholdStoppingCriterion(csa_peaks.gfa, .25)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will use an anatomically-based corpus callosum ROI as our seed mask to\ndemonstrate the method. In practice, this corpus callosum mask (labels == 2)\nshould be replaced with the desired ROI mask (e.g. gross tumor volume (GTV),\nlesion mask, or electrode mask).\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Make a corpus callosum seed mask for tracking\nseed_mask = labels == 2\nseeds = utils.seeds_from_mask(seed_mask, affine, density=[1, 1, 1])\n\n# Make a streamline bundle model of the corpus callosum ROI connectivity\nstreamlines = LocalTracking(csa_peaks, stopping_criterion, seeds, affine,\n step_size=2)\nstreamlines = Streamlines(streamlines)\n\n# Visualize the streamlines and the Path Length Map base ROI\n# (in this case also the seed ROI)\n\nstreamlines_actor = actor.line(streamlines, cmap.line_colors(streamlines))\nsurface_opacity = 0.5\nsurface_color = [0, 1, 1]\nseedroi_actor = actor.contour_from_roi(seed_mask, affine,\n surface_color, surface_opacity)\n\nscene = window.Scene()\nscene.add(streamlines_actor)\nscene.add(seedroi_actor)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If you set interactive to True (below), the scene will pop up in an\ninteractive window.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "interactive = False\nif interactive:\n window.show(scene)\n\nwindow.record(scene, n_frames=1, out_path='plm_roi_sls.png',\n size=(800, 800))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ ".. figure:: plm_roi_sls.png\n :align: center\n\n **A top view of corpus callosum streamlines with the blue transparent ROI in\n the center**.\n\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we calculate the Path Length Map using the corpus callosum streamline\nbundle and corpus callosum ROI.\n\nNOTE: the mask used to seed the tracking does not have to be the Path\nLength Map base ROI, as we do here, but it often makes sense for them to be the\nsame ROI if we want a map of the whole brain's distance back to our ROI.\n(e.g. we could test a hypothesis about the motor system by making a streamline\nbundle model of the cortico-spinal track (CST) and input a lesion mask as our\nPath Length Map base ROI to restrict the analysis to the CST)\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "path_length_map_base_roi = seed_mask\n\n# calculate the WMPL\nwmpl = path_length(streamlines, affine, path_length_map_base_roi)\n\n# save the WMPL as a nifti\nsave_nifti('example_cc_path_length_map.nii.gz', wmpl.astype(np.float32),\n affine)\n\n# get the T1 to show anatomical context of the WMPL\nt1_fname = get_fnames('stanford_t1')\nt1_data = load_nifti_data(t1_fname)\n\n\nfig = mpl.pyplot.figure()\nfig.subplots_adjust(left=0.05, right=0.95)\nax = AxesGrid(fig, 111,\n nrows_ncols=(1, 3),\n cbar_location=\"right\",\n cbar_mode=\"single\",\n cbar_size=\"10%\",\n cbar_pad=\"5%\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will mask our WMPL to ignore values less than zero because negative numbers\nindicate no path back to the ROI was found in the provided streamlines\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "wmpl_show = np.ma.masked_where(wmpl < 0, wmpl)\n\nslx, sly, slz = [60, 50, 35]\nax[0].matshow(np.rot90(t1_data[:, slx, :]), cmap=mpl.cm.bone)\nim = ax[0].matshow(np.rot90(wmpl_show[:, slx, :]),\n cmap=mpl.cm.cool, vmin=0, vmax=80)\n\nax[1].matshow(np.rot90(t1_data[:, sly, :]), cmap=mpl.cm.bone)\nim = ax[1].matshow(np.rot90(wmpl_show[:, sly, :]), cmap=mpl.cm.cool,\n vmin=0, vmax=80)\n\nax[2].matshow(np.rot90(t1_data[:, slz, :]), cmap=mpl.cm.bone)\nim = ax[2].matshow(np.rot90(wmpl_show[:, slz, :]),\n cmap=mpl.cm.cool, vmin=0, vmax=80)\n\nax.cbar_axes[0].colorbar(im)\nfor lax in ax:\n lax.set_xticks([])\n lax.set_yticks([])\nfig.savefig(\"Path_Length_Map.png\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ ".. figure:: Path_Length_Map.png\n :align: center\n\n **Path Length Map showing the shortest distance, along a streamline,\n from the corpus callosum ROI with the background set to -1**.\n\n## References\n\n.. [Jordan_2018_plm] Jordan K. et al., \"An Open-Source Tool for Anisotropic\nRadiation Therapy Planning in Neuro-oncology Using DW-MRI Tractography\",\nPREPRINT (biorxiv), 2018.\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 }