Empirical data: evaluate analyses

We have previously analyzed empirical datasets with CCA and PLS. Here, we illustrate the outcomes.

Note: These empirical datasets are not included in this release.

Setup

import os

import numpy as np
import xarray as xr
import pandas as pd

import pickle

from gemmr.sample_analysis.macros import *

import matplotlib
import holoviews as hv
from holoviews import opts
hv.extension('matplotlib')
hv.renderer('matplotlib').set_param(dpi=120)

from my_config import *

import warnings
from matplotlib import MatplotlibDeprecationWarning
warnings.simplefilter('ignore', MatplotlibDeprecationWarning)
warnings.filterwarnings(
    'ignore', 'aspect is not supported for Axes with xscale=log, yscale=linear', category=UserWarning
)  # holoviews emits this for log-linear plots

clr_cca2 = 'firebrick'
clr_pls2 = 'steelblue'

clr_cca3 = 'coral'
clr_pls3 = 'slateblue'

clr_cca4 = 'goldenrod'
clr_pls4 = 'mediumseagreen'

clr_perm1 = 'dimgrey'
clr_perm2 = 'black'

hcp_fMRI_cca_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/hcp_fMRI_cca_analyses.pickle'), 'rb'))
hcp_fMRI_pls_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/hcp_fMRI_pls_analyses.pickle'), 'rb'))

hcp_dMRI_cca_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/hcp_dMRI_cca_analyses.pickle'), 'rb'))
hcp_dMRI_pls_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/hcp_dMRI_pls_analyses.pickle'), 'rb'))

ukbb_fMRI_cca_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/ukbb_fMRI_cca_analyses.pickle'), 'rb'))
ukbb_fMRI_pls_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/ukbb_fMRI_pls_analyses.pickle'), 'rb'))

hcp_fMRI_optiNpcs_cca_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/hcp_fMRI_optiNpcs_cca_analyses.pickle'), 'rb'))
hcp_fMRI_optiNpcs_pls_analyses = pickle.load(open(os.path.expanduser('~/gemmr_data/empirical/hcp_fMRI_optiNpcs_pls_analyses.pickle'), 'rb'))

Plotting functions

class Artist_edgecolor:
    def __init__(self, color):
        self.color = color
    def __call__(self, plot, element):
        artist = plot.handles['artist']
        artist.set_edgecolor(self.color)
        artist.set_linewidth(1.5)


def plot_assocs_w_perm(analysis_results, qs=(.025, .975), y_label='y', color='steelblue',
                       color_perm=None, x_max=None, y_scale=1, cv=None, show_pvalue=False,
                       pvalue_deltay=0):

    n = analysis_results['whole_sample'].attrs['n_samples']
    p_tot = len(analysis_results['whole_sample'].x_feature) + len(analysis_results['whole_sample'].y_feature)

    ds = analysis_results['whole_sample']
    ds_subsampled = analysis_results['subsampled']
    p_value = analysis_results['p_value']

    between_assocs_perm_q = ds.between_assocs_perm.quantile(qs, 'perm')

    subsampled_between_assocs_mean = ds_subsampled.between_assocs.mean('rep')

    subsampled_between_assocs_q = ds_subsampled.between_assocs.quantile(qs, 'rep')
    subsampled_between_assocs_perm_q = ds_subsampled.between_assocs_perm.stack(it=('rep', 'perm')).quantile(qs, 'it')

    if x_max is None:
        x_max = n
    y_max = max(
        float(between_assocs_perm_q.max()),
        float(subsampled_between_assocs_q.max()),
        float(subsampled_between_assocs_perm_q.max())
    )

    panel = (
        hv.Area(
            (subsampled_between_assocs_perm_q.n / p_tot,
             y_scale * subsampled_between_assocs_perm_q.sel(quantile=qs[0]),
             y_scale * subsampled_between_assocs_perm_q.sel(quantile=qs[-1])),
            vdims=['y', 'y2'],
            label='perm'
        ).opts(color=clr_perm1, linewidth=1.5, hooks=[Artist_edgecolor(color)])
        * hv.Area(
            (subsampled_between_assocs_q.n / p_tot,
             y_scale * subsampled_between_assocs_q.sel(quantile=qs[0]),
             y_scale * subsampled_between_assocs_q.sel(quantile=qs[-1])),
            vdims=['y', 'y2'],
            label='data'
        ).opts(color=color, linewidth=0)
        * hv.Curve(
            (subsampled_between_assocs_mean.n / p_tot, y_scale * subsampled_between_assocs_mean),
            vdims='y'
        ).opts(color=color)
        * hv.ErrorBars(
            (n / p_tot,
             y_scale * float(between_assocs_perm_q.mean()),
             y_scale * float(between_assocs_perm_q.sel(quantile=qs[-1]) - between_assocs_perm_q.mean())
            )
        ).opts(color=clr_perm1, hooks=[Artist_edgecolor(color)])
        * hv.Scatter(
            ([n / p_tot], [y_scale * ds.between_assocs])
        ).opts(color=color)
    )

    if show_pvalue:
        panel *= hv.Text(x_max / p_tot, y_scale * y_max + pvalue_deltay, r'$p=%.3f$' % p_value, halign='right', valign='top', fontsize=8)

    if cv is not None:

        if color_perm is None:
            color_perm = color

        assoc_label, cv_method = cv

        between_assocs_cv_perm_q = ds['between_{}_cv_perm'.format(assoc_label)].sel(mode=0, cv=cv_method).quantile(qs, 'perm')

        subsampled_between_assocs_cv = ds_subsampled['between_{}_cv'.format(assoc_label)].sel(mode=0, cv=cv_method)
        subsampled_between_assocs_cv_mean = subsampled_between_assocs_cv.mean('rep')
        subsampled_between_assocs_cv_q = subsampled_between_assocs_cv.quantile(qs, 'rep')

        subsampled_between_assocs_cv_perm = ds_subsampled['between_{}_cv_perm'.format(assoc_label)].sel(mode=0, cv=cv_method)
        subsampled_between_assocs_cv_perm_q = subsampled_between_assocs_cv_perm.stack(it=('rep', 'perm')).quantile(qs, 'it')

        panel *= (
            hv.Area(
                (subsampled_between_assocs_cv_perm_q.n / p_tot,
                 y_scale * subsampled_between_assocs_cv_perm_q.sel(quantile=qs[0]),
                 y_scale * subsampled_between_assocs_cv_perm_q.sel(quantile=qs[-1])),
                vdims=['y', 'y2'],
                label='perm'
            ).opts(color=clr_perm1, linewidth=1.5, hooks=[Artist_edgecolor(color_perm)])
            * hv.Area(
                (subsampled_between_assocs_cv_q.n / p_tot,
                 y_scale * subsampled_between_assocs_cv_q.sel(quantile=qs[0]),
                 y_scale * subsampled_between_assocs_cv_q.sel(quantile=qs[-1])),
                vdims=['y', 'y2'],
                label='cv'
            ).opts(color=color_perm, linewidth=0)
            * hv.Curve(
                (subsampled_between_assocs_cv_mean.n / p_tot, y_scale * subsampled_between_assocs_cv_mean),
                vdims='y'
            ).opts(color=color_perm)
            * hv.ErrorBars(
                (n / p_tot,
                 y_scale * float(between_assocs_cv_perm_q.mean()),
                 y_scale * float(between_assocs_cv_perm_q.sel(quantile=qs[-1]) - between_assocs_cv_perm_q.mean())
                )
            ).opts(color=clr_perm1, hooks=[Artist_edgecolor(color_perm)])
            * hv.Scatter(
                ([n / p_tot], [y_scale * ds['between_{}_cv'.format(assoc_label)].sel(mode=0, cv=cv_method)])
            ).opts(color=color_perm, marker='^')
        )

    return panel.redim(
        x='Samples per feature',
        y=y_label,
        y2=y_label + '2'
    ).opts(
        opts.Area(alpha=.33, show_legend=True),
        opts.ErrorBars(alpha=.33),
        opts.Overlay(padding=.02),
    )

def plot_weight_cossims_n(ds_subsampled, color, qs=(.025, .975)):

    p_tot = len(ds_subsampled.x_feature) + len(ds_subsampled.y_feature)

    #xy_weight_sims_mean, xy_weight_sims_q = pairwise_weight_cosine_similarity(ds_subsampled, qs=qs)
    x_cossim_weights_stats = ds_subsampled.x_weights_pairwise_cossim_stats
    y_cossim_weights_stats = ds_subsampled.y_weights_pairwise_cossim_stats
    xy_weight_sims_mean = np.minimum(
        x_cossim_weights_stats.sel(stat='mean'),
        y_cossim_weights_stats.sel(stat='mean'),
    )
    xy_weight_sims_q025 = np.minimum(
        x_cossim_weights_stats.sel(stat='q2.5%'),
        y_cossim_weights_stats.sel(stat='q2.5%'),
    )
    xy_weight_sims_q975 = np.maximum(
        x_cossim_weights_stats.sel(stat='q97.5%'),
        y_cossim_weights_stats.sel(stat='q97.5%'),
    )

    vdim = hv.Dimension('weight_cossim', label='Weight stability')
    return (
        hv.Area(
            (ds_subsampled.n / p_tot, xy_weight_sims_q025, xy_weight_sims_q975),
            kdims='Samples per feature',
            vdims=[vdim, 'y2']
        )
        * hv.Curve(
            (xy_weight_sims_mean.n / p_tot, xy_weight_sims_mean),
            kdims='Samples per feature',
            vdims=vdim
        )
    ).opts(
        opts.Area(linewidth=0, color=color, alpha=.33),
        opts.Curve(color=color),
        opts.Overlay(logx=True)
    )

def twinx_scores_spearmansims(plot, element):
    ax = plot.handles['axis']
    ax.spines['top'].set_visible(False)
    try:
        color = plot.style.options['color']
    except:
        print("atist not found")
        color='black'
    twinax = ax.twinx()
    twinax.spines['top'].set_visible(False)
    twinax.set_ylim(-.1, 1)
    twinax.set_yticks([0, 0.2, 0.4, 0.6, 0.8, 1.0])
    twinax.set_yticklabels(['0.0', '0.2', '0.4', '0.6', '0.8', '1.0'], fontdict=dict(family='Helvetica', size=7, weight='light'))
    twinax.set_ylabel('Score stability', fontdict=dict(family='Helvetica', size=8, weight='normal', color=color))
    plot.handles['axis'] = twinax


def plot_score_spearmansims_n(ds, color):
    p_tot = len(ds.x_feature) + len(ds.y_feature)
    panel = (
        hv.Curve(
            (
                ds.n / p_tot,
                np.minimum(
                    np.abs(ds.x_scores_pairwise_spearmansim_stats.sel(stat='mean', mode=0)),
                    np.abs(ds.x_scores_pairwise_spearmansim_stats.sel(stat='mean', mode=0))
                )
            ),
            kdims='Samples per feature'
        ),
    ).opts(
        opts.Curve(initial_hooks=[twinx_scores_spearmansims], color=color)
    )
    return panel

def twinx_pc_bias(plot, element):
    ax = plot.handles['axis']
    ax.spines['top'].set_visible(False)
    try:
        color = plot.style.options['color']
    except:
        print("atist not found")
        color='black'
    twinax = ax.twinx()
    twinax.spines['top'].set_visible(False)
    twinax.set_ylim(-.1, 1)
    twinax.set_yticks([0, 0.2, 0.4, 0.6, 0.8, 1.0])
    twinax.set_yticklabels(['0.0', '0.2', '0.4', '0.6', '0.8', '1.0'], fontdict=dict(family='Helvetica', size=7, weight='light'))
    twinax.set_ylabel('Weight PC1 bias', fontdict=dict(family='Helvetica', size=8, weight='normal', color=color))
    plot.handles['axis'] = twinax


def plot_pc_bias(ds, color):
    p_tot = len(ds.x_feature) + len(ds.y_feature)
    panel = (
        hv.Curve(
            (
                ds.n / p_tot,
                np.maximum(
                    np.abs(ds.x_weights_pc_cossim.sel(x_feature=0, mode=0)),
                    np.abs(ds.y_weights_pc_cossim.sel(y_feature=0, mode=0))
                ).mean('rep')
            ),
            kdims='Samples per feature',
            label='data'
        )
        * hv.Curve(
            (
                ds.n / p_tot,
                np.maximum(
                    np.abs(ds.x_weights_pc_cossim_perm.sel(x_feature=0, mode=0)),
                    np.abs(ds.y_weights_pc_cossim_perm.sel(y_feature=0, mode=0))
                ).mean('perm').mean('rep')
            ),
            kdims='Samples per feature',
            label='perm'
        ).opts(linestyle='--')
    ).opts(
        opts.Curve(initial_hooks=[twinx_pc_bias], color=color)
    )
    return panel

# holoviews hooks

class Sup_ylabel:
    def __init__(self, label):
        self.label = label
    def __call__(self, plot, element):
        ax = plot.handles['axis']
        ax.text(-.6, .5, self.label, transform=ax.transAxes, rotation=90, fontsize=9, fontweight='bold', ha='center', va='center')


def ylabel_clr_cca(plot, element):
    ax = plot.handles['axis']
    ax.yaxis.label.set_color(clr_cca4)


def ylabel_clr_pls(plot, element):
    ax = plot.handles['axis']
    ax.yaxis.label.set_color(clr_pls4)


class Labelpos:
    def __init__(self, dim, labelpad=-10):
        if dim not in ['x', 'y']:
            raise ValueError('invalid dim: {}'.format(dim))
        self.dim = dim
        self.labelpad = labelpad
    def __call__(self, plot, element):
        ax = plot.handles['axis']
        if self.dim == 'x':
            ax.set_xlabel(ax.get_xlabel(), labelpad=self.labelpad)
        if self.dim == 'y':
            ax.set_ylabel(ax.get_ylabel(), labelpad=self.labelpad)


def cca_assocStrength_legend(plot, element):
    ypos = -0.6
    ax = plot.handles['axis']
    ax.text(.275, ypos, 'in-sample', ha='right', color=clr_cca, fontdict=dict(size=8), transform=ax.transAxes)
    ax.text(.5, ypos, 'CV', ha='center', color=clr_cca2, fontdict=dict(size=8), transform=ax.transAxes)
    ax.text(1-.275, ypos, 'perm', ha='left', color=clr_perm1, fontdict=dict(size=8), transform=ax.transAxes)


def pls_assocStrength_legend(plot, element):
    ypos = -0.6
    ax = plot.handles['axis']
    ax.text(.275, ypos, 'in-sample', ha='right', color=clr_pls, fontdict=dict(size=8), transform=ax.transAxes)
    ax.text(.5, ypos, 'CV', ha='center', color=clr_pls2, fontdict=dict(size=8), transform=ax.transAxes)
    ax.text(1-.275, ypos, 'perm', ha='left', color=clr_perm1, fontdict=dict(size=8), transform=ax.transAxes)

Figures

literature_markers = np.array([
    (461, .87),  # Smith et al. (2015)
    (498, .91),  # Rahim et al. (2017)
    (819, .73),  # Bijsterbosch et al. (2018): Yeo parcellation(d=108) full correlation network-FC
    (1001, .67),  # Bijsterbosch et al. (2019): Tnet ICA thrs full correlation
    (1094, .667),  # Li et al. (2019)
    (818, .77),  # Han et al. (2020): FC
])

# convert to samples per feature
literature_markers[:, 0] /= 200

def set_axis_position(plot, element):
    ax = plot.handles['axis']
    bbox = ax.get_position()

    if bbox.x0 < .1:
        x0 = .05
    elif bbox.x0 < .35:
        x0 = .26
    elif bbox.x0 < .6:
        x0 = .55
    else:
        x0 = .77

    if bbox.y0 < .1:
        y0 = .05
    elif bbox.y0 < .5:
        y0 = .05 + .3
    else:
        y0 = .05 + .6

    ax.set_position((x0, y0, .14, .2))


fig = (
    # --- HCP fMRI ---
    (
        plot_assocs_w_perm(hcp_fMRI_cca_analyses, cv=('corrs', 'kfold5'), y_label='Observed corr.',
                           color=clr_cca, color_perm=clr_cca2, pvalue_deltay=.15).opts(opts.Area(show_legend=False))
        * hv.Scatter(literature_markers, label='literature').opts(marker='s', color='None', edgecolor='black', linewidth=.5, s=10)
        * hv.Text(550/200, .925, 'literature', halign='left', valign='center', fontsize=7)
        * hv.Text(210/200, .82, 'in-sample', fontsize=7, halign='left', valign='top').opts(color=clr_cca)
        * hv.Text(700/200, .7, 'perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
        * hv.Text(450/200, .325, 'CV', fontsize=7, halign='right', valign='top').opts(color=clr_cca2)
        * hv.Text(700/200, -.1, 'CV perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
    ).opts(
        xlim=(1, 5.75), xlabel='', ylim=(None, 1), yticks=7, legend_position='top_left', show_legend=False,
        hooks=[set_axis_position, Sup_ylabel('HCP fMRI'), legend_frame_off, suptitle_cca,
               Format_log_axis('x', (1, 4))]
    )
    + (
        plot_weight_cossims_n(hcp_fMRI_cca_analyses['subsampled'].sel(mode=0), clr_cca4)
        * plot_pc_bias(hcp_fMRI_cca_analyses['subsampled'], clr_cca3).opts(show_legend=True)
    ).opts(xlim=(1, 5.75), xlabel='', ylim=(-.01, 1), yticks=6,
           hooks=[set_axis_position, ylabel_clr_cca, legend_frame_off, suptitle_cca, Format_log_axis('x', (1, 4))])
    #
    + (
        plot_assocs_w_perm(hcp_fMRI_pls_analyses, cv=('covs', 'kfold5'), y_scale=1./9000, y_label='Observed cov. / a.u.', color=clr_pls, color_perm=clr_pls2)
        * hv.Text(350/200, .875, 'in-sample', fontsize=7, halign='left', valign='top').opts(color=clr_pls)
        # * hv.Text(700/200, .7, 'perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
        * hv.Text(250/200, .075, 'CV', fontsize=7, halign='left', valign='top').opts(color=clr_pls2)
        # * hv.Text(700/200, -.1, 'CV perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
    ).opts(
        xlim=(1, 5.75), xlabel='', ylim=(None, 1.0), ylabel='Observed cov. / a.u.', show_legend=False, legend_position='bottom_left',
        hooks=[set_axis_position, legend_frame_off, suptitle_pls, Format_log_axis('x', (1, 4))])
    + (
        plot_weight_cossims_n(hcp_fMRI_pls_analyses['subsampled'].sel(mode=0), clr_pls4)
        * plot_pc_bias(hcp_fMRI_pls_analyses['subsampled'], clr_pls3).opts(show_legend=True)
    ).opts(
        xlim=(1, 5.75), xlabel='', yticks=6, ylim=(0, 1), legend_position='bottom_right',
        hooks=[set_axis_position, ylabel_clr_pls, legend_frame_off, suptitle_pls, Format_log_axis('x', (1, 4))]
    )
    # ----------------
    # --- HCP dMRI ---
    # ----------------
    + plot_assocs_w_perm(hcp_dMRI_cca_analyses, cv=('corrs', 'kfold5'), y_label='Observed corr.',
                         color=clr_cca, color_perm=clr_cca2, pvalue_deltay=.1).opts(
        xlim=(1 , 5.5), xlabel='', ylim=(None, 1), yticks=7, show_legend=False,
        hooks=[Sup_ylabel('HCP dMRI'), set_axis_position, Format_log_axis('x', (1, 4)), ]
    )
    + (
        plot_weight_cossims_n(hcp_dMRI_cca_analyses['subsampled'].sel(mode=0), clr_cca4)
        * plot_pc_bias(hcp_dMRI_cca_analyses['subsampled'], clr_cca3).opts(show_legend=True)
    ).opts(
        xlim=(1, 5.5), xlabel='', ylim=(-.01, 1), yticks=6, show_legend=False,
        hooks=[set_axis_position, ylabel_clr_cca, Format_log_axis('x', (1, 4))]
    )
    #
    + plot_assocs_w_perm(hcp_dMRI_pls_analyses, cv=('covs', 'kfold5'), y_scale=1./2000, y_label='Observed cov. / a.u.', color=clr_pls, color_perm=clr_pls2).opts(
        xlim=(1, 5.5), xlabel='', ylim=(None, 1), ylabel='Observed cov. / a.u.', show_legend=False,
        hooks=[set_axis_position, Format_log_axis('x', (1, 4))]
    )
    + (
        plot_weight_cossims_n(hcp_dMRI_pls_analyses['subsampled'].sel(mode=0), clr_pls4)
        * plot_pc_bias(hcp_dMRI_pls_analyses['subsampled'], clr_pls3).opts(show_legend=True)
    ).opts(
        xlim=(1, 5.5), xlabel='', yticks=6, ylim=(0, 1), show_legend=False,
        hooks=[set_axis_position, ylabel_clr_pls, Format_log_axis('x', (1, 4))]
    )
    # -----------------
    # --- UKBB fMRI ---
    # -----------------
    + plot_assocs_w_perm(ukbb_fMRI_cca_analyses, cv=('corrs', 'kfold5'), y_label='Observed corr.', color=clr_cca, color_perm=clr_cca2).opts(
        xlim=(1, 125), ylim=(None, 1), yticks=7, show_legend=False,
        hooks=[Sup_ylabel('UKBB fMRI'), set_axis_position, Format_log_axis('x', major_numticks=4, minor_numticks=9)]
    )
    + (
        plot_weight_cossims_n(ukbb_fMRI_cca_analyses['subsampled'].sel(mode=0), clr_cca4)
        * plot_pc_bias(ukbb_fMRI_cca_analyses['subsampled'], clr_cca3).opts(show_legend=True)
    ).opts(
        xlim=(1, 125), ylim=(-.01, 1), yticks=6, show_legend=False,
        hooks=[set_axis_position, ylabel_clr_cca, Format_log_axis('x', major_numticks=4, minor_numticks=9)]
    )
    #
    + plot_assocs_w_perm(ukbb_fMRI_pls_analyses, cv=('covs', 'kfold5'), y_scale=1./12000,
                         y_label='Observed cov. / a.u.', color=clr_pls, color_perm=clr_pls2).opts(
        xlim=(1, 125), ylim=(None, 1), ylabel='Observed cov. / a.u.', show_legend=False,
        hooks=[set_axis_position, Format_log_axis('x', major_numticks=4, minor_numticks=9)]
    )
    + (
        plot_weight_cossims_n(ukbb_fMRI_pls_analyses['subsampled'].sel(mode=0), clr_pls4)
        * plot_pc_bias(ukbb_fMRI_pls_analyses['subsampled'], clr_pls3).opts(show_legend=True)
    ).opts(
        xlim=(1, 125), yticks=6, ylim=(0, 1), xaxis='bottom', show_legend=False,
        hooks=[set_axis_position, ylabel_clr_pls, Format_log_axis('x', major_numticks=4, minor_numticks=9)]
    )
).cols(
    4
).opts(*fig_opts).opts(
    opts.Overlay(logx=True, logy=False),
    opts.Layout(fig_inches=(7, 5), sublabel_position=(-.55, .9), vspace=1.25)
)

hv.save(fig, 'fig/fig5_empirical_data_analysis.pdf')

fig

print(hcp_fMRI_cca_analyses['p_value'], hcp_fMRI_pls_analyses['p_value'])
print(hcp_dMRI_cca_analyses['p_value'], hcp_dMRI_pls_analyses['p_value'])
print(ukbb_fMRI_cca_analyses['p_value'], ukbb_fMRI_pls_analyses['p_value'])

0.000999000999000999 0.005994005994005994
0.000999000999000999 0.000999000999000999
0.000999000999000999 0.000999000999000999

def set_axis_position(plot, element):
    ax = plot.handles['axis']
    bbox = ax.get_position()

    if bbox.x0 < .1:
        x0 = .05
    elif bbox.x0 < .35:
        x0 = .26
    elif bbox.x0 < .6:
        x0 = .55
    else:
        x0 = .77

    ax.set_position((x0, .25, .14, .35))


fig_optiNpcs = (
    (
        plot_assocs_w_perm(hcp_fMRI_optiNpcs_cca_analyses, cv=('corrs', 'kfold5'), y_label='Observed corr.', color=clr_cca,
                        color_perm=clr_cca2, show_pvalue=False)
        * hv.Text(150/88, .92, 'in-sample', fontsize=7, halign='left', valign='bottom').opts(color=clr_cca)
        * hv.Text(300/88, .7, 'perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
        * hv.Text(250/88, .265, 'CV', fontsize=7, halign='right', valign='bottom').opts(color=clr_cca2)
        * hv.Text(700/88, -.15, 'CV perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
    ).opts(
        xlim=(1, 1150/88), ylim=(None, 1), yticks=7, show_legend=False,
        hooks=[suptitle_cca, Sup_ylabel('HCP fMRI (optim. # PCs)'), set_axis_position, Labelpos('x', 4), Format_log_axis('x')]
    )
    + (
        plot_weight_cossims_n(hcp_fMRI_optiNpcs_cca_analyses['subsampled'].sel(mode=0), clr_cca4)
        * plot_pc_bias(hcp_fMRI_optiNpcs_cca_analyses['subsampled'], clr_cca3).opts(show_legend=True)
    ).opts(
        xlim=(1, 1150/88), ylim=(-.01, 1), yticks=6, show_legend=False,
        hooks=[suptitle_cca, set_axis_position, ylabel_clr_cca, Labelpos('x', 3), Format_log_axis('x')]
    )
    #
    + (
        plot_assocs_w_perm(hcp_fMRI_optiNpcs_pls_analyses, cv=('covs', 'kfold5'), y_scale=1./12000, y_label='Observed cov. / a.u.',
                            color=clr_pls, color_perm=clr_pls2, show_pvalue=False)
        * hv.Text(190/88, .875, 'in-sample', fontsize=7, halign='left', valign='top').opts(color=clr_pls)
        # * hv.Text(700, .7, 'perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
        * hv.Text(120/88, 0, 'CV', fontsize=7, halign='left', valign='top').opts(color=clr_pls2)
        # * hv.Text(700, -.1, 'CV perm', fontsize=7, halign='right', valign='top').opts(color=clr_perm1)
    ).opts(
        xlim=(1, 1150/88), ylim=(None, 1), ylabel='Observed cov. / a.u.', show_legend=False,
        hooks=[suptitle_pls, set_axis_position, Labelpos('x', 3), Format_log_axis('x')]
    )
    + (
        plot_weight_cossims_n(hcp_fMRI_optiNpcs_pls_analyses['subsampled'].sel(mode=0), clr_pls4)
        * plot_pc_bias(hcp_fMRI_optiNpcs_pls_analyses['subsampled'], clr_pls3).opts(show_legend=True)
    ).opts(
        xlim=(1, 1150/88), show_legend=False, yticks=6, ylim=(0, 1), xaxis='bottom',
        hooks=[suptitle_pls, set_axis_position, ylabel_clr_pls, Labelpos('x', 3), Format_log_axis('x')]
    )
).cols(
    4
).opts(*fig_opts).opts(
    opts.Overlay(logx=True, logy=False),
    opts.Layout(fig_inches=(7, 6), sublabel_position=(-.55, .95), vspace=1.25)
)

hv.save(fig_optiNpcs, 'fig/figS_empirical_data_analysis_hcp_fMRI_optiNpcs.pdf')

fig_optiNpcs

hcp_fMRI_optiNpcs_cca_analyses['p_value'], hcp_fMRI_optiNpcs_pls_analyses['p_value']

(0.000999000999000999, 0.006993006993006993)