I think we may have followed some red herrings due to your original description being interpreted on our part more specifically than would have been ideal.
If we defer for a moment what a “peak” is, how one might be either discovered or represented, and whether it comes from FOD data or something else, and just think about the phrase:
and check in the data of the patients for changes in FODs
, this becomes a much more general question about how one might isolate some subset of an FOD image for quantification, but where that isolation specifically needs to occur in both the spatial and orientation domains.
Expressing the question in this way avoids the ambiguity of what a “peak” is, which is I’m pretty sure the principal source of the confusion. In retrospect our documentation page on “fixels” (which admittedly still needs a good image or two to demonstrate the concepts) would benefit from additionally disambiguating between “fixel” and “peak”…
We explicitly define “peak” as a location at which the FOD amplitude is maximal1. It has a direction on the unit (half-)sphere, and the FOD has an amplitude along that direction. The way that we identify peaks in command
sh2peaks is best demonstrated in Figure 1 of this manuscript.
While a fixel also has a direction, and quite commonly some measure of “magnitude” with it, it differs from a “peak” in a number of ways:
For a “peak”, the direction is by definition that in which the amplitude of the FOD is maximal. For a “fixel”, the direction is an estimate of the principal orientation of the fibres within the voxel that are ascribed to that fixel. This could theoretically come from any estimation mechanism. In the context of MRtrix3 command
fod2fixel, we by default use the mean orientation of the FOD lobe that was segmented to produce that fixel; this may differ from the peak orientation by quite a bit in cases of complex FOD shapes.
For a “peak”, there is only one scalar value associated with it, being the amplitude of the FOD at that orientation where the amplitude is a local maximum. For “fixels”, we can ascribe any scalar value we wish to each fixel. In a typical processing pipeline, the first such quantitative vvalue that we would derive is the Fibre Density (FD). This is not the peak amplitude of the FOD, but the integral of the FOD lobe in the half-sphere. This is invariant to fibre orientation dispersion, whereas the peak amplitude is not.
The way that we identify fixels in command
fod2fixel is described in Appendix 2 of this manuscript.
So, back to the confusion.
When @jdtournier and I interpret very literally a question like this:
What we are planning to do this by creating an FOD-image of the streamlines that are predicted to be damaged, establish the peak directions in this FOD-image, and check in the data of the patients for changes in FODs along these peak directions.
, what we would respond with would be something like:
tckmap with the
-tod option to produce Track (not Fibre) Orientation Distributions from a set of streamlines;
sh2peaks to identify the orientations of the “peaks” (maxima) of these functions in each voxel;
sh2amp command to optionally accept as input an image containing a set of orientations in each voxel, and in each voxel, evaluate the amplitudes of the subject’s FOD image only in those orientations specified by the image produced by step 2.
So the FOD amplitudes measured as the output of the experiment are not the amplitudes of the subject’s own FOD function peaks; they are the amplitudes of the subject’s own FOD function in those directions defined by the Track Orientation Distributions.
However, what I think you actually want, I would propose expressing something like:
“To evaluate quantitative FOD properties where traversed by a set of streamlines.”
If this is indeed still consistent with your actual intent, then my response would instead be something like:
Map streamlines to fixels (that were produced from segmentation of the subject’s own FOD functions) using the
Apply some threshold to these per-fixel streamline counts using
mrthreshold to produce a fixel mask;
Evaluate per-fixel quantitative values (e.g. Fibre Density) within those voxels included in the mask.
A re-interpretation of words, but a completely different recommendation; also closer to the conventional approach for such an experiment when dealing with voxel-wise quantitative values (map streamlines to voxels; apply a density threshold; extract voxel-wise values from within the mask).
I hope that gets us all back on the same page…
1 Unfortunately what might make this even more confusing is the fact that the way the term “peaks” is used in the context of MRtrix3 command names is different again. We have our convention of frequently naming commands based on the form of the data input / output. There, we currently use “
peaks” for commands where the nature of the data is a 4D image, with 3xN volumes where N is the maximum number of discrete orientations in each voxel, and each triplet of volumes encodes an orientation in 3D Euclidean space, optionally with some scalar parameter encoded based on the norm of the vector. We’ve had multiple discussions about this ambiguity, because there’s a discrepancy between what a peak of an orientation distribution function is, and the image data format just described in which such data have conventionally been encoded, which is more like “3-vectors”, but we’ve not justified changing the relevant command names to something better. The fixel directory format can be thought of as a more generalised solution to representing such data, and hence developments around such are of greater priority than how we process—or indeed name—what we hope will over time become a “legacy” format.