Brain oscillatory correlates of visuomotor adaptive learning

Sensorimotor adaptation involves the recalibration of the mapping between motor command and sensory feedback in response to movement errors. Although adaptation operates within individual movements on a trial-to-trial basis, it can also undergo learning when adaptive responses improve over the course of many trials. Brain oscillatory activities related to these "adaptation" and "learning" processes remain unclear. The main reason for this is that previous studies principally focused on the beta band, which confined the outcome message to trial-to-trial adaptation. To provide a wider understanding of adaptive learning, we decoded visuomotor tasks with constant, random or no perturbation from EEG recordings in different bandwidths and brain regions using a multiple kernel learning approach. These different experimental tasks were intended to separate trial-to-trial adaptation from the formation of the new visuomotor mapping across trials. We found changes in EEG power in the post-movement period during the course of the visuomotor-constant rotation task, in particular an increased (i) theta power in prefrontal region, (ii) beta power in supplementary motor area, and (iii) gamma power in motor regions. Classifying the visuomotor task with constant rotation versus those with random or no rotation, we were able to relate power changes in beta band mainly to trial-to-trial adaptation to error while changes in theta band would relate rather to the learning of the new mapping. Altogether, this suggested that there is a tight relationship between modulation of the synchronization of low (theta) and higher (essentially beta) frequency oscillations in prefrontal and sensorimotor regions, respectively, and adaptive learning.


Introduction
Adaptation is an essential feature of motor control in which the motor command is adjusted on a trial basis to compensate for disturbances in the external environment or in the motor system itself. It involves the recalibration, or updating, of the brain's internal model that predicts the upcoming state of the motor system from the current state and the ongoing motor command. Any mismatch between the predicted motor state and the actual motor state as estimated through feedback, labelled a sensory prediction error, is used for the update of the model and the adjustment of the motor command ( Diedrichsen et al., 2005 ;Wolpert et al., 1995 ;Wolpert and Miall, 1996 ). Adaptation can therefore be viewed as a change in the mapping that relates sensory inputs to motor outputs. Although this remapping operates fundamentally on a trial-to-trial basis, it can also engage a learning phase when sensorimotor mapping evolves over the course of many trials and stabilizes so as to become optimally tuned to specific environments and tasks ( Braun et al., 2009a ). This raises a fundamental question as to whether adaptation occurring on a trial basis and its improvement across trials arise through the same mechanisms or not.
A number of studies have investigated the oscillatory brain activities associated with trial-to-trial reach adaptation using paradigms There is also evidence that trial-to-trial adaptation triggers modulation of beta synchronization during movement planning (pre-movement), therefore linking beta synchronization to the processing of previous error and somatosensory information, both critical to the updating of motor plans ( Alayrangues et al., 2019 ;Jahani et al., 2020 ;Torrecillos et al., 2015 ). Altogether, the above results suggest that beta synchronization prior to and following a movement could express some sort of functional polymorphism, evaluating movement error and mediating the subsequent adaptation from both motor and sensory information.
None of the above studies addressed the question of the learning occurring when a new sensorimotor mapping between feedback and motor command is acquired and becomes stable. A few studies provided some insights into this issue, examining changes in spectral power of the other oscillatory bands by comparing early and late phases of adaptation to a constant perturbation. They reported an increase of gamma power during movement execution as well as an increased power of the slower frequencies (theta and alpha) during pre-movement, especially in the parietal and frontal cortices, during the course of learning ( Gentili et al., 2011 ;Perfetti et al., 2011 ;Thürer et al., 2018 ). Thus, rhythms outside the beta band may play a pivotal role in motor learning. Unfortunately, both adaptation and learning were intertwined in these aforementioned studies. Finally, it is worth to mention that besides spectral power, phase information is also involved in the neuronal encoding of motor processes ( Combrisson et al., 2017 ;Hammer et al., 2013 ;Jerbi et al., 2007 ). However, motor adaptation and learning have not yet been studied under this perspective.
In the present study, we aimed to come up with a more complete understanding of motor adaptation and learning in terms of neural oscillations and synchrony. To this end, we implemented a data-driven multivariate approach -multiple kernel learning (MKL) -that explored the different spectral features of the EEG signals -i.e. power and phase in different bandwidths and regions of interest -prior to (pre-movement) and following (post-movement) the movement during a visuomotor rotation task with constant perturbation. MKL is a powerful method that can decode states of interest from a combination of kernels ( Gönen and Alapydin, 2011 ), such as features of the EEG signals ( Schrouff et al., 2016 ). We examined features of adaptive behavior through an MKL model discriminating early from late stages of a visuomotor constant rotation task. We expected to find a combination of features including but not limited to pre-movement and post-movement beta power to bring information regarding adaptive behavior. However, it would have remained unclear from that modeling whether EEG features that contributed to the decision boundary were related either to trial-to-trial adaptation from errors or formation of the new visuomotor mapping since both processes occur concomitantly in a visuomotor constant rotation task. Accordingly, we also considered two other experimental conditions, including a normal movement condition that relies on the identity mapping and does not involve any adaptation, and a condition with a random rotational perturbation centered on identity mapping. In the latter, there was trial-to-trial adaptation from errors ( Albert et al., 2021 ;Diedrichsen et al., 2005 ;Donchin et al., 2003 ), but this adaptation did not lead to the acquisition of a new sensorimotor mapping since the imposed average mapping is the identity policy as in normal movement condition ( Braun et al., 2009b ). MKL modeling of the constant rotation condition against these additional conditions as well as further univariate analyses provided information on the EEG features most related to the trial-to-trial adaptation from errors and those relating mostly to the acquisition of the new visuomotor mapping.

Participants
19 right-handed healthy volunteers (9 females and 10 males, age range: 20 -33 years; mean ± SD age: 23 ± 3 years) participated in the study. All subjects had no history of neurological or psychiatric dis-orders and presented normal or corrected-to-normal vision. The study was conducted with the approval of the local ethics committee from Grenoble-Alpes University (IRB00010290-2019-02-12-60). Written informed consent was obtained from all participants.

Experimental setup
Participants were seated in front of a 27 ′ ' computer screen, 1.5 m away, and held a joystick with their right hand, which rested on a padded arm support. An opaque panel prevented the view of the hand and forearm. Participants were equipped with a 128electrodes cap (Biosemi©). The task, in line with previous studies (e.g. Perfetti et al. 2011, Tan et al. 2014a, consisted in performing target aiming movements with a joystick that controlled a green ball. Each movement started from the center of the screen to one of eight possible and equally spaced targets around a virtual circle (radius, 15 cmcirca 15-20°of wrist flexion). Each trial started with a 1500 to 1900 ms pre-cue period during which the eight targets were presented as red circles with transparent background. It was followed by a 1500 ms cue period during which the background of two neighbor targets became red and allowed movement planning. Then, all targets except one of the two red targets disappeared to indicate the GO signal. Participants had to reach the remaining red target as fast and as accurately as possible. Once the target was reached with a 200 ms stop inside the target or when the allowed time elapsed ( > 5000 ms), a 2500 ms inter-trial interval preceded the following trial. During this interval, the green ball disappeared and only the reached target remained visible. Subjects were asked to passively let the joystick return to its initial position. The task was implemented using a custom C ++ software based on Qt and Measurement Computing© libraries. The software recorded behavioral data (cursor positions) at 2048 Hz. EEG electrodes (Biosemi system) were placed according to the international 10-20 EEG system, and acquired at 2048 Hz. Data synchronization was controlled through triggering from behavioral software.

Experimental conditions
After 50 familiarization trials, participants performed the normal movement condition ( Norm ) and then the random rotation condition ( RdmRot ), including 80 trials each. In Norm, there was a normal relationship between the display and joystick. In RdmRot, rotation angle between the cursor and the actual movement was selected randomly between − 60°to + 60°in step of 20°with the average angle over the 80 trials being equal to 0° ( Fig. 1 ). Each rotation angle was presented 10 times except for the 0°angle that was presented 20 times. Finally, participants performed the condition of 60°constant rotation, which was divided in two runs of 80 trials each separated by a 2 ′ break ( CnstRot-1 and CnstRot-2 ; Fig. 1 ). In all conditions, each target was visited 10 times in a random order.
Early and late stages of adaptive learning were defined as the first 30 trials of CnstRot-1 ( Early-CnstRot ) and the last 30 trials of CnstRot-2 ( Late-CnstRot ), respectively, akin to previous studies ( Perfetti et al., 2011 ;Tan et al., 2014a ). Norm, RdmRot and CnstRot-1 were used to tease apart neural correlates related either to trial-to-trial adaptation or learning of the visuomotor transformation.

Behavioral analysis
Behavioral analysis was performed with custom Matlab routines (R2018b). Cursor positions were down-sampled at 100 Hz and then filtered through a dual low pass Butterworth filter with a 20 Hz cutoff frequency. Reaction time (RT) was calculated as the delay between GO signal and the time when the velocity crossed the threshold of three times its SD at rest. Trials in which participants exhibited an anticipatory behavior (RT < 100 ms) or did not reach the target were excluded Fig. 1. Experimental design. Subjects started with 80 trials of Norm condition (no distortion), followed by 80 trials of RdmRot condition (random rotational perturbation), then performed CnstRot condition, divided in two runs of 80 trials each. The first 30 trials of CnstRot-1 defined the Early-CnstRot stage while the last 30 trials of CnstRot-2 defined the Late-CnstRot stage. Every trial started with a pre-cue period displaying the 8 possible target locations. Then, two targets lighted up during the preparation period and finally only one remained, towards which participants performed their aiming movement. In Norm condition visual feedback was normal, while in RdmRot and CnstRot conditions visual feedback was perturbed using random and constant perturbation, respectively. Distribution of perturbations over trials are presented below each column. from further behavioral and EEG analyses. For each trial, several parameters were also computed to quantify adaptive behavior ( Braun et al., 2009a ;Perfetti et al., 2011 ;Tan et al., 2014a ): (i) movement duration (MT), defined as the time elapsed between movement onset and movement termination (i.e. the moment when target was reached, after the 200 ms stop within it); (ii) path length (PL), or normalized covered distance, computed as the total distance traveled by the cursor during the trial divided by the optimal path length, i.e. the length of the line connecting the starting position and the target; (iii) absolute initial angular error (AIE), defined as the absolute angle between the line connecting the initial cursor position to its position 200 ms after movement onset -and before any corrective movement -and the line connecting the initial cursor position to the target; and (iv) trial-to-trial adaptation rate (AR) as estimated from state-space modeling (c.f. appendix and Tan et al., 2014a ). State-space model performance was assessed using Akaike information criterion (AIC). All parameters respected normality criterion as assessed through Shapiro-Wilk tests. Accordingly, paired ttests and within-subjects ANOVAs combined with Fisher LSD post-hoc tests have been used to investigate differences in mean parameters (averaged across trials). All statistical tests were performed using Statistica 10 (StatSoft©), with a level of significance set at p < 0.05.

EEG preprocessing
EEG data preprocessing was performed using Python and open source MNE software ( https://mne.tools ). First, EEG was down-sampled at 200 Hz and referenced to the average signal across all electrodes. Raw EEG signals were band-pass filtered between 1 and 80 Hz with a notch filter at 50 Hz, and detrended. All channels were visually inspected to identify bad channels which were interpolated, and blinking artifacts were rejected through ICA analysis ( Delorme et al., 2007 ). Then, signals were epoched in two periods: 1) pre-movement, locked on GO signal, starting 2000 ms before and ending 700 ms after (to include first part of motion, before any corrective movement); 2) post-movement, locked on movement termination, starting when target was reached (200 ms before movement termination) and ending after the 2500 ms inter-trial interval. All epochs were visually inspected to remove residual noisy epochs. "Bad epochs " based on anticipatory behavior or trials with target not reached were also excluded. To keep a large enough number of trials across participants in our main investigation between Early-and Late-CnstRot, the first 30 "good epochs " of CnstRot-1 and the last 30 "good epochs " of CnstRot-2 were considered for each participant. Regarding other conditions, 3.6 ± 3.3, 7.4 ± 7,1 and 9.4 ± 6.7 over 80 epochs were discarded in Norm, RdmRot and CnstRot-1 conditions, respectively (averaged values for pre-and post-movement epochs).

EEG time-frequency analysis
EEG time-frequency analysis was performed using Matlab (R2018b) and SPM12 toolbox ( https://www.fil.ion.ucl.ac.uk/spm/ ). For every recording sites, each epoch was decomposed in its time frequency representation through a Hilbert transform between 1 and 80 Hz, with 1 Hz Fig. 2. Illustration of the EEG features in Early-and Late-CnstRot (on central electrode Cz). First column represents averaged power changes across subjects during pre-movement window (0 ms corresponds to GO signal). Theta (4-8 Hz) and gamma (31-80 Hz) waves presented concomitantly two synchronization peaks, after visual stimuli (immediately after the CUE signal at − 1500 ms, and immediately after the GO signal at 0 ms). Alpha (9-12 Hz) and beta (13-30 Hz) waves presented a desynchronization during the cue period, which became even more consistent after the go signal. Second column represents averaged ITPC across subjects during pre-movement window. Peaks of ITPC were essentially present in theta and alpha bands, right after the CUE and GO visual stimuli. Third column represents averaged power changes across subjects during post-movement window (0 ms corresponds to movement termination). Theta and gamma powers presented a slight progressive synchronization after movement. Beta and alpha powers, which were largely below the baseline at the end of the movement, presented a strong synchronization after movement, peaking between 1,5 s and 2.5 s post-movement. Fourth column represents averaged ITPC across subjects during post-movement window. A peak of ITPC was present mainly in theta and alpha bands right after the end of the movement. All these modulations are consistent with current literature (see Kilavik et al. 2013, Palmer et al. 2019, Tan et al. 2014a, Tzagarakis et al. 2010 for power modulation and Combrisson et al. 2017, Popovych et al. 2016 for ITPC). non overlapping intervals and using a two-way, zero phase-lag FIR filter from EEGLAB toolbox (order 3r where r is the ratio of the sampling rate to the low-frequency cutoff of the filter, rounded down). Phase and amplitude (power) signals of the complex transform were then extracted for each frequency. This method has been shown previously to be as accurate as power estimation using Morlet wavelet ( Bruns, 2004 ) while preserving phase information ( Cohen et al., 2009 ;Combrisson et al., 2017 ;Voytek et al., 2013 ).
Relative EEG power changes were then calculated for Early-and Late-CnstRot as the percentage change relative to a stable baseline, by dividing the power signals at each frequency and each time point by the baseline, and then subtracting 100 from the normalized value (expressed in percent). Positive values indicated an EEG power higher than the baseline and will be reported as a synchronization while negative values will be reported as a desynchronization. The baseline was defined for each frequency as the average power during the last 500 ms of the pre-cue period, pooled across trials (i.e. one global baseline was used for the 80 trials of CnstRot-1 and CnstRot-2, before extracting Early and Late subsets). In this 500 ms period, participants were at rest, had stopped performing motion since more than 1000 ms and were presented with only the 8 targets with transparent background. Power changes were then averaged over epochs. Then, inter-trial coherence of phase signals over epochs (ITPC) was computed. For each subject, mean ITPC value during baseline period (same period than for power) was then subtracted from ITPC signals. Typical examples of time-frequency maps (powers and ITPC averaged across subjects on central electrode Cz in Early-and Late-CnstRot) are presented on Fig. 2 . Power changes and ITPC were also estimated for conditions Norm, RdmRot, CnstRot-1 the same way.

MKL
modeling was performed using Pronto library ( http://www.mlnl.cs.ucl.ac.uk/pronto/ ) which is based on the Sim-pleMKL package ( Rakotomamonjy et al., 2008 ). MKL is a supervised machine learning approach that simultaneously learns from different kernels ( Gönen and Alapydin, 2011 ;Schrouff et al., 2016 ) in order to classify two classes, or in the present study, conditions. MKL uses the set of inputs to build kernels representing pair-wise similarity between observations (dot product in Pronto). Then, a support vector machine (SVM; Cortes and Vapnik 1995 ) allows to define a decision boundary, discriminating between classes, for each kernel. To determine this boundary, model parameters w m , representing the contribution of each feature (unitary element of each input, e.g. time point) is optimized. Each decision boundary (one per kernel) is then weighted by a parameter d m to define a global decision boundary. These two steps are implemented recursively in an optimization procedure, with a L1-norm sparsity constraint on the d m vector, encouraging a sparse selection of non-null kernel contributions. Final kernel contributions d m and feature contribution w m are then retrieved ( Fig. 3 ). However, contrary to d m , w m are not sparse in MKL modeling, meaning that every features contributed to the model. As such, interpreting w m is complex and raises several issues ( Haynes, 2015 ;Schrouff et al., 2016 ), so that we preferred focusing our interpretation of the results on the kernel contributions d m . MKL algorithm is based on SVM models, which includes a soft-margin hyper-parameter C. This hyper-parameter allows more or less misclassifications during training, affecting SVM decision boundary. To optimize this hyper-parameter, we performed a nested cross-validation scheme in which an inner loop was used for hyper-parameter selection leading to the highest model performance whereas an outer loop used the selected hyper-parameter to assess the performance. This optimization procedure selected C = 1 for around 80% of the folds in all our MKL implementations and C = 0.1 or 10 for the remaining ones. Thus, we chose to set C = 1 in order to keep the same hyper-parameter across all folds.
In the present study, for every two-class classification, MKL was fed with input signals representing either temporal band-related power Fig. 3. Illustration of multiple kernel learning modeling procedure for Early-vs. Late-CnstRot. For each subject but one, and for each electrode, time-frequency power (or ITPC) maps were computed in both Early-and Late-CnstRot stages. These maps were then averaged across frequencies (in theta, alpha, beta and gamma bands) and over predefined regions of interest (ROIs) to define the 44 feature vectors in the two conditions (labels of the feature vectors are the conditions). ROIs were defined as showed on head map on the left side of the figure. For each region/band couple, a linear kernel Km, representing the pair-wise similarity between features vectors of both conditions across subjects was built ( m = 1 … 44). All kernels and their associated labels were then used to train the model. First, features contributions w m are estimated to define a decision function f m per kernel. The weight of each decision function, d m , is then estimated to provide a final decision function f(x). These two steps are implemented recursively in an optimization procedure, leading to a sparse selection of non-null kernel contributions. The model was then applied to test data (feature vectors of the 19th subject without labels) to obtain associated predicted conditions. This whole process was repeated as many times as there were subjects, excluding a different subject each time to assess accuracy of the model. changes or temporal band-related ITPC, in the four main bands of interest (theta: 4-8 Hz, alpha: 9-12 Hz, beta: 13-30 Hz and gamma: 31-80 Hz). A region-of-interest (ROI) approach was used to limit influence of highly correlated data, averaging power changes and ITPC signals within different ROIs. ROIs definitions and examples of averaged (across subjects and ROIs) temporal band-related power changes are depicted in the left panel of Fig. 3 . From a technical standpoint, input signals were converted into nifti files in order to fit the Pronto imaging toolbox. Thus, 44 kernels (11 ROIs x 4 bands) of size 38 × 38 (19 subjects x 2 conditions) were built for each classification. To ensure that the scale of each kernel did not interfere in modeling, all kernels were meancentered and normalized. Adaptive learning as a whole was decoded classifying Early-against Late-CnstRot for power changes and ITPC signals and during pre-and post-movement separately, leading to four main MKL implementations. We did not include pre-and post-movement periods as well as power changes and ITPC signals in the same model to avoid over-parametrization.
When a classification was significant (see Section 2.8 .), further MKL classifications between conditions were performed to investigate more accurately the role of the identified kernels (CnstRot-1 vs. Norm, CnstRot-1 vs. RdmRot, and Norm vs. RdmRot). CnstRot-2 was discarded from these further classifications due to its intermediate position between CnstRot-1 and Norm both in term of error size and acquisition of the new mapping.

Decoding performance
Performances of the MKL models were assessed using a "leave-onesubject-per-class-out " cross validation scheme. This scheme is the best choice when subjects in different classes are correlated (i.e. within subject design), keeping testing and training sets independent. Indeed, to ensure proper cross-validation, it is crucial that correlated information (here the same subject) is not present both in the train and test sets. Such a dependency would lead to "leakage " and over-optimistic model performance. In our design, the only way to keep testing and training sets completely independent was to split our data by subjects. Thus, each model was trained using all subjects with their corresponding labels (here conditions) except one, and used to predict the labels of the subject left out. This process was repeated as many times as there were subjects. Total accuracy (TA) and class accuracy (CA) were subsequently obtained as the total number of correctly classified test samples divided by the total number of test samples, irrespective of class for TA and for a given class for CA.
Since training sets were not independent across folds (overlapping of training data between folds), the use of any parametric tests was excluded to assess the statistical significance ( Noirhomme et al., 2014 ;Pereira et al., 2009 ). P-value associated to each performance measure was estimated using a 1000-permutation testing framework in which cross-validation and accuracy were recomputed after randomly shuffling training labels ( Ojala and Garriga, 2010 ). The level of significance was set at p < 0.05 for TA.

Univariate analyses
We conducted further analyses on the kernels that were most contributing to the MKL models in order to refine the distinction between neural processes most related to trial-to-trial adaptation from those most linked to the learning of the new visuomotor transformation. It is however important to mention that these selective analyses were run on a subset of the kernels and as such suffered double dipping ( Kriegeskorte et al., 2009 ). Accordingly, univariate analyses served to make interpretation of MKL models easier, but should not be considered as hypothesis testing per se. In this perspective, t -tests and one-way ANOVAs were used to evaluate differences of average power or PLV amplitude over time in the kernels that contributed the most to distinguish between Early-CnstRot and Late-CnstRot as well as between Norm, RdmRot and CnstRot-1. Also, we assessed the correlation between this average amplitude and error size in CnstRot-1 and RdmRot conditions. To this end, for each condition, kernel and subject, trials were grouped into 8 bins according to the size of the AIE, and Fisher r-to-z transform was computed between bin's order and the average amplitude signal in the post-movement period. Significance of the correlations were afterwards examined using one sample ttests. In kernels associated with trial-to-trial adaptation from errors, correlation between error size and average magnitude should be present in both RdmRot and CnstRot-1 conditions while in kernels associated with learning of the new visuomotor transformation, this correlation should be present only in CnstRot-1 condition. Indeed in the latter, decrease of errors should more or less match with time course of the experiment and the increased level of acquisition of the transformation.

Behavioral results
Paired t -test on RT between Early-and Late-CnstRot was not significant ( p = 0.67), indicating that the level of readiness remained similar during the task (Mean ± SD RT: 0.36 ± 0.02 s; Fig. 4 ). Other behavioral measures showed that subjects adapted to the constant rotational perturbation. Both the AIE (41.0°± 2.1 to 22.0°± 1.2; p < 0.001) and the AR (0.0230 ± 0.004 to 0.009 ± 0.002; p < 0.001) decreased from Early-to Late-CnstRot. MT and PL followed the same pattern, decreasing significantly from Early-to Late-CnstRot (Early-CnstRot: MT = 1.15 s ± 0.05 and PL = 1.70 ± 0.05; Late-CnstRot: MT = 0.89 s ± 0.04 and PL = 1.36 ± 0.03; p < 0.001). Nevertheless, although the level of adaptation was largely improved in Late-CnstRot, it was still not comparable to the level of Norm condition ( Fig. 4 ). Hence, subjects adapted to the new visuomotor transformation although room for improvement remained.

MKL modeling
The main objective of the study was to decode the process of adaptive learning based on EEG spectral features. To do so, we ran four MKL models aiming at classifying early (Early-CnstRot) and late (Late-CnstRot) stages of a constant perturbation condition based on power and ITPC from multiple frequency bands and cortical regions during pre-and postmovement. The MKL model whose kernels were built from power of the different band/regions couples in post-movement period was statistically significant (TA = 86.84%, p = 0.002; CA = 84.21% and 89.47%, p = 0.004 and p = 0.002 for Early-and Late-CnstRot, respectively). MKL classification based on kernels built from power in pre-movement did not distinguish the two classes (TA = 65.79%, p = 0.11). As for ITPCbased modeling, both pre-and post-movement classifications failed to dissociate Early-from Late-CnstRot (TA = 60.53% and 26.32%, p = 0.26 and p = 0.97 for pre-and post-movement, respectively). These results suggested that motor adaptation and learning (on 160 trials) was related to modulation of brain waves' power in post-movement, only.
Kernel contribution to the MKL model of Early-vs. Late-CnstRot based on power in post-movement was assessed in terms of frequency bands, ROIs, and ROI for each frequency band (a.k.a. ROI × band), as depicted in Fig. 5 . Three frequency bands actively participated to the classification, including the theta, beta and gamma bands with a contribution of 29.4%, 18.6% and 51.7%, respectively ( Fig. 5 A). ROIs that contributed the most to the MKL model were mainly distributed over frontal and central regions (e.g., d FP = 20.5%, d RPM = 11.6%, d SMA = 17.9%, d RM = 10.3%, d LM = 20.6%; Fig. 5 B). Interestingly, when looking at ROI × band kernels ( Fig. 5 C), it appeared that frontal and premotor cortices were related to the theta band (d FP_theta = 19.9%, d RPM_theta = 9.0%), supplementary motor area to the beta band (d SMA_beta = 17.9%), whereas motor and post-central cortices were associated with the gamma band (d RM_gamma = 10.3%, d LM_gamma = 19.6%, d LP_gamma = 9.8%).
Power changes (averaged across subjects) in post-movement window are represented in Fig. 5 C for the three main contributing ROI × bandspecific kernels (FP theta, SMA beta and LM gamma, which accounted together for 57.4% of the total contribution) to qualitatively examine their evolution between Early-and Late-CnstRot. For FP theta, power was roughly at a baseline level in Early-CnstRot and became positive in Late-CnstRot, which reflected occurrence of post-movement theta band synchronization with adaptation. A similar trend was found for LM gamma, also indicating the occurrence of a post-movement synchronization in this frequency band during the course of adaptive learning. Regarding SMA beta, there was an increased post-movement rebound, or equivalently post-movement synchronization, from Early-to Late-CnstRot. These trends were confirmed by t -tests on average postmovement synchronization between Early-and Late-CnstRot in each kernels ( p < 0.001, p = 0.003 and p = 0.03 for FP theta, SM beta and LM gamma, respectively).

Behavioral results
Within-subjects ANOVA did not reveal ( p = 0.94) any RT difference between Norm, RdmRot and CnstRot-1 (Mean ± SD RT: 0.36 ± 0.01 s; Fig. 6  Averaged power changes across subjects of the three main contributing kernels are represented for Early-(red) and Late-CnstRot (green), as well as mean ± SE values (over the time period -0 ms corresponds to movement termination), showing the increase of post-movement synchronization in these three kernels during adaptation. Averaged power changes across and mean values in Norm (gray) condition are provided as reference. * p < 0.05, * * p < 0.01 in t-tests.

MKL modeling
We implemented three MKL models to uncover which of the postmovement power-based kernels that decoded early from late stages of CnstRot condition were most associated with trial-to-trial adaptation from errors or cumulative learning of the new visuomotor transformation: Norm vs. CnstRot-1, Norm vs. RdmRot and CnstRot-1 vs. RdmRot.
The MKL model Norm vs. CnstRot-1 showed a significant discrimination between the two conditions (TA = 82%, p = 0.002; CA = 79% and 84%, p = 0.02 and 0.006 for Norm and CnstRot-1, respectively). The contribution vector d m revealed a set of kernels that were equivalent to those of the classification of Early-against Late-CnstRot, with post-movement theta power in frontal pole and premotor areas, beta power in supplementary motor area and gamma power in motor and parietal areas most contributing to the decision function. There was also a noticeable difference (compared to Early-vs. Late-CnstRot) in the balance of the kernel contribution to the decision boundary, with more contributory signal in the theta band than in the higher-frequency beta and gamma bands (d theta = 66.3%, d beta = 10.9% and d gamma = 14.6%; Fig. 7 A). Qualitatively, all contributing kernels showed weaker power in CnstRot-1 compared to Norm ( Fig. 7 A). Given that Norm and CnstRot differed both in terms of movement error and level of acquisition of the imposed mapping, these results supported once again the idea that postmovement theta synchronization in frontal pole and premotor regions, as well as beta synchronization in SMA, and gamma synchronization in motor-related and parietal regions were related to adaptive learning.
The MKL model Norm vs. RdmRot led to a significant classification of the conditions (TA = 74%, p = 0.016; CA = 63% and 84%, p = 0.2188 and 0.006 for Norm and RdmRot respectively). Post-movement beta power in the supplementary motor area and gamma power in motor to parietal regions mainly contributed to the model (d beta = 46.2% and d gamma = 40.6%; Fig. 7 B). Contribution of theta band was still present, but to a much lesser extent than the higher-frequency bands (d theta = 11.9%), and emerged from parietal areas and not from prefrontal areas as in the case of Early-vs Late-CnstRot MKL model. Accordingly, this result suggested that both supplementary motor area beta power and motor to parietal areas gamma power after movement completion principally reflected trial-to-trial adaptation from errors. The observed lower amplitude of post-movement beta and gamma power in RdmRot compared to Norm suggested an inverse relationship between error and post-movement high-frequency band synchronization ( Fig. 7 B). This finding also shed new light on frontal and premotor theta power strongly involved in the classification of both Early-vs. Late-CnstRot and Norm vs. CnstRot-1 which may be rather associated with the level of acquisition of the sensorimotor mapping. However, this conclusion on the relationship between post-movement frontal theta band synchronization and update of the mapping lacks direct demonstration as the MKL model RdmRot vs. CnstRot-1 did not show a significant result (TA = 50.00%; p = 0.55).

Univariate analyses
Further comparisons between kernels identified through above mentioned MKL models are presented in Fig. 7 C. For FP theta kernel, withinsubjects ANOVA revealed an effect of condition on mean post-movement power (F 3, 54 = 5.25; p = 0.003). Post-hoc testing revealed a lower power in CnstRot-1 compared with the other conditions ( p < 0.001 and p = 0.03 in comparison to Norm and RdmRot respectively), while Norm and Rdm-Rot conditions presented similar levels ( p = 0.08). These differences suggested that increase in post-movement FP theta mean power would re- Fig. 7. Results of extra MKL models. (A) Norm vs. CnstRot-1 model: Theta, beta and gamma bands mainly contributed to the model (as in Early-vs. Late-Adapt model), although with a dominance of the theta contribution. Main identified kernels qualitatively showed a larger power in Norm compared to CnstRot-1. This supported the implication of these three kernels in both movement error and acquisition of the new mapping. (B) Norm vs. RdmRot model: beta and gamma bands mainly contributed to the model. SMA beta and LM gamma kernels showed a larger power in Norm compared to RdmRot, suggesting that these kernels were principally related to trial-to-trial adaptation. Note that theta band was also identified as contributing to the model. However, further analysis revealed a parietal contribution on it, which differed from the theta contribution in frontal areas of the Early-vs. Late-CnstRot MKL model. (C) Mean postmovement power for the three main contributing kernels in all conditions. CnstRot-1 postmovement average power was lower than both Norm and RdmRot powers in FP theta kernel suggesting a relationship with level of acquisition of the new transformation. Furthermore, Norm postmovement average power was lower than both RdmRot and CnstRot-1 levels in SMA beta kernel, relating beta modulation to the adaptation process. These extra-analyses remained inconclusive regarding LM gamma kernel. * p < 0.05, * * p < 0.01 in post-hoc Fisher LSD tests. flect principally the level of acquisition of the new mapping. Indeed, although the acquisition of the new mapping was underway in CnstRot-1 it remained largely incomplete while it was complete in Norm and Rdm-Rot condition and corresponded to the identity mapping (i.e. ⟨∅⟩ = 0 • ).
Regarding SMA beta kernel, within-subjects ANOVA also demonstrated an effect of condition (F 3, 54 = 5.10; p = 0.004). Post-hoc testing revealed a larger mean power in Norm compared to RdmRot and CnstRot-1 ( p = 0.001 and 0.003 respectively). These differences, in addition to the absence of significant difference between CnstRot-1 and RdmRot, corroborated our guess that mean amplitude of post-movement SMA beta reflected trial-to-trial adaptation from errors, which was similar in CnstRot-1 and RdmRot. On the other hand, within-subjects ANOVA for LM gamma kernel did not reveal any difference between conditions ( p = 0.20). Hence, this extra analysis remained inconclusive on whether post-movement power in LM gamma relate rather to trialto-trial adaptation or cumulative learning of the new mapping.
Finally, we found significant correlations between bin's order (sorted according to decreasing error size) and mean post-movement power for SMA beta in both CnstRot-1 and RdmRot conditions ( R = 0.38 ± 0.03, p = 0.006 and R = 0.24 ± 0.03, p = 0.04 respectively). The smaller the error, the larger the mean post-movement power in SMA in the two conditions. Thus, post-movement power in SMA was modulated by error size. There was also a significant correlation for FP theta in CnstRot-1, but not in RdmRot ( R = 0.23 ± 0.02, p = 0.03 and R = 0.12 ± 0.03, p = 0.28). This absence of relationship in FP theta for RdmRot condition indicated that postmovement theta modulation during adaptive learning rather reflected the formation of the new mapping. Regarding LM gamma kernel, power did not correlate with bins of decreasing error size in both conditions.

Discussion
The aim of this study was to come up with a more complete characterization of adaptive learning in terms of neural oscillations and synchrony of distributed brain regions, and try to tease apart two subtending processes, namely trial-to-trial adaptation from errors and cumulative learning of a new sensorimotor mapping. Our main MKL model indicated that post-movement synchronization of beta oscillations in supplementary motor area, gamma oscillations in motor to parietal regions and low-frequency, theta, oscillations in (pre)frontal regions contributed to adaptive learning. Further discriminant models suggested that higherfrequency -beta and gamma -synchrony in motor-related regions was associated rather to trial-to-trial adaptation, while low-frequency synchrony in (pre)frontal regions rather informed on acquisition of the new mapping. Univariate analyses also pointed in the same direction, except for the gamma band. Modulation of SMA synchrony in the beta band was similar when facing both constant and random perturbations, in particular with a power that increased as a function of decreasing error size. This reinforced the idea that post-movement beta power relates to error processing that instantiates trial-to-trial adaptation. Inversely, modulation of frontal synchrony in the theta band differed between the conditions of constant and random perturbations, with the latter condition being similar to the condition with normal movements. Increased theta power during the course of adaptation seems therefore to be a hallmark of the acquisition of a new visuomotor mapping. Finally, univariate outcomes were inconclusive concerning the link between postmovement gamma synchronization and either trial-to-trial adaptation or cumulative learning.
An important result was the increase in post-movement beta synchronization during the constant rotation task to a level comparable to normal condition. Hence, there was an attenuation of post-movement beta power at the beginning of adaptive learning, that tended to disappear as learning progressed, confirming earlier results ( Palmer et al., 2019 ;Tan et al., 2014aTan et al., , 2014b. Furthermore, we were able to relate these brain oscillatory changes to the adaptation process that occurs independently of cumulative learning of the new sensorimotor mapping, thus linking them to the processing of a completed movement with respect to its predicted outcome. This complements previous results which specifically linked post-movement beta synchronization to movement error and confidence in internal model estimations in motor and somatosensory cortices, suggesting that post-movement beta synchronization decrease could signal a need for adaptation of the motor output ( Palmer et al., 2019 ;Tan et al., 2016 ). The topography of our result was different, being confined mainly to SMA. An explanation for this topographic discrepancy is that we let the MKL to select the most discriminating regions (i.e. data-driven approach) while the above mentioned studies did an a priori choice with respect to EEG channels of interest (i.e. model-based approach). Accordingly, adaptation from errors may be mainly mediated through SMA although other sensorimotor regions also likely contribute to it. Interestingly, it is acknowledged that synchronized beta oscillations bind multiple sensorimotor regions into a large-scale network during motor behavior ( Brovelli et al., 2004 ). This also remains true with subcortical areas such as the subthalamic nucleus whose beta power increases in coherence with that of the sensorimotor cortex after movement offset ( Tan et al., 2014b ). Consequently, beta oscillatory activity may be synchronized between different sensorimotor cortical (SMA, motor, somatosensory) and subcortical areas to process movement error and mediates subsequent adaptation. In addition, our result strengthens the indication that regions of the medial frontal cortex, especially the SMA and the dorsal anterior cingulate cortex (dACC), play a critical role in processing errors and evaluating the outcomes of action ( Amiez et al., 2012 ;Bonini et al., 2014 ;Botvinick et al., 2004 ;Ridderinkhof et al., 2004 ). It is clear that the oscillatory activity we reported here is not limited to SMA but also includes influence from the dACC which lies underneath. Finally, we did not find relationship between adaptive features and pre-movement beta synchronization as suggested before ( Torrecillos et al., 2015 ). This discrepancy may reside in the fact that premovement period could be more involved in the integration of somatosensory information ( Alayrangues et al., 2019 ) than in processing of sensory prediction error, that are both essential for adaptive update of the upcoming movement. Thus, we propose that postmovement modulation of beta power reflects sensory error processing and subsequent signal that engages internal model update.
Our MKL main result also revealed that synchronization in postcentral regions within the gamma band during the post-movement period was involved in overall adaptive learning. However, we failed to reach a clear conclusion about the particular process that was encoded through gamma modulation. At the most, our study suggested that postmovement gamma modulation would relate to error processing that triggers adaptation on a trial basis rather than to learning a new visuomotor transformation. This would suggest that gamma oscillations may not only serve the function of encoding afferent (proprioceptive) feedback and properties (e.g. velocity, effort, force level) of the movement reported in previous motor control studies ( Muthukumaraswamy, 2010 ); see also ( Nowak et al. 2018 for review and references therein). Besides, the involvement of gamma oscillations in processing of action outcome that triggers motor adaptation recalls previous reports that found a role of gamma oscillations in encoding reward outcomes for adapting to challenging tasks ( Quilodran et al., 2008;Rothé et al., 2011 ). It might be that gamma band is a vehicle to encode outcome expectation in a broad sense, with different regions encoding different dimensions (reward, sensory feedback) of outcome expectancy. Indeed, these previous studies reported modulation of gamma power in frontal regions during adaptation while we located it in post-central regions. This assumption is neurophysiologically plausible as power in gamma band is highly localized (e.g. Sirota et al. 2008 ), reflecting specific computations of local groups of neurons in the neocortex. This calls for a more elaborated study on the extent to which gamma-modulated oscillations contribute to recalibration of internal model and adaptive behavior.
Another important finding of the present study was the implication of post-movement synchrony of prefrontal, and to a lower extent premotor, theta oscillations in adaptive learning. This finding is in line with the overall implication of the prefrontal cortex in the coordination of adaptive goal-directed behavior ( Koechlin, 2016 ;Miller et al., 2010 ;Ridderinkhof et al., 2004 ), and also replicates the previous finding of a modulation of prefrontal theta power during visuomotor adaptation ( Gentili et al., 2011 ;Perfetti et al., 2011 ). Although these latter studies have speculated on the role of theta enhancement in the formation of internal models, they were not designed to tease out whether it relates to changes in the state of the internal model that responds to error on a trial basis or to the accumulation of these changes over time to form a new mapping. The inclusion of a condition in which errors were random and equaled zero on average allowed us to specifically link the increase of prefrontal theta oscillatory activity to cumulative learning of the new sensorimotor mapping. This interpretation is in line with a large consensus about the implication of the prefrontal cortex, especially the ventromedial and lateral components, in encoding and learning predictive models mapping stimulus-action onto expected outcomes ( Anguera et al., 2009 ;Contreras-Vidal and Kerick, 2004 ;Koechlin, 2016 ). There is also firm evidence that new acquired motor memories are stored and consolidated in prefrontal and secondary motor cortices ( Dandolo and Schwabe, 2019 ;Pinsard et al., 2019 ). In particular, the ability to store newly learned behavior as memory-based constructs requires implementing top-down control signals from the prefrontal cortex to motor regions ( Narayanan and Laubach, 2006 ), these latter regions being involved in selecting motor plans in response to stimuli ( Koechlin et al., 2003 ). Finally, our interpretation is also consonant with a recent study having demonstrated that exploitation of learned associations between stimuli and responses during spatial context learning are implemented in prefrontal theta band activity ( Spaak and de Lange, 2020 ).
Even though we have emphasized that our approach distinguished between oscillatory activity related to trial-to-trial changes in the internal model (i.e. adaptation) and oscillatory activity responsible for accumulating these changes across trials (i.e. learning), it is also possible that another process, namely cognitive strategy, has affected the study outcomes. It has been shown that individuals employ cognitive strategies to eliminate the error at the beginning of adaptive learning, and that implicit adaptation increasingly takes over from explicit cognitive strategies as skill learning proceeds ( Miyamoto et al., 2020 ;Taylor et al., 2014 ;Taylor and Ivry, 2011 ). Thus, our primary outcome on comparison between early and late stages of constant rotation condition (increase of power in theta, beta and gamma bands) could to some extent also reflect the change in the balance between explicit and implicit processes which are at work during adaptive learning. Individuals may have also employed explicit cognitive strategies to adapt in the random rotation condition, possibly affecting our secondary comparison between constant and random rotation conditions. There is evidence that explicit cognitive strategies influence sensorimotor adaptation even in conditions where perturbations are poorly predictable ( Albert et al., 2021 ;Miyamoto et al., 2020 ). We have very few clues about the neural correlates of implicit and explicit processes in motor learning. To our knowledge, it has been proposed only recently that explicit cognitive and implicit processes are reflected into beta-band activities of distinct regions, the medial frontal region and the lateral central region, respectively ( Jahani et al., 2020 ). Hence, SMA-related increase in beta power during the course of adaptation may not only reflect adaptive changes in the sensorimotor mapping but may also relate to a reduction of movement error guided by a cognitive strategy. It will be important for future studies to confirm that the increase in beta and theta band activities during adaptive learning are hallmarks of trial-to-trial adaptation error and cumulative learning of the new sensorimotor mapping and not a byproduct of the changes in the interactions between explicit strategy and implicit motor adaptation. This will require designing adaptive learning experiments that isolate as much as possible explicit and implicit processes of adaptation, for instance by delaying presentation of the feedback for the former ( Brudner et al., 2016 ;Schween and Hegele, 2017 ) and limiting reaction time for the latter ( Haith et al., 2015 ;Leow et al., 2017 ). An-other solution to isolate the explicit component of adaptation may consist in asking subjects to report where they plan to aim before each trial ( Miyamoto et al., 2020 ;Taylor et al., 2014 ). Finally, we cannot rule out the possibility that subjects exposed across trials to random variations in visuomotor transformation were able to learn some meta-parameters related to the multiple transformations; the so-called structural learning that fasters learning in a new situation ( Bond and Taylor, 2017 ;Braun et al., 2010 ). Hence, our conditions of constant and random rotations may not have been perfectly orthogonal with respect to learning of the new mapping.
Our results on the relationship between specific sensorimotor processes and particular EEG bands and regions would not have been possible without a multi-band exploration of the tasks through machine learning MKL modeling. This stresses how important it is to study brain oscillations in adaptive learning, and more generally motor control, beyond beta band. This complements previous observations that not only beta band but also other oscillatory activities are important to enhance motor performance ( Chung et al., 2017 ). In this perspective, looking at the coupling between theta and gamma or beta is an important challenge, such couplings coordinating communication across brain regions and contributing to the formation of memory traces ( Lisman and Jensen, 2013 ). In particular, it has been shown that theta-gamma and theta-beta coupling can predict individual working memory performance ( Axmacher et al., 2010 ) which is known to be engaged in visuomotor adaptation ( Anguera et al., 2010 ;Christou et al., 2016 ). Closely related is the need to foster whole-brain approaches that integrate not only the sensorimotor regions that have been traditionally explored but also the frontal regions that likely generate theta oscillations. Some studies suggested that dACC would process errors to output a control signal that specifies the need for adjustments toward downstream regions, including the prefrontal cortex, responsible for implementing corresponding adjustments ( Ridderinkhof et al., 2004 ;Shenhav et al., 2016 ). Given this assumption, investigating the coupling between prefrontal theta and SMA beta oscillations could be valuable in the context of motor adaptation.

Conclusion
The current study advances our understanding of adaptive learning in humans by demonstrating changes of oscillatory activity in multiple bands and regions and linking these changes to specific "adaptive " and "learning " motor processes. Consistent with previous studies, the results indicated that beta power, here of the supplementary motor area, underpins error-dependent changes in trial-to-trial performance and as such reflects neural processes that evaluate motor error and likely signal a need for adaptation of the motor output. The other main outcome is the contribution of frontal theta power to adaptive learning, yet without being related to errors. This suggests a role of theta oscillations in storing the newly learned internal model. Whether these spatially and band-wise distinct oscillatory activities constitute neural correlates of implicit motor adaptation, only, or also reflect cognitive strategies that may operate during motor adaptation to learn from errors remains to be tested.

Data and code availability statement
All raw data are stored on SUMMER, the data server of Grenoble-Alpes university and will be made available without restrictions upon request to lucas.struber@univ-grenoble-alpes.fr or fabien.cignetti@univgrenoble-alpes.fr. Analysis codes can be found on gitlab, at the following address https://gricad-gitlab.univ-grenoble-alpes.fr/lstruber-research .