Decoding cerebro-spinal signatures of human behavior: Application to motor sequence learning

Mapping the neural patterns that drive human behavior is a key challenge in neuroscience. Even the simplest of our everyday actions stem from the dynamic and complex interplay of multiple neural structures across the central nervous system (CNS). Yet


Introduction
Mapping the neural patterns that drive human behavior is a central and daunting question of neuroscience.Each of our actions, no matter how simple, involves a dynamic and complex interplay between multiple neural structures throughout the central nervous system (CNS), from the brain to the spinal cord.Neuroimaging techniques, such as functional magnetic resonance imaging (fMRI), are instrumental in deciphering these mechanisms.Despite the strong emphasis on cerebral function in past research, non-invasive imaging of the spinal cord can also provide a window into spinal processes and has grown as a field in its own right (see Kinany et al., 2022b ;Landelle et al., 2021 for reviews).While imaging the spinal cord presents inherent challenges (e.g., small crosssectional dimensions, physiological noise, and field inhomogeneities), technological advances have allowed researchers to investigate various facets of spinal activity, notably during motor tasks ( Barry et al., 2021 ;Kinany et al., 2022aKinany et al., , 2019 ; ;Weber et al., 2016 ).Importantly, simulta-neous fMRI of the brain and spinal cord is also feasible, though it poses increasing challenges due to differences in ideal MRI protocols for these two regions (see Tinnermann et al., 2021 for review).This approach can offer a comprehensive characterization of these mechanisms at multiple levels and at a larger scale of the CNS ( Tinnermann et al., 2021 ).
To date, the handful of studies that have exploited cerebro-spinal signals have primarily relied on inferential univariate techniques (e.g., general linear model -GLM), to identify regional activations related to specific tasks.However, such methods fall short of characterizing the integration and segregation of large-scale networks mediating human behavior ( Friston et al., 2002 ).As such, there is a need for flexible multivariate methods that can disentangle the richness of cerebro-spinal activity and reveal the neural states that shape our actions.To address this knowledge gap, the current study parts ways from traditional techniques and introduces a framework to achieve a concise view of the cerebro-spinal signatures of human behavior, by combining (i) a datadriven functional connectivity method to reveal spatially and temporally overlapping networks, the so-called innovation-driven coactivation patterns (iCAPs) ( Karahano ğlu and Van De Ville, 2015 ;Kinany et al., 2020 ) and (ii) multivariate classification.The term innovation , in this context, refers to transient activity that is recovered using robust hemodynamicinformed deconvolution ( Karahano ğlu et al., 2013 ) We propose to demonstrate the relevance of this approach by exploring a key facet of our daily lives -the acquisition and consolidation of new motor skills.Specifically, we focus on motor sequence learning (MSL), a complex process in which complex motor skills emerge from sequences of individual movements that come to be executed jointly and effortlessly through repeated practice ( Diedrichsen and Kornysheva, 2015 ).MSL encompasses a wide range of motor, sensory, and cognitive processes supported by structures distributed along the neural axis ( Doyon et al., 2018( Doyon et al., , 2009 ; ;Doyon and Benali, 2005 ).This intricate machinery operates over multiple temporal scales, characterized by rapid improvements in the early stage of skill acquisition (i.e., within a single training session), followed later by a slower development of performance over several days or weeks, as skills became consolidate and optimized ( Lohse et al., 2014 ;Penhune and Steele, 2012 ;Wymbs and Grafton, 2015 ).
Traditionally, studies assessing learning-related functional plasticity have focused on cortical and subcortical regions ( Dayan and Cohen, 2011 ;Doyon et al., 2018 ), while the spinal cord was merely viewed as a hardwired relay.Nonetheless, evidence challenging the concept of the spinal cord as a passive relay has accumulated over the years ( Christiansen et al., 2017 ;Grau, 2014 ;Nielsen, 2016 ;Wolpaw and Tennissen, 2001 ) and physiological experiments in animals and humans have provided support for learning-dependent plasticity in the spinal cord by probing changes in reflexes (e.g., Hoffmann reflex) ( Lungu et al., 2010 ;Meunier et al., 2007 ;Nielsen et al., 1993 ;Nielsen, 2016 ;Wolpaw, 2007 ).These observations hinted at a more sophisticated role of the spinal cord in motor learning than initially assumed, and prompted the use of neuroimaging studies extending beyond the brain.Yet, only two studies so far have taken advantage of cerebrospinal fMRI to pursue this question by using classic inferential methods ( Khatibi et al., 2022 ;Vahdat et al., 2015 ).
Given its multifaceted and distributed nature, unveiling the mechanisms underlying MSL stands out as a challenge, providing a prime opportunity to validate our iCAP-driven approach.To this end, we deployed our data-driven multivariate framework in a simultaneous brain-spinal cord fMRI dataset acquired during early and late MSL ( Khatibi et al., 2022 ).We harnessed the dynamic content of CNS signals to uncover remarkably well-delineated cerebral and spinal networks, representing building blocks that are distinctly modulated and interacting to sustain evolving task demands over the course of motor learning.These neural states could be used to decode learning stages with a high accuracy, underscoring the potential of our approach to reveal meaningful correlates of motor learning.Considering its versatility, this framework can readily be extended to unveil interpretable cerebro-spinal signatures of aspects of human behaviors beyond those investigated in this study (e.g., pain processing, neurological conditions, etc.).

Participants
Thirty young healthy participants were enrolled in the study ( Khatibi et al., 2022 ).All participants were right-handed and had no history of neurological or psychiatric disease, no motor-system complications, nor were they using neurological medications.Individuals who previously participated in motor-learning experiments, as well as those with previous training in playing a musical instrument for more than three consecutive years in the last five years, were excluded from the study.Of these thirty participants, five subjects were subsequently excluded (two dropped-out, one exhibited excessive motion, and two had incomplete functional scans).Consequently, 25 subjects were included in the analyses (13 females, 24.8 ± 3 years old).All participants gave their written informed consent to participate, and the study was approved by the Ethics Committee at the centre de Recherche de l'Institut Universitaire de Gériatrie de Montréal (CRIUGM).

General procedure
Each participant participated in six experimental sessions, distributed over six consecutive days (Figure S1A).These sessions comprised two MRI scans: at the beginning (Day 1) and at the end (Day 6) of the study.At the beginning of the first experimental session, subjects were familiarized with the set-up by practicing a motor sequence learning (MSL) task in a mock scanner, using a different sequence than the one used during the actual experiment.Participants were then installed in the MRI scanner, where functional acquisitions were obtained during the performance of two different motor tasks (Sequence or Random, see below for details).During Day 2 to Day 5, participants performed a daily training of the sequential task while lying supine in the mock scanner.

Motor sequence learning paradigm
The motor task was performed using a joystick that the subject manipulated with the dominant (right) hand.At the beginning of the task, the cursor representing the joystick position was positioned in the center of the screen (i.e., neutral position) and indicated by a black circle.Four static targets, represented as circles, were positioned at equal distance from the center (3, 6, 9 and 12 0-clock, Figure S1B).Participants were instructed to move the cursor towards a target as soon as the corresponding circle was filled.Once the target was reached, it disappeared and the circle corresponding to the following target was highlighted.The sequence was not known in advance by the participants, who were simply asked to reach targets as quickly as possible while being as accurate as possible.The sequence was repeated ten times in each block, with 15 blocks in total, interleaved with rest periods of 20 s.Two different conditions were used: (i) an 8-item Sequence Learning condition (SEQ) in which the order of the targets was sequential and repeated every 8 reaching movements (Figure S1B), and (ii) a Random condition (RND) where targets appeared in a pseudo-random, non-sequential order.During each scanning session (Day 1 and Day 6), functional MRI acquisitions were acquired during both experimental conditions on two separate runs (SEQ and RND) and the order of the two runs was counterbalanced between subjects.

Data acquisition
Imaging data were acquired with a 3 Tesla Siemens Trio scanner (Erlangen, Germany), equipped with a 12-channel head coil paired with a 4-channel neck coil (both receive only).Extra foam pads were installed around the subjects' head and neck in order to minimize motion.The isocenter of the magnet was aligned to the lower edge of the head coil (i.e., approximately at the vertebral level C2/C3).
In order to investigate the cerebro-spinal correlates of motor learning, functional images were acquired using a custom fMRI protocol enabling combined brain and spinal cord imaging ( Finsterbusch et al., 2013 ).Briefly, this approach allows simultaneous functional measurements in two sub-volumes (brain, spinal cord) with different geometric properties (Figure S1C).As such, each sub-volume was tailored to the region-of-interest (i.e., brain or cervical spinal cord).Optimal shim settings (linear shim terms and resonance frequencies) were determined for the two sub-volumes at the beginning of each scanning session.Then, shim parameters were dynamically updated during the functional acquisitions.Details regarding shimming procedure are discussed in ( Tinnermann et al., 2017 ).Based on a reference acquisition, a slice-specific z-shim approach was also used to maximize the signal intensity within the spinal cord ( Finsterbusch et al., 2012 ).Forty to 43 slices (depending on the participant size) were acquired for both sub-volumes.Specifically, the brain sub-volume contained 30 to 33 axial slices (FOV = 220 ×220 mm 2 , in-plane resolution = 2 ×2 mm 2 , slice thickness = 3.5 mm, no gap between slices, flip angle = 90º, GRAPPA factor 2, TE = 30 ms).Such parameters allowed to cover the wholebrain of most participants, except for a few subjects for whom the box was tilted upward by approximately 10°(around the x axis), hence cutting a small part of the upper motor cortex in the field of view.As for the spinal sub-volume, 8 to 10 slices were acquired (FOV = 132 ×132 mm 2 , in-plane resolution = 1.2 ×1.2 mm 2 , slice thickness = 5 mm, 4 mm gap between slices, flip angle = 90º, GRAPPA factor 2, TE = 33 ms).The imaged region extended from the C3 to T1 spinal cord segment levels.The Repetition Time (TR) was 3140 ms for S01, 3200 ms for S02 to S04 and 3050 ms for the other subjects, as the TR was slightly optimized in the first acquisitions.Importantly, the same TR was always used within the same subject.To limit noise, only signal from the head coil was used for the brain images, while only signal from the neck coil was used for the spinal cord images.A 3D-MPRAGE T1-weighted anatomical image was also acquired (175 sagittal slices, head to upper thoracic spine, TR = 2300 ms, TE = 3.45 ms, flip angle = 9º, TI = 1.1 s, FOV = 192 ×240 ×320 mm 3 , voxel resolution 1 ×1 ×1 mm 3 ).

Data preprocessing
All preprocessing steps were performed independently for the brain and the spinal cord using the Oxford Center for fMRI of the Brain's (FM-RIB) Software Library (FSL) ( Jenkinson et al., 2012 ) and the Spinal Cord Toolbox (SCT) (versions 3.1 to estimate warping fields and for motion correction (see Khatibi et al., 2022 ) and 4.3 for other processing steps) ( De Leener et al., 2017 ).
(1) Motion correction Brain: Following removal of the first two volumes to allow for T1 equilibration effects, non-brain regions were removed using BET ( Smith, 2002 ) and motion correction was applied using MCFLIRT ( Jenkinson et al., 2002 ).
Spinal cord: Similarly to the brain, the first two fMRI volumes were discarded and slice-wise motion correction was performed with the mean image as reference, using the SCT ( De Leener et al., 2017 ).The corrected time series were visually inspected and parameters were adapted if needed.
Motion parameters were extracted independently for the brain and spinal cord and used to compute the mean framewise displacement (FD BR and FD SC ) for each subject, day and condition.A repeated measures ANOVA with Day (1 vs 6) and Condition (SEQ vs RND) as withinsubject factors was conducted to assess changes in FD BR and FD SC .Performances across days and conditions were then evaluated using paired t-tests (Bonferroni corrected for multiple comparisons).
(2) Image denoising Brain: FSL's FEAT tool was used to regress out signals of no interest from fMRI time series.The following noise regressors were included: six motion correction parameters (translations and rotations), CSF and white matter signals.Residuals were spatially smoothed using an isotropic Gaussian kernel with a FWHM of 4 mm Spinal cord: Noise regression was also performed in the spinal cord using the following confounds: physiological noise modeling (PNM) regressors ( Brooks et al., 2008 ), slice-wise motion correlation parameters (x and y), CSF signal and mean global signal outside the spinal cord.The resulting residuals were then spatially smoothed using a 3D Gaussian kernel with a full width half maximum (FWHM) of 2 ×2 ×4 mm 3 , along the spinal cord.
(3) Estimating warping fields for normalization Brain: Functional-to-anatomical and anatomical-to-template transformations were carried out using FLIRT ( Jenkinson et al., 2002 ).Registration from high resolution structural to standard space (3 mm resolution) was then further refined using FNIRT nonlinear registration.
Spinal cord: Spinal cord segmentation was performed automatically on the T1 anatomical images using a two-step process: the spinal centerline was extracted using a first segmentation and subsequently used to yield a smoothed image on which a second segmentation was obtained.Segmentations were visually inspected and manual adjustments were performed when necessary.Landmarks were placed at the C4 and C7 vertebrae and used to normalize the structural scans to the MNI-Poly-AMU template ( Fonov et al., 2014 ).The spinal cord was semiautomatically segmented (i.e., with manual corrections) on the mean motion-corrected fMRI volume.All segmentations were performed by one rater and inspected by a second examiner, in order to optimize segmentation quality.This segmentation was used to support the registration of the mean motion-corrected to the T1 anatomical scan.Functional-to-template warping fields were obtained by concatenating the functional-to-anatomical and anatomical-to-template transformations.

Task performance
Subject's performance during each block of practice of the two experimental conditions (Sequence and Random) was evaluated by measuring the movement's jerkiness and the duration of each block.We expected that motor skill learning would be reflected by a reduction in movement jerkiness; an indicator that the sequence of movements was being learned and that participants were able to predict the next stimulus and adjust the trajectory of their movements accordingly.Briefly, jerkiness was computed as the third-time derivative of the pointer displacement, for each time point.Values were then integrated to obtain an Integrated Absolute Jerk (IAJ) index ( Goldvasser et al., 2001 ) for each task block.Likewise, we hypothesized that block durations would decrease over the course of the experiment, reflecting participants' improved ability to perform the movements more quickly with practice.A repeated measures ANOVA with Blocks (1:7 vs 9:15) x Day (1 vs 6) x Condition (SEQ vs RND) as within-subject factors was conducted to assess improvement in performance, for both IAJ and block duration.Performances across days and conditions were then evaluated using paired t-tests (Bonferroni corrected for multiple comparisons).As for within-session learning, a similar analysis was conducted to compare the level of jerkiness and average block duration in the first seven blocks against the last seven blocks, for both days and conditions.

Extracting innovation-driven coactivation patterns (iCAPs)
We deployed the (SPiC)iCAP frameworks to extract spinal and brain iCAPs ( Karahano ğlu and Van De Ville, 2015 ; Kinany et al., 2020 ).Briefly, the raw BOLD signals were denoised, and a regularized HRF deconvolution was employed to recover the activity-inducing signals, thanks to the Total Activation (TA) framework ( Karahano ğlu et al., 2013 ).Transient activity (i.e., the so-called innovation signal) was obtained as the temporal derivative of the activity-inducing time courses.A two-step thresholding process was then used to select significant innovation frames: (i) temporal thresholding: a surrogate distribution was generated by applying TA to phase randomized data and a 5% confidence interval was used to pick significant voxels; (ii) spatial thresholding: only the frames containing more that 5% of active voxels were kept.This procedure was applied to each run independently, for both the spinal cord and the brain.Significant innovation frames were normalized to the MNI-Poly-AMU and MNI template, respectively.Of note, a mask, shown in Figure S2, was used to constrain brain analyses to regions relevant to the study of motor control and learning, as defined based on previous literature ( Dayan and Cohen, 2011 ;Doyon and Benali, 2005 ;Hardwick et al., 2013 ;Vahdat et al., 2015 ).As for the spinal cord, the common imaged region was used and extended from C5 to T1 spinal cord segment levels.The spatially normalized maps were masked and concatenated (i.e., both days and conditions, for all subjects, independently for the brain and the spinal cord) and temporal K-means clustering was conducted to obtain group level iCAP maps.Similarly to Zöller et al. (2019aZöller et al. ( , 2019b) ) , the optimal number of clusters was determined by means of consensus clustering ( Monti et al., 2003 ), which uses a subsampling scheme to estimate clustering stability (Figure S3).Specifically, K-means clustering was applied to multiple subsets of the data to generate a consensus matrix indicating the fraction of subsets for which two frames were clustered together.Ideally, two frames should either (i) always be clustered together (i.e., fraction value of 1) or (ii) never be clustered together (i.e., fraction value of 0).This procedure was repeated for different values of K, yielding a consensus matrix for each K of the selected range.The average consensus (i.e., mean consensus per cluster) and the cumulative distribution function of values in the consensus matrices were then probed to select the best clustering solution.At the brain level, 17 iCAPs were recently reported to describe restingstate activity ( Tarun et al., 2021 ;Zöller et al., 2019b ).In order to obtain a similar level of detail, we applied consensus clustering for K BR ∈ [14, 20] (i.e., range of ± 3 compared to these previous studies).Results are presented in Figure S3 and prompted the selection of K BR = 15 for further analyses.As for the spinal cord, we relied on a priori knowledge on anatomy and explored a range corresponding to the potential number of dorsal and ventral components, expected to be around 10, as 5 spinal cord segment levels were imaged.Specifically, consensus clustering was evaluated for K SC ∈ [7, 13] (i.e., range of ± 3 compared to anatomical expectations) and we opted for K SC = 11 based on the consensus clustering quality measures.In order to ensure that the spatial organization of iCAPs was stable across conditions (i.e., for both days and both tasks), we computed condition-specific iCAPs maps (i.e., mean of the frames of one condition, assigned to a particular iCAP).Stability was evaluated using cosine similarity and Dice coefficients (maps binarized with a threshold at Z = 2).Finally, subject-specific activity-inducing time courses were obtained for all iCAPs (brain and spinal cord) using spatio-temporal transient-informed regression ( Zöller et al., 2019a ).As both the RND and SEQ tasks were self-paced, time courses were averaged per task / rest block to allow comparisons between subjects.

Classification of cerebro-spinal correlates of MSL
A multivariate data-driven approach was employed to highlight brain and spinal correlates of MSL.Specifically, we used Linear Discriminant Analysis (LDA) to investigate whether distinct neural activity patterns were associated with the execution of RND or SEQ tasks, in Day 1 and Day 6. Features were defined as the mean difference of activity during task and rest blocks for each iCAP (i.e., 11 spinal cord features and 15 brain features).Values were Z-scored and fed to a two-class LDA classifier, with leave-one-subject-out cross validation.Confusion matrices were computed using all cross-validation folds and statistical significance was assessed using non-parametric permutation tests.Specifically, 1000 classifications were performed with randomly assigned labels and used to build a null distribution on which significance thresholds were calculated.The significance of the LDA weights was evaluated with the same procedure.A similar procedure was employed to compare MSL (i.e., SEQ > RND) between sessions.

Investigating iCAP-informed correlations
In order to explore whether distinct cerebro-spinal interactions were associated with MSL, we performed a series of correlation analyses between the spinal cord and the brain.In particular, we examined task-dependent functional connectivity by leveraging the ability of the iCAP framework to retrieve activity-inducing signals unblurred from the hemodynamic lag.Functional connectivity estimates were performed using iCAPs as seeds.Specifically, time courses from each spinal iCAP were correlated with whole-brain activity-inducing signals from the brain.In order to rule out that connectivity merely reflected shared task input (i.e., correlation with task timings) ( Cole et al., 2019 ), we focused our analysis only on the portions of the time courses corresponding to task blocks.To minimize the impact of state transition, one TR was removed at both ends of each block, and all task blocks were concatenated.Pearson's correlation coefficients were computed between these task-only time courses, for all subjects, days and conditions.Values were then Fisher Z-transformed and MSL-related correlation maps were obtained by contrasting the SEQ and RND conditions within each subject, for each day.Group-level maps were generated using a mixed-effects general linear model ( Beckmann et al., 2003 ), with FSL's FLAME (FM-RIB's Local Analysis of Mixed Effects).Correction for multiple comparisons was conducted using Gaussian Random Field (GRF) theory ( Z > 3, p < 0.001), with FSL's smoothest and cluster functions.

Results
We assessed motion using framewise displacement, separately for the spinal cord (FD SC ) and the brain (FD BR ), across the different days and conditions (Figure S4).A two-way repeated-measures ANOVA revealed a significant interaction effect between day and condition for both the spinal cord (F(1,24) = 248.1,p < 0.001) and the brain (F(1,24) = 133.09,p < 0.001).Paired t-tests indicated that on Day 6, motion estimates were significantly smaller for RND than for SEQ, both in the spinal cord (t(24) = 4, p < 0.01, difference of 0.08 ± 0.02 mm, mean ± SE across subjects) and in the brain (t(24) = 4.76, p < 0.001, difference of 0.03 ± 0.01 mm).

Behavioral performance: motor sequence learning leads to smoother task execution
To probe early and late MSL, subjects were scanned at the beginning (Day 1) and end (Day 6) of a MSL training regimen (Figure S1A).Scans were acquired during sequential (SEQ) or random (RND) movement, performed with the right hand (15 blocks of 80 trials interleaved with rest periods).Participants' performance over time and conditions was evaluated using a measure of movement jerkiness (i.e., integrated absolute jerk [IAJ]) ( Fig. 1 A).Distinct learning curves were observed, as emphasized by a significant three-way interaction between block, day and condition (repeated measures ANOVA, F(1,24) = 859.06,p < 0.001).In particular, learning effects were present across days (i.e., significant interaction between day and condition, F(1,24) = 888.15,p < 0.001) as well as within the same day (i.e., significant interaction between block and condition (F(1,24) = 986.96,p < 0.001).Paired t-tests were used to further analyze these learning dynamics and confirmed an improvement between sessions (t(24) = 2.85, p < 0.05, Fig. 1 B) as well as within session (on Day 1: t(24) = 4.44, on Day 6: t(24) = 4.40, p < 0.01, Fig. 1 C), for the SEQ, but not for the RND condition.As expected, these observations suggest that repeated practice of the SEQ task, as opposed to the non-specific effect of random motor practice, led to smoother movements, hence indicating that subjects had learned the repeating motor sequence.In order to disentangle learning effects from motor practice, MSL was defined as SEQ > RND in further analyses.Block duration (i.e., reflecting speed) was also analyzed as a complementary measure of task execution (Fig. S5).We observed a significant three-way interaction between block, day and condition (Fig. S5A,  F(1,24) = 3458.7,p < 0.001), with significant interaction between day and condition (F(1,24) = 3556.5,p < 0.001) and between block and condition (F(1,24) = 3624.8,p < 0.001).As for speed, SEQ movements were consistently performed faster than RND ones ( p < 0.001) (Fig. S5B-C).An improvement between sessions was observed for SEQ (t(24) = 8.46, p < 0.001) and, to a lesser extent, for RND (t(24) = 3.30, p < 0.05) (Fig. S5B).Over the course of Day 1, there was a significant increase in movement speed for the SEQ condition (t(24) = 6.19, p < 0.001), but not for the RND condition (Fig. S5C).During Day 6, both SEQ (t(24) = 6.87, p < 0.001) and RND (t(24) = 5.70, p < 0.001) movements showed a decrease in block duration.

Brain and spinal activity can be decomposed into meaningful networks
We applied the iCAP framework in the brain and spinal cord independently, so as to retrieve regional iCAPs across all subjects and conditions.Fifteen brain iCAPs were extracted and corresponded to well delineated cortical and subcortical regions known to be involved in various aspects of motor learning and practice ( Fig. 2 A, see Table S1 for details of the structures involved in each network).Notably, a task-specific network was identified (iCAP 6) and included contralateral sensorimotor cortices, the premotor areas bilaterally and the ipsilateral cerebellum.Of note, there were several additional iCAPs comprising parts of the cerebellum (iCAPs BR 3, 7, 8, 13 and 14).Moreover, multiple areas known to be involved in various aspects of motor control and learning were present, such as the striatal regions of the basal ganglia (iCAPs BR 1 and 11), the premotor cortex (iCAPs BR 9 and 10) and areas supporting visuospatial perception (iCAPs BR 5 and 10).The anterior cingulate cortex (ACC), which notably plays a role in movement coordination ( Wenderoth et al., 2005 ) and error detection ( Seidler et al., 2013 ) was also well represented (iCAP BR 2, 11 and 12).Other regions comprised the somatosensory areas (iCAPs BR 2, 9 and 15), fronto-parietal networks (iCAP BR 4) and hippocampus (iCAPs BR 1 and 13).At the level of the spinal cord, eleven iCAPs were identified ( Fig. 2 B), covering a region extending between C5 to T1 spinal cord segment levels.As previously reported ( Kinany et al., 2020 ;Kong et al., 2014 ), all maps spanned a limited rostro-caudal extent.ICAPs appeared to be divided into dorsal and ventral components, likely reflecting sensory and motor mechanisms.Although all patterns appeared to spread bilaterally when visualized at a threshold of Z > 2, there was a noticeable trend towards right-sided lateralization (i.e., ipsilateral to the task).This was supported by significant differences in both the number of voxels ( p < 0.05) and amplitude ( p < 0.001) between the left (1437 ± 101 voxels, mean number of voxels ± SD over iCAPs; 3.03 ± 0.17, mean Z score ± SD over iCAPs) and right (1654 ± 188 voxels; 3.14 ± 0.14) hemicords.
To ensure that these networks could reliably be used to characterize cerebro-spinal activity under different conditions, we assessed their spatial similarity across days and tasks.Spatial patterns were highly stable, both in the brain (cosine similarity = 0.97 ± 0.01, Dice coefficient = 0.72 ± 0.20, mean over iCAPs ± SD) and in the spinal cord (cosine similarity = 0.99 ± 0.002, Dice coefficient = 0.86 ± 0.03).Given the spatial robustness of the functional organization as described using iCAPs, we then probed the temporal expression of these building blocks in each subject.Task-related changes in activity that were differentially modulated by SEQ and RND task practice were observed in the brain and the spinal cord ( Fig. 2 C), hence suggesting the ability of the iCAP framework to capture learning dynamics.
To rule out the possibility that iCAP signals were confounded by motion, we investigated whether mean iCAP task activities were correlated with the corresponding framewise displacement across subjects (i.e., correlation of FD SC with the spinal iCAPs and FD BR with the brain iCAPs).Our analysis revealed no significant correlation between motion and iCAP activities, for all iCAPs, regions, days, and conditions.

Cerebro-spinal dynamics are predictive of early and late motor learning
In order to further explore the relevance of these dynamics, we conducted a multivariate classification analysis to determine whether specific activity profiles were associated with different learning stages ( Fig. 3 ).Using mean iCAP task activities as features, we first aimed to quantify changes in cerebro-spinal activity pertaining to within-session MSL, as reflected by the difference between SEQ and RND task practice.On both days, the task could be decoded with a high accuracy, as emphasized by diagonal confusion matrices ( Fig. 3 A).On Day 1, the average accuracy was 64% ( p < 0.05 against chance level, non-parametric permutation testing), while it reached 86% on Day 6 ( p < 0.001).Linear Discriminant Analysis (LDA) weights were stable between cross-validation folds (cosine similarity = 0.98 ± 0.008 on Day 1 and 0.99 ± 0.02 on Day 6, mean ± SD over folds) and they revealed distinct activity patterns linked to early (Day 1,) and late (Day 6) MSL ( Figs. 3 A, S7A-B).On Day 1, iCAP SC 10 and iCAP BR 2 were consistently engaged in classifying MSL (i.e., SEQ > RND), as highlighted by a high discriminative power.These iCAPs correspond, respectively, to the ventral side of C8 (ipsilateral to the task, as highlighted in Figure S6) and to ACC and somatosensory regions.On Day 6, numerous features carried discriminant power, with significant weights found for six iCAPs.In particular, two spinal iCAPs (iCAPs SC 4 and 6, corresponding to ventral networks at the C6 and C7 spinal cord segment levels, with iCAP SC 4 spreading bilaterally and iCAP SC 6 focused on the ipsilateral side, Fig. S6) and two brain iCAPs (iCAPs BR 1 and 15, corresponding to striatal regions of the basal ganglia and to sensorimotor regions) were positively associated with MSL.On the other hand, iCAPs BR 6 and 10, which mainly include motor and visuospatial regions, appeared to be less prominent during the late phase of MSL.
To better characterize the temporal evolution of learning-specific neural correlates, we also classified MSL activity patterns between sessions ( Fig. 3 B).For that purpose, iCAP activity during the SEQ and RND Fig. 2. iCAP spatio-temporal properties -The iCAP framework was employed to identify cerebral (A.) and spinal (B.) networks for all subjects (for both days, 1 and 6, and both conditions, SEQ and RND).The average consensus value of the frames forming each iCAP is indicated in italic.A. Spatial patterns for the 15 iCAPs that were extracted at the level of the brain.Networks correspond to regions (see Table S1 for details) with distinct functional and anatomical relevance.Maps are overlaid on the MNI template.ACC = Anterior Cingulate Cortex.B. Spatial patterns for the 11 spinal iCAPs.Networks were clearly segregated in the rostrocaudal direction, reflecting the segmental anatomy of the spinal cord.Dorsal (i.e., sensory) and ventral (i.e., motor) iCAPs are present.Maps are overlaid on the MNI-Poly-AMU template ( Fonov et al., 2014 ) and spinal cord segment levels are provided as reference ( Cadotte et al., 2015 ).L = left, R = right, D = dorsal, V = ventral.C. Examples of taskrelated timecourses for one brain iCAP (top panels) and one spinal cord iCAP (bottom panels), corresponding to practice of the SEQ (left panels) and RND (right panels) tasks, on Day 1.The brain iCAP corresponds to a motor task-specific iCAP (iCAP BR 6), while the spinal cord iCAP (iCAP SC 10) is related to the ventral side of spinal cord segment level C8.Timecourses are normalized by the mean value during rest for this particular iCAP.Each bar corresponds to the mean activity in one block (task in orange for SEQ, blue for RND, and rest in gray).Values are presented as mean over subjects ± SE. task practice, on both days, were taken as features.Again, the high classification accuracy (on average 70%, p < 0.01 against chance level, nonparametric permutation testing) underlined that cerebro-spinal iCAPs were differentially activated during early (fast) and late (slow) MSL ( Figs. 3 B, S7C).LDA weights were stable across cross-validation folds (cosine similarity = 0.98 ± 0.005) and indicated that iCAP SC 6 (ventral side of spinal cord segment level C7) along with iCAPs BR 1 and 15 (basal ganglia and sensorimotor) were involved in late motor learning (MSL D6 > MSL D1), in line with within-session observations.In contrast, the fast motor learning phase (MSL D1 > MSL D6) displayed a functional pattern including iCAP SC 10 (ventral side of C8) as well as iCAPs BR 6 and 10 (motor and visuospatial).
In summary, these findings underscored the functional relevance of the iCAP-derived organization, demonstrating that it can be leveraged to decode learning phases with a high accuracy.Furthermore, we could delineate the neural signatures underlying learning across cortical, subcortical, and spinal levels.This provides a clear and concise system-level description of the complex neuro-functional changes underlying MSL.

Late motor learning changes brain-spinal cord interactions
To better understand the mechanisms by which the brain and the spinal cord support MSL, we next investigated whether their relationship exhibited learning-dependent characteristics.To this end, we used the spinal iCAP networks as seeds and examined their functional connectivity with supraspinal regions across days and conditions.Specifically, we compared connectivity patterns related to performance of both sequential and random conditions on both days.No significant difference was observed on Day 1.In contrast, the same analysis revealed significant learning-dependent connectivity changes on Day 6 ( Fig. 4 A).In particular, we observed an increased interaction between iCAP SC 3 (i.e., dorsal side of C6) and the left cerebellum (Crus I), as well as between iCAP SC 8 (i.e., ventral side of C8) and the left sensorimotor cortex (hand area).Both spinal iCAPs exhibited a lateralized focus on the right side, ipsilateral to the task (Fig. S6).In order to gain additional insights into these changes, we finally assessed the mean correlation between the reported cerebral clusters and the associated spinal iCAPs, for the different days Fig. 3. Decoding cerebro-spinal correlates of early and late MSL. A. Independently within each session (i.e., Day 1 and Day 6), two-class LDA classifiers (leave-onesubject-out cross validation) were employed to discriminate between SEQ and RND task practice using iCAP activities as features.The accuracy of these classifications is highlighted by the confusion matrices, for both days (left panel).On Day 1, the average accuracy was 64%, while it reached 86% on Day 6. LDA weights (mean over cross-validation folds ± SD) are displayed in the middle panel.Significant values are indicated and the corresponding iCAPs are summarized in the right panel.B. A two-class LDA classifier (leave-one-subject-out cross validation) was used to discriminate between activity patterns associated with MSL (defined as the contrast SEQ > RND) between the two sessions (i.e., Day 1 and Day 6).The resulting confusion matrix is displayed in the left panel and shows an average accuracy of 70%.LDA weights (mean over cross-validation folds ± SD) and the iCAPs corresponding to significant weights are presented in the middle and right panels, respectively.* corresponds to p < 0.05, * * to p < 0.01 and * * * to p < 0.001 (non-parametric permutation testing).ACC = anterior cingulate cortex, S1 = primary somatosensory cortex, M1 = primary motor cortex.and conditions ( Fig. 4 B).This underscored that, in the late phase of motor learning, a positive synchronization with the cerebellum emerged in SEQ task practice, while RND task performance led to a negative correlation between the caudal cervical cord and the sensorimotor cortex.
Altogether, these results emphasized the potential of iCAP-driven functional connectivity to shed light on the specific brain-spinal cord interactions that arose during the early and late learning phases of motor skill acquisition.

Discussion
In this study, we aimed to evaluate the potential of a data-driven framework that combines large-scale cerebral and spinal networks derived from iCAPs with multivariate classification to decipher how CNS structures dynamically operate to mediate human behavior.Drawing on a recent dataset ( Khatibi et al., 2022 ), in which fMRI time series were acquired simultaneously from the brain and spinal cord during early (Day 1) and late (Day 6) stages of MSL, we illustrate the ability of this approach to provide a strikingly concise view of the neural states linked to learning-induced functional plasticity at multiple CNS levels and timescales.A compelling observation was a shift from the caudal spinal cord segment level C8 at the onset of learning to the more rostral C6 and C7 spinal cord segment levels following long-term practice of the sequence.In the brain, extended learning generated a transition from motor and visuospatial networks to striatal and sensorimotor networks.Using iCAP-informed functional connectivity analyses, we also uncovered changes in neural synchronization in late learning, with the interaction between the spinal cord and cerebellum becoming stronger during sequential practice while the spinal cord became less synchronized with cortical sensorimotor regions during random move- ments.Taken together, our findings not only corroborate observations derived from inferential univariate approaches, but support the use of our data-driven framework to disentangle cerebro-spinal signals into well-delineated and interpretable signatures, here in the context of motor learning, and more generally to probe the large-scale neural correlates of human behavior in healthy and impaired populations.

Modular description of task-relevant brain and spinal functional networks
To parse out the building blocks of cerebral and spinal activities, we deployed a dynamic connectivity approach based on hemodynamicinformed deconvolution of BOLD time courses ( Karahano ğlu and Van De Ville, 2015 ; Kinany et al., 2020 ).Even though this methodology has previously been used in resting-state studies ( Kinany et al., 2022a( Kinany et al., , 2020 in the spinal cord; Piguet et al., 2021 ;Pirondini et al., 2022 ;Tarun et al., 2021 ;Zöller et al., 2019b in the brain), this is the first demonstration, to our knowledge, that it can be extended to probe task-related signals.
Here, we have highlighted its ability to isolate the constituent elements of neural activity, which are dynamically combined to meet the demands of the tasks.
Unsurprisingly, the 15 brain networks did not correspond to commonly observed resting-state patterns, but rather exhibited motor-and learning-relevant structures ( Fig. 2 A, Table S1).Notably, iCAP BR 6 included contralateral sensorimotor cortices, bilateral premotor cortices and superior parietal lobules, as well as ipsilateral cerebellum.These regions correspond to well-known structures involved in movement generation ( Miall, 2013 ); their lateralization advocating for the task-specific nature of this network.In addition, the iCAP-derived organization delineated networks that comprised striatal regions of the basal ganglia, as well as distinct anatomically and functionally meaningful subregions of the cerebellum.This is also in line with knowledge on the neural substrates mediating motor skill learning, as both the cortico-striatal and cortico-cerebellar circuits have been proposed to play a critical role in this form of procedural memory ( Dayan and Cohen, 2011 ;Doyon et al., 2018 ).Several networks also included parietal regions involved in representation of the hands and in visuospatial guidance ( Hardwick et al., 2013 ), an expected result owing to the nature of the motor task (visual presentation of targets to be reached using a joystick).The latter regions were combined with the premotor cortex, also implicated in spatial planning of movement ( Tanji, 2001 ).Finally, the anterior cingulate cortex (ACC), known to participate in movement coordination ( Wenderoth et al., 2005 ) and error processing ( Seidler et al., 2013 ), was also observed in several networks.
With regards to the 11 spinal iCAPs ( Fig. 2 B), we isolated networks corresponding to spinal cord segment levels and exhibiting both dorsal (i.e., sensory) and ventral (i.e., motor) components, similar to restingstate studies ( Barry et al., 2014 ;Kinany et al., 2020 ;Kong et al., 2014 ;Vahdat et al., 2020 ).Interestingly, all networks were bilateral, despite the unilaterality of the task.Although unilateral networks could emerge with higher numbers of components ( Kinany et al., 2020 ;Kong et al., 2014 ), the current study instead resorted to a limited level of granularity, with the aim of accessing simple large-scale building blocks of spinal cord's functional architecture.
Importantly, we also harnessed the time-varying nature of the iCAP method to confirm that the temporal dynamics of both cerebral and spinal modules displayed fluctuations in activity related to task performance ( Fig. 2 C), further underscoring their functional relevance.

Data-driven signatures of early and late motor learning
Capitalizing on this modular description, we implemented a multivariate classification approach to identify patterns reflecting MSL at different timescales.To this end, we relied on the availability of fMRI time series acquired during RND and SEQ tasks in the early (Day 1) and late (Day 6) stages of learning.Importantly, the behavioral results ( Fig. 1 ) confirmed that within-and between-session improvements were specific to the SEQ task, hence suggesting that we could distinguish skill learning from the subjects' motoric performance.In line with these observations, we showed that activity patterns associated with the SEQ and RND tasks in each session could be classified with high accuracy ( Fig. 3 A).This entailed an increase in discriminability from Day 1 (64%) to Day 6 (86%), suggesting a stabilization of the underlying neural structures as the skills became consolidated and increasingly automatic ( Diedrichsen and Kornysheva, 2015 ).Further advocating for the specificities of the neural activity patterns associated with early and late MSL, we also noted that they could be classified between days with a mean accuracy of 70% ( Fig. 3 B).Nonetheless, the observation that late MSL could occasionally be misclassified as early MSL suggests that there may also be partial overlap in the learning mechanisms across sessions.
Critically, distinct cerebro-spinal structures were involved in discriminating between tasks.In the brain, we found a number of iCAPS showing significant discriminative power.On Day 1, one cortical iCAP (iCAP BR 2) was found to be associated with SEQ practice.Composed of regions in the ACC and S1, we posit that this network may contribute to the initial learning of the sequence through processes supporting sensory feedback and error correction ( Hardwick et al., 2013 ).Likewise, we observed that networks encompassing premotor areas (iCAPs BR 6 and 10) were tied to early MSL, in line with evidence pointing at a dominant contribution of these regions during the formation of sequence-specific knowledge ( Lohse et al., 2014 ;Penhune and Steele, 2012 ;Wymbs and Grafton, 2015 ).Following training, a global decrease was observed in task-specific sensorimotor areas (iCAP BR 6), along with visuospatial regions (iCAP BR 10), potentially indicating that sequential movements were done effortlessly after several days of practice, with minimal reliance on visual information.These diminutions occurring with late learning were accompanied by an increase of activity in other cortical and subcortical domains, in particular in a network comprising striatal and hippocampal regions (iCAP BR 1), and in an iCAP including bilateral sensorimotor regions (iCAP BR 15).Importantly, cortico-striatal systems are presumed to mediate MSL consolidation ( Doyon et al., 2018( Doyon et al., , 2003 ; ;Doyon and Benali, 2005 ;Lohse et al., 2014 ), notably by means of interactions with the hippocampus ( Albouy et al., 2013 ).Furthermore, these observations are in agreement with a reported shift from a task-specific cortical network in short-term learning, to a bi-hemispheric cortical and subcortical network in long-term learning ( Floyer-Lea and Matthews, 2005 ).
In the spinal cord, our findings emphasized several networks positively associated with MSL, with disparate regions dominating during the early and late learning intervals.Specifically, a shift from caudal (C8, iCAP sc 10) to rostral (C6-C7, iCAPs sc 4 and 6) spinal cord segment levels occurred from Day 1 to Day 6.These networks were located ventrally, pointing to their motor essence, and coincided with innervation sites of fingers and wrist muscles, respectively ( Kendall et al., 2005 ;Kinany et al., 2019 ).We hypothesized that the shift in rostro-caudal activity arising throughout learning may indicate a temporal evolution of the motor strategies employed by the participants to achieve smooth sequential motion.Initially, accurate movements would rely on distal muscles, as finger grip enables precise control of the joystick and allows for prompt trajectory correction in the event of errors.In contrast, a more optimal muscle control scheme would develop following practice, leading to a decrease in finger activations, coupled with the emergence of synergistic strategies involving more proximal muscles, such as the wrist actuators.As a matter of fact, synergies have been proposed as a way to simplify complex motor control, by simultaneously recruiting muscles (i.e., muscle synergies, D 'Avella et al., 2003 ;Mussa-Ivaldi et al., 1994 ) or joints (i.e., postural synergies, Santello et al., 1998 ) in a coordinated fashion.Here, such an hypothesis is also supported by changes in muscular activations in the different phases of learning ( Khatibi et al., 2022 ).
In conclusion, we have exposed distributed networks of cortical, subcortical and spinal structures that selectively contribute to the acquisition and consolidation of novel skills.While this study corroborates prior GLM-based investigations of task-related activity during MSL ( Khatibi et al., 2022 ), it also extends them by revealing a concise lowdimensional representation of large-scale learning-induced plasticity in a data-driven manner.From a methodological viewpoint, these results lend support to the use of this framework to delineate clear and interpretable signatures of MSL and, more generally, of other aspects of human behavior.

Integrated brain-spinal cord networks in late learning
Finally, we leveraged the ability of the iCAP framework to retrieve time series unblurred from the hemodynamic lag, and thus conducted an iCAP-informed correlation analysis to further characterize the cerebrospinal interplay underlying learning.This revealed that two brain-spinal cord networks emerged following long-term practice ( Fig. 4 A), involving an increased correlation between the spinal cord and the cerebellum, accompanied by an increased anti-correlation with the sensorimotor cortex during performance of the RND condition.Interestingly, cerebrospinal networks have been formerly reported, during task ( Khatibi et al., 2022 ;Vahdat et al., 2015 ) and at rest ( Khatibi et al., 2022 ;Vahdat et al., 2020Vahdat et al., , 2015 ) ).In particular, connectivity changes have been highlighted in networks comprising the spinal cord, the cerebellum and the sensorimotor cortex within a session of sequential finger tapping practice ( Vahdat et al., 2015 ), hinting at the importance of these interactions in skill acquisition.
Yet, our results indicate the involvement of the same integrated networks in the later stage of learning when using a different, wristcontrolled task.Accordingly, the regions involved in the reported spinocerebellar network are congruent with the nature of the task.Indeed, experiments in primates and humans have linked Crus I with forelimb movements ( Lu et al., 2007 ;Mottolese et al., 2013 ) and the dorsal side of the C6 spinal cord segment level with sensory signals from muscles, joints and skin from the same region.We speculate that the enhanced synchronization between these regions may reflect the integration of muscular synergies relying on wrist muscles, in line with previous research suggesting a role of spino-cerebellar interactions in synergy recruitment ( Jörntell, 2016 ).According to this hypothesis, the spinal cord adapts to the optimal muscle coordination patterns based on the task's kinematic and biomechanical constraints, while the cerebellum times the use of these synergistic components.This notion is supported by evidence that cerebellar damage impacts the spatiotemporal structure of muscle synergies ( Berger et al., 2020 ).While these conclusions remain conjectural, future studies could explicitly compare synergistic and non-synergistic hand movements, similar to previous work in the brain ( Ehrsson et al., 2002 ).
In contrast to the spino-cerebellar network, the increased negative correlation between activity in the ventral part of the C8 spinal cord segment level and the sensorimotor cortex may rely on corticospinal inhibitory processes.Indeed, changes after long-term practice were in that case observed for the random condition and not the sequential one, which may suggest that participants had to inhibit the sequence program during the random performance.Although there is no direct inhibitory corticomotoneuronal connections, neurons originating from the sensorimotor cortex also project to inhibitory spinal interneurons, which in turn synapse directly with motoneurons to control their excitability ( Goulding et al., 2014 ;Jankowska, 2001 ;Nielsen, 2004 ).This sensorimotor integration at the spinal level may have functional significance in regulating ongoing motor activity during visuomotor training through the gating of spinal motoneurons via inhibitory mechanisms ( Floeter et al., 2013 ;Perez et al., 2005 ).Given the change in speed reported for the random sequence from Day 1 to Day 6, it is also plausible that the changes observed during the random practice may have been influenced by a confounding effect of sequence practice (during Day 2 to Day 5) on the random task performance on Day 6.

Limitations
Although a thorough discussion and investigation of learning versus performance is beyond the scope of the present methodological study, it is important to stress that these two constructs are often intertwined and difficult to disentangle ( Doyon et al., 2018 ;Orban et al., 2010 ).Here, we set out to use a serial reaction time (SRT) task, a well-established and widely-used approach for studying learning (see Janacsek et al., 2020 for review) that allows for the control of several confounding factors (e.g., perceptuo-motor) through the comparison of SEQ and RND conditions.However, we observed a slight improvement in speed for the RND condition on Day 6, which could be due to general task execution improvement linked to the SEQ task practice during Days 2 to 5. Despite this limitation, our approach revealed marked improvements in movement smoothness that were specific to the SEQ task, providing evidence that we were able to isolate learning-specific effects.
From a methodological standpoint, denoising spinal cord fMRI time courses is notoriously challenging, and physiological noise such as breathing and cardiac signals is a prominent issue ( Brooks et al., 2008 ;Fratini et al., 2014 ).Although we followed recommended practices for physiological noise removal ( Eippert et al., 2017 ;Kong et al., 2012 ), no study to date has investigated the impact of different denoising approaches on the spatiotemporal properties of iCAPs.Likewise, we carefully removed motion-related noise and established that motion was not correlated with iCAP activity.However, given the observed differences in motion across conditions, it is possible that residual motion effects could have influenced our results to some extent, despite our efforts to minimize this potential impact.
While our data-driven multivariate framework using iCAPs allowed us to uncover cortical, subcortical, and spinal functional networks during motor learning, it is important to note that our approach is sensitive to the number of iCAPs selected.In our study, we chose the number of iCAPs based on stability metrics and our goal of achieving a low-dimensional space.This may explain why spinal iCAPs exhibited a stronger ipsilateral focus, yet still spread bilaterally.While this choice allowed us to reveal large-scale networks and concisely summarize learning-related mechanisms, it cannot provide detailed information on the underlying sensorimotor processes.Future studies could explore the impact of varying the number of iCAPs on the results, to assess the balance between model complexity and capturing finer underlying neural states.Moreover, our choice of applying a mask on regions related to motor control and motor learning at the brain level, while it also helped limit the dimensionality of our analyses, may have masked other unexpected -but potentially interesting -regions, something that could also be explored in future work.
The iCAP analysis was conducted independently on the brain and spinal cord, drawing on prior work in each respective region ( Kinany et al., 2022a( Kinany et al., , 2020b ; ;Zöller et al., 2019a ).However, an interesting avenue would be to explore the potential of our approach to extract joint brain and spinal cord components.Yet, conducting such analyses poses unique challenges due to differences between the two regions (acquisition parameters, signal-to-noise ratio, sensitivity to noise confounds, number of voxels), and to the high computational cost associated with a large number of voxels.To date, only one study has attempted to investigate joint cerebro-spinal cord networks, using ICA ( Vahdat et al., 2020 ) and future work could focus on establishing a robust pipeline to perform such analyses, beginning with resting-state data and eventually addressing more complex questions, such as those explored here.

Broader potential for fundamental and clinical applications
A key strength of this framework is its versatility, which opens up numerous possibilities for future research.By combining simultaneous cerebro-spinal fMRI with data-driven multivariate techniques, it offers the prospect of uncovering the complexities of large-scale CNS activity and translating them into relevant signatures.We have demonstrated the potential of this approach in the context of MSL, but it can be readily applied to other facets of human behavior, such as pain processing or proprioceptive perception.
Interestingly, its data-driven nature can also allow researchers to bypass certain limitations of traditional inferential analyses (e.g., GLM).While the latter demands knowledge of the underlying task paradigm, our framework can be extended beyond the study of task execution itself.Thus, it offers the flexibility to harness resting-state recordings to delve into the intrinsic functional organization of the CNS.In the context of MSL, for instance, this could be deployed to reveal the long-term plastic reorganization that takes place in the brain and spinal cord following learning (i.e., by decoding changes in the organization at rest before and after training).
Altogether, this can have profound implications for clinical applications, enabling the study of cerebro-spinal signatures of neurological conditions in patients with various levels of disability.By rendering a concise view of the impact of these conditions on large-scale cerebrospinal networks and connectivity, relative to control populations, this could not only play a crucial role in advancing our understanding of the impaired CNS but also help delineate biomarkers of its functional integrity.

Conclusion
Using a data-driven multivariate framework, we uncovered cortical, subcortical and spinal networks and demonstrated that their dynamics could be leveraged to accurately decode the cerebro-spinal signatures of early versus late MSL.Our approach offers a complement to traditional univariate analyses by providing a concise and modular view of neural function and plasticity throughout the neural axis.The potential of such approaches to disentangle CNS activations paves the way to future investigations of the large-scale cerebro-spinal networks underlying healthy and impaired human sensorimotor behavior.

Declaration of Competing Interest
None.

Fig. 1 .
Fig. 1.Behavioral performance.Performance was assessed using an Integrated Absolute Jerk (IAJ) measure, which characterizes the smoothness of movement trajectories (i.e., a decrease in IAJ values reflects an increase in smoothness).A. IAJ averaged across all subjects for each block, illustrating learning curves for both days and conditions.B. Comparison of subjects' performance (i.e., mean IAJ values) between days for both conditions.C. Comparison of the performance between the early (first seven blocks: B1-B7) and late (last seven blocks: B9-B15) phases of learning, for both days and both conditions.Values are presented as mean ± SE. n.s.= non significant, * corresponds to p < 0.05, * * to p < 0.01, and * * * to p < 0.001 (paired t-tests corrected for multiple comparisons).

Fig. 4 .
Fig. 4. Brain-spinal cord interactions during motor learning.An iCAP-informed correlation analysis (i.e., correlating spinal iCAPs time courses with brain activity inducing signals) was used to inspect learning-induced brain-spine connectivity patterns (i.e., SEQ vs RND) during the two sessions.A. While no significant difference was observed during Day 1, two brain-spinal cord networks emerged in late MSL (i.e., stronger during SEQ task practice than RND task practice on Day 6).The first network (left) encompasses iCAP SC 3 and a region in the left cerebellum (Crus I).The second one (right) corresponds to iCAP SC 8 and a brain region in the left sensorimotor cortex (hand area).Spinal iCAPs are presented along with the corresponding brain statistical maps (corrected for multiple comparison using Gaussian Random Field (GRF) theory, p < 0.001).B. For the two brain-spinal cord networks, we present bar plots showing the mean Fisher Z-transformed correlations values independently for the different days and conditions.Values are provided as mean ± SE.SMC = sensorimotor cortex.