I assume it is ok and recommended to use the two commands together when using the 5ttgen hsvs defaults, and needing a lut table ouput, please confirm
By the way, which are the defaults for the 5ttgen hsvs for hippocamal and thalamic segmentation ( or -hippocampi and -thalami options)
When it is correct/indicated to use labelsgmfix after 5ttgen hsvs, when specific -hippocampi and -thalami options are used? I assume after using ‘first’ it is ok, and after using ‘aseg’ is not, how about after ‘subfields’ and ‘nuclei’, respectively.
If within 5ttgen hsvs, FSL’s FIRST is being used to segment the sub-cortical grey matter (as is the default as long as it is installed), then most likely you would want the delineations of those structures within your parcellation image to be based on the same method; otherwise you could have streamlines terminating within such a structure as determined by the 5TT image but not being assigned to that structure in the connectome.
However I’m not sure what you mean here by “needing a lut table output”?
I think I didn’t add the “defaults” to the command help page because it’s not quite as simple as that; there’s more complex logic in place based on both the availability of data in the input dataset and the availability of installed software that determines what methods are used for what structures. I could possibly add more verbose details into the Description field I suppose.
But in the simplest case:
The FreeSurfer hippocampal subfields (and amygdalae if included) module will be used (if present);
FIRST will be used for the thalami, even if the FreeSurfer thalamic nuclei segmentation is present (I wasn’t happy with the results I saw from the latter, but am happy to be proven wrong in my decision here).
I’m not sure I follow your question here fully. labelsgmfix and 5ttgen are separate operations that deal with different types of data in parallel; one does not come “after” the other. It would be reasonable to say that whatever technique is used to segment a particular structure in one method, it would be preferable to use the same technique in the other method. Which I suppose exposes a feature-incompleteness, which is maybe what you were trying to get at: labelsgmfix currently just runs FIRST and substitutes its parcels, whereas if the hippocampi / amygdalae / thalami are segmented using FreeSurfer modules in 5ttgen hsvs, it would theoretically be preferable for their delineations in a parcellation image to be derived from the same source. So maybe I need to add some options to labelsgmfix to import data from those modules… Added to the ever-expanding issue list.
Thank you, Rob.
Yes, it would be nice to have labelsgmfix options to adjust freesurfer labels to (exactly) the 5ttgen hsvs output. I am referring not only to the basal ganglia/thalami, or mesial temporal structures, but also the cortex (gm-wm boundary). The OBHM abstract suggests (from the figure) that with 5ttgen hsvs, the cortex is segmented by Freesurfer. If true, I would expect the cortical ribbon in the act.mif image to be similar to that in aparc.mgz, but on close inspection, I see small differences in some places, which can affect, as you said, the connections of tracts to the parcels, if one uses Freesurfer-derived parcels. I here attach such image
Fig 1 shows the aparc-aseg (Freesurfer) derived parcelation image, fig 2 shows a gm-wm mask derived from the act.mif vol obtained using 5ttgen hsvs command, and fig 3 shows a decimated tdi vol overlapping the parcellation image, with tracts cropped at the gm-wm boundary (seed dynamic, 10 mil fibers, sift2 processed). The cursor indicates a region where the tracts do not abut the cortical ribbon, likely due to a local incongruity between the act.mif (fig2) and aparc (fig 1). Thank you,
Octavian
it would be nice to have labelsgmfix options to adjust freesurfer labels to (exactly) the 5ttgen hsvs output.
To nitpick (as I do regularly), actually requires an additional layer of development in order to be “exact”. Currently, both 5ttgen hsvs and labelsgmfix will internally execute FSL’s run_first_all command in order to obtain sub-cortical grey matter segmentations. If there is any imprecision / stochastic influence within that command, or indeed if the input images to that command are different in any way (manual T1w input vs. orig.mgz), the two may not match exactly. What would be preferable would be for both scripts to have a command-line option whereby a pre-calculated run_first_all output can be provided as input.
The OBHM abstract suggests (from the figure) that with 5ttgen hsvs, the cortex is segmented by Freesurfer. If true, I would expect the cortical ribbon in the act.mif image to be similar to that in aparc.mgz, but on close inspection, I see small differences in some places
There’s a number of potential sources of discrepancy here, most of which aren’t actually problematic.
aparc.mgz is generated by FreeSurfer, not MRtrix3; hence I can’t provide any guarantees on its accuracy nor can I even describe with confidence exactly how it is generated. I invested a decent amount of effort in getting the algorithm for going from a closed surface representation to a partial volume fraction image as accurate as I could, but it’s not a trivial exercise (especially pial surfaces on adjacent gyri where geodesically distant vertices can just about touch ), so I can’t guarantee that there’s no issue; but it would require stronger evidence than what you have to prove that this algorithm is giving an incorrect result given the other points below.
aparc.mgz is an integer label image, meaning that each voxel can only be classified as cortex or not cortex. Conversely, the output of 5ttgen hsvs has the capacity to encode partial volume fractions in each voxel. As such, even if there is no other underlying source of discrepancy between the two, there must by definition be voxels with a non-zero GM fraction in the 5TT image but absent in aparc.mgz, and voxels present in aparc.mgz that in the 5TT image have a GM fraction less than 1.0.
If aparc.mgz were to include within it all voxels that have any intersection with the reconstructed cortical surface, then the cortex would actually look thicker in that image, because where a voxel has a tiny intersection with the surface in one corner, that entire voxel gets included, and you then perceive the extremity of the cortical surface as being in the opposing corner of the voxel. Hence, there is some heuristic within FreeSurfer that determines whether or not a voxel should or should not be included in aparc.mgz. As mentioned previous this is out of my control and I don’t actually know how it’s done, but I expect it to be something like voxels requiring at least a 0.5 partial volume fraction.
In all of these cases, there’s the issue of looking at inherently 3D data within a 2D slice, which can be misleading depending on the shape of the geometry through-plane.
The cursor indicates a region where the tracts do not abut the cortical ribbon, likely due to a local incongruity between the act.mif (fig2) and aparc (fig 1).
This is precisely why tck2connectome by default does not just look at the contents of the parcellation image within the two voxles in which the streamline endpoints reside. A more advanced heuristic is necessary because even in this ideal case, where the 5TT image and the parcellation image are derived from exactly the same source, the fact that partial volume fractions are encoded & interpolation can be performed on one image but not the other means that it’s still possible for streamlines to terminate in a voxel that’s not present in the parcellation image. Historically what others have done to address this is dilate the labels within the parcellation image so that all voxels contain labels. The default tck2connectome approach is the “radial search”, which just expands a sphere at the location of the streamline termination point looking for the nearest label. It’s far from ideal, but as long as one is not natively using the surface representation for both streamline termination and assignment to parcels, some heuristic is required.
For the sake of the argument / better understanding: an alternative (not fantastic) way that this mismatch can be mitigated is by using the 5ttgen freesurfer algorithm. That directly uses aparc.mgz to generate the 5TT image, and hence there would be no mismatch between the two by definition. However it would also mean that the tissue segmentations utilised by ACT would have a similar jagged-voxel effect to that present in aparc.mgz, and hence those non-smooth structures would manifest in the tractogram also.
Thank you, Rob, the answer is very detailed, as always, and the issue more complex than I thought, as always. I have been using -assignment_radial_search (although I am not sure one should choose -assignment_forward_search max_dist instead, to preserve streamline direction, any thoughts?). The gap shown in the pictures may also be due to the decimated .tck file used for visualization purposes, it may be the whole track vol would fill the gap better. One option, for the future, would be to make labelsgmfix dependent on 5ttgen output, and add an option combining the freesurfer/aparc.mgz cortical segmentation (as is) with the FSL subcortical segmentation. I agree this may help mostly testing of the incongruencies between outside parcelation, and segmentation, and the run-to-run variability, but would not necessarily add to the ground truth. Thank you, Octavian
Even though I left the implementation there, as I described in another thread, I don’t think I’d actually advise using -assignment_forward_search; it’s there for conceptual demonstration more than anything else.
The gap shown in the pictures may also be due to the decimated .tck file used for visualization purposes, it may be the whole track vol would fill the gap better.
This depends on what you mean by “decimated”. If you’re just reduced streamline count, then that’s unlikely to “close the gaps”, given that you are showing enough streamlines such that gaps do not appear between them elsewhere in the image. If in reference to vertex downsampling, whether explicit via tckresample or implicitly within mrview (which happens seamlessly), the endpoint vertices are always preserved.
One option, for the future, would be to make labelsgmfix dependent on 5ttgen output
There isn’t enough information in the 5ttgen output for labelsgmfix to do its job. What’s needed (and would be more faithful to the software tools being utilised) would be for FSL’s run_first_all to be executed first, and then its output utilised by both 5ttgen and labelsgmfix.
I’ve also noticed a mismatch between the labelsgmfix output and the 5TT hsvs output, particularly in the SGM structures.
What is the current recommendation on the input T1 image to labelsgmfix to get optimal overlap with the 5TT hsvs image? I’m currently using the 001.mgz (from FreeSurfer’s recon-all /mri/orig/ directory). It seems to do an okay job, but doesn’t provide perfect overlap - just wondering if this is the best input to get optimal matching
The two can never have perfect correspondence if only because the 5TT image can represent partial volume fractions whereas a hard parcellation image necessitates a binary representation per parcel; I’m presuming that you’re observing a difference of greater magnitude than this?
When executed from within 5ttgen hsvs, FIRST is provided with image norm.mgz. I’m not intimately familiar with the data provenance of FreeSurfer so can’t comment on what difference they may have with 001.mgz, but the MGH format can in fact embed some provenance tracking that MRtrix3 can pull into its own reserved command_history tag and that may provide a clue. You could also check the transforms and voxel sizes of those images. But if you use norm.mgz and get different results, that would suggest some stochasticity in the FIRST processing, which could be evaluated separately to those scripts.
Better again would be to use FIRST once as a standalone tool and then have both 5ttgen and labelsgmfix pull data from the same source; I’ve just never gotten around to implementing it…
Thanks Rob - yes, the differences I was seeing were fairly large. In some instances, about half of the segmentation from the parcellation image did not overlap with the corresponding structure in the 5TT image.
Using the norm.mgz image as an input for parcellation generation helped fix this problem - there are still (small) differences between the two segmentations, but as pointed out, these differences are probably due to partial volume effects
When using the -sgm_amyg_hipp option, the overlap in SGM structures between the parcellation image and the 5TTvis image is ideal. However, this leads to a void in the cortical ribbon (shown at the position of the cursor, and is particularly obvious in the coronal view). Conversely, when the -sgm_amyg_hipp option is not used, the overlap in SGM structures between the parcellation image and the 5TT image is not particularly good (which is to be expected) and there is no void in the parcellation image.
I think the void occurs because FastSurfer assigns the “missing” region as hippocampus or amygdala, and these SGM structures are then replaced with FSL’s SGM segmentations (which do not extend into the cortical ribbon, thereby creating the void). The potential issue here is that, while the -sgm_amyg_hipp option works well for SGM structures, it does result in GM regions of the 5TT image that are unassigned in the parcellation image - which might have implications for tractography/connectomics…
I have received interest in filling the parcellation image void, but am not sure about the best way to do this…extending adjacent nodes along the cortical ribbon to meet in the middle of the void might be one option (which is neat in that retain the 84 node atlas, but it also jeopardises the accuracy of the adjacent nodes)…another simpler option might be to assign the missing areas as dummy nodes (which won’t change the accuracy of existing nodes, but will increase the number of nodes in the atlas)…
), so I can’t guarantee that there’s no issue; but it would require stronger evidence than what you have to prove that this algorithm is giving an incorrect result given the other points below.
