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.
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:
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.
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:
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:
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?
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.