rmldj / hcp-utils

Utilities to use HCP and HCP-like data with nilearn and other Python tools
MIT License
41 stars 7 forks source link

Generating brain figure of subset of parcellation #2

Closed amsinha9 closed 4 years ago

amsinha9 commented 4 years ago

Is there a way to specify what functional networks to color/highlight in a brain figure using hcp.view_parcellation(mesh_sub.inflated, hcp.ca_network)? For example, if analyses identify specific networks within the Cole-Anticevic Network partition that we want to illustrate (say the first 5 out of the 12 networks), is there currently a way to do this using hcp-utils?

Similarly, is there a way to specify parcels from a given parcellation that we want to highlight in generating brain figures?

amsinha9 commented 4 years ago

More explicitly stated, is there a way to use hcp.view_parcellation() and specify parcel/network indices to only highlight a subset of parcels/regions instead of showing all 379 parcels/12 networks?

I would like to use this toolbox for generating figures in a manuscript, so any help would be greatly appreciated.

rmldj commented 4 years ago

If you want to highlight/show only some signal from a specific region you can do it using the following snippet from the docs:

plotting.view_surf(mesh_sub.inflated, 
    hcp.cortex_data(hcp.mask(Xn[29], hcp.yeo7.map_all==2)), 
    threshold=0.1, bg_map=mesh_sub.sulc)

For the moment unfortunately hcp.view_parcellation() does not allow for showing only a subset of regions/networks. I will think about how to do it...

amsinha9 commented 4 years ago

Thanks for the response. Given that the field is moving away from volume-based ROIs to surface-based parcellations, it would be extremely helpful if there was a way we could give a list of network/parcel indices of a parcellation scheme to highlight/color networks/parcels found to be significant, perturbed, etc. as opposed to providing ROI coordinates and a corresponding adjacency matrix to generate brain figures of brain connections as in the Matlab toolbox, BrainNet Viewer, and most other brain visualization tools used today.

rmldj commented 4 years ago

For the moment you can try something like this:

import numpy as np
import matplotlib
import nilearn.plotting as plotting

def view_parcellation(meshLR, parcellation, subset=None, color=(0.9,0.9,0.9,1)):
    """
    View the given parcellation on an a whole brain surface mesh.
    """
    # for some parcellations the numerical ids need not be consecutive
    cortex_map = hcp.cortex_data(parcellation.map_all)
    ids = np.unique(cortex_map)
    normalized_cortex_map = np.zeros_like(cortex_map)
    rgba = np.zeros((len(ids), 4))
    for i in range(len(ids)):
        ind = cortex_map==ids[i]
        normalized_cortex_map[ind] = i
        if subset is None or i in subset:
            rgba[i,:] = parcellation.rgba[ids[i]]
        else:
            rgba[i,:] = color

    cmap = matplotlib.colors.ListedColormap(rgba)
    return plotting.view_surf(meshLR, normalized_cortex_map, symmetric_cmap=False, cmap=cmap)

Then you can write e.g.

view_parcellation(hcp.mesh.inflated, hcp.yeo7, subset=[2,4])

which gives you view_parcellation_with_subset

amsinha9 commented 4 years ago

Yes that works really well! Thanks for coming up with a workaround solution, I appreciate it.

amsinha9 commented 4 years ago

Quick follow-up, with:

hcp.view_parcellation(hcp.mesh.inflated, hcp.mmp)

we can generate a mesh brain figure that depicts all 379 parcels from Glasser parcellation. However, to visualize only parcels 100, 200 and 300 within hcp.mmp using the above function, view_parcellation:

view_parcellation(hcp.mesh.inflated, hcp.mmp, subset=[100,200, 300])

the legend to the right of the brain figure only shows one of the corresponding parcel label colors.

newplot

Is there a way to resolve this?

rmldj commented 4 years ago

This is really an issue within the surface plotting code from nilearn. But this is not really a bug there - just as for mmp there are 379 distinct colors, and due to the limited pixel size of the colorbar (which I guess is beyond control) only a subset of colors can appear - and which ones do is governed by rounding.

So I guess that the only viable solution would be to construct a new "subset" parcellation which would have just the three regions and the rest of the brain unassigned. However in order to cure the above issue one would have to change the numerical id's so that in the subset parcellation 100->1, 200->2, 300->3.

rmldj commented 4 years ago

I guess including such a subset parcellation would make sense in an update of the package...

For the moment you can use the following function

from sklearn.utils import Bunch

def make_subset_parcellation(parcellation, subset):
    """
    Takes the given parcellation and produces a new one where only parcels in the subset are retained.
    Remaining part of the brain is set to 0 (unassigned).
    """
    ids = np.arange(len(subset)+1)
    nontrivial_ids = np.arange(1, len(subset)+1)
    # labels
    labels = dict()
    labels[0] = 'None'
    for i, c in enumerate(subset):
        labels[i+1] = parcellation.labels[c]
    # colors
    rgba = dict()
    rgba[0] = np.array([0.9,0.9,0.9,1])
    for i, c in enumerate(subset):
        rgba[i+1] = parcellation.rgba[c]
    # mapping
    map_all = np.zeros_like(parcellation.map_all)
    for i, c in enumerate(subset):
        map_all[parcellation.map_all==c] = i+1

    new_parcellation = Bunch()
    new_parcellation.map_all = map_all
    new_parcellation.labels = labels
    new_parcellation.ids = ids
    new_parcellation.nontrivial_ids = nontrivial_ids
    new_parcellation.rgba = rgba

    return new_parcellation

Now once you do

submmp = make_subset_parcellation(hcp.mmp, [100, 200, 300])

you can use the regular visualization

hcp.view_parcellation(hcp.mesh.inflated, submmp)

You can also display the labels through

hcp.parcellation_labels(submmp)
amsinha9 commented 4 years ago

The make_subset_parcellation function works well for this purpose, thank you! Yes, given that we often report findings on network connectivity changes/differences between subsets of regions in a given parcellation, it would be a helpful tool to incorporate into hcp-utils, however, make_subset_parcellation does the job well.

Thanks again.