I have been running into issues trying to obtain the centroid coordinates of an ROI that spans a large curved region. The center of mass is outside the ROI because of its shape. See attached for an example in volumetric space, the crosshairs point to the calculated centroid coordinate and the ROI is colored in yellow.
My approach has been to use
mri_surfcluster in freesurfer to obtain centroid coordinates in MNI from the cortical surface as follows:
mri_surfcluster --in $SUBJECTS_DIR/fsaverage6/surf/rh.curv \
--clabel $SUBJECTS_DIR/Schaefer2018_100Parcels_7Networks_order_fsaverage6_labels/rh.7Networks_RH_Default_PFCm_2.label \
--sum $SUBJECTS_DIR/Schaefer2018_100Parcels_7Networks_order_fsaverage6_labels/summary.txt \
--thmin 0.000001 \
--hemi rh \
--subject fsaverage6 \
The reason why I am deriving centroid coordinates from the cortical surface is because the ROI was originally derived on the cortical surface.
I have also tried deriving centroid coordinates after transforming the ROI into volumetric space. I used
fslstats on the volume with the
-C flag to obtain the coordinates. I pretty much receive the same results where the centroid is outside the ROI.
The ROI has been obtained from Schaefer et al. 2018 Cerebral Cortex 100 area parcellation (https://github.com/ThomasYeoLab/CBIG/tree/master/stable_projects/brain_parcellation/Schaefer2018_LocalGlobal).
Any advice or suggestions would be helpful, thank you!
This appears to be a case where the geometric center (centroid) of your ROI exists outside of the region actually segmented – this is due to the non-convex morphology of the cortical surface label. Thus the centroid measurements you are obtaining from the two different tools are actually correct. Further answers to your question may depend on what you hope to accomplish with the centroid coordinate.
All the best,
Thanks for your reply @jclauneuro, I intend to have the centroid inside the segmented region, not outside. From the centroid coordinate, I need to create uniform 4-6mm spheres. I will be doing this across multiple such ROIs across the brain.
One idea that was suggested to me was to first make the surface “flat” and then obtain the center of mass, while this makes sense I don’t know how to execute it. Any suggestions or tips would be appreciated.
The centroid is outside of the intended region because your cortical parcellation captures both mesial and lateral aspects of the superior frontal gyrus – so your centroid ends up in the WM in between.
This could be resolved by subparcellating the ROI. One option would be to do this manually, but probably the best option would be to use one of the more densely parcellated ROIs available at the link you provided (I see they have a 400 parcel version and up to 1000 parcels).
Unless you’re stuck using the 100 parcel version, you should consider using either the 400 or 1000 parcel one.
I’m not familiar with how to do this on a flattened cortical surface but I’m sure someone on this forum has experience with it.
Sounds to me like you would want to extract the medoid, i.e. the voxel in the parcel that minimizes the average distance with other voxels inside the parcel. That voxel belongs by definition to the parcel, unlike the centroid.
Thanks @pbellec, that sounds exactly right. However, I haven’t found an algorithm or tool to perform this appropriately. I looked at the kmedoids matlab function (https://www.mathworks.com/help/stats/kmedoids.html) which takes n-by-p 2D matrices as input. Voxel coordinates are 3D, so I could potentially iterate this several times (let’s say I do it along each z slice for x-by-y matrices, and then re-iterate it over another plane, and find the coordinate that converges?)…not sure if this would work. Does freesurfer or any neuroimaging tool have a medoid search tool?
Advice appreciated, thank you!
I don’t know of any tools currently written for neuroimages that will give you the medioid, but there is recent work that includes what appears to be a fairly easily implementable algorithm that beats O(n^2).
If you want to find the centroid by geodesic distance along the mesh, as opposed to Euclidean distances in the volume, you could find an all-to-all shortest-path-length matrix (
networkx has options), take the mean along one axis, and then find the minimum.