M Couraud

D Cattaert

F Paclet

P Y Oudeyer

A De Rugy

Model and experiments to optimize co-adaptation in a simplified myoelectric control system

Keywords: Myoelectric control, Co-adaptation, Visuomotor rotation, Sensorimotor adaptation, Electromyography (EMG)

Objective. To compensate for a limb lost in an amputation, myoelectric prostheses use surface electromyography (EMG) from the remaining muscles to control the prosthesis. Despite considerable progress, myoelectric controls remain markedly di↵erent from the way we normally control movements, and require intense user adaptation. To overcome this, our goal is to explore concurrent machine co-adaptation techniques that are developed in the field of brain-machine interface, and that are beginning to be used in myoelectric controls. Approach. We combined a simplified myoelectric control with a perturbation for which human adaptation is well characterized and modeled, in order to explore coadaptation settings in a principled manner. Results. First, we reproduced results obtained in a classical visuomotor rotation paradigm in our simplified myoelectric context, where we rotate the muscle pulling vectors used to reconstruct wrist force from EMG. Then, a model of human adaptation in response to directional error was used to simulate various co-adaptation settings, where perturbations and machine co-adaptation are both applied on muscle pulling vectors. These simulations established that a relatively low gain of machine co-adaptation that minimizes final errors generates slow and incomplete adaptation, while higher gains increase adaptation rate but also errors by amplifying noise. After experimental verification on real subjects, we tested a variable gain that cumulates the advantages of both, and implemented it with directionally tuned neurons similar to those used to model human adaptation. This enables machine co-adaptation to locally improve myoelectric control, and to absorb more challenging perturbations. Significance. The simplified context used here enabled to explore co-adaptation settings in both simulations and experiments, and to raise important considerations such as the need for a variable gain encoded locally. The benefits and limits of extending this approach to more complex and functional myoelectric contexts are discussed.

Introduction

Myoelectric control systems are primarily designed to enable amputees to control the actuators of a prosthesis with the activity of their remaining muscles, recorded from surface electromyograms (EMG). To overcome the unnatural sequential aspect associated with successive movements recognized by conventional pattern recognitionbased controls, the field is now quickly moving toward simultaneous proportional control of multiple degrees of freedom [START_REF] Farina | The extraction of neural information from the surface emg for the control of upper-limb prostheses: emerging avenues and challenges[END_REF][START_REF] Jiang | Myoelectric control of artificial limbs-is there a need to change focus[END_REF][START_REF] Jiang | Extracting simultaneous and proportional neural control information for multiple-dof prostheses from the surface electromyographic signal[END_REF]. This proportional control can be achieved using various regression techniques that capture the natural relationship between EMG and movements [START_REF] Hahne | Linear and nonlinear regression techniques for simultaneous and proportional myoelectric control[END_REF]. However, getting relevant data to establish this relationship is not always possible, as for instance with amputees who can no longer produce movements with their missing limb. The activity of remaining muscles might still be regressed against movements intended by the amputees while producing muscle contractions, but injuries could be such that the resulting EMG-movement relationship might be too deteriorated to be used to control movement. And, if the activity of the remaining muscles gathered while amputees are simulating movements, are su cient to elicit a decent initial myoelectric control, signal nonstationarities due for instance to electrode shifts, sweat or varying postures, remain a major di culty to deal with [START_REF] Vidovic | Covariate shift adaptation in emg pattern recognition for prosthetic device control[END_REF].

While a strategy could be to let the sensorimotor system adapt to these imperfections and changes in signal quality, and despite some degree of flexibility of muscle patterns [START_REF] He | User adaptation in long-term, open-loop myoelectric training: implications for emg pattern recognition in prosthesis control[END_REF][START_REF] Ison | Proportional myoelectric control of robots: muscle synergy development drives performance enhancement, retainment, and generalization[END_REF][START_REF] Nazarpour | Flexible cortical control of taskspecific muscle synergies[END_REF], there is also evidence for limited adaptation capabilities. For instance, perturbations of the natural relationship between EMG and movement are mostly resolved within the existing repertoire of muscle activities, and the system has great di culties in learning perturbations for which no such solution exist [START_REF] Berger | Di↵erences in adaptation rates after virtual surgeries provide direct evidence for modularity[END_REF]de Rugy et al., 2012a). Relying solely on subjects' adaptation when faced to poor or deteriorating myoelectric controls appears therefore unsatisfactory.

Similar issues occur in brain-machine interfaces, where the usual strategy consists in designing a decoder that is as intuitive as possible, and rely on brain plasticity to adapt to the decoder's imperfections when using it to control movements. As for myoelectric controls, however, it is not always practical to establish an e cient initial decoder, and signal quality can deteriorate. Brain adaptation capabilities are also known to be limited, and recent studies illustrate that the structure of existing networks constitutes important constrains on adaptation to novel or perturbed decoders. For instance, a perturbation that requires novel relationships between the activities of the di↵erent motor cortical neurons used by the decoder is more di cult to learn than a perturbation of similar magnitude that maintain the original relationships [START_REF] Sadtler | Neural constraints on learning[END_REF]. When only a portion of the recorded neurons are used for the decoder, adaptation was also found that produced natural correlation with the activities of other recorded neurons, such that the solution generated by the nervous system was essentially already present in the repertoire of cortical activities [START_REF] Chase | Behavioral and neural correlates of visuomotor adaptation observed through a brain-computer interface in primary motor cortex[END_REF][START_REF] Hwang | Volitional control of neural activity relies on the natural motor repertoire[END_REF].

To overcome these limitations in terms of adaptation capabilities, a strategy is to complement human adaptation by an adaptation of the machine, or decoder. This strategy is currently being explored with increasing success in the context of brain machine interface [START_REF] Iturrate | Exploiting task constraints for self-calibrated brain-machine interface control using error-related potentials[END_REF][START_REF] Dangi | Design and analysis of closed-loop decoder adaptation algorithms for brain-machine interfaces[END_REF][START_REF] Gilja | A high-performance neural prosthesis enabled by control algorithm design[END_REF][START_REF] Orsborn | Closed-loop decoder adaptation on intermediate time-scales facilitates rapid bmi performance improvements independent of decoder initialization conditions[END_REF][START_REF] Shanechi | Robust brain-machine interface design using optimal feedback control modeling and adaptive point process filtering[END_REF][START_REF] Shenoy | Combining decoder design and neural adaptation in brain-machine interfaces[END_REF][START_REF] Vidaurre | Machine learning based co-adaptive learning: Towards a cure for bci illiteracy[END_REF], and starts also being developed in the context of myoelectric controls [START_REF] Hahne | Concurrent adaptation of human and machine improves simultaneous and proportional myoelectric control[END_REF][START_REF] Rezazadeh | Co-adaptive and a↵ective human-machine interface for improving training performances of virtual myoelectric forearm prosthesis[END_REF][START_REF] Tuga | Concept of a co-adaptive training environment for human-machine interfaces based on emgcontrol[END_REF]. However, many questions remain as to the particular form and setting of the machine co-adaptation that would best complement the human adaptation.

Here, we used a simplified myoelectric context, combined with a perturbation for which human adaptation is well characterized and modeled, to explore co-adaptation strategies in both simulated and real experiments. The myoelectric task consists in aiming for two dimensional targets on a computer screen with a cursor controlled with isometric force developed at the wrist joint, and reconstructed online from EMGs combined to a simplified representation of the muscles' biomechanics obtained with linear regressions. This virtual biomechanics technique (de Rugy et al., 2012b) represents the intuitive decoder that will receive experimental perturbations, and that will also be the site of machine co-adaptation. First we show that a 45 degree rotation of the muscle pulling vectors that form the virtual biomechanics, elicit a pattern of errors that is consistent with that obtained in a typical visuomotor rotation paradigm [START_REF] De Rugy | Changes in muscle directional tuning parallel feedforward adaptation to a visuomotor rotation[END_REF][START_REF] Krakauer | Learning of visuomotor transformations for vectorial planning of reaching trajectories[END_REF]. Then we used a model of this adaptation, via connection weights updating in a simple sensorimotor network of directionally tuned neurons (de Rugy, 2010; [START_REF] Thoroughman | Learning of action through adaptive combination of motor primitives[END_REF], to simulate various perturbations and co-adaptation settings. These simulations established that an optimal gain of machine co-adaptation that minimizes final error generates slow and incomplete adaptation, while higher gains can produce faster adaptation but also amplify noise. After verifying these predictions experimentally, we set up and demonstrate the e↵ectiveness of a variable and local gain of machine co-adaptation for absorbing more challenging perturbations that, by design, cannot be dealt with by the sole subjects' adaptation. Finally, we discuss the benefits and limits of extending this approach to more complex and functional myoelectric contexts.

Materials and Methods

Our task involved participants to hit 1 of 16 di↵erent possible targets by moving a cursor on the 2D space of a computer screen (figure 1(a)). The position of this cursor was determined either by the real force generated by the subjects at the wrist joint (i.e., 'force-driven' condition), or by the force reconstructed from muscle activities recorded while generating wrist forces (i.e., 'EMG-driven' condition). Participants were required to perform a fast and direct cursor motion aiming on the target while they received continuous visual feedback of the cursor position, and a feedback on movement time was delivered after each trial. This movement time corresponded to the time taken by the cursor to go from 10% to 90% percent of target distance from the center of the screen. Moreover, angular error was calculated for each trial as the di↵erence between the target angle and the initial direction of the movement, measured by movement direction at 80ms from movement initiation (de [START_REF] De Rugy | Changes in muscle directional tuning parallel feedforward adaptation to a visuomotor rotation[END_REF]. This was designed to prevent the influence of online feedback correction that might intervene as soon as 100ms from movement onset, and concentrate on the feedforward component of the movement for this measure and the (co-)adaptation it will elicit. Indeed, the model of human adaptation used here considers only the feedforward component of the sensorimotor adaptation system, and a more comprehensive model that includes adaptation in response to both feedback and feedforward error correction mechanisms is beyond the scope of the present contribution.

In the following, section 2.1 presents our simplified myoelectric context, the virtual biomechanics technique used to reconstruct forces from EMGs, and introduce the di↵erent perturbations applied on this biomechanics during simulations and experiments. Section 2.2 explains the model of human adaptation used to simulate various co-adaptation settings. Section 2.3 describes the algorithm used to perform machine co-adaptation, and the di↵erent gain values chosen to study co-adaptation e↵ects. Section 2.4 concentrates on the experiments conducted on real subjects, including the apparatus, recordings, and experiments details. Finally, section 2.5 indicates the analysis conducted on simulated and real data.

Virtual biomechanics

Pulling vectors. Cursor position corresponded to wrist forces reconstructed using a previously designed virtual biomechanics technique (de Rugy et al., 2012b), whereby muscles activities are multiplied by a set of pulling vectors that represent the directions and extents to which muscles are pulling per unit of activation. First, muscles activities were recorded while subjects produced forces toward targets evenly distributed in 16 directions at the same distance from the center. Muscles tuning curves illustrated figure 1(b) were then obtained by averaging, for each muscle and target direction, muscle activity when the force-driven cursor was maintained within the target (i.e. within 10% of target distance from the center of the target). Then, the virtual biomechanics was determined as the set of five pulling vectors (one per recorded wrist muscle, figure 1(c)) that minimizes the sum quadratic errors between target positions and the reconstructed reaches. In practice, we used linear regression to obtain this set of vectors.

Force reconstruction. In the EMG-driven condition, the cursor position corresponded to the wrist forces reconstructed online. Cartesian coordinates of this position were obtained by multiplying the rectified and filtered EMGs signals by the virtual biomechanics (i.e. the cartesian coordinates of associated vectors):

 X pos Y pos =  X vec1 • • • X vec5 Y vec1 • • • Y vec5 ⇥ 2 6 4 EM G 1 . . . EM G 5 3 7 5 (1)
Perturbations of the pulling vectors. Perturbations of the virtual biomechanics were introduced in order to simulate and test the e↵ects of di↵erent co-adaptation processes on task performance (figure 2). We used two types of perturbations: a +45 counterclockwise rotation of all vectors that compose the virtual biomechanics (figure 2(b)), and a vertical alignment of all these vectors (figure 2(c)). The overall rotation by 45 of the virtual biomechanics produces an e↵ect similar to the visuomotor rotation paradigm that is classically used to study sensorimotor adaptation, and for which human adaptation is well characterized and modelled (e.g. [START_REF] Krakauer | Learning of visuomotor transformations for vectorial planning of reaching trajectories[END_REF][START_REF] De Rugy | Changes in muscle directional tuning parallel feedforward adaptation to a visuomotor rotation[END_REF]. In contrast, the vertical alignment of all vectors provided a situation that is impossible to resolve by human adaptation only, as the cursor would be constrained on the vertical axis irrespective of the subjects' muscle activities. This situation is a good test because it requires a co-adaptation context to be resolved, whereby the decoder complements the participants' actions to reach targets that are not situated on that axis.

Model of human adaptation

We used a model of human adaptation (see Algorithm 1) in order to simulate the co-adaptation with di↵erent settings, and determine the parameters of machine co-adaptation that should optimize performance. The first aim of the model was to reproduce a typical learning curve obtained in response to visuomotor rotation, with rapid, exponential reduction of error, which stabilizes at an asymptotic level of residual error (figure 3(a)). This asymptotic level can be explained by the natural tendency of the system to return to baseline behavior in the absence of errors [START_REF] Galea | Dissociating the roles of the cerebellum and motor cortex during adaptive learning: the motor cortex retains what the cerebellum learns[END_REF], a behavior usually modeled by the forgetting factor, such as used below in equation (4). These general features of sensorimotor adaptation in response to perturbations are well described by the state space model, which reproduces a typical learning curve describing how motor intention is adapted on each trial as a function of the error on the preceding trial [START_REF] Donchin | Quantifying generalization from trial-by-trial behavior of adaptive systems that learn with basis functions: theory and experiments in human motor control[END_REF][START_REF] Smith | Interacting adaptive processes with di↵erent timescales underlie short-term motor learning[END_REF][START_REF] Thoroughman | Learning of action through adaptive combination of motor primitives[END_REF][START_REF] Haith | Theoretical models of motor control and motor learning[END_REF]. In order to account for the local property of the adaptation, whereby generalization to a learned rotation downgrades with distance from the practiced target [START_REF] Krakauer | Learning of visuomotor transformations for vectorial planning of reaching trajectories[END_REF], we implemented the adaptation in a network of directionally tuned neurons similar to that proposed in de Rugy (2010). The network was made of 20 sensory neurons (figure 4), each receiving target direction as input signal, and connected to a single output neuron that defines the motor intention. As our control space was a set of target that described a circle, the network had to manage circularity (figure 4 left). This was done using the following circular tuning function (cf. [START_REF] Amirikian | Directional tuning profiles of motor cortical cells[END_REF])

r i = exp k⇥cos(✓targ ✓i) exp k (2)
where r i is the firing rate of neuron i, k its receptive field, ✓ i its preferred direction and ✓ targ the input target direction. We set k = 15 in order to reproduce the width of the tuning functions involved in rotation adaptation in human (±23 , [START_REF] Tanaka | Adaptation to visuomotor rotation through interaction between posterior parietal and motor cortical areas[END_REF]. Each sensory neuron was then connected to the central motor neuron with di↵erent weights (w i ). The output of the network, the motor intention illustrated figure 4(b), is given by

✓ output = arg( 20 X i=1 r i ⇥ ! v i ) (3) with arg( ! v i ) = w i .
To simulate human adaptation, the direction that corresponds to this motor intention (figure 4(b)) was used to select muscle activity from the muscle tuning curves of the simulated subject (figure 4(c)-(d)), which was then coupled to the virtual biomechanics (figure 4(e)) to provide the movement output (figure 4(f)). In order to simulate subject-specific variability, a Gaussian white noise term with variance similar to that of directional error recorded for this subject at baseline was added to the movement output. The resulting position reached was then subtracted from the target position to obtain angular error, which was used to adapt both the human model (figure 4(h)) and the machine (figure 4(g); section 2.3). The network adaptation illustrated figure 4(h) was elicited by updating weights after each trial according to

w i.t+1 = w i.t + ↵ ⇥ ✏ t ⇥ r i.t ⇥ (w i.t w i.0 ) (4)
where coe cients ↵ (error sensibility) and (forgetting factor) are determined empirically (cf. section 3 'Results'), w i.0 is the initial weight of neuron i, w i its weight at trial t, r i.t its firing rate and ✏ is the angular error between target direction and movement direction at trial t. x mean(r ⇥ cos(✓))

6:

y mean(r ⇥ sin(✓))

7:

NetworkOutput arctan(x, y) . Model Adaptation : weights update 28:

w w + ↵ ⇥ ✏ ⇥ r ⇥ (w w0) . (4) 29:
listError[thisTrial] ✏

machine co-adaptation

Co-adaptation algorithm. We used an iterative co-adaptation algorithm to modify the virtual biomechanics according to the error produced in the preceding trial (figure 4(g)). The direction of each virtual biomechanics vector was changed ( ✓ i ) as a function of the angular error produced ✏ t , the associated muscle implication (i.e., the relative muscle activity M i.t involved in the trial producing that error), and a global gain value (g):

✓ i = g ⇥ ✏ t ⇥ M i.t (5)
where M i.t is the proportion of muscle activity, such that

P 5 i=1 M i.t = 1. Gain values.
The gain parameter g determines the speed and stability of adaptation, such that high values correspond to fast adaptation. For instance, a gain g of 1 with a muscle activity M i.t of 1 would elicit a rotation of the pulling vector of that muscle that would totally correct the error in only one trial. In practice, however, muscle sharing is such that a gain of 1 would not completely correct for the error in one trial, and the correction would be distributed over the contributing muscles. In our study, we tested two fixed values (a low or optimal one, and a high one), and a variable gain. The two fixed values were defined using simulations, such that the low (optimal) gain is the one that should minimize mean squared error calculated on the adapted state, and the high gain is the one that should allow maximal amount of angular error correction in the minimal number of trials (ideally, correcting all the error in one trial). The variable gain was inspired by the Resilient propagation algorithm (Rprop algorithm, [START_REF] Igel | Empirical evaluation of the improved rprop learning algorithms[END_REF], whereby the gain was updated by comparing the sign of the error at trial t-1 to the sign of the error at trial t: if the two signs were identical, the gain value was increased because the correction was beneficial but not su cient

g t = g t 1 ⇥ ⌘ + (6)
with ⌘ + = 1.2. In contrast, if the signs were di↵erent, the correction was considered excessive, and the gain value was therefore decreased

g t = g t 1 ⇥ ⌘ (7)
with ⌘ = 0.5. However, an adjustment of a global gain g based on the Rprop algorithm would not take into account the local requirement of error correction in the context of adaptation to a visuomotor rotation. Indeed, the nature of the perturbation, as well as that of the human adaptation, could be such that a requirement to change the gain of the adaptation in an area of the workspace would be detrimental if applied to another area of the workspace. To enable a local adjustment of gain based on the Rprop, we designed a model based on two 20-neurons networks comparable to that used to model human adaptation: one to locally encode the error (local errors network), and a second to locally encode the gain (local gains network; figure 5 and Algorithm 2). As for the human network, these networks were composed of directionally tuned neurons, and the networks' output was the angle of the resultant vector (equation ( 3)). However, the width of the tuning curves for the local errors network were adjusted to increase local resolution (in equation ( 2), k was replaced by k er = 45). The tuning curves of the local gains network was identical to the human network (k g = 15). Furthermore, while the output of the human network corresponds directly to the response of the human model used in the co-adaptation process, the di↵erence between the networks input and the output was used in these two local networks. In the local errors network, this di↵erence is used to retrieve the sign of the previous error encountered for a similar target direction. In the local gains network, this di↵erence is the new gain value associated with the target direction. Note that before starting the co-adaptation process, the weights of the local gains network were tuned so that the input-output di↵erence matched the initial gain. In the local errors network, the connection weights were updated after each trial as a function of the error obtained on that trial, according to (see equation ( 4) for variable significations):

w i.t+1 = w i.t + sign(✏ t ) ⇥ r i.t (8)
These weights were limited by a maximal and minimal value (limit = ±0.5) so that in most cases, w t+1 = w 0 0.5 when sign(✏ t ) is negative, and w t+1 = w 0 + 0.5 when sign(✏ t ) is positive. Therefore for a given input direction, the di↵erence between the input and output indicates the sign of the error encountered for that movement direction. This procedure was used to ensure that sign of successive errors would alternate between negative and positive values. The width of the tuning curves of this network were reduced using k er = 45 to ensure that the sign of the error encoded for a target direction is specific to that direction, and do not generalize to other target directions like in human model (where k = 15). When probed for a di↵erent direction on another trial, this network will return the sign of the error encountered for that direction. This sign will then be compared to the sign of the actual error produced on that trial, in order to update the weights of the local gains network according to:

w i.t = ( w i.t 1 + ⌘ + ⇥ g t 1 ⇥ r i.t
same signs

w i.t 1 + ⌘ ⇥ g t 1 ⇥ r i.t else (9)
where ⌘ + et ⌘ , the coe cients applied when the sign of the current error is similar (figure 5, case 1) or di↵erent (figure 5, case 2) from the sign of the previous error for a similar direction, were determined so that the dynamic of this network mimics that of a single Rprop. This was achieved with ⌘ + = 0.28 and ⌘ = 0.7. The output is calculated from equation (3), and the gain is the di↵erence between output and input. Moreover, we used an inferior limit of 0.04 for the gain in order to avoid negative or null gains to occur. Tuning curves with similar width to that used for the model of human adaptation (i.e. k g = 15) were used for this second network to benefit from similar generalization order properties. Critically, the weights of the local errors network were updated after its output was used by the local gains network, such that the sign of the current error is compared to that of the previous error for the same direction. The local gain coded by the local gains network was then used to adapt the virtual biomechanics. Figure 5 illustrates the two issues: (case 1) sign of current error and sign of previous error were similar (figure 5(a,b)), in this case, weights of the local gains network were updated aiming to increase gain value for that direction (figure 5(c,d)), and the local errors network recorded the current error according to its sign (figure 5(e)); (case 2) signs of current and previous error were di↵erent (figure 5(f,g)), weights of the local gains network were updated aiming to decrease gain value for that direction (figure 5(h,i)) and the local errors network recorded the current error according to its sign (figure 5(j)).

Experimental plan

Subjects. Twenty three subjects participated in the experiments. All were right-handed (mean ±sd score at handedness Edinburgh test were 87% ±13%; Oldfield, 1971) and had normal or corrected to normal vision. Five subjects constituted the "control" group, in which they had to adapt to a 45 rotation in the absence of machine coadaptation."Low Gain", "High Gain" and "Global Variable Gain" groups were formed by four subjects each, who had to adapt to the 45 rotation in the presence of coadaptation, with a low, a high and a global variable gain, respectively. The last group formed by six subjects was the "Local Variable Gain" group, which involved subjects to adapt to a vertical alignment perturbation with co-adaptation using a variable gain coded locally.

Apparatus. Subjects sat 80 cm from a computer display positioned at eye level and their right hand was kept in the custom-made manipulandum described previously (de Rugy et al., 2012b). The elbow was approximately at 110 with the forearm parallel to the table, and the hand was placed in a semi-pronated position. The wrist was fixed by an array of adjustable supports to fit the hand at the metacarpal-phalangeal gain thisT rial gain joints (12 contacts) and the wrist just proximal to the radial head (10 contacts). This apparatus allowed subjects to contract their muscles without limb movements. Moreover, flexion/extension and abduction/adduction forces were recorded using a 6-degrees-of-freedom force/torque transducer (FT DELTA SI-660-60, ATI Industrial Automation) coupled with the wrist manipulandum. These forces determined cursor position on the x and y axes respectively. EMG. Electromyograms were recorded from Extensor Carpi Radialis brevis (ECRb) and longus (ECRl), Extensor Carpi Ulnaris (ECU), Flexor Capri Ulnaris (FCU) and Flexor Carpi Radialis (FCR), with self-adhesive surface electrodes using the Delsys BagnoliTM EMG System. Wrist muscles were identified using palpation during resisted isometric contraction according to the SENIAM recommendations (www.seniam.org). Signals were, pre-amplified 10 times, amplified 1000 times, bandpass filtered from 20Hz to 450Hz, and sampled at 2kHz. To obtain muscle activity envelopes, EMG signals were centered with zero mean, rectified, and filtered using a 2 nd order Butterworth filter with 1.5 Hz cuto↵ frequency.

Experiments. The experiments consisted of two phases: a calibration phase and an experimental phase (see figure 3(b)). During the calibration phase, EMGs were recorded while subjects controlled the cursor with wrist forces (force-driven condition). This phase provided the data necessary to establish the virtual biomechanics as described section 2.1. During the experimental phase, subjects controlled the cursor with their EMGs (EMG-driven condition). This phase is composed of three stages illustrated in figure 3(b): baseline (5 blocks of 16 targets), learning (20 blocks of 16 targets) and after e↵ect (3 blocks of 16 targets). During the baseline stage, subjects got used to control the cursor with their EMGs, and produced their baseline behavior in that condition. During the learning stage, subjects adapted to the perturbation of the virtual biomechanics with or without co-adaptation of the machine. During the after e↵ect, the perturbation was removed and there was no machine co-adaptation. This last stage allowed to indirectly assess the e↵ectiveness of the machine co-adaptation by the end of the learning phase. Indeed, e cient co-adaptation should absorb the perturbation and revert the virtual biomechanics to that of the initial EMG to force mapping, thereby eliminating after-e↵ects. This condition was also particularly useful to distinguish the subjects' adaptation from that of the machine, which is not easy when both are concurrently operating.

Analysis

To evaluate the e↵ects of co-adaptation, di↵erent analyzes were conducted: (i) the speed of adaptation was assessed as the averaged error on the first five blocks of the learning stage (figure 3(a), point 2); (ii) the final performance was assessed as the averaged error on the last ten blocks of the learning stage (figure 3(a), point 3); and (iii) the after e↵ect was assessed as the averaged error on the first block of the after e↵ect stage (figure 3(a), point 4). To compare angular errors obtained with di↵erent gain values we performed a Wilcoxon rank-sum test. Noise amplification was assessed by performing a Wilcoxon rank-sum test on the standard deviation of directional error over the last 10 blocks of the learning stages. To determine if an after e↵ect was present, we compared data against an error of zero by using a one-sample t-test. To compare the initial versus final virtual biomechanics vectors directions, we calculated the angular di↵erence between the initial angle and the final angle for each pulling vector of each subject (figure 6). Then, we performed di↵erent statistical tests on simulated and real data to accommodate for di↵erences in sample size between these two groups. For simulated data we performed a Wilcoxon signed rank test and for real data we compared data against zero by using a one-sample t-test. The significance level was fixed at ↵ = 0.05 and median values (Med) were given throughout.

Results

Before doing the experiments on real subjects, we used the model of human adaptation to simulate various co-adaptation settings. At first, we set the coe cients ↵ and used to update the connection weights in the model of human adaptation so as to fit the time course of adaptation obtained in de [START_REF] De Rugy | Changes in muscle directional tuning parallel feedforward adaptation to a visuomotor rotation[END_REF]. The resulting values were ↵ = 0.170 and = 0.004, and the corresponding fit is illustrated in black figure 7. Because of the slight di↵erence in method and results (i.e., in the present study we applied the rotation perturbation on the pulling vectors used to reconstruct force from EMGs rather than directly on force), we subsequently adjusted the coe cients ↵ and to fit the data obtained here from our control group. This resulted in ↵ = 0.172 and = 0.003 (cf. red data and fit figure 7), the two parameter values that we used subsequently to simulate co-adaptation settings. It is important to note that despite small di↵erences, the fits obtained for a standard visuomotor rotation and for the present rotation of the muscle pulling vectors are very close, which justifies our protocol and model of human adaptation to explore co-adaptation.

Without and with fixed gains of co-adaptation

For the first set of simulations, we used a simple perturbation that consisted in a 45 counterclockwise rotation of all virtual pulling vectors. In order to simulate subjectspecific co-adaptation, we used the virtual biomechanics and the muscle patterns of real subjects, as well as their specific level of noise. This was done by setting the standard deviation of the white Gaussian noise term added to the simulated reached direction (figure 4(f)) to that of the directional error obtained on the baseline blocks of the corresponding subjects (see section 2.2). Three di↵erent simulations were conducted with the noise level of each of the five control subjects, thereby giving results for 15 simulated subjects. Figure 8(a)-(c) displays raw data simulated for a representative subject in the di↵erent co-adaptation conditions, and provides a good illustration of the level of noise recorded in the experiments.

Simulated and real results conducted without co-adaptation are qualitatively similar.

Results of simulation conducted in the absence of co-adaptation (blue points and lines figure 8(d)-(f); n = 15) are qualitatively similar to those of control group (figure 8(g)-(i); n = 5). These results show a quick initial adaptation, which stabilizes at a relatively low level or residual error (figure 8 (d),(g); points 2 and 3). An after e↵ect was present when the perturbation was removed, both in the simulated and in the real subjects (figure 8(d),(g); point 4; t (14) = 33.3, p < 0.001 and t (4) = 21.7, p < 0.001 respectively). Moreover, the level of variability of directional error remains relatively low for both the simulated and the real subjects during adaptation (figure 8

(e)-(h)).

Optimal fixed co-adaptation results in slow and incomplete adaptation.

The optimal gain of machine co-adaptation was defined as the one that elicits the lowest final mean squared error (MSE), which accounts for systematic bias as well as the variability around it, thereby promoting both accuracy and precision. To determine this optimal low gain, we evaluated 100 simulations per gain, for gain values between 0 and 1.5, conducted with the muscle patterns, the virtual biomechanics, and the noise level of our control subjects. The resulting averaged MSE calculated on the last 10 blocks are displayed figure 9 . As can be seen on that figure, minimal errors were obtained for gain values around 0.05, but there were also substantial variations in this optimal value due to subject-specific characteristics and noise level. To account for this inter-subject variability, each subject of the low gain group was tested with his/her specific optimal gain value, determined with his/her pulling vectors, muscle pattern, and noise level.

Results of simulated subjects (figure 8; n = 12) conducted in the presence of co-adaptation with the low gain (red points and lines) show a faster adaptation than without co-adaptation (figure 8(d); point 2; U = 0, p < 0.001 with Med without = 25.5 and Med lowGain = 14.7) and a final averaged angular error closer to zero (figure 8(d); point 3; U = 0, p < 0.001 with Med w ithout = 6.5 and Med lowGain = 1.6). However, a residual after e↵ect was predicted by the simulations (figure 8(d); point 4). Indeed, we found a significant di↵erence between the initial error on the first block after the perturbation has been removed and an error of zero (t (11) = 7.8, p < 0.001). This after e↵ect was significantly lower than the after e↵ect observed for the simulated "control" group (figure 8(d); point 4; U = 0, p < 0.001 with Med without = 34.9 and Med lowGain = 10.5). Figure 10(b) illustrates that during simulation of co-adaptation with the low gain, the angles of the virtual biomechanics vectors gradually changes toward their initial position. However, as the gain of the machine co-adaptation was low, the simulated human adaptation is correcting for part of the error, such that the behavior stabilizes before the virtual biomechanics has reverted to its initial, unperturbed values. Moreover, we found a significant di↵erence between initial and final vector angles (figure 8(k); W = 1830, p < 0.001) and the median di↵erence angle was equal to 13.6 .

Results of real subjects (n = 4) show a similar pattern of results. Indeed adaptation with a low, "optimal" co-adaptation gain was faster than without co-adaptation (figure 8(g); point 2; U = 0, p = 0.016 with Med without = 24.2 and Med lowGain = 13.9), and the final averaged angular error was closer to zero (figure 8(g); point 3; U = 0, p = 0.016 with Med without = 6.1 and Med lowGain = 0.4). However, although results show that average performance of "low gain" group were lower than performance observed for the "control" group (figure 8(g); point 4; U = 0, p = 0.016 with Med without = 35.1 and Med lowGain = 10.6), the after e↵ect predicted by simulation was not significantly di↵erent to the null error (t

(3) = 1.7, p = 0.183). Nevertheless, angular di↵erence between initial and final biomechanics vectors direction was significantly di↵erent to zero (figure 8(l); t (19) = 3.7, p = 0.001)) and was comparable to that obtained in simulations (figure 8(k)). Therefore, although optimal low gain allows better final performance, simulated results suggest that a residual human adaptation still exists due to slowness of machine co-adaptation.

High fixed co-adaptation results in faster but unstable adaptation.

The high gain of co-adaptation corresponded to the value that allowed correcting the maximal amount of angular error in the minimal number of trials, ideally, correcting all the error in one trial. However, because the machine co-adaptation equation ( 5) is contingent upon muscle sharing, this gain is di↵erent for each target and for each subject. For each target direction, we therefore determined the specific gain value that corrects for the entire error in that direction in one trial, before averaging these target-specific gains over all directions for each subject.

Results of simulated subjects (n = 12) with the high gain of co-adaptation (green points and lines figure 8) show an adaptation that is faster than without co-adaptation (figure 8(d); point 2; U = 0, p < 0.001 with Med without = 25.5 and Med highGain = 1.9), and than with the low gain of co-adaptation (U = 0, p < 0.001 with Med lowGain = 14.7 and Med highGain = 1.9). Final averaged angular errors were significantly lower than that of the "control" group (figure 8(d); point 3; U = 0, p < 0.001 with Med without = 6.5 and Med highGain = 0.02) and that of the "low gain" group (U = 0, p < 0.001 with Med lowGain = 1.6 and Med highGain = 0.02). There was also no after e↵ect (t (11) = 0.85, p = 0.415), and di↵erence between initial and final vector angles was not significant (figure 8(k); W = 398, p = 0.144). Moreover, the directional errors in the first after e↵ect block were significantly lower than the after e↵ect observed for the "control" group (figure 8(d); point 4; U = 0, p < 0.001 with Med without = 34.9 and Med highGain = 0.2), and than that observed for the "low gain" group (U = 4, p < 0.001 with Med lowGain = 10.6 and Med highGain = 0.2). However, standard deviation (figure 8(e)) showed that this high gain of co-adaptation induced a level of variability that is significantly higher than for the "control" and the "low gain" groups (U = 14 and U = 15, respectively, ps < 0.001 with Med without = 18.3, Med lowGain = 15.0 and Med highGain = 28.9). Indeed, as the gain of machine co-adaptation is high, small errors that might be due to noise rather than to a genuine error in the motor output, might induce unnecessary large corrections in the virtual biomechanics. This excessive sensitivity of the machine coadaptation to error is responsible for the unstable direction of the virtual biomechanics vectors displayed figure 10(b), and explained the noise amplification observed in the simulations (figure 8(e)-(f); point 3).

Results of real subjects (n = 4) show a similar pattern of results. Indeed, adaptation with the high co-adaptation gain was faster than without co-adaptation (figure 8(g); point 2; U = 0, p = 0.016 with Med without = 24.2 and Med highGain = 0.1), and than with the low gain of co-adaptation (U = 0, p = 0.029 with Med lowGain = 13.9 and Med highGain = 0.1). The final averaged angular error was closer to zero than the "control" group (figure 8(g); point 3; U = 0, p = 0.016 with Med without = 6.1 and Med highGain = 1.1), but not di↵erent from the "low" gain coadaptation group (U = 3, p = 0.2 with Med lowGain = 0.4 and Med highGain = 1.1). Noise amplification was higher than for the "control" group, as indicated by standard deviations of directional error that were significantly higher for the high co-adaptation gain than for the "control" group (figure 8(h); point 3; U = 0, p = 0.016 with Med without = 14.0 and Med highGain = 26.4). Moreover, there was no after e↵ect (figure 8(g); point 4) as directional errors on the first block after perturbation was removed were not significantly di↵erent from zero (t

(3) = 0.5, p = 0.654), and as on simulated data, di↵erence between initial and final vector angles was not significantly di↵erent to zero (t (19) = 0.388, p = 0.702). In sum, although the high gain enables faster adaptation and better final performance, it also amplifies noise and results in more variables errors.

With variable gain values

Global variable gain can deal e ciently with global rotation, but not with more complex perturbations.

We tested the e↵ects of co-adaptation with global variable gain values (Rprop) to a 45 rotation on simulated and real subjects (n = 12 and n = 4, respectively, figure 11). Results show a rapid adaptation (figure 11(a,d); point 2) with a final error level near to null error (figure 11(a,d); point 3) and without any noise amplification (figure 11(b,e)). In sum, we cumulate low and high gain advantages: fast adaptation and final performance close to null error without noise amplification nor after e↵ect. Indeed, Rprop simulated results show a faster adaptation (figure 11(a); point 2) than "control" and "low gain" groups (U = 0 and U = 3, respectively, ps < 0.001 with Med without = 25.5; Med lowGain = 14.7 and Med R propGain = 2.9) and was no significantly di↵erent from "high gain" group (U = 44, p = 0.112 with Med highGain = 1.9 and Med R propGain = 2.9). In the same way, variability on the last 10 blocks (figure 11(b); point 3) was lower than high gain variability (U = 7, p < 0.001 with Med highGain = 28.9 and Med R propGain = 17.7). Finally, even if the angle di↵erence between initial and final virtual biomechanics vectors was significant (figure 11(h); W = 694, p = 0.011), an after e↵ect was no present(figure 11 (a); point 4; t (11) = 1.9, p = 0.081). Results obtained on real subjects confirm these simulated results, showing a faster adaptation than for the "control" and "low gain" groups (figure 11(d); point 2; Us = 0, p = 0.016 and p = 0.029, respectively, with Med without = 24.2; Med lowGain = 13.9 and Med RpropGain = 3.4), and that was not significantly di↵erent from the "high gain" group (U = 4, p = 0.343 with Med highGain = 0.1 and Med RpropGain = 3.4). Noise amplification at the end of learning was lower than with a high gain (figure 11 (e); point 3; U = 0, p < 0.029 with Med highGain = 26.4 and Med RpropGain = 15.3). Finally, there was no after e↵ect (figure 11 (d); point 4; t

(3) = 1.8, p = 0.163), and angular di↵erence between initial and final vectors was not significantly di↵erent from zero (figure 11 (i); t (19) = 0.245, p = 0.809).This co-adaptation with a global variable gain was therefore easily able to deal with a global 45 perturbation applied to all muscle pulling vectors. Next, we chose to focus on a more challenging perturbation that elicits di↵erent errors for di↵erent areas of the workspace, and that, by design, cannot be compensated for by the sole human adaptation. This perturbation consisted in aligning all muscle pulling vectors to a single, vertical axis, such that irrespective of the recorded muscle activities, the cursor could only move on that axis. As with fixed gain value, we began by testing co-adaptation through simulations. We first simulated adaptation to the vertical alignment of biomechanics vectors without co-adaptation and with global gain algorithm (Rprop). Raw simulated data confirm that without co-adaptation (blue data), subjects could not adapt their behavior to perform the task (figure 12(a)). Simulations conducted with the Rprop algorithm implemented with a single, global gain (figure 12(b)) result in an averaged error that is close to zero from early on, at the beginning the adaptation (purple data figure 12(d)). However, instead of an indication of an e cient adaptation, this small initial averaged error is due to the fact that large consecutive errors of opposing signs are cancelling each other out. Indeed, high variable errors and MSE are obtained for the initial part of the adaptation, and maintained for about half of the learning blocks (figure 12(e),(f)). This is because early in the adaptation, consecutive errors of opposite signs elicit a rapid decrease of the variable gain to its minimal value (see red dot figure 10(c)), which resulted in only minimal changes in the muscle pulling vectors angles (figure 10(d) and figure 12 (k)).

From early on in this condition, the only adaptation that operates is therefore that of the subjects, which is applied on a set of pulling vectors that are just slightly deviated from the perturbed, vertical position (figure 12 (j), (k)).

Local variable gain can absorb complex perturbations.

In contrast, simulations of co-adaptation with a gain that is coded locally (n = 15; figure 12(c) green data) show a fast adaptation with reduced variable errors in response to the more complex vertical alignment perturbation (figure 12(d); points 2 and 3). Results of real subjects (n = 6; figure 12(g)-(i), yellow data) were similar to those of simulations (figure 12(d)-(i), black data). Indeed, there were no after e↵ect (figure 12(g); point 4), and the noise was not amplified in the learning blocks, nor in the after e↵ect blocks, as indicated by low variability (figure 12(h)-(i)). The network used to implement this local variable gains allow to keep high gain values for target directions that required it, while simultaneously enabling low gain for others (figure 10(c)). These gains induced a better convergence of the direction of the muscle pulling vectors toward their initial positions (figure 10(d) and figure 12(l)).

Discussion

To overcome di culties associated with establishing and maintaining an e cient decoder despite limited neural recordings and changes in signal quality, closed-loop adaptation of the decoder are increasingly studied and used to complement subjects' adaptation in brain machine interface [START_REF] Iturrate | Exploiting task constraints for self-calibrated brain-machine interface control using error-related potentials[END_REF][START_REF] Dangi | Design and analysis of closed-loop decoder adaptation algorithms for brain-machine interfaces[END_REF][START_REF] Orsborn | Closed-loop decoder adaptation on intermediate time-scales facilitates rapid bmi performance improvements independent of decoder initialization conditions[END_REF][START_REF] Shenoy | Combining decoder design and neural adaptation in brain-machine interfaces[END_REF]. As similar issues occur for myoelectric controls, co-adaptation strategies are starting to be explored there as well [START_REF] Hahne | Concurrent adaptation of human and machine improves simultaneous and proportional myoelectric control[END_REF][START_REF] Rezazadeh | Co-adaptive and a↵ective human-machine interface for improving training performances of virtual myoelectric forearm prosthesis[END_REF][START_REF] Tuga | Concept of a co-adaptive training environment for human-machine interfaces based on emgcontrol[END_REF]. In both situations, however, the dynamics of the concurrent adaptation of the human and machine needs to be carefully considered since it is likely to be complex and to greatly impact the e↵ectiveness of particular co-adaptation settings. Here, we combined a simplified myoelectric control with a perturbation protocol for which human adaptation is well characterized and modeled, in order to explore co-adaptation strategies in a principled manner, through both simulations and real experiments.

First, we showed that adaptation to a 45 degree rotation of all muscle pulling vectors used to reconstruct wrist forces in our study was similar to that observed in response to a visuomotor rotation applied to the task controlled with real forces (cf. figure 7, de [START_REF] De Rugy | Changes in muscle directional tuning parallel feedforward adaptation to a visuomotor rotation[END_REF]. This is critical because visuomotor rotation is a well-studied protocol, for which human adaptation is well characterized and modelled [START_REF] Krakauer | Motor learning and consolidation: the case of visuomotor rotation[END_REF][START_REF] Krakauer | Learning of visuomotor transformations for vectorial planning of reaching trajectories[END_REF]. For instance, the typical learning curve obtained in response to sensorimotor perturbation is reproduced by a state space model that describes how motor intention is adjusted on a trial-by-trial basis as a function on errors [START_REF] Donchin | Quantifying generalization from trial-by-trial behavior of adaptive systems that learn with basis functions: theory and experiments in human motor control[END_REF][START_REF] Smith | Interacting adaptive processes with di↵erent timescales underlie short-term motor learning[END_REF][START_REF] Thoroughman | Learning of action through adaptive combination of motor primitives[END_REF]. To additionally account for the local aspect associated with adaptation of movement direction, whereby generalization reduces away from the trained direction [START_REF] Krakauer | Learning of visuomotor transformations for vectorial planning of reaching trajectories[END_REF], the trial-by-trial adaptation can be implemented by updating connection weights in a network of directionally tuned neurons [START_REF] Georgopoulos | Mental rotation of the neuronal population vector[END_REF][START_REF] Tanaka | Adaptation to visuomotor rotation through interaction between posterior parietal and motor cortical areas[END_REF][START_REF] Taylor | Feedback-dependent generalization[END_REF]. Here, we used a previous implementation of such network (de Rugy, 2010), in which the width of the directionally tuned neurons and the gain of the adaptation were set so as to reproduce the adaptation observed here and in previous study (cf figure 7). With this model of human adaptation reproducing results observed with perturbation of the muscle pulling vectors in our simplified context, we are well equipped to explore various co-adaptation settings. Indeed, we can simulate concurrent adaptation with whatever algorithm and settings of the machine co-adaptation applied as a function of directional error, and test the most interesting settings on real subjects.

Simulations conducted to determine the optimal gain of co-adaptation that elicits the lowest mean squared error at the final stage of the adaptation revealed that this gain improves performance but remains perfectible. Indeed, if applying this relatively low gain of machine co-adaptation elicits a more complete adaptation, by reducing the residual error that is associated with the forgetting factor in the model of human adaptation, the rate of adaptation is only moderately improved. Four or five cycles of 16 trials are still necessary to reduce averaged error to approximately nil, and the persistence of an after e↵ect indicates that even after the 320 trials of practice, the modelled adaptation of the subjects is still responsible for part of this result. Together with experimental results that were in good agreement with these predictions (figure 8), these features are calling for a more e cient machine co-adaptation. Subsequent simulations conducted with a higher gain, set to provide the highest possible machine co-adaptation in the lowest number of trials, were able to elicit much faster adaptation, with averaged directional error approximately nil within the first cycle of 16 trials. However, this dramatic improvement in adaptation rate came at a great expense in term of movement variability, which translates into higher mean squared errors. This is because the high gain of machine co-adaptation is amplifying the natural noise of the system, applying substantial corrections in response to errors that might pertain more to noise in the system than to genuine errors in motor intention. A comparable result was noted in [START_REF] Hahne | Concurrent adaptation of human and machine improves simultaneous and proportional myoelectric control[END_REF], when the recursive least square algorithm used to compute linear regressions from incoming data was set in a way that gives too much emphases on new data, resulting in faster but also unstable adaptation patterns.

The previous results obtained with fixed gains of machine co-adaptation are therefore calling for a variable gain, capable of producing high corrections when needed early in the adaptation, followed by smaller corrections to stabilize performance. This can be e ciently achieved with the Resilient propagation algorithm (Rprop algorithm), where the rate of learning is adjusted according to the sign of successive errors [START_REF] Igel | Empirical evaluation of the improved rprop learning algorithms[END_REF][START_REF] Riedmiller | A direct adaptive method for faster backpropagation learning: The rprop algorithm, Neural Networks[END_REF]. Interestingly, human sensitivity to errors has recently been shown to also depend on the history of errors, and a memory of errors appears critical to explain the observed faster relearning of previously encountered perturbations [START_REF] Herzfeld | A memory of errors in sensorimotor learning[END_REF]. Simulations conducted with a variable gain set according to the Rprop confirmed that the advantages of high and low gains can be cumulated (figure 11), but also revealed that it should be implemented in a local manner to be able to deal with more complex perturbations. Indeed, if a single gain modified with Rprop can elicit faster initial adaptation with little noise amplification for a global perturbation of all pulling vectors by 45 degree in the same direction, this method is unable to deal e ciently with a perturbation that aligned all pulling vectors on the vertical axis (figure 12). This is because such perturbation elicits directional errors of di↵erent sign for di↵erent areas of the workspace, such that the randomized presentation of targets can elicit large errors of di↵erent consecutive signs that prematurely reduce the learning rate. Simulations conducted in this context showed that the aligned perturbation can nevertheless be partly absorbed, but this is only due to subjects adaptation working on a set of pulling vectors that have been slightly deviated from vertical by the initial machine co-adaptation, before consecutive errors of changing signs had reduced machine co-adaptation to zero. As a result, the pattern of directional error is quite variable (i) for a large part of the adaptation, where only the subjects adaptation operates and takes time to absorb the pattern of di↵erent local errors, and (ii) in the aftere↵ect, where errors of various signs and directions initially occur when the vectors are reverted to their original directions, as an evidence of subjects adaptation that took place in the preceding phase.

To respond to these problems, we implemented a local version of the variable gain with two networks of directionally tuned neurons similar to that used to model subject adaptation: a first one that locally encodes the error, and a second one that locally encodes the gain of machine co-adaptation. As the first network enables comparing the sign of the current error with that produced earlier in the same area of the workspace, the second network can use this information to update and memorize the gain of machine co-adaptation according Rprop for that specific area. Simulations and experiments revealed that this mechanism was capable of correcting e ciently for the aligned perturbation that, by design, cannot be dealt with by the sole adaptation of subjects (figure 12). The fact that the adjusted pulling vectors were close to the original ones by the end of the adaptation indicates that this mechanism should converge toward a sensible solution in cases where the original mapping cannot be determined, as for instance with amputees that can no longer produce force of movement with their missing limb.

Although promising, several limitations would have to be carefully addressed in order to extend this approach to more complex and functional myoelectric controls. In particular, our simplified myoelectric context concentrates on the directional error produced by the feedforward component of movement production, and rests upon an initial calibration of the muscle pulling vectors. Each of those aspects raises its own concerns. While our linear decoder was biologically interpretable in terms of muscle pulling vectors, less interpretable non-linear decoders could be used to map EMG to movement in more complex situations. For instance, support vector regression techniques enables myoelectric control in situations where electrodes are positioned randomly, and where recordings are performed on moving limbs, which implies relative displacement between surface electrodes and recorded muscles, and a more complex EMG-movement relationship [START_REF] Ameri | Support vector regression for improved real-time, simultaneous myoelectric control[END_REF]. As already mentioned, initial calibration is not always practical. Although it could be removed, or roughly guessed and adapted as new data are coming in, this usually requires an accurate estimate of movement intentions to be mapped onto recorded control signals. This is for instance the case in [START_REF] Hahne | Concurrent adaptation of human and machine improves simultaneous and proportional myoelectric control[END_REF], which uses a recursive least square algorithm to update, with incoming data, a linear regression model between EMG recordings and target positions. This strategy seems to work well for slow target movements, for which we can reasonably assume that the ongoing movement intention correspond to the target. However, movement intention is more di cult to assess for faster natural movements, which typically combine feedforward mechanisms with visuomotor feedback corrections of varying magnitude and delay [START_REF] Dimitriou | The temporal evolution of feedback gains rapidly update to task demands[END_REF][START_REF] Scott | Feedback control during voluntary motor actions[END_REF]. This has been partly addressed in brain machine interface, where before being used for decoder adaptation, the estimation of movement intention is improved either with movement likelihood from history, using Kalman filter [START_REF] Gilja | A high-performance neural prosthesis enabled by control algorithm design[END_REF], or with an optimal feedback control model [START_REF] Shanechi | Robust brain-machine interface design using optimal feedback control modeling and adaptive point process filtering[END_REF]. Although similar approaches could be employed for myoelectric controls, an adaptation of the decoder dealing simultaneously with movement direction and magnitude would require, to enable exploration of various co-adaptation settings through simulations as was done here, a model of human adaptation also dealing with both. Indeed, adaptation to perturbation of movement direction and extent are known to rely on di↵erent neural substrates, to generalize di↵erently, and even to be represented in di↵erent coordinate systems [START_REF] Davare | Dissociable contribution of the parietal and frontal cortex to coding movement direction and amplitude[END_REF][START_REF] Krakauer | Learning of visuomotor transformations for vectorial planning of reaching trajectories[END_REF][START_REF] Vindras | Error parsing in visuomotor pointing reveals independent processing of amplitude and direction[END_REF][START_REF] Poh | Distinct coordinate systems for adaptations of movement direction and extent[END_REF]. The model used here for direction only would therefore be insu cient. Biomimetic models of sensorimotor cortex are being developed for the specific purpose of informing the development of rehabilitation strategies and brain machine interfaces [START_REF] Dura-Bernal | Towards a real-time interface between a biomimetic model of sensorimotor cortex and a robotic arm[END_REF], 2016). However, the important contribution from spinal networks [START_REF] Raphael | Spinal-like regulator facilitates control of a two-degree-of-freedom wrist[END_REF][START_REF] Tsianos | Useful properties of spinal circuits for learning and performing planar reaches[END_REF], occurring downstream from the brain but clearly implicated in muscle activation and movement control, would also need to be included for the case of myoelectric control. Finally, the two dimensional workspace used here and in most brain machine interfaces would have to be extended to the three dimensional space of our environment, and to the multiple degrees of freedom of a prosthesis [START_REF] Wodlinger | Ten-dimensional anthropomorphic arm control in a human brain-machine interface: di culties, solutions, and limitations[END_REF]. In sum, while the simplified context used here enabled to explore co-adaptation in both simulations and experiments, and to raise important considerations, application to more complex and functional settings will require intense research to solve some of the issues mentioned here. a, Representation of the set of virtual biomechanics vectors at three di↵erent stages: initial calibrated stage (full arrows), perturbed stage (dashed arrows) and final stage (dotted arrows). (b, Angular di↵erences between initial and final states (shaded areas) and between final and perturbed states (filled areas) are projected on a referent vector ( vertical grey arrow used as a reference for the initial direction). Resultant vector (large black arrow) is calculated and illustrates the di↵erence between initial and final virtual biomechanics. A persistence of this di↵erence, such as the one illustrated here, could explain the presence of an after e↵ect after the perturbation has been removed.

Real results

Simulated results

Visuomotor rotation :

Real results

Simulated results

Virtual biomechanics rotation :

Figure 1 .

 1 Figure1. Online force reconstruction from muscle virtual biomechanics and realtime electromyograms (EMGs). Subjects produced isometric force at the wrist joint to reach for 1 of 16 equally distributed target directions (a). In force driven condition, the cursor moves according to the force developed at the wrist joint (wrist force). EMGs recorded during this task are used to generate the muscle tuning curves (b) that are used to extract the muscle pulling vectors (i.e. virtual biomechanics, c). Online force reconstruction is then obtained by multiplying the rectified filtered EMG signals to the pulling vector of each muscle (d). Note that EMG signals are shown from 1 trial only, whereas EMGs from 3 consecutive trials to each of the 16 targets are used to compute muscle tuning curves. Also note that for the purpose of illustration, the same EMG signals are used for o✏ine averaging and online force reconstruction, whereas in the experiments, EMG for those 2 processes come from di↵erent acquisition blocks (i.e., in force-driven and EMG-driven conditions, respectively). Adapted from deRugy et al. (2012b)

Figure 2 .Figure 4 .Figure 5 .Figure 6 .

 2456 Figure 2. Virtual biomechanics perturbations. a. Initial calibrated virtual biomechanics. b. Perturbation 1: counterclockwise 45 rotation of all muscles pulling vectors. c. Perturbation 2: vertical alignment of the pulling vectors that form the virtual biomechanics.

Figure 7 .Figure 9 .

 79 Figure7. Mean angular error according to block number in two di↵erent experiments. Red data points (±sd : shaded areas) correspond to real data obtained by de[START_REF] De Rugy | Changes in muscle directional tuning parallel feedforward adaptation to a visuomotor rotation[END_REF] where subjects control cursor via wrist forces and have to adapt to a 45 visuomotor rotation between the direction of wrist forces and that of the cursor used to reach targets. Black data points (±sd : shaded areas) correspond to real data obtained in our experiments with control subjects. The two lines represent simulated results, with the model parameters tuned to fit the real data of the corresponding experimental situation.

Acknowledgments

This work has been funded by a PEPS Idex Bordeaux /CNRS.

Typical learning curve obtained in response to a visuomotor rotation by 45 , and after-e↵ect observed when the perturbation is removed. During baseline block A, a human subject performs the task with a mean error near to 0 (point 1). Then, the subject produces a large initial error in response to the perturbation (block B), quickly adapts his behavior (point 2), and stabilizes at an asymptotic level of error di↵erent from zero (point 3). When the perturbation is removed, an after-e↵ect occurs, and the subject re-adapts his behavior toward initial performance (point 4).

b. Experimental plan. In force-driven condition, subject's wrist forces are used to control the cursor movement. In EMG maximal, subject exerts maximal forces in each direction in order to normalize EMG during EMG condition. Then EMG/Force relationships are used to calibrate the EMGs in three consecutive blocks (calibration). In EMG-driven condition, subject's EMGs are used to control the cursor movements. Subject exerts to reach all targets during 5 blocks (baseline). Then the biomechanics is rotated (perturbation) for the next 20 blocks, and the subject adapts to this perturbation (learning). Finally the initial biomechanics is applied again and re-adapt to this original condition for 3 blocks (after-e↵ect)