Note
Click here to download the full example code
Parallel Transport Tractography (PTT) [Aydogan2021]
from dipy.direction import peaks_from_model
from dipy.data import default_sphere
from dipy.io.streamline import save_trk
from dipy.io.stateful_tractogram import Space, StatefulTractogram
from dipy.data import get_sphere
from dipy.direction import PTTDirectionGetter
from dipy.reconst.shm import CsaOdfModel
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, load_nifti_data
from dipy.reconst.csdeconv import (ConstrainedSphericalDeconvModel,
auto_response_ssst)
from dipy.tracking import utils
from dipy.tracking.local_tracking import LocalTracking
from dipy.tracking.streamline import Streamlines
from dipy.tracking.stopping_criterion import ThresholdStoppingCriterion
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')
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)
seed_mask = (labels == 2)
white_matter = (labels == 1) | (labels == 2)
seeds = utils.seeds_from_mask(seed_mask, affine, density=1)
response, ratio = auto_response_ssst(gtab, data, roi_radii=10, fa_thr=0.7)
csd_model = ConstrainedSphericalDeconvModel(gtab, response, sh_order=6)
csd_fit = csd_model.fit(data, mask=white_matter)
0%| | 0/58788 [00:00<?, ?it/s]
1%|1 | 643/58788 [00:00<00:09, 6422.73it/s]
3%|2 | 1530/58788 [00:00<00:07, 7853.65it/s]
4%|4 | 2443/58788 [00:00<00:06, 8431.34it/s]
6%|5 | 3348/58788 [00:00<00:06, 8671.79it/s]
7%|7 | 4268/58788 [00:00<00:06, 8860.71it/s]
9%|8 | 5193/58788 [00:00<00:05, 8992.74it/s]
10%|# | 6124/58788 [00:00<00:05, 9096.39it/s]
12%|#2 | 7064/58788 [00:00<00:05, 9191.11it/s]
14%|#3 | 8012/58788 [00:00<00:05, 9279.94it/s]
15%|#5 | 8940/58788 [00:01<00:05, 9257.52it/s]
17%|#6 | 9878/58788 [00:01<00:05, 9294.45it/s]
18%|#8 | 10823/58788 [00:01<00:05, 9340.29it/s]
20%|## | 11758/58788 [00:01<00:05, 9335.06it/s]
22%|##1 | 12702/58788 [00:01<00:04, 9364.68it/s]
23%|##3 | 13666/58788 [00:01<00:04, 9445.47it/s]
25%|##4 | 14627/58788 [00:01<00:04, 9491.27it/s]
26%|##6 | 15577/58788 [00:01<00:04, 9467.71it/s]
28%|##8 | 16543/58788 [00:01<00:04, 9522.29it/s]
30%|##9 | 17502/58788 [00:01<00:04, 9542.36it/s]
31%|###1 | 18457/58788 [00:02<00:04, 9489.56it/s]
33%|###3 | 19407/58788 [00:02<00:04, 9452.42it/s]
35%|###4 | 20353/58788 [00:02<00:04, 9440.66it/s]
36%|###6 | 21298/58788 [00:02<00:03, 9436.38it/s]
38%|###7 | 22266/58788 [00:02<00:03, 9507.31it/s]
39%|###9 | 23217/58788 [00:02<00:03, 9437.38it/s]
41%|####1 | 24161/58788 [00:02<00:03, 9411.71it/s]
43%|####2 | 25114/58788 [00:02<00:03, 9446.75it/s]
44%|####4 | 26059/58788 [00:02<00:03, 9440.33it/s]
46%|####5 | 27004/58788 [00:02<00:03, 9382.84it/s]
48%|####7 | 27943/58788 [00:03<00:03, 9349.28it/s]
49%|####9 | 28879/58788 [00:03<00:03, 9234.36it/s]
51%|##### | 29803/58788 [00:03<00:03, 9104.07it/s]
52%|#####2 | 30714/58788 [00:03<00:03, 9071.47it/s]
54%|#####3 | 31626/58788 [00:03<00:02, 9085.19it/s]
55%|#####5 | 32551/58788 [00:03<00:02, 9133.51it/s]
57%|#####6 | 33504/58788 [00:03<00:02, 9248.82it/s]
59%|#####8 | 34449/58788 [00:03<00:02, 9308.42it/s]
60%|###### | 35394/58788 [00:03<00:02, 9350.24it/s]
62%|######1 | 36350/58788 [00:03<00:02, 9412.81it/s]
63%|######3 | 37304/58788 [00:04<00:02, 9449.86it/s]
65%|######5 | 38257/58788 [00:04<00:02, 9472.59it/s]
67%|######6 | 39225/58788 [00:04<00:02, 9533.16it/s]
68%|######8 | 40182/58788 [00:04<00:01, 9541.69it/s]
70%|######9 | 41137/58788 [00:04<00:01, 9532.98it/s]
72%|#######1 | 42119/58788 [00:04<00:01, 9616.45it/s]
73%|#######3 | 43081/58788 [00:04<00:01, 9565.64it/s]
75%|#######4 | 44039/58788 [00:04<00:01, 9568.41it/s]
77%|#######6 | 44996/58788 [00:04<00:01, 9545.14it/s]
78%|#######8 | 45951/58788 [00:04<00:01, 9491.55it/s]
80%|#######9 | 46901/58788 [00:05<00:01, 9449.36it/s]
81%|########1 | 47847/58788 [00:05<00:01, 9446.39it/s]
83%|########2 | 48792/58788 [00:05<00:01, 9432.49it/s]
85%|########4 | 49739/58788 [00:05<00:00, 9441.22it/s]
86%|########6 | 50708/58788 [00:05<00:00, 9513.17it/s]
88%|########7 | 51660/58788 [00:05<00:00, 9483.05it/s]
89%|########9 | 52609/58788 [00:05<00:00, 9338.64it/s]
91%|#########1| 53544/58788 [00:05<00:00, 9258.99it/s]
93%|#########2| 54471/58788 [00:05<00:00, 9164.29it/s]
94%|#########4| 55388/58788 [00:05<00:00, 9145.18it/s]
96%|#########5| 56303/58788 [00:06<00:00, 9055.36it/s]
97%|#########7| 57209/58788 [00:06<00:00, 8953.44it/s]
99%|#########8| 58105/58788 [00:06<00:00, 8835.14it/s]
100%|##########| 58788/58788 [00:06<00:00, 9243.01it/s]
We use the GFA of the CSA model to build a stopping criterion.
csa_model = CsaOdfModel(gtab, sh_order=6)
gfa = csa_model.fit(data, mask=white_matter).gfa
stopping_criterion = ThresholdStoppingCriterion(gfa, .25)
Prepare the PTT direction getter using the fiber ODF (FOD) obtain with CSD. Start the local tractography using PTT direction getter.
sphere = get_sphere(name='repulsion724')
fod = csd_fit.odf(sphere)
pmf = fod.clip(min=0)
ptt_dg = PTTDirectionGetter.from_pmf(pmf, max_angle=15, probe_length=0.5,
sphere=sphere)
# Parallel Transport Tractography
streamline_generator = LocalTracking(direction_getter=ptt_dg,
stopping_criterion=stopping_criterion,
seeds=seeds,
affine=affine,
step_size=0.2)
streamlines = Streamlines(streamline_generator)
sft = StatefulTractogram(streamlines, hardi_img, Space.RASMM)
save_trk(sft, "tractogram_ptt_dg_pmf.trk")
if has_fury:
scene = window.Scene()
scene.add(actor.line(streamlines, colormap.line_colors(streamlines)))
window.record(scene, out_path='tractogram_ptt_dg_pmf.png',
size=(800, 800))
if interactive:
window.show(scene)
Aydogan DB, Shi Y. Parallel Transport Tractography. IEEE Trans Med Imaging. 2021 Feb;40(2):635-647. doi: 10.1109/TMI.2020.3034038. Epub 2021 Feb 2. PMID: 33104507; PMCID: PMC7931442.
Total running time of the script: ( 1 minutes 30.201 seconds)