How to use the brainspace.plotting.plot_hemispheres function in brainspace

# PCA, Laplacian eigenmaps and diffusion mapping
embeddings = ['pca', 'le', 'dm']

gradients_embedding = [None] * len(embeddings)
for i, emb in enumerate(embeddings):
    gm = GradientMaps(kernel='normalized_angle', approach=emb, random_state=0)
    gm.fit(conn_matrix)

    gradients_embedding[i] = map_to_labels(gm.gradients_[:, 0], labeling, mask=mask,
                                           fill=np.nan)


# sphinx_gallery_thumbnail_number = 2
label_text = ['PCA', 'LE', 'DM']
plot_hemispheres(surf_lh, surf_rh, array_name=gradients_embedding, size=(1200, 800),
                 cmap='viridis_r', color_bar=True, label_text=label_text)


###############################################################################
# Gradient alignment
# +++++++++++++++++++
#
# A more principled way of increasing comparability across gradients are
# alignment techniques. BrainSpace provides two alignment techniques:
# Procrustes analysis, and joint alignment. For this example we will load
# functional connectivity data of a second subject group and align it with the
# first group.

conn_matrix2 = load_group_fc('schaefer', scale=400, group='holdout')
gp = GradientMaps(kernel='normalized_angle', alignment='procrustes')
gj = GradientMaps(kernel='normalized_angle', alignment='joint')

# Compute parcellation centroids and append to spheres
aop.get_parcellation_centroids(sphere_lh, parcellation_lh, mask=mask_lh,
                               non_centroid=0, append=True,
                               array_name='centroids')
aop.get_parcellation_centroids(sphere_rh, parcellation_rh, mask=mask_rh,
                               non_centroid=0, append=True,
                               array_name='centroids')

mask_centroids_lh = sphere_lh.get_array('centroids') > 0
mask_centroids_rh = sphere_rh.get_array('centroids') > 0

centroids_lh = sphere_lh.Points[mask_centroids_lh]
centroids_rh = sphere_lh.Points[mask_centroids_rh]

# We can see the centroids on the sphere surfaces
plot_hemispheres(sphere_lh, sphere_rh, array_name='centroids',
                 interactive=False, embed_nb=True, size=(800, 200),
                 cmap_name='binary')


###############################################################################
# Now, let's generate 2000 random samples using spin permutations.

from brainspace.null_models import SpinRandomization

n_spins = 2000
sp = SpinRandomization(n_rep=n_spins, random_state=0)
sp.fit(centroids_lh, points_rh=centroids_rh)
gradient_spins_lh, gradient_spins_rh = sp.randomize(gradient[:200],
                                                    x_rh=gradient[200:])

# Let's create some rotations
n_permutations = 1000

sp = SpinPermutations(n_rep=n_permutations, random_state=0)
sp.fit(sphere_lh, points_rh=sphere_rh)

t1wt2w_rotated = np.hstack(sp.randomize(t1wt2w_lh, t1wt2w_rh))
thickness_rotated = np.hstack(sp.randomize(thickness_lh, thickness_rh))


###############################################################################
# As an illustration of the rotation, let’s plot the original t1w/t2w data

# Plot original data
plot_hemispheres(surf_lh, surf_rh, array_name=t1wt2w, size=(1200, 300), cmap='viridis',
                 nan_color=(0.5, 0.5, 0.5, 1), color_bar=True)


###############################################################################
# as well as a few rotated versions.

# sphinx_gallery_thumbnail_number = 2
# Plot some rotations
plot_hemispheres(surf_lh, surf_rh, array_name=t1wt2w_rotated[:3], size=(1200, 800),
                 cmap='viridis', nan_color=(0.5, 0.5, 0.5, 1), color_bar=True,
                 label_text=['Rot0', 'Rot1', 'Rot2'])


###############################################################################
#
# .. warning::

from brainspace.datasets import load_group_fc, load_parcellation, load_conte69

# First load mean connectivity matrix and Schaefer parcellation
conn_matrix = load_group_fc('schaefer', scale=400)
labeling = load_parcellation('schaefer', scale=400, join=True)

# and load the conte69 surfaces
surf_lh, surf_rh = load_conte69()


###############################################################################
# Let’s first look at the parcellation scheme we’re using.

from brainspace.plotting import plot_hemispheres

plot_hemispheres(surf_lh, surf_rh, array_name=labeling, size=(1200, 200),
                 cmap='tab20', zoom=1.85)


###############################################################################
# and let’s construct our gradients.

from brainspace.gradient import GradientMaps

# Ask for 10 gradients (default)
gm = GradientMaps(n_components=10, random_state=0)
gm.fit(conn_matrix)


###############################################################################
# Note that the default parameters are diffusion embedding approach, 10
# components, and no kernel (use raw data). Once you have your gradients, a

from brainspace.datasets import load_group_fc, load_parcellation, load_conte69

# First load mean connectivity matrix and Schaefer parcellation
conn_matrix = load_group_fc('schaefer', scale=400)
labeling = load_parcellation('schaefer', scale=400, join=True)

# and load the conte69 surfaces
surf_lh, surf_rh = load_conte69()


###############################################################################
# Let’s first look at the parcellation scheme we’re using.

from brainspace.plotting import plot_hemispheres

plot_hemispheres(surf_lh, surf_rh, array_name=labeling, size=(1200, 200),
                 cmap='tab20', zoom=1.85)


###############################################################################
# and let’s construct our gradients.

from brainspace.gradient import GradientMaps

# Ask for 10 gradients (default)
gm = GradientMaps(n_components=10, random_state=0)
gm.fit(conn_matrix)


###############################################################################
# Note that the default parameters are normalized angle kernel, diffusion
# embedding approach, 10 components. Once you have your gradients, a good first

# Visualize results
from brainspace.datasets import load_fsa5
from brainspace.plotting import plot_hemispheres
from brainspace.utils.parcellation import map_to_labels

# Map gradients to original parcels
grad = [None] * 2
for i, g in enumerate(gm.gradients_.T):
    grad[i] = map_to_labels(g, labeling, mask=mask, fill=np.nan)


# Load fsaverage5 surfaces
surf_lh, surf_rh = load_fsa5()

# sphinx_gallery_thumbnail_number = 2
plot_hemispheres(surf_lh, surf_rh, array_name=grad, size=(1200, 400), cmap='viridis_r',
                 color_bar=True, label_text=['Grad1', 'Grad2'], zoom=1.5)

# components, and no kernel (use raw data). Once you have your gradients, a
# good first step is to simply inspect what they look like. Let’s have a look
# at the first two gradients.

import numpy as np

from brainspace.utils.parcellation import map_to_labels

mask = labeling != 0

grad = [None] * 2
for i in range(2):
    # map the gradient to the parcels
    grad[i] = map_to_labels(gm.gradients_[:, i], labeling, mask=mask, fill=np.nan)

plot_hemispheres(surf_lh, surf_rh, array_name=grad, size=(1200, 400), cmap='viridis_r',
                 color_bar=True, label_text=['Grad1', 'Grad2'], zoom=1.55)


###############################################################################
# But which gradients should you keep for your analysis? In some cases you may
# have an a priori interest in some previously defined set of gradients. When
# you do not have a pre-defined set, you can instead look at the lambdas
# (eigenvalues) of each component in a scree plot. Higher eigenvalues (or lower
# in Laplacian eigenmaps) are more important, so one can choose a cut-off based
# on a scree plot.

import matplotlib.pyplot as plt

fig, ax = plt.subplots(1, figsize=(5, 4))
ax.scatter(range(gm.lambdas_.size), gm.lambdas_)
ax.set_xlabel('Component Nb')

# embedding approach, 10 components. Once you have your gradients, a good first
# step is to simply inspect what they look like. Let’s have a look at the first
# two gradients.

import numpy as np

from brainspace.utils.parcellation import map_to_labels

mask = labeling != 0

grad = [None] * 2
for i in range(2):
    # map the gradient to the parcels
    grad[i] = map_to_labels(gm.gradients_[:, i], labeling, mask=mask, fill=np.nan)

plot_hemispheres(surf_lh, surf_rh, array_name=grad, size=(1200, 400), cmap='viridis_r',
                 color_bar=True, label_text=['Grad1', 'Grad2'], zoom=1.55)


###############################################################################
# But which gradients should you keep for your analysis? In some cases you may
# have an a priori interest in some previously defined set of gradients. When
# you do not have a pre-defined set, you can instead look at the lambdas
# (eigenvalues) of each component in a scree plot. Higher eigenvalues (or lower
# in Laplacian eigenmaps) are more important, so one can choose a cut-off based
# on a scree plot.

import matplotlib.pyplot as plt

fig, ax = plt.subplots(1, figsize=(5, 4))
ax.scatter(range(gm.lambdas_.size), gm.lambdas_)
ax.set_xlabel('Component Nb')

n_pts_lh = surf_lh.n_points
n_parcellations = len(list_parcels)

# first we are going to append the parcellation to the hemispheres
name_parcels = [None] * n_parcellations
for i, np in enumerate(list_parcels):
    array_name = 'Schaefer{0}'.format(np)
    name_parcels[i] = array_name

    surf_lh.append_array(labelings[i][:n_pts_lh], name=array_name, at='p')
    surf_rh.append_array(labelings[i][n_pts_lh:], name=array_name, at='p')

# Then plot the data on the surface
plot_hemispheres(surf_lh, surf_rh, array_name=name_parcels, interactive=False,
                 embed_nb=True, size=(800, 800), cmap_name='tab20')


###############################################################################
# We have 4 mean connectivity matrices built from each parcellation.

import matplotlib.pyplot as plt

fig, ax = plt.subplots(1, n_parcellations, figsize=(10, 25))

# The mean connectivity matrix built from the HCP data for each parcellation
for i in range(n_parcellations):
    ax[i].imshow(conn_matrices[i], cmap='hot', interpolation='bilinear')

    ax[i].set_title(name_parcels[i])
    ax[i].xaxis.set_visible(False)

# For each embedding approach, we compute 2 gradients and append them to the
# left and right surfaces. Note that we have used 'normalized_angle' to build
# the affinity matrix.
#
# Next, for all embedding approaches, we display the first gradient

array_names = ['pca_grad0', 'le_grad0', 'dm_grad0']
plot_hemispheres(surf_lh, surf_rh, array_name=array_names, interactive=False,
                 embed_nb=True, size=(800, 600), cmap_name='viridis')


###############################################################################
# And the second gradient

array_names = ['pca_grad1', 'le_grad1', 'dm_grad1']
plot_hemispheres(surf_lh, surf_rh, array_name=array_names, interactive=False,
                 embed_nb=True, size=(800, 600), cmap_name='viridis')