Toward an Autonomous Brain Machine Interface: Integrating Sensorimotor Reward Modulation and Reinforcement Learning

Surgery.

Three bonnet macaques (Macaca radiata) were implanted bilaterally in the primary motor cortex with chronic 96-channel platinum microelectrode arrays (10 × 10 array separated by ∼400 μm, 1.5 mm electrode length, 400 kOhm impedance, ICS-96 connectors, Blackrock Microsystems). The implantation of large numbers of electrodes was made possible because of our previously described techniques (Chhatbar et al., 2010). We briefly summarize here. All surgical procedures were conducted in compliance with guidelines set forth by the National Institutes of Health Guide for the Care and Use of Laboratory Animals and were further approved by the State University of New York Downstate Institutional Animal Care and Use Committee. All surgical procedures were performed under general anesthesia, and aseptic conditions were maintained throughout. Anesthesia and animal preparation were performed directly or were supervised by members of the State University of New York Division of Comparative Medicine veterinary staff. Ketamine was used to induce anesthesia; isofluorane and fentanyl were used in maintenance of anesthesia. Dexamethasone was used to prevent inflammation during the procedures, and diuretics, such as mannitol and furosemide, were available to further reduce cerebral swelling if needed. All subjects were observed hourly for the first 12 hours after implantation and were provided with a course of antibiotics (baytril and bicilin) and analgesics (buprenorphine and rimadyl) commensurate with the recommendations of the Division of Comparative Medicine veterinary staff.

Three to 6 months before this electrode array implantation, an initial implantation with a footed titanium headpost (Crist Instrument) to allow head fixation during training was conducted, again under aseptic conditions. Head restraint of the animal is required for our experiments to ensure minimization of movement artifacts on our neural recording system, as well as to track the movement of the eyes. Implantation was performed following training to a sustained performance level of at least 90% correctly completed trials per training session for manual tasks.

Extracellular unit recordings.

Subjects were allowed to recover for 2–3 weeks after implantation surgery, after which recordings of single-unit, multiunit, and local field potential (LFP) activity were recorded while the subjects performed the below tasks. Recordings were made using multichannel acquisition processor systems (Plexon). Neural signals were amplified, bandpass filtered (170 Hz to 8 kHz for single and multiunit activity and 0.7–300 Hz for LFPs), sampled at 40 kHz for single-unit/multiunit activity activity and 2 KHz for LFP, thresholded and single-unit/multiunits were sorted based on their waveforms using principal component (PC)-based methods in Sort-Client software (Plexon). For data analysis, we used bilateral M1 units from Monkeys A and C, and ipsilateral M1 units from Monkey Z for our manual experiments. In our observation experiments, we analyzed data from the contralateral M1 (with respect to the right arm) of Monkey A and the ipsilateral M1 (with respect to the right arm) of Monkey Z. For the purposes of this work, we did not specifically segregate between single units and multiunits. An average of 180 M1 units were recorded per session/d.

Electromyography.

Surface gold disc electrodes (Grass Technologies) were sewn onto elastic bands and placed on the skin overlying muscle groups. EMG was recorded from the following muscle groups: latissimus dorsi, biceps, deltoid, triceps, forearm extensors, and forearm flexors. EMG signals were acquired through the Plexon system at a sampling rate of 2 kHz.

Experimental setup and behavioral training.

Three Bonnet macaques (M. radiata, 2 females, 1 male) were trained to perform a center-out reaching task while their right arm rested inside the Kinarm exoskeletal robotic manipulandum (BKIN Technologies). There were two main types of experiments: manual and observational tasks. Visual feedback of the current hand position was provided by a hand feedback cursor on the monitor that precisely colocated with the tip of the monkey's middle finger during manual trials. The manual task (see Fig. 1a), consisted of right hand movements from a center target to a peripheral target located 5 cm to the right of a start position. The target radius was 0.8 cm. Trials were initiated by entering the center target (with neutral color, see Fig. 1, green) and holding for 325 ms (center hold). The center hold was followed by the color cue period (100–300 ms depending on the animal's temperament). The color-cued peripheral target was displayed, and the color of the center target changed from the neutral color to the same color as the cued peripheral target, informing the monkey whether the trial would be rewarding or nonrewarding. The monkey was required to maintain its hold on the color-cued center for 325–400 ms, again depending on the animal's temperament. The implicit GO cue was when the center target disappeared after the 300 ms color cue period, at which time the monkey could move to the peripheral target, where it had to hold for 325 ms before receiving a liquid reward or no reward. A liquid reward was provided only after a successful reach for a rewarding trial. If the monkey failed to complete a trial correctly, the same trial was repeated, giving incentive to perform nonrewarding trials correctly the first time. Trials types were randomized otherwise.

For the first observational task (observational task 1, OT1), the rewarding and nonrewarding trials were color-coded as in the manual task; however, in the observational tasks, the hand feedback cursor would automatically move on its own to the peripheral target to the right of the center target while the KINARM was locked into place, so that the monkey could not make active reaching movements (see Fig. 2a). The left arm was restrained with a padded cuff for all experiment types. In this task, the cursor moved at a constant speed toward the peripheral target. The monkeys knew when the KINARM was locked and did not attempt to make reaching movements during these trials. For the second observational task (OT2), the color of the targets was maintained to be the same for rewarding and nonrewarding trials (i.e., there was no color cueing in OT2). The cue of reward versus no reward was the cursor moving toward or away from a peripheral target, respectively. In this version of the task, there were two possible peripheral targets, as can be seen in Figure 3a.

During every trial, eye-tracking was conducted using an IR-sensitive camera. A trial was aborted if the monkey failed to look at the projection screen where the tasks was displayed during either the color cue period (manual task and OT1) or during the initial feedback cursor movement (OT2).

Data analysis.

Multivariate linear regression was performed on the neural firing rates (100 ms bins) to fit and predict shoulder and elbow positions acquired during the manual task (for fits and predictions, see Table 1). The mean of each 100 ms of position data was fit by 10 bins of neural data, corresponding to 1 s of causal information (Francis and Chapin, 2006) (for fits and predictions, see Table 1). Multivariate linear regression was also performed on the neural data (100 ms bins) to fit and predict EMGs of the right latissimus dorsi and right biceps brachii acquired during the manual task and OT1 (for fits and predictions, see Table 2).

The total number of units acquired in each task per brain area ranged from 31 to 214. For the manual task: Monkey A had 172 contralateral M1 units and 126 ipsilateral M1 units (total of 298 units). Monkey C had 89 contralateral M1 units and 100 ipsilateral M1 units (total 189 units). Monkey Z had 52 ipsilateral M1 units. Hence, the total number of M1 units was 539. For observation task 1: Monkey A had 214 contralateral M1 units and Monkey Z had 51 ipsilateral M1 units. For observation task 2: Monkey A had 54 contralateral M1 units and Monkey Z had 51 ipsilateral M1 units. The amount of units available slowly decreased over time after implantation.

For the manual task, we pruned the data in the following manner to be sure that the differences between rewarding and nonrewarding trials were not due to differences in kinematics. Nonrewarded trials were pruned, so that only trials with maximum velocity, path length, and time to reward within one SD of the average rewarding trials were selected. All the trials whose maximum velocity peak occurred at or after 1200 ms (qualitatively/visually selected) after the initiation of the trial were eliminated to remove trials with delayed reach time. Trials with double peaks in the velocity profile were also removed. Only neural data from pruned trials were selected for analysis. The separability between rewarding and nonrewarding trials was evident without pruning the data (data not shown). However, the data were pruned to show that the separability was not purely due to kinematic differences between the trials.

We analyzed the neural data comparing rewarding to nonrewarding trials. The observational tasks lasted longer compared with the manual tasks due to our choice of the cursor speed. In the manual task, we considered data (binned at 100 ms) starting 200 ms before the color cue and ending 1500 ms (includes 300 ms after average reach time to the target) after the presentation of the color cue for each trial, whereas in the observation task, we considered data (binned at 100 ms) starting 200 ms before the color cue and ending 2700 ms (includes 300 ms after reach time to the target) after the presentation of the color cue for OT1, and movement onset for OT2 for each trial. There was no statistical difference (two-sample t test, p < 0.05) between the neural data considered 200 ms before the color cue (manual and OT1), or the start of movement (OT2) for rewarding versus nonrewarding trials.

The square root transform was performed on all units' binned data to bring the values closer to a Gaussian distribution (Yu et al., 2009). Reward-modulated units (units with a significantly different firing rate between rewarding and nonrewarding trials for a state in the trial: two-sample t test, p < 0.05) were further separated based on whether their average firing rate was higher for rewarding or nonrewarding trials. Units from these groups were selected as samples for the figures.

For Figure 5a–c, the following states were considered: before the color cue (200 ms), color cue present (300 ms), movement of cursor (task-dependent), reward period (300 ms), after reward period (300 ms). For each state, a two-sample t test (p < 0.05) was used to determine whether a significant difference existed between rewarding and nonrewarding trials.

Classifiers.

PC analysis was performed on all the units' data from 100 ms bins during the times taken from the pruned data for the manual files and all the completed trials in the observational tasks. The neural data were z-scored before running the princomp function in MATLAB (MathWorks). As an example, if we had 100 units for a given experiment and there were 100 good trials that were 2000 ms long, we would then have a data matrix that was 2000 rows (nins) by 100 (units). This data matrix was then passed into the MATLAB function princomp. PC score values were separated into rewarding and (pruned, for the manual task) nonrewarding trials. Support vector machines, logistic regression, and linear discriminant analysis (linear classify function of MATLAB) were tested to obtain the best prediction of rewarding vs nonrewarding trials by using PC scores as inputs. The best results were obtained from the linear classify function in MATLAB. After training, the function classifies each bin of sampled data PC scores into one of the training groups: rewarding or nonrewarding. The smallest number of PC scores that gave the best prediction values was selected for use in the algorithm, which was 10. For each monkey, we applied leave-one-out cross validation on its data to quantify the classification performance.

LFP analysis.

LFP signals were collected from 32 channels of the M1 array through the plexon system at a sampling rate of 2k (filtered between 0.7 Hz and 300 Hz). Event-related time frequency (ERTF) analysis was performed using methods adapted from Canolty et al. (2012). Briefly, the average of the 32 channel LFP signal was filtered at frequencies centered from 1 to 128 Hz (in log space) using Gabor filters with a proportional filter bandwidth of 0.25 (bandwidth = 0.25 × the center frequency). For each frequency channel, the power of the filtered signal was realigned to the start of each trial and then averaged across trials. The trial-averaged power at each frequency was then normalized with respect to the value for the 500-0 ms pretrial period and expressed in dB unit.

Reinforcement learning (RL)-based BMI.

Briefly, the theory of RL indicates that an agent, such as an animal, or in our case the RL-BMI system, should act in a manner that leads to the most rewards while interacting with its environment. The term “environment” in our case includes the neural activation patterns from M1, and the type of reinforcement learning architecture we are considering here is termed actor critic, where the actor is the motor BMI and the critic is the evaluative feedback. The logic used by the actor to perform an action given a state, neural firing pattern, is called the policy. An action performed by the actor under a given policy leads it to a new state in the environment, and the consequence of such an action is used as feedback to modify its behavior/policy, which is learning.

Temporal difference learning is a branch of reinforcement learning that allows moment-to-moment updating given a simple evaluative feedback signal, such as the one we are deriving from our classifier. Specifically we used Q learning. The state-action value, Q^π(s, a), is the expected return starting from state “s” given that the RL agent executes the action “a” in state “s” under a policy π (Sutton and Barto, 1998). Specifically, we used an ε-greedy policy as the actor and the Q learning paradigm augmented with Eligibility Trace Q(λ), as the actor's update rule. An eligibility trace is extremely useful in dealing with the credit assignment problem (Sutton and Barto, 1998). The action with the highest Q value is selected 1 − ε percentage of the time (exploitation), whereas a random action is performed ε percentage of the time (exploration) under the ε-greedy policy. There are also ways to change ε given the systems performance, but such full descriptions are outside the scope of this work.

In Q learning, the TD error equation is as follows: Where r = {− 1,1} is the immediate reward, γ = the discount rate and its allowable range is [0,1], (s, a) = the previous state and the action performed in state s under an ε-greedy policy π, respectively, and (s', a') = the current state and a' is the ε-greedy action in state s', respectively.

The TD error is used as feedback to update the estimates of the state-action values (Q values) as follows: In our architecture, r is the class label predicted by a reward classifier (critic) whose input is the M1 neural activity. Specifically, when population firing is classified as rewarding, r is set to 1, whereas when the neural activity is classified as nonrewarding, r is set to −0.1. As such, a classifier outputs a binary evaluative measure by decoding the neural signal, which critiques the executed action. The architecture suggested here conforms to a broader definition of the actor-critic architecture as it has a separate memory structure to explicitly represent the policy independent of the entity providing the evaluative signal. The scalar evaluative signal is the sole output of the critic and drives all learning in the actor. The suggested architecture can easily be modified to conform to the stricter definition of actor-critic wherein the critic represents the estimated value function and the evaluative feedback provided by the critic is used to update itself along with the actor (Sutton and Barto, 1998). One can also envision versions where the user gives feedback on the critic's performance as a perfect source of feedback to update the critic and subsequently the actor when necessary.

Simulations.

One of our future biomedical engineering goals is to use neural activity from M1 simultaneously to control movement, via a motor BMI, as well as to update this BMI via an evaluative feedback signal, also derived from M1. One architecture that is well suited for this type of updating would be a BMI that works via reinforcement learning (DiGiovanna et al., 2009; Sanchez et al., 2011; Tarigoppula et al., 2012; Mahmoudi et al., 2013; Pohlmeyer et al., 2014), as this would only necessitate an evaluative signal, such as rewarding or nonrewarding, rather than a full error signal, such as the difference on a moment-to-moment basis between the desired movement trajectory and the actual one made by the BMI. This later full error signal is what most BMIs to date use. This decrease in the amount of information necessary for updating the BMI system makes it more plausible that such a system could autonomously correct itself in real-world changing environments. One can easily imagine combining the best of the supervised learning world with the best of the reinforcement-learning world, but we leave these architectures as future work at this point (see discussion below on methods of the RL agent).

Here we simply demonstrate that the amount of evaluative feedback obtained in our experiments from M1 would be sufficient for a reinforcement learning BMI in theory to work. Toward this goal, we used a simulation of the motor cortical output that we have previously used (Tarigoppula et al., 2012) in testing RL-BMI systems. We have previously used M1 information for BMI purposes (Chhatbar and Francis, 2013), including RL-based systems (Sanchez et al., 2011), and thus know that the movement-related activation is present, a well-known fact at present. We therefore needed to test whether the evaluative feedback from M1 would be sufficient as well. The RL system we used for our proof of concept was a simple one-step system. It used only the neural output from M1 at the start of the trial during target presentation. From that neural activation pattern, it decided which target to approach. This type of one step system has been shown previously in real-time (Pohlmeyer et al., 2012). Our rationale for using the simulation rather than the actual monkeys for this RL-BMI proof of concept is that the monkeys used for these reward-based experiments had lost their chronic recordings to a large extent by the time this work would have been conducted.

We used the classification rates obtained in our OT2 for this simulation, as it is reasonable to expect no difference, from the animal's point of view, between OT2 and an online one step RL-BMI. Our group has previously described the use of a RL paradigm in which an RL agent performed a 4 target-8 action center out reaching task by decoding the firing rate of a simulated M1 neuronal ensemble (Tarigoppula et al., 2012). We used this same neural model in our current work, and thus only briefly describe it here. In this M1 model, a group of neurons was simulated using the Izhikevich model neuron. The neural ensemble consisted of 80 neurons; 60% of the neurons had unimodal tuning curves, 15% were had bimodal tuning curves, and 25% had assigned asymmetric tuning curves, as observed by Amirikian et al. (2000). A tuning curve directed the degree of neural modulation, given the direction of the target with respect to the cursor's position. Preferred directions of these neurons were assigned randomly. A spike was detected every time the membrane potential of a neuron surpassed 30 mV (Tarigoppula et al., 2012). The task was identical to OT2 in spatial arrangement and cursor motion; however, the cursor was controlled by the RL system.

The target direction in a given trial changed each neuron's firing rate with respect to its baseline activity based on their respective tuning curves. That is, given a target in the left direction, the neurons that had their preferred direction to the left fired at their maximum firing rate, whereas the remaining neurons modulated their firing based on their tuning curve. Using the output of the simulated neural ensemble as the input to an artificial neural network, the Q value for each potential action was determined. Specifically, a multilayer perceptron (MLP) with a single hidden layer consisting of 120 units was used to calculate the Q value, given an input from the neural ensemble. Here, 99% of the time, the action with the highest Q value was executed (the “greedy” part of the ε-greedy policy). The other 1% of the time, a random action was taken (the exploratory rate, the “ε” part of the ε-greedy policy). The exploratory rate, the percentage of steps in which an action is executed randomly regardless of its optimality, was set at 1% (“ε” part of ε-greedy policy). The random exploration allows for discovery of new solutions by the RL agent, useful especially in an altering environment. Update of the weights of the MLP was performed by back-propagation (Bryson and Ho, 1969) of a qualitative error signal “TD error × eligibility trace,” calculated using the immediate reward it received based on the correct or incorrect action performed. A correct action resulted in evaluative feedback to the RL agent of either 1 (rewarding) or −0.1 (nonrewarding) with a probability determined by the success rate of our M1 classifier for OT2, equivalent to 70% correct feedback. This means that, 70% of the time in our simulation, the RL agent was given the correct evaluative feedback of rewarding or nonrewarding and 30% of the time it was given false information.