Four monkeys participated in the experiments. They were chronically implanted with multielectrode arrays in motor (M1) and somatosensory cortices (S1)1,35. Extracellular electrical activity from 570-783 neurons was sampled chronically from monkeys M (410-501 neurons), C (196-214), O (156-229) and K (140-156) up to five years after implantation36. B2s operated using at least 550 neurons, while the B3s operated with at least 775 neurons. During experimental sessions, monkeys sat in separate rooms where each of them viewed a realistic monkey avatar arm on a computer screen (Fig. 1A). Movements of the avatar arm were generated using either a hand-held joystick (hand control) or from modulations in VLSBA in a mode called brain control4,5. In brain control, a neural decoding algorithm called the unscented Kalman filter (UKF)37 extracted reach kinematics from neuronal ensemble activity in real-time.

Experiment 1: Shared control

In the first experiment (shared control task, Fig. 1B), two monkeys worked together to move the avatar arm in 2D space. Each monkey contributed equally (50%) to both the X and Y coordinates. The monkeys moved the avatar arm from the screen centre towards a circular target (Movies S1, S2). These movements were enacted either using the joystick (dyad M&O) or the Brainet (monkey dyads M&O and M&C). To obtain training data for brain control, we started each session with a passive observation epoch during which monkeys watched the computer-generated avatar arm movements along centre-out trajectories38. One UKF was fit per monkey and the UKF outputs (X and Y coordinates) were combined across monkeys to produce the conjoint brain-controlled movements of the avatar arm.

The behaviour of the dyad MO improved significantly with training (Fig. 2A–D). Over several weeks of training (4 weeks of hand control followed by 3 weeks of brain control), we observed a significant reduction in target acquisition time between early and late sessions (Figs. 2A,B) (p < 0.001 and p < 0.01; KS test for hand control and brain control, respectively), trial duration (p < 0.001 and p < 0.02) and inter-reward interval (p < 0.004 and p < 0.02). The mean trial duration across weeks of training significantly decreased for both hand control (∆t = 0.13 s, 1-way ANOVA: p < 0.05) (Fig. 2C) and brain control training (∆t = 0.7 s, 1-way ANOVA: p < 0.05) (Fig. 2D). To assess the performance of the B2 in manual control and brain control modes, we also calculated performance metrics reported in the previous literature39,40. Table 1 shows the overall mean target acquisition time, target hit rate and the movement precision required, which was quantified as the ratio of target size to the size of the workspace40.

Table 1 Full size table

Figure 2 B2 shared control task. (A,B) Boxplots comparing target acquisition time (left panel), trial duration (centre panel), inter-reward interval (right panel) in hand control (A) and in brain control (B). (C,D) Changes in trial duration with conjoint training across weeks. Mean ± SEM trial duration for each of the four weeks of hand control experiments (C) and three weeks of brain control (D) for dyad M&O. Dashed line shows linear regression fit. Trial duration reduced significantly with improvement in coordination between monkeys. (E,F) Trial-averaged movement profiles from an early session (dashed) and a late session (solid) for monkey M (red) and O (blue) aligned to the time of target onset (time = 0 ms.). Plots are shown up to 1.5 s for hand control (with mean target acquisition time of 0.8 s) and 3 s for brain control (with mean target acquisition time of 2.2 s). The target was located 9 cm from the centre (y-axis). As the dyad trained together over a span of 7 weeks the reaction time lag, the time of movement onset of monkey O relative to monkey M, decreased during both hand control (E) and brain control (F) (G,H) Neural activity PETH from an early session and a late session for hand control (G) and brain control (H) Each row is a single neuron, colour denotes normalized firing amplitude. Monkey M and O neurons marked by black (M) and grey (O) vertical bands. Notice that neuronal modulations are more intense and synchronized across the monkey pair in the later session as compared to the early session. (I) Changes in reaction time lag between monkey M and O across experiments. Lag was derived from peak of cross-correlation of two monkeys’ behavioural traces (E,F) Trends fit with linear regression. As the dyad trained together the reaction time lag reduced to 10 ± 27 ms. (J) Lag in neural activity between two monkeys over the same experiments as in (E) computed again by finding the peak of cross correlogram on each session. As the dyad trained together the reaction time lag reduced to 4.2 ± 5 ms. Full size image

As monkeys M and O learned the shared control task, their coordination improved in both hand control and brain control modes (Fig. 2E–J). In the initial session, with monkey M initiating movements earlier than monkey O (Fig. 2E,I; hand control), the lag between the monkeys’ reaction times was 200 ± 12 ms (mean ± SEM). With conjoint training, reaction time lag between the two monkeys steadily decreased to 10 ± 27 ms over 21 sessions (13 sessions of hand control and 8 sessions of brain control, Fig. 2I). This gradual reduction in lag was statistically significant (hand control: linear regression: r2 = 0.57; p < 0.01, brain control: linear regression: r2 = 0.53; P < 0.05). As the monkey dyad became more coordinated, the cortical neural activity lag between the two brains decreased and stabilized near zero (Fig. 2G,H,J). More specifically, cortical activity lag decreased during hand control (linear regression: r2 = 0.41; P < 0.05) and remained approximately zero (4.2 ± 5.1 ms) throughout brain control (purple data in Fig. 2J, linear regression: R2 = 0.01; P > 0.05).

For the other dyad (M&C), the behavioural responses were synchronous and lag remained close to zero (3 ± 35 ms; linear regression: P = 0.44) throughout all sessions (8 sessions) (Fig. S1, Movie S2). Dyad M&C performed brain control after passive observations of the avatar arm movements and without a requirement to coordinate their behaviours in the manual task. This training sequence probably aided better coordination between the subjects. Table 1 shows the mean target acquisition time, target hit rate and the movement precision required for the dyad M&C.

Concurrent changes in neural activity during Brainet operation could be related to multiple factors: common visual feedback of the avatar arm movements, common representation of the target (spatial location and temporal onset) and common BMI outputs (movements to the same location). To understand if the interplay between the common sensory input and common motor behaviour could account for the genesis of concurrent neural activity we carried out further analyses. First, we investigated the relationship between the behavioural outcomes in each trial and concurrent neuronal modulations between the monkeys in the same trials. We found correct task outcomes were associated with greater neural correlations between monkeys, both in hand control (Fig. 3A; P < 0.01; unpaired t-test) and brain control (Fig. 3B; P < 0.01). To clarify the role of concurrent neuronal activity further, we conducted an analysis where a k-NN classifier predicted the trial outcome, based on the neural correlations between monkeys, both before (700 ms before target appearance until target appearance) and after target onset (from 100 to 800 ms after target appearance). The neuronal correlations were significant predictors of trial outcome (P < 0.05, 1-proportion z-test) not only after target appearance but also before (Fig. 3C,D). Notably, the effect of neural synchrony on the trial outcome was stronger in brain control mode (Fig. 3D). We suggest that the presence of neural correlation prior to target onset indicated that the monkeys attended to the avatar arm since the common visual input from the avatar arm was the only source of neural correlation during that period. This interpretation in terms of an increased attention agrees well with the better performance on the trials with high neural correlation: the more the monkeys attend to the avatar arm and the other screen events, the better they perform.

Figure 3 (A-B) Cross correlogram of Monkey M vs. Monkey O neural activity on correct trials (upper left panel) and error trials (lower left panel) for both a representative hand control (A) and brain control session (B) Mean cross correlogram for the correct and error trial group shown in panel on right of each panel. The correlation is higher in correct trials indicative of increased synchronous neuronal modulations between monkeys. (C,D) k-NN prediction of trial outcome (reward or error) using the mean neuronal cross correlogram on a single trial (i.e. left panels in A-B) between the two monkeys either prior to target appearance (“Pre-target”) or after target appearance (“Post-target”). Chance level prediction (95% confidence interval) shown in yellow band. * denotes P < 0.01, unpaired t-test. Neuronal synchrony between monkeys before and after target onset was predictive of trial outcome. (E,F) Extra correlation analysis. Velocity profile in a trial for dyad M&O were cross-correlated. The average cross correlation (trial specific correlation, red trace) was estimated for all trials: hand control (E) and brain control (F) Extra correlation was the excess correlation (in the red trace) that cannot be accounted for by the distribution of across-trial correlation (enclosed by the grey dashed lines). The vertical line shows the time lag at peak correlation. (G-H) Partitioned control task: (G) Panels show the average trajectory to the 4 target locations, denoted by colour, from the first (left column) and the last (right column) session. Trajectories derived from the monkeys’ controlled axes predictions (controlled) are shown in the upper panels and non-controlled axes predictions (non-controlled) are shown in the lower panels. Non-controlled traces were initially more similar to the controlled traces (left panels) but over time became shrunk and convoluted (right panels). (H) The fraction of rewarded trials was computed based on trials when the avatar arm is moving according to the predictions of the controlled axis (blue) or the non-controlled axes (red). Mean fraction correct trials shown for each of the three weeks of experiments. Shown separately is fraction correct amongst all trials (upper panel) and only among trials where a reward was achieved (lower panel). The percentage of trials in which the complementary trace reached the target (red bars) decreased significantly over training. Full size image

Next, we examined whether the increase in concurrent neuronal activity between the two monkeys could emerge due to the fact that the two animals move simultaneously, without necessarily paying attention to each other, following visual stimulus onset. If this scenario were true, there would not be any trial specific neuronal correlation beyond the level that could be observed if the behavioural trials of one of the monkeys were randomly shuffled (provided that target locations still matched for two monkeys). To investigate this possibility, we used movement velocities as a proxy to test the correlation between monkeys. To determine the additional component to the coordination between the monkeys, over and above the coordination that comes about as a result of jointly performing reaching movements to a specific target location, we calculated a parameter called “extra-correlation”. Extra-correlation is the correlation between the two monkeys’ movement velocities in a trial that cannot be explained by the across-trial correlation. Across-trial correlation was computed by correlating velocity of movement of a monkey from one trial with that of the partner monkey on another randomly chosen trial in which the monkeys moved to the same target location. This was performed for several such combinations of trials (40-50 combinations) to obtain a distribution of the across-trial correlation between the velocities of the dyad. We observed that extra correlation in movement velocities, which was significant in 3/13 sessions during hand control (Fig. 3E), was however, significant during all brain control sessions for both dyad MO (Fig. 3F) and MC (Fig. S2). These results suggest that increases in concurrent activity is not only due to motor outputs triggered by the same visual stimulus, but also due to a significant trial specific coordination in movement velocities, especially during brain control.

Experiment 2: Partitioned control

In the next experiment, a B2 (dyad M&C) performed a partitioned control task where each monkey performed a subtask of a 2D movement (Fig. 1C). One monkey controlled only the X position of the avatar arm (X-monkey) and the other monkey controlled only the Y position (Y-monkey). The targets appeared on the diagonals at 45°, 135°, 225°, or 315°, which meant that both the X-monkey and Y-monkey had to simultaneously and accurately enact movements along their individual axes to achieve a correct trial. Each session started with a period of passive observations of the avatar arm movements. Based on the data recorded during this period, one UKF was fit per monkey. During brain control, the UKF outputs (X for monkey M and Y for monkey C) were combined. The unused outputs (Y for monkey M and X for monkey C) were computed, but not shown to the monkeys as any kind of feedback. We compared the used and unused X and Y outputs in an offline analysis to estimate cortical adaptation to preferentially represent the coordinate being controlled through the BMI.

During brain control mode, despite being given control of only one axis, a monkey could in principle generate motor plans at the neural level that encoded movements in both the X and Y dimensions. Alternatively, each monkey could parse the 2D trajectory into a controlled dimension (which it had to attend to) and non-controlled dimension (which it could disregard). For example, if monkey M controlled movements along X-axis, we asked how this animal’s cortical neurons changed their contribution to movements along Y-axis, which they did not control. That provided us with a way to quantify how much functional plasticity was occurring in the monkeys’ cortex as a result of partitioned Brainet operation. To measure that, we first plotted the actual avatar movement traces generated by the B2 (Fig. 3G; upper panels) which corresponds to the axes the monkeys controlled (X for monkey M and Y for monkey C). We then compared this plot with the traces for the complementary pair of dimensions (Y for monkey M and X for monkey C; Fig. 3G: lower panels). We found that the members of the B2 contributed more along the dimension they controlled and less along the dimension they did not control. This result is evident from the comparison of the actual avatar movement traces generated by the B2 (Fig. 3G; upper panels) to the complementary traces (Fig. 3G; lower panels). These traces were similar in the early sessions (Fig. 3G; left panels) but diverged over time (Fig. 3G, right panels). By the 3rd week of training, average complementary traces were shrunk and convoluted compared to the average actual traces. As a result of this cortical adaptation, the complementary trajectories seldom reached the target. The percentage of trials in which the complementary trace reached the target (red bars in Fig. 3H) decreased from 26% to 12% over training (regression, r2 = 0.42 upper panel of Fig. 3H, P < 0.05). Similar reduction was observed when only the rewarded trials were included (r2 = 0.38, P < 0.05, lower panel Fig. 3H).

Experiment 3: Brainet triad control of an avatar arm in 3-D space

In the third experiment, we built a 3D super-task out of individual 2D BMI tasks. Each of the three monkeys forming a B3 viewed the 2D projection of a spherical target in 3D space, from an X-Y, Y-Z, or X-Z reference frame (Fig. 1D). The monkeys moved the avatar arms on their displays in their 2D space to the projected target location. As they reached the projected target, a 3D avatar arm that received inputs from all the three UKF decoders reached the true target location in 3D space (Movie S3). This design reduced task difficulty for any one individual monkey by breaking down the 3D task into simpler 2D subtasks.

Each session started with a period of passive observation during which each monkey watched the computer-generated avatar arm perform the corresponding subtask. A UKF model was fit per monkey and its subtask and the individual monkey outputs (X, Y for monkey M; Y, Z for monkey C; and X, Z for monkey K) were combined to produce 3D movement of the avatar arm. Put another way, during brain control, pairs of monkeys shared equal control over one of the axes: Monkeys M and K controlled the X-axis, monkeys M and C controlled the Y-axis and monkeys C and K controlled the Z-axis.

Clear behavioural improvements occurred during the brain control epochs over a span of three weeks of training. We observed a significant reduction in target acquisition time (p < 0.02; KS test), trial duration (p < 0.03; KS test) and inter-reward interval (p < 0.002; KS test) between early and late sessions (Fig. 4A). The mean trial duration significantly decreased (4.25 to 3.65 s, 1-way ANOVA: p < 0.05, Fig. 4B panel on the left) over the span of three weeks and the mean reward rate increased from 6 to 10 trials per minute (Fig. 4B, panel on the right). Across 11 sessions, the B3 significantly improved its performance from 20% correct trials in the first session to 78% correct trials in the last session (P < 0.01, 1-way ANOVA with bootstrapping, green+black bars in Fig. 4D).

Figure 4 Triad (B3) control of 3D movements. (A) Boxplots comparing target acquisition time (left panel), trial duration (centre panel), inter-reward interval (right panel). (B) Reduction in trial duration (left panel) and concurrent increase in the reward rate (right panel) with conjoint training across weeks. (C) Normalized contribution of each of the three monkeys across a representative subset of 30 trials. The relative contribution of each monkey varied from trial to trial. (D) Fraction of trials that were correctly performed by a dyad (black) as a triad (green), or incorrectly (purple) shifted across the 11 triad experiments. The fraction of total trials with a rewarded outcome in which all three monkeys contributed (green), or those in which two monkeys contributed (black) increased significantly within each week and across sessions whereas the fraction of erroneous/unattempted trials reduced significantly. (E,F) Decoded trajectories and neural data from the triad experiment. (E) Mean X,Y,Z traces produced by individual monkeys shown separately by colour among trials where all monkeys contributed (left column) or when only monkeys M and C contributed (right column). Mean X,Y,Z (the value used to move avatar) shown in black. Distance to target in each axis is 5 screen-cm. When one monkey (monkey K) opted out, the working dyad generated higher-amplitude trajectories (Right column, X axis and Z axis) as opposed to when all the members contributed (left column). (F) PETHs aligned on target onset for same trial subsets as in (E) Rows represent individual neurons and colour indicates normalized firing rate (z-score). Neurons from different monkeys marked by colour along right edge (same colours as in (E)). Increased effort by the working dyad also resulted in stronger cortical modulations between the members (right panel) as compared to when all the members contributed (left panel). Full size image

Table 1 shows the overall mean target acquisition time, target hit rate and the movement precision required for the B3. Due to the addition of the third dimension, movement precision required of B3 was scaled by a factor of 5.5 as compared to B2. However, the target acquisition times increased only by a factor of 1.9 (Dyad MC vs. Triad MCK). Target hit rate remained nearly the same for both paradigms. The increase in target acquisition times indicates the cost of synchronizing an additional monkey into a B3. However, as a result of the design, the B3 performed reach movements in conditions that required greater movement precision, without compromising the hit rate.

Even though the relative contribution of each monkey varied from trial to trial (Fig. 4C), the highest rate of success was attained when all three monkeys contributed (Fig. 4D). Furthermore, the percent of total trials with a rewarded outcome in which all three monkeys contributed to the final movement outcome (green bars in Fig. 3D) grew from 7% to 50%, the largest increase of any possible dyad or triad tested in this study. The increase was significant within each week of training (post-hoc Tukey test: P < 0.05) and across weeks (1-way ANOVA with bootstrapping: P < 0.001).

The B3 design was resilient to any one individual monkey underperforming in a given trial because any monkey dyad could successfully move the avatar arm to the true target location in 3D space. We observed that the number of rewarded trials for the dyad M&C was lower as compared to the triad but improved with sessions as well (from 12% to 27%, P < 0.01, 1-way ANOVA with bootstrapping). This improvement suggests that the B3 could benefit from the redundancy that was built into the design in the form of shared control. However, when one monkey opted out, the working dyad would need to produce higher-amplitude trajectories (Fig. 4E, right column) as opposed to when all the members contributed (Fig. 4E, left column). Increased effort by the working dyad also resulted in stronger cortical modulations between the members (Fig. 4F, right panel) as compared to when all the members contributed (Fig. 4F, left panel).

Neuron-dropping curves: decoding improves with number of neurons

During B2 and B3 operations, cortical ensembles in each monkey exhibited clear task-related activity during both passive observation (Fig. 5A,C) and brain control modes (Fig. 5B,D). Analysis of these neural signals confirmed that the accuracy of arm movement decoding would improve when VLSBA was recorded and combined from multiple brains (Fig. 5C). This finding extends our previous results, using neuron dropping curve (NDC) analyses, where we showed that decoding accuracy consistently improved when larger neuronal samples were recorded from a single brain4,5. Here, this analysis has been extended to visualize the effect of different relative quantities of two and three brains performing together as part of a Brainet. NDCs were constructed for the decoding of avatar arm position during passive observations and brain control mode for the B2 (Fig. 5E,F) and B3 (Fig. 5G-I). NDCs were plotted as families of curves, where each curve represented decoding accuracy for a neuronal sample composed of a variable-size ensemble from one monkey (Mk1) and a fraction of the full ensemble from the other monkey (Mk2). The NDCs indicated that decoding accuracy benefited from mixing the contributions from different brains as well as the overall neuronal mass. The best accuracy was typically achieved when all neurons from all monkeys were combined (Fig. 5E-I), with the only exception being the Z-axis prediction in B3 (Fig. 5I). When only a small number of Mk1 cells (fewer than 10) were added to Mk2 ensembles during passive observation, decoding accuracy stayed the same or slightly decreased (Fig. 5E-I). However, by scaling to 102 or 103 neurons, this trend shifted such that additional Mk1 neurons added to Mk2 ensembles yielded more accurate predictions.

Figure 5 Neuronal representations during B2 and B3 shared control experiments. Neuronal modulations from monkey M (orange) and C (green) during a 10 second window of passive observation. (A) and brain control (B) Centre and peripheral target onset times are denoted by grey vertical lines. Cortical ensembles in each monkey exhibited clear task-related activity during both passive observation and brain control modes. (C) Passively observed trajectory (black) compared with predicted trajectory using only Monkey M (green), only Monkey C (red), or both (blue) neuronal ensembles. Grey vertical lines from (A) again denote relevant task events. Accuracy of arm movement decoding improved when VLSBA was recorded and combined from multiple brains. (D) Decoded X and Y trajectories from monkey M (black, solid) and monkey C (red), as well as the average of the two (dotted) during 8 second window of brain control experiment. (E,F) Neuron dropping curves (NDC) showing effect of ensemble size of each of the two monkeys on prediction accuracy during B2 passive observation. The number of neurons used from Mk1 marked by x-axis. The percent of Mk2 population used for predictions denoted by colour (see Legend). Accuracy of predictions measured as correlation coefficient r. The decoding accuracy benefited from mixing the contributions from different brains as well as the overall neuronal mass. (G-I) Same NDC analysis as (E,F) except for prediction of X, Y and Z position during B3 passive observation. Prediction of X (G), Y (H) and Z (I) shown separately. Again (as seen in E-F), the decoding accuracy benefited from mixing the contributions from different brains as well as the overall neuronal mass. Full size image

The NDC analysis indicated that Brainet decoding accuracy improved with an increase in the number of neurons. However, a larger total number of neurons was not the sole factor that contributed to the improvement, particularly, the improvement in performance over time. The total number of recorded neurons per monkey did not change considerably throughout the study (Monkey M: 410-501; Monkey C: 196-214; Monkey O: 156-229; Monkey K: 140-156). Yet, we observed a significant improvement in performance (Fig. 2A,B, Fig. 4A,B), as well as improvement in coordinated behaviour (Fig. 2I) and occurrence of concurrent neuronal activity (Fig. 2J). Coordinated behavior was observed even during hand control when the recorded neurons did not influence behavior.

Additionally, improvements in Brainet performance were related to conjoint behavior of the monkeys (as in team work) rather than improvements in individual skills or individual decoding accuracy. All monkeys involved in these experiments had been previously overtrained on the single-BMI task, which makes individual improvements in BMI control an unlikely explanation for the overall improvement of the Brainet performance. Furthermore, we assessed the decoding accuracy in individual monkeys across sessions of manual control (Fig. 6A) and passive observations (Fig. 6B–D). The decoding accuracy showed some fluctuations, which is typical as the motivation levels vary between days of the week. However, there was no systematic trend overall, which contrasts to a steady improvement in Brainet performance over time. These observations suggest that coordinated behavior rather than individual improvements the major factor, which led to better performance.