I came across this reply by Rob, where he implies that it’s possible to “collapse” the measures (FD, FC, FDC) of a specified pathway and end up with a single measurement for the specified tract for each subject. I’m interested in doing exactly that, but I couldn’t find any more detailed information about that specific idea in the rest of the thread.
So anyway, I produced a fixel mask for a specific region of my FOD template, and ran fixelcfestats successfully using the mask. My question is, is it somehow possible to derive a single (average) measurement of FD/FC/FDC for the masked region for each subject separately? If so, how would I go around doing that?
I hope my question makes sense (I know it’s a bit convoluted ), and thanks for the help!!
If you already have a fixel mask of the “fixel region” you’re after, you’re very close to what you want to get. The key thing to realise is that many MRtrix commands designed to work on “voxel images” can work on fixel data files just as well… as long as those data files come with (a copy of) the same index.mif file. That would be the case already for your fixel mask as well as the data files (FD, FC, FDC), as they also share that requirement of a common index.mif for the fixel stats (and those you’ve been able to run successfully).
Your fixel mask, for instance, would be a binary “fixel data file” sitting in one of your fixel folders. Lets say it’s called “fixelmask.mif”. You’d also have fixel folder with all your FD fixel data files, if you’ve been running FBA analysis. Let’s say a subject’s data file in your FD fixel folder is named “subject001.mif” or something. You can then simply run mrstats (name doesn’t suggest being able to act on fixels, but it works nonetheless!) on subject001.mif and add the -mask option, while providing the latter with fixelmask.mif. The output would give you your mean FD within that fixel mask, for that subject (and a range of other stats as well, like the median, or min/max, etc… even the fixel count in your mask).
To already answer another potential question that might be coming up in this context: if you’re going to be averaging stuff within a (fixel) region, and wondering whether to use FC or log(FC), the answer is log(FC).
For a further example of this in an actual study: see e.g. https://academic.oup.com/brain/article/141/3/888/4788771 , where we used this for further comparison of MCI vs controls, with regions (“fixel masks”) actually being obtained from an FBA of AD vs controls. See for instance figures 6 and 7: even though the regions are shown as streamlines there, these are actually derived from the FBA (fixel) results, and were used also as fixel masks to extract FDC values specifically.
Hi,
to deepen in this issue: is it normal to obtain slightly different mean FD values (using a fixelmask) depending on whether we calculate FD in the subject fixel directory or the template fixel directory?
I calculate FD in the subject fixel directory (sub01/fixel_in_template_space_reoriented) where I have directions.mif, index.mif, fd.mif and fixelmask_sub01.mif and I get the following:
If I calculate the same in the template fd fixel directory (fd) where I have directions.mif, index.mif, sub01.mif,sub02.mif…, the values differ:
I first transform the fixelmask to the template:
It’s entirely technically valid both to produce a fixel mask in subject space and compute mean FD within it, and to produce a fixel mask in template space and compute mean FD within it. The reason the latter is more conventionally employed is because the delineation of the pathway of interest in derivation of that mask only needs to be performed once, and can then be applied across all subjects. The former could be done if you could specify an automated & robust mechanism for consistently segmenting the same pathway across subjects within the space of each individual subject, but even doing so there would still likely be residual doubts regarding whether any effect observed was due to a difference in FD across subjects, or a difference in the behaviour of the segmentation method across subjects.
In this specific case the difference is arising through use of the fixelcorrespondence command, which is really a quite basic process in its current form, and may or may not behave as one may subjectively expect it to in any regions where the subject FOD and template FOD are not highly similar.
This made me wonder whether I’m using the template appropriately…
I’m creating seed and target ROIs for each subject separately and then tracking (in template space). In order to create the individual ROIs, I have to compute several transformations: MNI space=> Individual anatomical space=> Individual diffusion space => Template diffusion space…So I have different tracts for each subject even if I’m tracking in template space.
Are there more straightforward ways of transforming a ROI from MNI space to template space and end with a single ROI in template space (and thus, a single tract)?
I’m creating seed and target ROIs for each subject separately and then tracking (in template space).
This is generally ill-advised. The generally-accepted preferred approach is to reconstruct a pathway entirely in subject space (potentially warping atlas information from a template into the individual subject space to do so), and then transform the streamlines into template space. You might then e.g. combine streamlines from across subjects in order to define a single ROI in template space corresponding to that pathway, such that the position and spatial extent of the ROI is then defined only once in template space rather than varying across subjects.
Hi,
Thank you for the response. I’m a bit confused with this information. I will try to summarise to figure out which step I’m misinterpreting:
Tracking:
Should be done in native diffusion space.This means that (for subject=sub01) the FOD input for tckgen is sub01/wmfod.mif… as defined in step 14 before warping them to template space, right?
Extracting FD, FC,FDC values from a tck:
1.fod2fixel: fixel directory has to be created from sub01/wmfod.mif etc since we tracked in native diffusion space. This will create sub01/fixel_directory/
2.tck2fixel: From here, I interpret that fixel_folder_in should be sub01/fixel_directory/ and we could create a fixelmask_native.mif as fixel_data_out.
Since I want to extract the FC, FDC values (which only exist in template space), I need to transform my fixelmask (which is in native diffusion space) to template space. Here comes my question:
3.fixelcorrespondence:
From here, I interpret that subject_data could be fixelmask_native.mif and that template_directory should be Template/Fixel_mask/.
However, the description says “It is assumed that the subject image has already been spatially normalised and is aligned with the template”. I don’t think fixelmask_native.mif meets this assumption, right? It’s only after
I think perhaps the issue is that you are considering the individual steps as being entirely independent; whereas what needs to happen is to decide the overall structure of the experiment, and then determine the requisite steps to complete that experiment.
Given you are now speaking about also extracting tract-wise FC and FDC values, which are only defined in template space voxels, I would personally focus on deriving one mask in template space per tract, from which you can then extract mean FD / FC / FDC per subject. Going to the effort of defining an individual fixel mask per tract per subject, only to then try to map those into template space, just sounds like a source of unnecessary work and unwanted variance to me. Is there any strong motivation to perform any of these calculations in individual native spaces?
Hello Thijs,
I am very much interested in performing a similar Tract of interest analysis. I understand that we can extract the FD, FC or FDC within a tract using a fixel mask.
Could you please guide me on how to generate a tract specific fixel mask as shown in figures 6 & 7 from the Atlas?
Are these binary fixel masks? If so, how do i explicitly calculate FD in a fixel with crossing fibers (for example posterior cingulum and splenium of corpus callosum)?
Yep, that is correct. So this is conceptually the same as a simple region-of-interest analysis in the world of “voxels”. But here, rather than a binary (voxel) ROI, you’ll need a binary fixel ROI (or “mask”) indeed. We typically create this directly on the fixel-analysis mask in template space you’d already have sitting there when you’re running the main steps of an FBA or similar pipeline.
There’re multiple ways to get these.
In some previous works, we’ve used the actual significant fixels from one (FBA) analysis to “become” the binary fixel ROI(s) of another analysis or statistical test. Key here is realising that you can use most commands that work on images, also on the fixel “data files” in each fixel folder. If you want to get a binary fixel mask from your significant fixels for example, you can simply run mrthreshold with the -abs 0.95 option on the fixel data file (.mif, in the fixel folder) that has the p-values. This will yield another fixel data file, that you can store either in the same fixel folder, or alternatively copy in a new fixel folder wherein you also make a copy of index.mif and directions.mif of the original fixel folder.
The other scenario to get a binary fixel mask is when you’ve got some a priori hypothesis or interest for a tract of interest (e.g. “the cingulum”). Again multiple ways to go at this, but the simplest is:
Perform tractography (tckgen) on the FOD template of your study, using inclusion/exclusion regions and parameters of your choice to obtain a good set of streamlines representing that structure. It’s important the structure of interest itself is really densely filled with streamlines. If there’s some spurious streamlines that go beyond or outside the structure of interest, that’s not a worry (see next steps!)
Map this tractogram of the bundle of interest to your fixel analysis mask you have in the FBA pipeline. This is done using tck2fixel. This is essentially a “track density image” (TDI), but then on fixels, rather than on voxels. It’ll generate another fixel data file (that once more also belongs together with index.mif and directions.mif in the fixel folder). For each fixel, it’ll have the count of streamlines going “through” that fixel (which means: which go through the relevant voxel, and along an angle that is close enough to the fixel). Notably, this is indeed fixel-specific: streamlines may add to the count of one fixel in a voxel, but not to the other fixel(s) in the same voxel.
Threshold this fixel-wise TDI image using mrthreshold with -abs ... and some number of minimal streamlines that need to have passed through the fixel for the fixel to be “in” the binary fixel mask. This is the step where you can play with this threshold to get a nice continuous, clean and solid “bundle” in your final binary fixel mask, even when there were a few spurious false positive streamlines before. This will then output the binary fixel mask (similarly as above in the example of thresholding the p-value fixel image).
So at this point, you’ve got a binary fixel mask for your bundle/tract of interest (or multiple, if you did this for multiple bundles). That means the FD fixel files match the binary fixel masks in ordering of storing all the numbers: this is per individual fixel. So essentially, there’s no worry about voxels with multiple fixels here: each entry is an individual fixel, both in the FD data files as well as the ROI binary files. That means you can provide these just to e.g. mrstats: i.e. run mrstats on a subject’s FD fixel data file, and add -mask ... where you provide the binary fixel mask. With mrstats for example, you can then compute the mean of all FD values in the binary fixel ROI. You can run this for example in a foreach loop, and extract specifically the mean values; and then copy-paste those values in whichever software you prefer for stats (e.g. R or Matlab or similar…).
Those are the basics; but it takes some experimenting to get used to the process of doing this. Hope it helps!
I had some hypothesis specific ROIs in NIFTI format from a previous study. I converted them into fixel masks using voxel2fixel. As you suggested i ran mrstats with -mask option to -output mean and it works!!
Now, i managed to get the mean FD within the fixel mask. I see that there are quite some fixels within this fixelmask that have crossing fibers. Its not clear to me as to how the mean FD within a fixel with multiple fibers/crossing fibers is calculated. For example: if there are two crossing fibers within a fixel, will the mean FD in this fixel correspond to average(FD_fiber1,FD_fiber2)?
Aha, yes, that would work for sure. But now I also better understand your context, and I realise where part of your initial question came from, I think.
So, there’s something to consider here, I think. Having your hypothesised ROIs initially as voxel images, basically means you’ve got no “fixel specificity” in your hypothesis in a way. The question maybe is: do those ROIs represent actual “tracts”? Or are they just more “blobs” (for lack of a better word) somewhere in the WM? If they’re supposed to represent tracts, and you’ve got an idea which ones, you might still pursue something like the tractography-to-fixel approach I mentioned above. You could even limit the tractography itself to be masked in those ROIs. Or alternatively, it’s even possible to just track those tracts, convert them (as above, via a fixel TDI) to a binary fixel mask, and then: compute the intersection of your a priori ROI(s)-to-fixel with the tract-to-fixel binary mask. Similarly to how you can “simply” use mrstats, you can also use mrcalc to e.g. compute the intersection of 2 binary fixel masks!
So a number of options here, if there’s some definition of a known “tract” involved. The purpose then would be to end up with a fixel mask that has some degree of fixel specificity in areas: i.e. voxels where the fixel mask itself doesn’t per se include all fixels in the voxel.
Otherwise, you’ve got indeed this situation:
In that case, to answer what happens for this one:
Well, it just averages the FD over all fixels in the complete binary fixel mask, regardless of how they are grouped (or not) in individual voxels. So for example, if you mask had just 4 fixels, where 1 sits in a voxel alone and the 3 other ones sit together in a single voxel, the result of mrstats -output mean is the mean of the FD of all 4 fixels, equally weighted. There’s no “separate” voxel-wise averaging that happens first.
This might be important depending on what you actually want to compute, especially for FD very specifically (because it has some properties of a density, even though it’s far from a proper density)! If your initial ROI(s) is/are just WM “blobs” without a specific tract hypothesis, then in lots of scenarios you may not want this behaviour. You might want to e.g. average over all voxels in an ROI, but after first summing (not averaging) all FDs of fixels within voxels. This would give you the average FD of that blob/chunk of WM tissue. Due to averaging over all voxels, it normalises for size of the ROI, where size comes at 1 (same) unit of volume per voxel. If that is what you’re after, then actually you’re taking an unnecessary detour going via fixels! And one that introduces imprecisions at that…
In this case, you can rather apply mrstats -output meandirectly to the first volume of the FOD image, using your voxel ROI with -mask.
I hope that all makes sense. It’s a bit complicated at times, and it gets more complicated due to FD not being an actual density. It’s sensitised to intra-axonal signal, but not specific to intra-axonal volume. This is generally still ok, if you deal with it reasonably, e.g. most of the time in FBA. It gets more problematic when things are integrated over different regions, e.g. in ROI analyses (or certain kinds) and connectomics.
I’m new to Mrtrix-could you please explain how you obtained the fixel masks from FBA, if these fixel masks are tract-specific regions? (i.e. fixel masks covering each white matter tract, as I am looking to calculate FBA metrics for each WM tract)
In particular, I saw the response about inclusion and exclusion ROIs. How can you make these? Is there a specific guide/atlas one can use? I would like to get the tracts identified in Mito. et al, 2018 (MCI and PD versus control)
No worries, but it’s quite challenging if you’re not yet used to some of the basic (and slightly more advanced) usage of some of these steps… It might be good to first get yourself somewhat familiarised with a basic FBA pipeline; in particular the use of the fixel “format” (aka fixel folders). Otherwise, I can imagine some of the content in the posts above might be puzzling. Just a few quick remarks maybe, to makes sure you’re on the right track generally:
So it’s important here to distinguish between 2 very different scenarios. I’m mentioning this, because your wording at least seems to slightly (maybe!) mix them up. Ignore this if you were absolutely sure what you were after, but in any case, these are the 2 different kinds of scenarios:
“obtain fixel masks from FBA”; this reads a bit as the following scenario: you’ve performed an FBA, and thus you’ve already obtained a bunch of significant fixels. You then want to extract metrics from (subgroups of) these fixels. This is what was done in Mito et al. 2018. The FBA there had itself pulled out significant fixels for a controls versus AD comparison. These fixels themselves were then split up in groups of fixels that roughly belonged to different “identified” structures in that result. Those are the visualisations in the glass brains you see in that paper. For those “identified” structures, we computed then average metrics for the controls and MCI population, to then compare those. We did this because a “full blown” FBA of controls versus MCI wasn’t sensitive enough to pick things up. So the controls versus AD result itself was in a way then used as a more specific hypothesis for controls versus MCRI. But in practice, this was then an exercise to “segment” the FBA result into sub-structures. @rmito did a lot of manual work combining the FBA result with anatomical knowledge. And all of that using a range of commands to get this done in practice.
“fixel masks covering each white matter tract, as I am looking to calculate FBA metrics for each WM tract”; this on the other hand reads (to me) more like: you have some tracts in mind a priori. So the tracts or fixels done come “from” an FBA. You just have (as in the FBA pipeline) at some point a mask of all the fixels in the brain, and then you have a set of tracts with names in mind (e.g. “the corticospinal tract” or “the optical radiations” or… etc…). This is a different exercise: you have to perform targeted tractography of those tracts, and then use that outcome to map it onto fixels (and via some threshold turn it into a binary fixel mask).
Scenario 2 would indeed use inclusion and exclusion ROIs: that’s part of performing targeted tractography. You essentially define rules for tractography so you obtain your tract of interest. Inclusions regions are typically located at both ends of the tract, and possibly some “check points” in between. Exclusion regions are typically defined out of experience of tracking without them: you’ll see tractography getting so creative that it finds alternate (false positive, nonsense) paths between your inclusion regions. Once you know where this happens, you define an exclusion region along those false pathways to forbid them to happen next time. There’s many different review papers on tractography of specific structures that list some of these with tips. But generally, all it takes is an anatomy book and common sense: if you know roughly where a tract starts/ends, you can just use the ROI tool to draw regions there, and provide them via the -include option to tckgen. If you’ve never done this before, allow yourself a day or two to get experience with this.
The tracts identified in @rmito’s work were defined in more hybrid ways, since those already had the FBA output up front. This is more scenario 1, whereas (I think) you’re after scenario 2.
I hope that helps somewhat! It’s hard to teach this in writing beyond this though. Here’s an ISMRM educational I did for ISMRM a few years ago, which brings up targetted tractography at some point, if you weren’t familiar with it yet. To turn this into fixels you’ll need tck2fixel to make a fixel-wise tract-density image, followed by mrthreshold to turn this into a binary fixel mask. That can then be used with mrstats to pull out e.g. the average FD or similar.
I ran my fb analysis (on some tracts of interest, in this case the corticospinal tract) and get some nice results. Now I wanted to extract the fd values (only in the significant fixels) and used the approach you described here, using mrthreshold on my fwe_…mif file. However when I’m using the resulting image as the mask in mrstats, it’s telling me that the dimensions between my mask and the directions.mif + index.mif files of my input file mismatch (between dimensions 0 and 2). I already did the same approach but used the binarised fixel mask from the whole tract, and there it worked fine, so I can’t figure out what’s different with my new thresholded mask.
I’m running my analyses on data from a longitudinal study, so my stats input file is the difference image of timepoints 1 and 2. Could this be the problem?
I hope you can help me with this.
Thank you and best,
Kristin
If utilising a fixel mask to perform statistics on values from a fixel data file only using values from within specific fixels, then there is no reason why the index or directions files should be involved at all. If you have a fixel mask that is Nx1x1 (with N being the number of fixels), and a fixel data file (e.g. containing Fibre Density (FD) values) that is also Nx1x1, then mrstats will have no issue: the binary mask just controls which fixel values contribute to the statistics and which do not. The index and directions images serve a fundamentally purpose: they essentially specify the spatial location and fibre orientation for each of the N fixels in those data files. But as long as fixel index i in the mask image is the same fixel as index i in the FD data file, then there’s actually no need to know the spatial location of each fixel in order to calculate those statistics.
If you’re still uncertain, maybe include here the exact commands that you are attempting to execute, along with the output of mrinfo for each of the files involved.
Thanks Rob for your clear answer!
Turned out I should have written my bash script a bit clearer. This way it tried to run mrstats on the direction and index file as well, which of course had to result in an error, but everything else worked fine. So just a stupid mistake that cost me a few hours staring at mrinfo outputs trying to figure out what went wrong
And one that introduces imprecisions at that…