This example is meant to be an introduction to some of the streamline tools
available in DIPY. Some of the functions covered in this example are
To get started we’ll need to have a set of streamlines to work with. We’ll use
EuDX along with the CsaOdfModel to make some streamlines. Let’s import the
modules and download the data we’ll be using. We’ve loaded an image called Now we can use EuDX to track all of the white matter. To keep things reasonably
fast we use The first of the tracking utilities we’ll cover here is We can use some of DIPY’s visualization tools to display the ROI we targeted
above and all the streamlines that pass through that ROI. The ROI is the yellow
region near the center of the axial image. Once we’ve targeted the corpus callosum ROI, we might want to find out which
regions of the brain are connected by these streamlines. To do this we can use
the We’ve set Because we’re typically only interested in connections between gray matter
regions, and because the label 0 represents background and the labels 1 and 2
represent white matter, we discard the first three rows and columns of the
connectivity matrix. We can now display this matrix using matplotlib. We display it using a log
scale to make small values in the matrix easier to see. In our example track there are more streamlines connecting regions 11 and
54 than any other pair of regions. These labels represent the left and right
superior frontal gyrus respectively. These two regions are large, close
together, have lots of corpus callosum fibers and are easy to track so this
result should not be a surprise to anyone. However, the interpretation of streamline counts can be tricky. The
relationship between the underlying biology and the streamline counts will
depend on several factors, including how the tracking was done, and the correct
way to interpret these kinds of connectivity matrices is still an open question
in the diffusion imaging literature. The next function we’ll demonstrate is Let’s save this density map and the streamlines so that they can be
visualized together. In order to save the streamlines in a “.trk” file we’ll
need to move them to “trackvis space”, or the representation of streamlines
specified by the trackvis Track File format. Footnotes Example source code You can download Connectivity Matrices, ROI Intersections and Density Maps
target
, connectivity_matrix
and density_map
. target
allows one
to filter streamlines that either pass through or do not pass through some
region of the brain, connectivity_matrix
groups and counts streamlines
based on where in the brain they begin and end, and finally, density map counts
the number of streamlines that pass through every voxel of some image.import numpy as np
from scipy.ndimage import binary_dilation
from dipy.core.gradients import gradient_table
from dipy.data import get_fnames
from dipy.io.gradients import read_bvals_bvecs
from dipy.io.image import load_nifti_data, load_nifti, save_nifti
from dipy.direction import peaks
from dipy.reconst import shm
from dipy.tracking import utils
from dipy.tracking.local_tracking import LocalTracking
from dipy.tracking.stopping_criterion import BinaryStoppingCriterion
from dipy.tracking.streamline import Streamlines
hardi_fname, hardi_bval_fname, hardi_bvec_fname = get_fnames('stanford_hardi')
label_fname = get_fnames('stanford_labels')
t1_fname = get_fnames('stanford_t1')
data, affine, hardi_img = load_nifti(hardi_fname, return_img=True)
labels = load_nifti_data(label_fname)
t1_data = load_nifti_data(t1_fname)
bvals, bvecs = read_bvals_bvecs(hardi_bval_fname, hardi_bvec_fname)
gtab = gradient_table(bvals, bvecs)
labels_img
which is a map of tissue types such
that every integer value in the array labels
represents an anatomical
structure or tissue type [1]. For this example, the image was created so that
white matter voxels have values of either 1 or 2. We’ll use
peaks_from_model
to apply the CsaOdfModel
to each white matter voxel
and estimate fiber orientations which we can use for tracking. We will also
dilate this mask by 1 voxel to ensure streamlines reach the grey matter.white_matter = binary_dilation((labels == 1) | (labels == 2))
csamodel = shm.CsaOdfModel(gtab, 6)
csapeaks = peaks.peaks_from_model(model=csamodel,
data=data,
sphere=peaks.default_sphere,
relative_peak_threshold=.8,
min_separation_angle=45,
mask=white_matter)
density=1
which will result in 1 seeds per voxel. The stopping
criterion, determining when the tracking stops, is set to stop when the
tracking exits the white matter.affine = np.eye(4)
seeds = utils.seeds_from_mask(white_matter, affine, density=1)
stopping_criterion = BinaryStoppingCriterion(white_matter)
streamline_generator = LocalTracking(csapeaks, stopping_criterion, seeds,
affine=affine, step_size=0.5)
streamlines = Streamlines(streamline_generator)
target
. This
function takes a set of streamlines and a region of interest (ROI) and returns
only those streamlines that pass through the ROI. The ROI should be an array
such that the voxels that belong to the ROI are True
and all other voxels
are False
(this type of binary array is sometimes called a mask). This
function can also exclude all the streamlines that pass through an ROI by
setting the include
flag to False
. In this example we’ll target the
streamlines of the corpus callosum. Our labels
array has a sagittal slice
of the corpus callosum identified by the label value 2. We’ll create an ROI
mask from that label and create two sets of streamlines, those that intersect
with the ROI and those that don’t.cc_slice = labels == 2
cc_streamlines = utils.target(streamlines, affine, cc_slice)
cc_streamlines = Streamlines(cc_streamlines)
other_streamlines = utils.target(streamlines, affine, cc_slice,
include=False)
other_streamlines = Streamlines(other_streamlines)
assert len(other_streamlines) + len(cc_streamlines) == len(streamlines)
from dipy.viz import window, actor, colormap as cmap
# Enables/disables interactive visualization
interactive = False
# Make display objects
color = cmap.line_colors(cc_streamlines)
cc_streamlines_actor = actor.line(cc_streamlines,
cmap.line_colors(cc_streamlines))
cc_ROI_actor = actor.contour_from_roi(cc_slice, color=(1., 1., 0.),
opacity=0.5)
vol_actor = actor.slicer(t1_data)
vol_actor.display(x=40)
vol_actor2 = vol_actor.copy()
vol_actor2.display(z=35)
# Add display objects to canvas
scene = window.Scene()
scene.add(vol_actor)
scene.add(vol_actor2)
scene.add(cc_streamlines_actor)
scene.add(cc_ROI_actor)
# Save figures
window.record(scene, n_frames=1, out_path='corpuscallosum_axial.png',
size=(800, 800))
if interactive:
window.show(scene)
scene.set_camera(position=[-1, 0, 0], focal_point=[0, 0, 0], view_up=[0, 0, 1])
window.record(scene, n_frames=1, out_path='corpuscallosum_sagittal.png',
size=(800, 800))
if interactive:
window.show(scene)
connectivity_matrix
function. This function takes a set of streamlines
and an array of labels as arguments. It returns the number of streamlines that
start and end at each pair of labels and it can return the streamlines grouped
by their endpoints. Notice that this function only considers the endpoints of
each streamline.M, grouping = utils.connectivity_matrix(cc_streamlines, affine,
labels.astype(np.uint8),
return_mapping=True,
mapping_as_streamlines=True)
M[:3, :] = 0
M[:, :3] = 0
return_mapping
and mapping_as_streamlines
to True
so that
connectivity_matrix
returns all the streamlines in cc_streamlines
grouped by their endpoint.import numpy as np
import matplotlib.pyplot as plt
plt.imshow(np.log1p(M), interpolation='nearest')
plt.savefig("connectivity.png")
density_map
. This function allows
one to represent the spatial distribution of a track by counting the density of
streamlines in each voxel. For example, let’s take the track connecting the
left and right superior frontal gyrus.lr_superiorfrontal_track = grouping[11, 54]
shape = labels.shape
dm = utils.density_map(lr_superiorfrontal_track, affine, shape)
from dipy.io.stateful_tractogram import Space, StatefulTractogram
from dipy.io.streamline import save_trk
# Save density map
save_nifti("lr-superiorfrontal-dm.nii.gz", dm.astype("int16"), affine)
lr_sf_trk = Streamlines(lr_superiorfrontal_track)
# Save streamlines
sft = StatefulTractogram(lr_sf_trk, hardi_img, Space.VOX)
save_trk(sft, "lr-superiorfrontal.trk")
the full source code of this example
. This same script is also included in the dipy source distribution under the doc/examples/
directory.