#!/usr/bin/env python # coding: utf-8 # # Hyperalignment # In[ ]: """ Whole-brain searchlight hyperalignment This script performs searchlight-based hyperalignment using Neuroboros and Procrustes alignment on surface data. For each subject, data from all other subjects are used to build the common space, and the left-out subject is aligned to this common space (leave-one-subject-out scheme). Outputs: - Aligned training data - Aligned test data - Corresponding labels (sequences) Author: Theo Marins """ # ============================================================================ # Imports # ============================================================================ import os import pickle import warnings from time import time import numpy as np import neuroboros as nb import nibabel as nib from tqdm import tqdm from sklearn.exceptions import ConvergenceWarning from hyperalignment.searchlight import ( searchlight_procrustes, searchlight_template ) # Silence convergence warnings (expected in ridge/procrustes) warnings.filterwarnings("ignore", category=ConvergenceWarning) # ============================================================================ # Paths # ============================================================================ subject_data_file = os.path.join(PKL_DIR, "darrays.pkl") run_label_file = os.path.join(PKL_DIR, "run_label.pkl") subject_label_file = os.path.join(PKL_DIR, "subject_label.pkl") session_label_file = os.path.join(PKL_DIR, "session_label.pkl") sequence_label_file = os.path.join(PKL_DIR, "sequence_label.pkl") # ============================================================================ # Load data # ============================================================================ with open(subject_data_file, "rb") as f: subject_data = pickle.load(f) with open(run_label_file, "rb") as f: run_label = pickle.load(f) with open(subject_label_file, "rb") as f: subject_label = pickle.load(f) with open(session_label_file, "rb") as f: session_label = pickle.load(f) with open(sequence_label_file, "rb") as f: sequence_label = pickle.load(f) print("Loaded data for sessions:", subject_data.keys()) # ============================================================================ # Parameters # ============================================================================ sessions = ["01"] subjects = list(subject_data["01"].keys()) runs = list(subject_data["01"]["01"].keys()) radius = 20 space = "fsavg-ico32" # ============================================================================ # Masks # ============================================================================ mask_l = nb.mask("l", space=space, mask=True) mask_r = nb.mask("r", space=space, mask=True) # ============================================================================ # Output containers # ============================================================================ train_data_dict = {} test_data_dict = {} train_label_dict = {} test_label_dict = {} # ============================================================================ # Hyperalignment helper # ============================================================================ def hyperalign_train_data(train_data, sls, dists, radius=20): """ Perform searchlight hyperalignment on training subjects. Parameters ---------- train_data : array Shape (n_subjects, n_samples, n_vertices) sls, dists : searchlight definitions radius : int Searchlight radius Returns ------- template : array Group template aligned_train : array Hyperaligned training data """ print(" Building template...") template = searchlight_template( dss=train_data, sls=sls, dists=dists, radius=radius, tmpl_kind="pca", n_jobs=1 ) aligned_train = [] for i, subj_data in enumerate(train_data): print(f" Aligning training subject {i}") W = searchlight_procrustes( subj_data, template, sls, dists, radius=radius ) aligned_train.append(subj_data @ W) aligned_train = np.stack(aligned_train, axis=0) return template, aligned_train # ============================================================================ # Main loop # ============================================================================ start_time = time() for hemi_name, hemi_key in zip(["left", "right"], ["l", "r"]): print("\n" + "=" * 60) print(f"Processing {hemi_name.upper()} hemisphere") print("=" * 60) mask = mask_l if hemi_key == "l" else mask_r # Define searchlights sls, dists = nb.sls( hemi_key, radius, space=space, mask=True, return_dists=True ) train_data_dict[hemi_key] = {} test_data_dict[hemi_key] = {} train_label_dict[hemi_key] = {} test_label_dict[hemi_key] = {} # Leave-one-subject-out for test_subject in subjects: print(f"\nTest subject: {test_subject}") train_data = [] train_labels = [] test_labels = [] # ------------------------------------------------------------ # Assemble data # ------------------------------------------------------------ for ses in sessions: for subject in subjects: subj_runs = [] subj_labels = [] for run in runs: X = subject_data[ses][subject][run][hemi_name][:, mask] y = sequence_label[ses][subject][run] subj_runs.append(X) # subj_labels.extend(y) subj_data = np.concatenate(subj_runs, axis=0) if subject == test_subject: test_data = subj_data test_labels = subj_labels else: train_data.append(subj_data) train_labels.extend(subj_labels) train_data = np.stack(train_data, axis=0) # ------------------------------------------------------------ # Hyperalignment (train subjects) # ------------------------------------------------------------ print("Running hyperalignment...") template, train_data_ha = hyperalign_train_data( train_data, sls, dists, radius ) # ------------------------------------------------------------ # Align test subject # ------------------------------------------------------------ W_test = searchlight_procrustes( test_data, template, sls, dists, radius=radius ) test_data_ha = test_data @ W_test # ------------------------------------------------------------ # Store # ------------------------------------------------------------ train_data_dict[hemi_key][test_subject] = train_data_ha test_data_dict[hemi_key][test_subject] = test_data_ha train_label_dict[hemi_key][test_subject] = train_labels test_label_dict[hemi_key][test_subject] = test_labels print(f"\nCompleted in {time() - start_time:.2f} seconds.") # ============================================================================ # Save results # ============================================================================ with open(os.path.join(PKL_DIR, "train_data.pkl"), "wb") as f: pickle.dump(train_data_dict, f) with open(os.path.join(PKL_DIR, "test_data.pkl"), "wb") as f: pickle.dump(test_data_dict, f) with open(os.path.join(PKL_DIR, "train_labels.pkl"), "wb") as f: pickle.dump(train_label_dict, f) with open(os.path.join(PKL_DIR, "test_labels.pkl"), "wb") as f: pickle.dump(test_label_dict, f) print(f"\nSaved hyperaligned data to: {PKL_DIR}")