Beginner: DWI preprocessing for tracking (TrackVis)

Hi MRtrix experts,

I am working with a DTI study with the following parameters:

  • b=3000 s/mm2
  • 32 directions acquired twice (each time a b0 scan was obtained at the begining)
  • No extra b0s (therefore no chance of applying TOPUP)
  • Halfscan: 61%

I trying to use MRtrix3 software to solve crossing fibres problem in my tractography (particulary in Aslant tracks), but I am a beginner in this software. So, I have followed the “DWI denoising” and “Basic DWI processing” tutorials to get a final .tck which I have converted to .trk using TrackConverter tool (https://github.com/MarcCote/tractconverter) for viewing on TrackVis tool. However, I didn’t get good enough results.

Particulary, I ran the following pipeline:

  1. dwidenoise dti_native.nii.gz dti_denoise.mif
  2. dwipreproc -rpe_none -fslgrad dti.bvec dti.bval AP dti_denoise.mif dti_preproc.nii
  3. dwi2mask -fslgrad dti.bvec dti.bval dti_preproc.mif dti_mask_mrtrix.mif
  4. dwi2response tournier -fslgrad dti.bvec dti.bval dti_preproc.mif RF.txt
  5. dwi2fod csd dti_preproc.mif RF.txt dti_preproc_FOD.mif -mask dti_mask_mrtrix.mif -fslgrad dti.bvec dti.bval –force
  6. tckgen dti_preproc_FOD.mif tracks_80000_iFOD2.tck -seed_image dti_preproc.mif -mask dti_mask_mrtrix.mif -fslgrad dti.bvec dti.bval -force -algorithm iFOD2 -number 80000
  7. TractConverter.py –i tracks_80000_iFOD2.tck -o T80000_iFOD2.trk -a corregistered_T1_image.nii

Note: In the step 3, I modified the dti_mask_mrtrix with mrview to get a more precise mask. In the step 6 I put 80000 streamlines to have sufficient fibers for tractography.

The results of each step are shown below:


Figure 1 (Response Function from step 4)


Figure 2 (FOD from step 5)


Figure 3 (2D tractography from step 6)


Figure 4 (3D Aslant tracks in Trackvis).

After all the steps I obtained a complete reconstrution of Aslant tracks (Figure 4). But It has a poor number of tracks and noise in comparison with this:


Figure 5 (image preprocessed from another software without CSD algorithms).

Is there something which I can change to improve my final .tck (.trk)?

Best regards,

Jose

Hi Jose,

I think the problem you have is that for whole brain tractography, 80000 fibers is a low number. I think that you have two options:

a) Do the same but with 5000000 tracts (aproximate number)

or

b) Generate tractography for only the tract that you are interested in. In this case, your seed image should be an ROI containing a part of the tract (note that with this approach you can use more ROIs for exclusion criteria, other inclusion criterie, etc…)

I hope this helps.

Regards,

Manuel

Thank!

I did some tractography tests and I obtained good results generating tractography for only the left aslant tract. I used tckgen command with a seed image containing part of the aslant tract. I obtained the following results on right aslant tract depending on the number of streamlines I generated (with iFOD2):
500.000 streamlines: 19,2 ml/3120 voxels/690 tracts
1.000.000 streamlines: 19,46 ml/3175 voxels/ 732 tracts
2.000.000 streamlines: 18,78 ml/3066 voxels/ 729 tracts
3.000.000 streamlines: 21,25 ml/3325 voxels/ 1150 tracts

What could be the correct number os streamlines? Would I have to try with more streamlines although I use only a aslant seed image (no whole brain)?

Thank for all.

Best regards,

Jose.

Hi Jose,

I can’t answer to this, I’m not expert in the brain anatomy, but I have a doubt, why if you use seeding in the aslant tract with 500000 streamlines, the final result has 690 tracts? from my understanding supose it to have 500000.

Regards,

Manuel

Because I am seeding in a region with contain part of the aslant tract, not completely and not only this. So, the 500.000 tracts are used in generation of all tract with a part of them in the ROI.

In any case, these results are obtained converting the .tck file into .trk file and opening it with Trackvis. So, I think the the number of streamlines of MRTrix3 are not equivalent to the number of tracts of Trackvis.

Could J. D. Tournier have an idea more clear about this topic?

Thank you, Manuel.

Best regards,

Jose.

OK, there’s a lot to discuss here, I’ll try to keep it concise…

Why not just acquire 64 unique directions? 32 directions limits you to lmax=6, when you could go to lmax=8 with 64 - might be a good idea, especially given you’re operating at b=3000, this might start to make a difference. I’m guessing this is a Philips system (given your use of the term Halfscan)? There’s nothing to stop you from inputting your own directions - although I’m told it’s a pretty tedious process…

Also, beware the overplus option on Phillips scanners: do not ever turn that on. It doesn’t look like you have thankfully, otherwise I’d expect your response to look terrible - but it looks fine.

And in general, the tracrography looks sensible, so I don’t think changing the acquisition will solve any of the issues you’re asking about. It’s just that if you have the opportunity to change things, I would recommend going to a true 64 direction acquisition.

OK, from this point on, what you’re asking is dependent on exactly what it is you want to do. I can see you want to delineate the Aslant tract, but what do you want to do with the result? Do you need to extract quantitative information about it, in particular some measure of its ‘connectivity’? This makes quite a big difference to the advice you might get…

So first off: why focus on the Aslant tract specifically? As I understand it (and I may very well be wrong), there’s actually quite a bit of controversy as to whether it even exists… So it seems to me to be very much the wrong tract to practice on and optimise your analysis pipeline. Also, I remember Maxime Descoteaux showing results from their recent ‘Tractometer’ analysis: the FiberFox phantom they used for this was constructed with a limited subset of known in vivo tracts, not including the Aslant tract, yet just about every tracking algorithm tested produced a very convincing depiction of it…

All I’m saying is the evidence for its existence is in my opinion pretty weak, if you consider the diffusion MRI literature only. There may well be convincing independent evidence for it, I’m just not aware of it - but then I haven’t looked for it either, so please don’t take my word for it…


But assuming this isn’t an issue, you also need to be aware of the glaring differences between deterministic DTI and probabilitic CSD tractography: with deterministic approaches, you get the same result for the same seed. In fact, it’s even more restrictive than that: you get the same result for any seed along that streamline. So if you do whole-brain tracking, there are a finite number of possible unique streamlines you can get.

When you go multi-fibre (with CSD or other higher-order models), you now have the potential to choose from multiple directions at each step along the streamline. So the number of possibilities increases substantially. If on top of that you introduce probabilistic propagation, then no two streamlines will ever come out the same. You can keep generating streamlines, and observe new unique trajectories being produced. @rsmith may want to comment here, but if I recall correctly, using a simple Haussdorf-based clustering approach, new clusters (i.e. with at least 5 members) not yet seen in the output were still being produced at the 300,000,000 streamlines mark (or maybe he made it to a billion, I can’t remember). The point being, there’s a lot of possibilities in the data if you allow for the uncertainty you know is there.

So to characterise that uncertainty requires a lot more samples. Note this isn’t unique to tractography: any probabilistic sampling approach needs to produce a big enough sample for inferences to be stable. And the number of samples needed will grow with the complexity of the probability space you’re trying to characterise. It’s just that in tractography, the probability space is ridiculously/wonderfully complex…


So to come back to your original question:

If all you’re after is a clean delineation of the tract of interest, there is no substitute for using well-informed and anatomically justifiable regions of interest. This might consist of a seed ROI in the middle of the tract, an initial direction of tracking to avoid wasting time on unrelated structures, one or two inclusion ROIs near the expected end-points, and maybe some exclusion ROI to remove clearly spurious streamlines if they become a problem. This is basically the classical tract-editing technique, and it requires imposing a lot of prior anatomical knowledge to help constrain the algorithm.

If you want to perform quantification on these results, it depends what kind. If you’re after the mean FA within the tract or some other scalar measure, then the above will probably be fine. If on the other hand you want to quantify some measure of connectivity, then you can’t just do the above and expect the streamline count, or the proportion of streamlines to seeds, or any other metric related to streamline count to be meaningful (see this article on the topic) - unless you process your data differently. In your case, the most appropriate would probably be the SIFT2 technique - @rsmith mentioned this in another post recently. Essentially this would allow you to combine a whole-brain tractogram with a more focused, dense, tract-specific tractogram, and analyse them both jointly to derive a measure of connectivity that also matches the dMRI data over the whole brain - see the SIFT2 paper for details.


Finally, onto the specific questions that came up:

The problem here is that very few of the streamlines generated will end up meeting the criteria for inclusion in your Aslant tract. Increasing the size of the ROI doesn’t necessarily help, if the extra volume isn’t part of the Aslant tract to begin with - you’re just wasting time tracking unrelated structures. As I mentioned above, you can help the algorithm out by providing an initial direction of tracking, which would help increase the proportion of ‘successful’ streamlines. But really, looking at numbers like these is simply not all that informative, for the reason mentioned above (as per Derek Jones’ aforementioned article on streamline counts).

There is no such thing as ‘the correct number of streamlines’. This is a bit like asking how long is a piece of string… If you generate more streamlines, you will get more streamlines. If you use a smaller seed ROI that is located deeper in the core of the tract you’re interested in, you will most likely get a higher proportion of successful streamlines. It all depends on so many factors that it’s essentially meaningless. The only thing that really matters is that the streamlines sample the tract(s) you’re interested in with sufficient density.

If you use something like SIFT/SIFT2 on the other hand, your numbers will be meaningful, but only in relation to the total number generated over the whole brain. But I can’t comment on what the right ‘percentage’ of streamlines might be in any given tract.

I have no idea about this, but I would be surprised if the TrackConverter was somehow messing with the data without any warning. There may be issues if converting the full 500,000 (although I doubt it), but if you’re converting just those streamlines you’ve selected, I can’t see that being an issue. But then I’ve no experience with either TrackConverter or TrackViz, so I can’t really say much more than that. Other users may have something sensible to say here.

Wow :open_mouth:, thank you for your awesome answer, Donald! (and Manuel, of course). All is explained very well. It will help me a lot.
I’ll try to work with this advices.

Best regards,

Jose.

Hi, everyone.

After seeing SIFT2 paper, FAQ section, and other parts of web docs about SIFT, I am not sure what do in my situation (to get volume and number of voxels of aslant track with some reliability).

On the one hand, in SIFT workflow on web, it says: “Generate a whole-brain tractogram; * Apply SIFT; * Extract the pathway(s) of interest using tckedit. * Get the streamline count using tckinfo. The SIFT algorithm is not directly applicable to targeted tracking data”. But If I obtain the aslant track from a whole-brain tractogram generated with 10 millions, for instance, the track does not sufficient density. And with more millions, it does not improve a lot. Also, tckinfo doesn’t show the volume or the number of voxels.

On the other hand, in SIFT2 paper, it combine a whole-brain tdi with a targeted seeding tdi to get a combined reconstruction and apply SIFT/SIFT2. In my case, I have tried to generate a whole-brain .tck, with -seed_dynamic (with FOD image) option and 10 millions of streamlines, and a aslant-specific .tck with -seed_image (Aslant ROI) option and 1 million os streamlines. In both cases I have generated the corresponding track density image .mif (with -vox 0.2) just in case.
The problem is, how should I combine these tractograms? I tried with mrcalc but it does not accept .tck format.
I would like use the combined .tck to use it with tcksift2 command, apply the generated weights with tckedit and convert the final .tck into .trk to view it on TrackVis (as TrackVis show me the volume, number of voxels, etc. of a track).

I hope the post is not too confusing and you can help me.

Thank and best regards,

Jose.

I am just in preparation to do the same thing. According to the some post here and help page of tckedit, the .tck file combination can be achieved by tckedit.

Regards,

Antonin

Yes, tckedit is a bit of a ‘workhorse’: there’s a lot of different tricks and manipulations you can do with it, one of which is combining track files by simply providing multiple input files.

I would like use the combined .tck to use it with tcksift2 command, apply the generated weights with tckedit and convert the final .tck into .trk to view it on TrackVis (as TrackVis show me the volume, number of voxels, etc. of a track).

I suspect the end goal of your experiment needs to be thought out a little more. Firstly, the weights calculated by tcksift2 could only be ‘applied’ to a tractogram if a track file format had the capacity to represent and store such information, which the .tck format does not; hence this information must be handled explicitly. Secondly, altering the weights of individual streamlines will most likely not affect the volume / number of voxels of a tract of interest, because there’s no way to incorporate such information into calculation of such a measure. If you’re looking for a ‘quantitative measure of connectivity’, it’s most likely the streamline weights themselves that you’re looking for…

Thanks to both of you!. I have been a bit busy and I haven’t could reply before.

After reading your posts and testing somethings, I have some new doubts about streamlines weights:

  1. SIFT2 paper says the weights are a “cross-sectional area of the white matter axons connecting regions”. How is it calculated internally from a tractogram and a FOD images? What is the measurement unit?

  2. For instance, in my case, using tcksift2 command I obtain a txt file with the weights of one million of streamlines (streamlines which constitute some pathways). How could I use this information to get a quantitative measure of connectivity in a particular pathway? Should I isolate the pathway tracts, use tcksift2 and calculate the mean or the sum of the corresponding weights?

  3. I saw you can use the weights txt file as a input parameter in tckmap to obtain a image with a sum of streamline weights or a sum of streamline volumes per voxel (even a sum of lenghts using -precise option). Is it posible to save these new data into a txt file? Could they be used as a reliable quantitative measure of connectivity? What are the measurement units used in these cases?

Regards,

Jose

I will get this publication out eventually, I promise… :confounded:
Grant writing season… :rage:

  1. The unit is “AFD per mm”. This is captured within the SIFT proportionality coefficient, which is:
    (Sum of fibre density) / (Sum of track density)
    It seems a little obscure, but AFD is fundamentally a measure of volume (despite the T2 confound & lack of physical units), and hence “AFD per mm” is a measure of cross-sectional area.

  2. Almost: You need to run tcksift2 first (it only works on a whole-brain tractogram; never select a pathway of interest prior to running tcksift / tcksift2), then isolate the pathway streamlines, then take the sum of the weights within your pathway of interest (since the quantitative measure is the total fibre cross-sectional area of the pathway).

  3. Again here it helps to think in physical units. If, for each streamline, the contribution to a particular voxel is multiplied by both the per-streamline weight from SIFT2 (a cross-sectional area) using the -tck_weights_in option, and the length of that streamline within that particular voxel by using the -precise option, what you get is a volume; specifically a fibre volume in units of AFD. Summed across all streamlines traversing a voxel, this gives the fibre volume within that voxel corresponding to your pathway of interest.
    Note however that this is not a measure of “connectivity” in and of itself, but is better thought of as a per-voxel fibre volume. For instance: Summing these across all voxels would give you the fibre volume of the entire pathway; but this quantity is also dependent on the length of the pathway, which shouldn’t affect a measure of ‘connectivity’ as much as the cross-sectional area does.

Rob

Thank you, Rob!

I am trying to obtain the weights within my pathway of interest, but I am not sure If I am doing it correctly.
For instance, after generate a whole brain tractogram (tckgen with 1 million of streamlines) and the tractogram with the pathway isolate (tckgen with 100.000 streamlines), I used tckedit to join both tractograms in one (as SIFT2 paper, I think). Then, I used tcksift2 in this combined tractogram to get the .txt file with the weights of all streamlines. Finally, I used tckedit again with all weights, the combined tractogram, and some including spheres (to isolate my pathway) as inputs to generate the weights of my pathway of interest (-tck_weights_out option). Is it correct?

Moreover, after reading some papers, I am little confused about the differences between “streamline”, “fiber”, “tract”, “axons” concepts. Especially with respect to “density” and “bundles” of these structures (for instance, Is a streamline the representation/modeling of a axon bundle?). Could you tell me a little about the differences?

I add another question to the previous two: If, for instance, I would want to report the sum of cross-sectional area of my pathway as a quantitative measure in a paper, how should I write this? For instance, the IFOF has a connectivity of 3000 AFD per mm? (writing corresponding references).

Best regards,

Jose.

Sounds about right to me. You just need to sum the weights in the file to get the measure you’re after.

You’re not alone, I think there is a lot confusion about this in the field, much of it due to slightly loose use of language.

In my mind, ‘fibres’ refer to single physical fibres, such as a single axon - but it could also refer to e.g. a capillary tube in a hardware phantom. An ‘axon’ is very clearly a real, biological entity. A ‘tract’ is a collection of white matter axons, and should also be reserved for actual biological bundles of white matter axons, along with words like ‘fascicles’ or ‘pathways’. The word ‘bundle’ is a bit more vague, and can probably be applied to describe ‘tracts’ or fascicles, but also their abstract representations - but in this case the context should make it clear (i.e. bundle of fibres vs. bundle of streamlines).

On the other hand, ‘streamlines’ are clearly just abstract 3D paths, which are used as representations of white matter pathways. Another word that I find problematic is ‘track’, since that should also refer to a streamline (i.e. an abstract concept), but sounds sufficiently similar to ‘tract’ (a real biological entity) to lead to confusion; for this reason I tend to avoid it if I can.

What I think is really important is to make sure the distinction between streamlines and axons is clear: one is an abstract approximation of the trajectory that the other might follow. But crucially, there is no one-to-one relationship between axons and streamlines. A single streamline might be used to represent a full tract, for instance. As part of making this distinction clear, I think it’s really important to use consistent terminology: ‘streamlines’ for abstract paths (‘track’ is too problematic), while ‘fibre’, ‘axon’, ‘tract’, ‘pathway’, ‘fascicle’ should be reserved exclusively to describe biological entities - things that can be seen under a microscope. Unfortunately, the word ‘fibre’ is often used as a synonym for ‘streamline’, which I think can lead to misinterpretations, particularly from readers unfamiliar with the field.

The word ‘density’ is also very tricky - we’ve had many discussions between ourselves on this very issue. When looking at the voxel-wise apparent fibre density (AFD), what we’re really reporting is a measure proportional to the intra-axonal volume within that voxel. The word ‘density’ is justified at the single-voxel level by assuming that fibres traverse the voxel, so the effective cross-sectional area of the fibres should be proportional to their volume within each voxel. Things get a bit messy when summing these values across voxels, since then it really does become a measure of the volume occupied by the fibres. However, when using SIFT/SIFT2, the length of the tract is inherently factored out since we’re measuring the (weighted) sum of the streamlines, ensuring their density matches with the voxel-wise AFD, and in this case we end up with a measure of cross-sectional area for the tract - not really density either…

So it’s clearly not a simple issue to reason about, and the terminology is far from perfect. As I mentioned earlier, we had many discussions about whether the term ‘apparent fibre density’ was appropriate, or whether we should switch to something like ‘apparent fibre volume’ - but I think in practice, the interpretation ascribed to these measures will vary depending on what was done and how they were obtained. Most of the time we are genuinely more interested in the density or cross-sectional area than the volume, so many of our tools are geared towards obtaining robust measures of that, rather than volume per se. But to be fair, this is an ongoing debate…

Interesting question… The measures you’d get would be in arbitrary units, and probably best described as proportional to the effective cross-sectional area of the tract you’ve tried to delineate. But there are no units as such… So it would be ‘the IFOF has total apparent fibre cross-section of 3000 a.u.’. I wouldn’t necessarily call this ‘connectivity’ as such, certainly not when reporting the Results: although it is a very tempting interpretation, and cross-sectional area probably is the most closely related measure, I’d leave these issues of interpretation to the Discussion.

But this is just my opinion, others may very well have other opinions on the matter…

However, when using SIFT/SIFT2, the length of the tract is inherently factored out since we’re measuring the (weighted) sum of the streamlines, ensuring their density matches with the voxel-wise AFD, and in this case we end up with a measure of cross-sectional area for the tract - not really density either…

I’d say that referring to these quantities as “connection density” is comparably problematic to referring to AFD itself as “density”. But it does introduce more ambiguity since “density” in the context of discussions such as these can now refer to something with L2 or L3 units.

I add another question to the previous two: If, for instance, I would want to report the sum of cross-sectional area of my pathway as a quantitative measure in a paper, how should I write this? For instance, the IFOF has a connectivity of 3000 AFD per mm? (writing corresponding references).

Interesting question… The measures you’d get would be in arbitrary units, and probably best described as proportional to the effective cross-sectional area of the tract you’ve tried to delineate. But there are no units as such… So it would be ‘the IFOF has total apparent fibre cross-section of 3000 a.u.’

Technically it does actually have units of AFD per mm. But I wouldn’t go writing as such in an article; at the very least not until I’ve published explaining why that is in fact the case and why it’s a useful measure… I’d avoid ‘units’ as best as possible, and indeed would avoid quoting values if possible; focus on relative values. If absolutely necessary, just say it’s a parameter that’s proportional to the (hopefully intra-cellular due to high b-value) cross-sectional area of the pathway.

Thank to both of you! Now, I understand better the process, the issue of the units, and the concepts under the theory. I will be looking out for your future papers! :wink:

Cheers,

José

I’m a bit late to the party (can you tell I’m catching up with a month worth of posts on the forum? :grin:), but to add to this: AFD itself is technically of course not a unit of its own… we (e.g. @Dave) typically add an “a.u.” to graphs or colour scales’ numbers that depict AFD. So the unit of “AFD per mm” becomes “a.u. per mm”, which is again quite an a.u. … :poop: