Particle Filtering Tractography (PFT) [Girard2014] uses tissue partial
volume estimation (PVE) to reconstruct trajectories connecting the gray matter,
and not incorrectly stopping in the white matter or in the corticospinal fluid.
It relies on a stopping criterion that identifies the tissue where the
streamline stopped. If the streamline correctly stopped in the gray matter, the
trajectory is kept. If the streamline incorrectly stopped in the white matter
or in the corticospinal fluid, PFT uses anatomical information to find an
alternative streamline segment to extend the trajectory. When this segment is
found, the tractography continues until the streamline correctly stops in the
gray matter. PFT finds an alternative streamline segment whenever the stopping criterion
returns a position classified as ‘INVALIDPOINT’. This example is an extension of An introduction to the Probabilistic Direction Getter and
Using Various Stopping Criterion for Tractography examples. We begin by loading the
data, fitting a Constrained Spherical Deconvolution (CSD) reconstruction
model, creating the probabilistic direction getter and defining the seeds. Continuous map criterion (CMC) [Girard2014] and Anatomically-constrained
tractography (ACT) [Smith2012] both uses PVEs information from
anatomical images to determine when the tractography stops.
Both stopping criterion use a trilinear interpolation
at the tracking position. CMC stopping criterion uses a probability derived
from the PVE maps to determine if the streamline reaches a ‘valid’ or ‘invalid’
region. ACT uses a fixed threshold on the PVE maps. Both stopping criterion can
be used in conjunction with PFT. In this example, we used CMC. Girard, G., Whittingstall, K., Deriche, R., & Descoteaux, M.
Towards quantitative connectivity analysis: reducing tractography biases.
NeuroImage, 98, 266-278, 2014. Smith, R. E., Tournier, J.-D., Calamante, F., & Connelly, A.
Anatomically-constrained tractography: Improved diffusion MRI
streamlines tractography through effective use of anatomical
information. NeuroImage, 63(3), 1924-1938, 2012. Example source code You can download Particle Filtering Tractography
import numpy as np
from dipy.core.gradients import gradient_table
from dipy.data import get_fnames, default_sphere
from dipy.direction import ProbabilisticDirectionGetter
from dipy.io.gradients import read_bvals_bvecs
from dipy.io.image import load_nifti, load_nifti_data
from dipy.io.stateful_tractogram import Space, StatefulTractogram
from dipy.io.streamline import save_trk
from dipy.reconst.csdeconv import (ConstrainedSphericalDeconvModel,
auto_response_ssst)
from dipy.tracking.local_tracking import (LocalTracking,
ParticleFilteringTracking)
from dipy.tracking.streamline import Streamlines
from dipy.tracking import utils
from dipy.viz import window, actor, colormap, has_fury
# Enables/disables interactive visualization
interactive = False
hardi_fname, hardi_bval_fname, hardi_bvec_fname = get_fnames('stanford_hardi')
label_fname = get_fnames('stanford_labels')
f_pve_csf, f_pve_gm, f_pve_wm = get_fnames('stanford_pve_maps')
data, affine, hardi_img = load_nifti(hardi_fname, return_img=True)
labels = load_nifti_data(label_fname)
bvals, bvecs = read_bvals_bvecs(hardi_bval_fname, hardi_bvec_fname)
gtab = gradient_table(bvals, bvecs)
pve_csf_data = load_nifti_data(f_pve_csf)
pve_gm_data = load_nifti_data(f_pve_gm)
pve_wm_data, _, voxel_size = load_nifti(f_pve_wm, return_voxsize=True)
shape = labels.shape
response, ratio = auto_response_ssst(gtab, data, roi_radii=10, fa_thr=0.7)
csd_model = ConstrainedSphericalDeconvModel(gtab, response)
csd_fit = csd_model.fit(data, mask=pve_wm_data)
dg = ProbabilisticDirectionGetter.from_shcoeff(csd_fit.shm_coeff,
max_angle=20.,
sphere=default_sphere)
seed_mask = (labels == 2)
seed_mask[pve_wm_data < 0.5] = 0
seeds = utils.seeds_from_mask(seed_mask, affine, density=2)
CMC/ACT Stopping Criterion
from dipy.tracking.stopping_criterion import CmcStoppingCriterion
voxel_size = np.average(voxel_size[1:4])
step_size = 0.2
cmc_criterion = CmcStoppingCriterion.from_pve(pve_wm_data,
pve_gm_data,
pve_csf_data,
step_size=step_size,
average_voxel_size=voxel_size)
# Particle Filtering Tractography
pft_streamline_generator = ParticleFilteringTracking(dg,
cmc_criterion,
seeds,
affine,
max_cross=1,
step_size=step_size,
maxlen=1000,
pft_back_tracking_dist=2,
pft_front_tracking_dist=1,
particle_count=15,
return_all=False)
streamlines = Streamlines(pft_streamline_generator)
sft = StatefulTractogram(streamlines, hardi_img, Space.RASMM)
save_trk(sft, "tractogram_pft.trk")
if has_fury:
scene = window.Scene()
scene.add(actor.line(streamlines, colormap.line_colors(streamlines)))
window.record(scene, out_path='tractogram_pft.png',
size=(800, 800))
if interactive:
window.show(scene)
# Local Probabilistic Tractography
prob_streamline_generator = LocalTracking(dg,
cmc_criterion,
seeds,
affine,
max_cross=1,
step_size=step_size,
maxlen=1000,
return_all=False)
streamlines = Streamlines(prob_streamline_generator)
sft = StatefulTractogram(streamlines, hardi_img, Space.RASMM)
save_trk(sft, "tractogram_probabilistic_cmc.trk")
if has_fury:
scene = window.Scene()
scene.add(actor.line(streamlines, colormap.line_colors(streamlines)))
window.record(scene, out_path='tractogram_probabilistic_cmc.png',
size=(800, 800))
if interactive:
window.show(scene)
References
the full source code of this example
. This same script is also included in the dipy source distribution under the doc/examples/
directory.