ThomasYeoLab / CBIG

MIT License
579 stars 382 forks source link

Network ordering of Kong2018_MSHBM #4

Open rubykong opened 5 years ago

rubykong commented 5 years ago

We are trying to release a few scripts to provide the network label with the network name. The only problem might be that the network structure could be very different if you are training from scratch on a different dataset. Here are two cases:

Use our estimated group priors to estimate individual parcellation

The group atlas will be utilized to initialize the algorithm, therefore, the network order in the individual parcellation will be the same as the group parcellation. We provided two group priors which are generated by a group atlas estimated by GSP dataset in fsaverage5 space, and a group atlas estimated by HCP dataset in fsLR_32k space. I will provide these group atlases with corresponding network labels.

Use your own data to train the model and estimate group priors

If you train on your own data, then your group atlas will be different from ours, the network ordering will also be different. In this way, you will have to run a Hungarian match algorithm to match the networks in your estimated group atlas with our atlases. To do that, you can use the following code:

1) If it's fsaverage5, you can directly use the following command to reorder the color table and visualize it:

my_par = load('<output_dir>/group/group.mat');
ref_par = load('/Kong2019_MSHBM/examples/results/estimate_group_priors/group/group.mat');
[output, assign, cost, dice_overlap] = CBIG_HungarianClusterMatch([ref_par.lh_labels; ref_par.rh_labels], [my_par.lh_labels; my_par.rh_labels], 1);
colors_old=ref_par.colors(2:end,:);
colors_new=colors_old(assign,:);
colors_new=[[0,0,0];color_new];
CBIG_DrawSurfaceMaps(my_par.lh_labels, my_par.rh_labels, 'fsaverage5', 'inflated',-Inf,Inf,colors_new);

2) If your data is in fslr, sorry, I didn't release the HCP group atlas in the current release (will do it asap), so contact me (roo.cone@gmail.com) to get the group atlas file HCP_40sub_1000iter_17nets_cen_sm4.mat.

my_par = load(fullfile(output_dir,'group','group.mat'));
ref_par = load('HCP_40sub_1000iter_17nets_cen_sm4.mat');
[output, assign, cost, dice_overlap] = CBIG_HungarianClusterMatch([ref_par.lh_labels; ref_par.rh_labels], [my_par.lh_labels; my_par.rh_labels], 1);
colors_old=ref_par.colors(2:end,:);
colors_new=colors_old(assign,:);
CBIG_DrawSurfaceMaps_fslr(my_par.lh_labels, my_par.rh_labels, 'fs_LR_32k', 'inflated',-Inf,Inf,colors_new);

I will try to release all above commands and other useful scripts later.

mwaskom commented 5 years ago

Thanks for making the code for this method open and so easy to use! I have a question that is related to this issue, although a bit distinct. While the code you have posted here is useful in the case of matching an experiment-specific set of priors to the GSP parcellation, I believe that the order of the clusters in the GSP prior set itself doesn't match those in the original Yeo2011 publication.

I refer to the priors in stable_projects/brain_parcellation/Kong2019_MSHBM/lib/group_priors/GSP_37/Params_Final.mat

And likewise, the provided group parcellation in stable_projects/brain_parcellation/Kong2019_MSHBM/examples/results/estimate_group_priors/group/group.mat.

The latter mat file does contain a colormap that properly plots the parcellation with the expected colors. But the labels don't match what I expect.

If I want an individual-specific parcellation with labels that match Yeo2011, should I manually implement a matching algorithm based on the colors / eye / some other system? Or am I missing something in the CBIG code that will transform the labels for me?

rubykong commented 5 years ago

@mwaskom The color table in stable_projects/brain_parcellation/Kong2019_MSHBM/examples/results/estimate_group_priors/group/group.mat matches with the label name in Yeo2011 parcellation. For your reference, the network name of GSP priors is:

network_name = {'DorsalAttentionA','DefaultB','SomatomotorA','DefaultA','TemporalParietal','VisualB','VentralAttentionA','DorsalAttentionB','LimbicB','VentralAttentionB','ControlA','VisualA','SomatomotorB','ControlC','ControlB','LimbicA','DefaultC'};
peter3200 commented 2 years ago

Thank you for providing these excellent resources! We would like to perform the Hungarian matching algorithm on individual parcellations in fsaverage6 resolution. Are there group priors available in fsaverage6 space? Or alternatively, how would I go about resampling the fsaverage5 GSP priors to fsaverage6 space?

Thanks again!

rubykong commented 2 years ago

@peter3200 We use nearest neighbor algorithm to upsample parcellation from fsaverage5 to fsaverage6. However, this will cause some medial wall issues. Specifically, if we upsample parcellation from fs5 to fs6, we will get a upsampled medial wall area in fs6, which differs from medial wall area of fs6 template defined by Freesurfer. The difference is not big, but you might be careful with it.

Here is the code for upsampling:

lh_FS_data = your_input_fs5_parcellation_left;
rh_FS_data = your_input_fs5_parcellation_right;
FS_mesh = "fsaverage5";

if(size(lh_FS_data,1) ~= 1);
    lh_FS_data = lh_FS_data';
end
if(size(rh_FS_data,1) ~= 1);
    rh_FS_data = rh_FS_data';
end 
lh_mesh6 = CBIG_ReadNCAvgMesh('lh', 'fsaverage6', 'sphere', 'cortex');
rh_mesh6 = CBIG_ReadNCAvgMesh('rh', 'fsaverage6', 'sphere', 'cortex');
lh_FS_mesh = CBIG_ReadNCAvgMesh('lh', FS_mesh, 'sphere', 'cortex');
rh_FS_mesh = CBIG_ReadNCAvgMesh('rh', FS_mesh, 'sphere', 'cortex');

lh_data6 = MARS_NNInterpolate_kdTree(lh_mesh6.vertices, lh_FS_mesh, lh_FS_data);
rh_data6 = MARS_NNInterpolate_kdTree(rh_mesh6.vertices, rh_FS_mesh, rh_FS_data);