# 17 Schaeffer et al .pdf

Nom original:

**17_Schaeffer et al.pdf**Titre:

**Switching Markov decoders for asynchronous trajectory reconstruction from ECoG signals in monkeys for BCI applications**Auteur:

**Marie-Caroline Schaeffer**

Ce document au format PDF 1.7 a été généré par Elsevier / Acrobat Distiller 8.0.0 (Windows), et a été envoyé sur fichier-pdf.fr le 16/07/2018 à 14:10, depuis l'adresse IP 132.168.x.x.
La présente page de téléchargement du fichier a été vue 685 fois.

Taille du document: 1.1 Mo (13 pages).

Confidentialité: fichier public

### Aperçu du document

Journal of Physiology - Paris 110 (2016) 348–360

Contents lists available at ScienceDirect

Journal of Physiology - Paris

journal homepage: www.elsevier.com/locate/jphysparis

Original Research Paper

Switching Markov decoders for asynchronous trajectory reconstruction

from ECoG signals in monkeys for BCI applications

Marie-Caroline Schaeffer, Tetiana Aksenova ⇑

Univ. Grenoble Alpes, CEA, LETI, CLINATEC, MINATEC Campus, 38000 Grenoble, France

a r t i c l e

i n f o

Article history:

Received 30 June 2016

Received in revised form 7 December 2016

Accepted 6 March 2017

Available online 10 March 2017

Keywords:

Asynchronous brain-computer-interface

Hidden Markov model

Mixture of experts

Piece-wise linear models

Switching Kalman filter

a b s t r a c t

Brain-Computer Interfaces (BCIs) are systems which translate brain neural activity into commands for

external devices. BCI users generally alternate between No-Control (NC) and Intentional Control (IC) periods. NC/IC discrimination is crucial for clinical BCIs, particularly when they provide neural control over

complex effectors such as exoskeletons. Numerous BCI decoders focus on the estimation of

continuously-valued limb trajectories from neural signals. The integration of NC support into continuous

decoders is investigated in the present article. Most discrete/continuous BCI hybrid decoders rely on

static state models which don’t exploit the dynamic of NC/IC state succession. A hybrid decoder, referred

to as Markov Switching Linear Model (MSLM), is proposed in the present article. The MSLM assumes that

the NC/IC state sequence is generated by a first-order Markov chain, and performs dynamic NC/IC

state detection. Linear continuous movement models are probabilistically combined using the NC and

IC state posterior probabilities yielded by the state decoder. The proposed decoder is evaluated

for the task of asynchronous wrist position decoding from high dimensional space-time-frequency

ElectroCorticoGraphic (ECoG) features in monkeys. The MSLM is compared with another dynamic hybrid

decoder proposed in the literature, namely a Switching Kalman Filter (SKF). A comparison is additionally

drawn with a Wiener filter decoder which infers NC states by thresholding trajectory estimates. The

MSLM decoder is found to outperform both the SKF and the thresholded Wiener filter decoder in terms

of False Positive Ratio and NC/IC state detection error. It additionally surpasses the SKF with respect to

the Pearson Correlation Coefficient and Root Mean Squared Error between true and estimated continuous

trajectories.

Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction

Brain-Computer Interfaces (BCI) are systems that allow severely

motor-impaired patients to use their brain activity to control

external devices, and thereby to interact with their environment

(Leuthardt et al., 2006). Among BCI systems, movement-control

BCIs aim at providing users with control over upper- or lowerlimb orthoses or prostheses (Mak and Wolpaw, 2009). Several

steps are usually necessary to translate the user’s neuronal activity

into commands for effector control. The first stage consists in

acquiring an observation of the user’s brain activity (Mak and

Wolpaw, 2009). Microelectrode arrays (MEA), which are directly

implanted in the cortex, permit the observation of neuronal Single

or Multi-Unit Activity (SUA/MUA) and of Local Field Potentials in

localized areas. The high spatial resolution of MEAs comes at the

price of high invasiveness (Mak and Wolpaw, 2009). Electroen⇑ Corresponding author.

E-mail address: tetiana.aksenova@cea.fr (T. Aksenova).

http://dx.doi.org/10.1016/j.jphysparis.2017.03.002

0928-4257/Ó 2017 Elsevier Ltd. All rights reserved.

cephalography (EEG) is a non-invasive acquisition method, but it

exhibits a poor spatial resolution along with a high sensitivity to

artifacts and noise (Leuthardt et al., 2006). Finally, minimally invasive electrocorticographic (ECoG) arrays measure neural signals

from the cortex surface. They are implanted under or over the dura

mater (Schalk and Leuthardt, 2011). Several pre-clinical (Chao

et al., 2010; Shimoda et al., 2012) and clinical studies (Schalk

et al., 2008; W. Wang et al., 2013) have confirmed the potential

of ECoG for long-term stable and accurate BCI systems.

Movement-control ECoG-based BCIs are considered in the present

paper.

Brain signal acquisition is followed by the extraction of features

specifically related to the user’s intentions (Mak and Wolpaw,

2009). ECoG recordings are generally described by potentially

high-dimensional combinations of time, frequency and space characteristics (Bundy et al., 2016; Chao et al., 2010; Eliseyev and

Aksenova, 2014). A decoder is then applied to estimate the user’s

intention from the brain features. After being optionally enhanced

by post-processing methods (Bashashati et al., 2007a), the estimate

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

of the user’s intention is converted into commands sent to control

the effector (Leuthardt et al., 2006).

MEA-based control of a limb prosthesis or orthosis generally

relies on the direct extraction of trajectories’ kinematic parameters

from the neural features (Li, 2014). By contrast, EEG-based control

of limb prostheses or orthoses is typically based on the discrimination between several motor imageries (e.g., Baxter et al., 2013;

Bhattacharyya et al., 2015; Hortal et al., 2015). Best strategy for

ECoG decoding is still unclear (Chin et al., 2007; Ashmore et al.,

2012). The reported trajectory reconstructions obtained with

ECoG-based models of limbs’ position or velocity (Pistohl et al.,

2008; Bundy et al., 2016) suggest that kinematic decoding may

be applicable to ECoG-based BCIs. Continuously-valued kinematic

parameters are generally estimated from the neural features by

means of regression techniques. Although several teams have

investigated the use of nonlinear models (e.g., Eliseyev and

Aksenova, 2014; Li et al., 2009), estimation of limbs’ position or

velocity is generally performed using static or dynamic linear models. Wiener filtering has been utilized for both SUA/MUA (Velliste

et al., 2008; Wodlinger et al., 2015) and ECoG static decoding (W.

Wang et al., 2013; Williams et al., 2013). Similarly, Kalman filters

and variants have been applied for dynamic decoding of SUA/

MUA (Wu et al., 2002; Hochberg et al., 2012; Aggarwal et al.,

2013; Velliste et al., 2014) and ECoG signals (Pistohl et al., 2008;

P.T. Wang et al., 2013).

One challenge of BCI clinical applications is the ability to provide users with asynchronous control over the effector

(Leuthardt et al., 2006). Unlike synchronous BCIs which are periodically controllable by users, asynchronous (or self-paced) BCI decoders are continuously available (Graimann et al., 2010). In

asynchronous settings, users generally alternate between Intentional Control (IC) and No-Control (NC) periods (Mason et al.,

2006). Kinematic decoders usually output disturbing non-neutral

values when applied to NC states (Chao et al., 2010; Shimoda

et al., 2012; Velliste et al., 2014). NC support is thus highly desirable (Leeb et al., 2007).

NC support has initially been studied for brain-switches, i.e. BCI

systems which integrate NC detection into discrete decoders (classifiers) (Bashashati et al., 2007b; Pfurtscheller et al., 2010; SolisEscalante et al., 2010). Reported methods for NC/IC state discrimination included Linear Discriminant Analysis (LDA) (Pfurtscheller

et al., 2010), Support Vector Machine (Solis-Escalante et al.,

2010), logistic regression (logit) (Blokland et al., 2014), etc.

Specific decoding strategies are required to integrate NC support into continuous decoders. Several approaches have been proposed to ensure that the outputs of continuous kinematic decoders

are neutral during NC states. A first strategy consists in postprocessing the estimated continuous trajectories, e.g. by applying

dynamic link functions (Wang et al., 2011). Another approach

relies on hybrid discrete/continuous decoders. Such hybrid decoders exploit a state detector to bring out NC and IC periods, and

a continuous movement model is applied when appropriate. State

and movement decoders are mostly combined using a WinnerTakes-All strategy (Velliste et al., 2014; Flamary and

Rakotomamonjy, 2012; P.T. Wang et al., 2013). The movement

model which corresponds to the most likely class at a given time

is selected and applied on the neural features. The use of hybrid

decoders for the integration of NC support in continuous decoding

approaches was reported for both SUA/MUA (Velliste et al., 2014;

Wood et al., 2005; Aggarwal et al., 2013; Suway et al., 2013) and

ECoG signals (e.g., Bundy et al., 2016; Flamary and

Rakotomamonjy, 2012; P.T. Wang et al., 2013). State decoders were

used to combine static or dynamic linear IC movement models.

Wiener filters were coupled with a LDA-based class detector in

SUA/MUA signals (Suway et al., 2013), and were mixed with logistic regression (Bundy et al., 2016) or linear regression of state

349

labels on neural signals (Flamary and Rakotomamonjy, 2012) in

ECoG signals. Kalman filters and variants were cascaded with

LDA in SUA/MUA signals (Velliste et al., 2014; Aggarwal et al.,

2013) or with a Bayes classifier in ECoG signals (P.T. Wang et al.,

2013).

Most hybrid decoders relied on the assumption that successive

states are temporally independent (Velliste et al., 2014; Flamary

and Rakotomamonjy, 2012; Bundy et al., 2016). The interest of

making more realistic hypotheses about temporal dependencies

in the state sequence was investigated for classification-based

BCIs. Dynamic classifiers such as Hidden Markov Models (HMM)

(Kemere et al., 2008; Hotson et al., 2016) or Conditional Random

Fields (CRF) (Hasan and Gan, 2011; Saa et al., 2016) were specifically used for the estimation of NC/IC states, possibly with several

NC- or IC-related sub-states, in SUA/MUA (Kemere et al., 2008),

EEG (Hasan and Gan, 2011) and ECoG signals (Saa et al., 2016).

When compared with generic static state decoders, they are

expected to be less liable to output spurious detections (e.g., Saa

et al., 2016).

By contrast, few hybrid decoders made use of a sequential

assumption for NC/IC state estimation, e.g. the Switching Particle

Filter (SPF) (Wood et al., 2005) and the Switching Kalman Filter

(SKF) (Srinivasan et al., 2007). The SKF-based hybrid decoder proposed in Srinivasan et al. (2007) assumed that the NC/IC state

sequence was generated by a first order Markov chain. SKFs extend

Kalman Filters by conditioning the emission and transition models

on a Markovian hidden switching state variable (Murphy, 1998).

The SKF (Srinivasan et al., 2007) was tested for EEG-based asynchronous wheelchair control from simulated data. Validation on

real data was not carried out. Although latent classes were not

associated to NC and IC periods, the use of a SKF was additionally

reported for the task of 2D wrist position decoding from SUA/MUA

signals in monkeys (Wu et al., 2004). Similarly, the SPF relied on a

Markov assumption to model NC/IC state succession (Wood et al.,

2005). The state posterior probabilities issued by a LDA-based

dynamic state decoder were used to weight a Particle Filter. The

SPF-based decoder was applied on monkeys’ SUA/MUA signals

for offline real movement decoding. Despite its decoding performance, the high computational cost of the SPF (Kotecha and

Djuric´, 2003) limits its applicability for real-time decoding.

The present article addresses the issue of NC support in continuous ECoG decoders, and specifically investigates the integration of

dynamic state detection. Previously reported hybrid decoders with

sequential state assumptions (Srinivasan et al., 2007; Wood et al.,

2005) were tested on MUA/SUA signals or simulated data. Significant differences between MEA and ECoG feature representations

limit the portability of SUA/MUA decoding methods to ECoG.

Strategies for the integration of sequential state decoding into

ECoG-based continuous kinematic decoders are explored in the

present article. In particular, a Markov Switching Linear Model

(MSLM) is proposed as an ECoG hybrid decoder with dynamic NC

support. The sequence of NC and IC states is assumed to be generated by a first-order Markov chain. A linear model between neural

features and trajectory is conditioned on the current state. The trajectory decoder applied during IC states is exclusively trained on IC

samples. Contrary to the traditional Winner-Takes-All approach,

the MSLM relies on a probabilistic combination rule to integrate

state discrimination into continuous trajectory decoding. As the

MSLM is a piecewise linear model, it can be seen as an approach

to introduce non-linearities into brain signal decoders.

The MSLM decoding performance was evaluated for asynchronous wrist position reconstruction from both epidural and

subdural ECoG data in monkeys (2 monkeys for each data set).

The MSLM was compared to the SKF. Similarly to the MSLM, the

SKF continuous model is conditioned on a first-order hidden Markov chain. Its relevance for ECoG data asynchronous decoding had

350

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

not yet been investigated. The decoding performance of a Wienerbased decoder which trivially supports NC periods by means of trajectory thresholding was additionally analyzed for comparison

purposes. The MSLM outperforms both the SKF and thresholded

Wiener Filter (WF) approaches in terms of state detection error

and of False Positive Ratios. It additionally surpasses the SKF with

respect to the Pearson Correlation Coefficient and Root Mean

Squared Error between true and estimated continuous trajectories.

2. Methods

Notations: matrices are denoted by bold upper-case romans

(e.g., X), vectors by bold lower-case romans (e.g., x) and scalar by

lower-cases italics (e.g., x). xij and xi refer to the ði; jÞth entry of

matrix X and to the ith element of vector x, respectively.

The MSLM relies on the state posterior probabilities yielded by a

dynamic state decoder to combine continuous linear statedependent models. Combination strategy is built on the Mixture

of Experts (ME) framework (Waterhouse, 1997) extended to time

series.

2.1. Mixture of experts framework

First introduced by Jacobs et al. (1991), MEs combine several

functions (‘‘experts”) to model the dependence between input

and output variables (Jacobs et al., 1991; Waterhouse, 1997).

Experts are weighted according to the input variable (Bishop,

2006). Optimal weights are computed for each sample according

to a given criterion.

Let the vector xt be the explanatory (input) variable,

t

x 2 X Rm , and the vector yt be the continuous response (output)

variable, yt 2 Y Rn . Samples are indexed by t 2 N. MEs assume

that the input space X is partitioned into K (possibly overlapped)

S

regions: X ¼ Kk¼1 X k , and that a different sub-process generates

the output vector from the input vector in each region

(Waterhouse, 1997). Mixtures of Experts model the dependence

between xt and yt by means of a set of K local functions

K

ff k g1 ; f k : X ! Y, referred to as ‘‘experts”. Expert k is relevant for

the region X k :

yt ¼

K

X

dzt ;k f k ðxt Þ þ t

k¼1

where ðzt Þt2N is the expert sequence, zt ¼ k if yt has been generated

by expert k. dzt ;k ¼ 1 if zt ¼ k and dzt ;k ¼ 0 otherwise. t is the observation noise.

^t ¼ Eðyt jxt Þ (Bishop, 2006) of the response

The Bayes estimate y

variable yt is computed via the decomposition of the conditional

expectation (Waterhouse, 1997):

Eðyt jxt Þ ¼

¼

XK

k¼1

Eðyt ; zt ¼ kjxt Þ ¼

k¼1

^tk :

g tk y

XK

XK

k¼1

Pðzt ¼ kjxt ÞEðyt jxt ; zt ¼ kÞ

^tk ¼ Eðyt jxt ; zt ¼ kÞ is the estimate issued by expert k and Eð:Þ

Here, y

indicates the expected value.

The input-dependent mixing coefficients g tk ¼ Pðzt ¼ kjxt Þ com^ tk . Mixing coefficients can be

bine (‘‘gate”) the experts’ estimations y

interpreted as the conditional probability of expert f k having generated the output value yt , given the value of xt (Waterhouse,

P

1997). They satisfy the constraints g tk 2 ½0; 1 ; k g tk ¼ 1 (Bishop,

2006). The structure which computes the mixing coefficients is

often referred to as ‘‘gating network”.

Linear or nonlinear experts and gating networks can be integrated into the ME flexible structure (Yuksel et al., 2012). ME-

based approaches were used in several BCI studies which didn’t

support NC periods, e.g. for multi-model movement decoding from

monkeys’ SUA/MUA signals (Kim et al., 2003).

Several extensions have been specifically proposed to adapt

MEs for sequential data modeling (Yuksel et al., 2012). Recurrent

gating networks (Cacciatore and Nowlan, 1994; Meila and Jordan,

1996) were integrated into ME and for example applied for movement modeling (Meila and Jordan, 1996) and efficient control of

switching systems (Cacciatore and Nowlan, 1994). Particularly,

Markov Mixtures of Experts (MME) (Meila and Jordan, 1996)

model a temporal sequence of input-output data by dynamic

switching between local static models. The expert sequence

ðzt Þt2N is assumed to be generated by a first-order Markov chain.

The inputs and outputs corrupted by noise are measured, but

states are hidden and must be estimated from measurements.

2.2. Hidden Markov model-based gating network

The MSLM proposed herein is a variant of MME which relies on

Hidden Markov Models (HMMs) for dynamic gating of the experts’

estimates. It relies on the underlying hypothesis that neural signals

are generated by cognitive states. Thus, the distribution of the

extracted neural features xt is conditioned on unobserved discrete

states zt ; Pðxt jzt Þ. By contrast, the previously reported MMEs

(Bengio and Frasconi, 1996; Meila and Jordan, 1996) assumed that

the input variable xt conditions the probability of switching from

one state to another, Pðztþ1 jzt ; xt Þ.

Hidden Markov Models (HMMs) are a powerful tool for the

modeling of stochastic time-varying processes (Rabiner and

Juang, 1986). Two discrete-time stochastic processes are involved

in discrete-time HMMs. The first process is an unobservable (‘‘hidden”) dynamic sequence ðzt Þt2N . By contrast, the second sequence

ðxt Þt2N is observable, and is referred to as the ‘‘observation

sequence”. The MSLM state decoder models the state ðzt Þt2N and

feature sequences ðxt Þt2N by a HMM, where zt 2 N is the hidden

variable and xt 2 Rm is the observed variable.

The HMM hidden variable zt can take K distinct discrete values

(‘‘states”), zt 2 Z with Z ¼ fz1 ; z2 ; . . . ; zK g Z. In the typical case of a

first-order HMM, the value of ztþ1 depends exclusively on the current value zt : Pðztþ1 ¼ kjz1:t Þ ¼ Pðztþ1 ¼ kjzt Þ; t 2 N. The hidden

state sequence is then fully characterized by the transition matrix

A ¼ ðaij Þ; aij ¼ Pðztþ1 ¼ jjzt ¼ iÞ, and by the probability p associated

to the value of the first hidden state: pi ¼ Pðz1 ¼ zi Þ. The probability distribution of the observed variable xt 2 Rm is exclusively conditioned on the current hidden state zt : Pðxt jz1:t Þ ¼ Pðxt jzt Þ; t 2 N,

where the conditional probability Pðxt jzt Þ is referred to as ‘‘state

emission density”. Various parametric (Gaussian, Student) and

non-parametric state emission distributions have been integrated

into HMMs (Rabiner, 1989).

Let B gather the parameters necessary to characterize the distributions Pðxt jzt ¼ jÞ; j ¼ 1; . . . ; K. Then a HMM is fully described by

the parameter set H g ¼ fA; B; pg. HMM training typically relies

on Maximum-Likelihood (ML) estimation (Ghahramani, 2001).

While state labels are hidden during HMM application, they are

either hidden or available during training. When the state labels

T

ðzt Þt¼1 are observed in the training data set (complete training data

set), HMM training is supervised with respect to the state

sequence. While the application of specific algorithms (e.g. the

Expectation-Maximization algorithm) is often required for unsupervised HMM training, direct ML estimation can be performed

when a complete training data set is available. Once the HMM

has been trained, an inference algorithm is applied to estimate

T

the hidden state sequence ðzt Þt¼1 associated to a new observation

sequence

T

ðxt Þt¼1 .

T

If the full sequence ðxt Þt¼1 has been observed,

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

i.e. if past, current and future observations are known, estimation

of states zt ; t ¼ 1; . . . ; T is generally performed using the Viterbi

algorithm or the forward-backward algorithm. Both Viterbi and

forward-backward algorithms cannot be applied for online processing, because they rely on a back-tracking procedure from the

last estimated state ^zT ¼ argmaxk Pðzt ¼ kjx1:T Þ (Rabiner, 1989). By

contrast, the probability Pðzt ¼ kjx1:t Þ yielded by the forward algorithm (Rabiner, 1989) is only conditioned on past and current

observations and can be applied online.

As online estimation of users’ intentions is required for BCI systems, the MSLM input-dependent mixing coefficients are exclusively conditioned on past and current neural observations:

g tk ¼ Pðzt ¼ kjx1:t Þ. Mixing coefficients are first decomposed using

Bayes’ rule:

Pðzt ¼ kjx1:t Þ ¼

Pðzt ¼ k; x1:t Þ

:

Pðx1:t Þ

Both the denominator Pðx1:t Þ and numerator Pðzt ¼ k; x1:t Þ are

then computed using the forward algorithm.

2.3. Linear experts

The MSLM decoder applies linear experts to model the dependence between the explanatory and response variables xt 2 Rm

and yt 2 Rn :

yt ¼

XK

d t ðBk xt

k¼1 z ;k

þ tk Þ k ¼ 1; . . . ; K

351

subset fXk ; Yk g, where Xk and Yk gather training samples observed

at times t such that zt ¼ k. Training of the HMM-based gating network is performed on the data set fX; zg. Transition and emission

probabilities are separately estimated (Dietterich, 2009). The transition matrix A is estimated by counting transition frequencies in

the sequence z (Dietterich, 2009).

3. Application

The MSLM decoder was tested for a task of asynchronous 3D

trajectory decoding from ECoG signals. It was compared to the

SKF. The SKF relies on the assumption that its observation and

transition models are conditioned on a Markovian hidden state

(Murphy, 1998). Trajectory estimation is based on a probabilistic

combination of state-specific models (Murphy, 1998). The SKF

can therefore be considered as an analogue of the MSLM in the

frame of dynamic continuous models. SKF state estimation is based

on the consistency of the state-specific dynamic models with the

observed neural signals. Thus, unlike the MSLM, it doesn’t directly

exploit differences between the distributions of neural features

generated during each state. A Wiener filter with postprocessing-based NC support was additionally implemented for

comparison purposes. In contrast with the MSLM and the SKF, state

values are not taken into account when training the WF continuous

decoder. As state detection is performed at the post-processing

stage, its accuracy depends on the quality of trajectory estimates.

where Bk 2 Rn m and k 2 Rn is the observation noise.

3.1. Data sets

2.4. MSLM training

Model identification and test were carried out on publicly available subdural and epidural ECoG datasets acquired and distributed

by the Laboratory for Adaptive Intelligence, RIKEN Brain Science

Institut, Saitama, Japan (http://neurotycho.org/).

The cortical activity of 4 Non-Human Primates (Monkeys A, K, B

and C) was recorded by chronic ECoG arrays during a food reaching

task. Monkeys A and K were respectively implanted with a 32- and

64-channel subdural ECoG array. Both Monkeys B and C were

epidurally implanted with a 64-channel ECoG array. Fig. 1 shows

the position of the implanted ECoG arrays. The monkeys were sitting in front of an experimenter during the data acquisition sessions. The experimenter presented food to the monkeys at

random intervals. The monkeys performed reaching movements

to grab the food and to bring it to their mouth. ECoG signals were

acquired at a sampling rate of 1 kHz. A motion tracking system

tracked subjects’ wrist coordinates at a sampling frequency of

Let H ¼ He ; H g gather the experts’ and HMM-based gating

network’s parameters He ¼ fBk gKk¼1 and H g ¼ fA; B; pg. Whereas

state labels are hidden in test data sets, the case where they are

observed in the training data set is here considered. Fully supervised training of the MSLM relies on a complete training data set

fX; Y; zg, where X 2 RT m ; Y 2 RT n , and z 2 NT K gather the

T

T

T

observed sequences ðxt Þt¼1 ; ðyt Þt¼1 and ðzt Þt¼1 respectively. The

Maximum-Likelihood (ML) criterion is routinely used for both

EM (Yuksel et al., 2012) and HMM training (Ghahramani, 2001).

When a complete training data set fX; Y; zg is available, computa^ ML ðX; Y; zÞ ¼ argmax ln PðX; Y; zjHÞ is

tion of the ML estimate H

H

equivalent to separate ML estimations of experts and gate parameters (see Appendix A). Training of expert k is based on the training

Fig. 1. ECoG array implantation for Monkeys A, K, B and C (reproduced from Chao et al., 2010; Shimoda et al., 2012). Solid gray circles indicate reference electrodes. Ps:

principal sulcus; As: arcuate sulcus; Cs: central sulcus; IPs: intraparietal sulcus. A: subdural data set, Monkeys A (a) and K (b). B: epidural data set, Monkeys B (a) and C (b).

352

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

120 Hz. Body-centered 3D wrist trajectories downsampled at

20 Hz were extracted from the motion tracking system outputs.

Full descriptions of the experimental set-ups can be found in

Chao et al. (2010) and Shimoda et al. (2012).

One session was discarded from the subdural data set because of

very short idle states. Four epidural sessions presented multi-second

long null velocity periods at essentially different altitude levels. As

the experimental description didn’t permit to label these periods

as NC or IC states with certainty, these sessions were also excluded

from analysis. The final data set was thus composed of 7 subdural

sessions and 16 epidural sessions. The data set was analyzed session

per session (23 sessions, 4 monkeys). Each recording was split into a

training and a test set: the first 70% of each session were used to

identify decoder parameters, and the last 30% were exclusively used

to assess the decoders’ generalization performance.

3.2. Feature extraction

3.2.1. Neural signal

Time-frequency features were extracted for each channel from

Dt-long ECoG sliding epochs (Dt ¼ 1 s, sliding step: 100 ms) following (Eliseyev and Aksenova, 2014). A Complex Continuous

Wavelet Transform (CCWT) was applied (Morlet wavelet)

(Eliseyev et al., 2012). The frequency content was analyzed

between 1 and 250 Hz. A redundant sampling of this frequency

domain was achieved via 38 daughter wavelets chosen with a

logarithmic scale. A logarithmic transform was applied on the

CCWT-transformed signals’ absolute value. Average values were

computed in 100 ms sliding windows (100 ms step), resulting in

a 10-point description of ECoG 1s-long time epochs for each

frequency band and each channel. Thus an ECoG epoch was

described by the temporal-frequency-spatial feature vector

xt 2 Rm ; m ¼ 32 10 38 and m ¼ 64 10 38 for Monkey A

and Monkeys K, B and C respectively.

Trajectories. Wrist positions yt 2 R3 were issued by the motion

tracking system. Velocity was derived from position using a

forward-difference approximation (Eberly, 2014). NC and IC states

were manually labeled.

3.3. Decoder implementation

3.3.1. MSLM

A MSLM decoder with K ¼ 2 was trained using the complete

training data set fX; Y; zg, where zt ¼ 0 for NC samples and zt ¼ 1

for IC samples.

The IC expert was dedicated to wrist position decoding during

IC periods. The NC expert yielded the NC neutral position value

NC , estimated as the average value of yt computed over NC states.

y

A Gaussian noise was associated to the IC expert. A few NC transition samples were added to the IC samples for training. Partial

Least Squares (PLS) (Höskuldsson, 1988) regression was chosen

as an approximate ML estimator for the IC expert because of the

high dimension of input variable xt , which makes Ordinary Least

Squares (OLS) solutions unstable. PLS regression is a training

method able to handle high dimensional and/or correlated

explanatory variables, especially when the explanatory variable’s

dimension is higher than the number of training samples. PLS

regression has been performed on high-dimensional data sets in

several BCI studies (Bundy et al., 2016; Chao et al., 2010;

Eliseyev and Aksenova, 2014). It is based on the projection of xt

and yt onto low-dimensional subspaces. The subspaces are found

jointly on the criterion that they maximize the covariance between

the respective projections of xt and yt . The optimal subspace

dimension was here found by 6-fold cross-validation on the training data set.

The movement models were gated by a hybrid logit-HMM

model. Conventional HMMs directly model the emission distributions Pðxt jzt ¼ jÞ; j; . . . ; K. Variants based on discriminative modeling were proposed by several teams (e.g., using Neural Networks

(Renals et al., 1994) or Support Vector Machines (Valstar and

Pantic, 2007)). A discriminative approach consists in modeling

the class conditional probabilities Pðzt ¼ jjxt Þ; j; . . . ; K (e.g., Valstar

and Pantic, 2007). By contrast, ‘‘generative” classifiers infer

equiprobable decision boundaries from the probabilities

Pðxt jzt ¼ jÞ; j; . . . ; K. To clarify the relevance of both approaches

for the considered NC/IC classification task, pilot tests were completed on a reduced data set (7 acquisition sessions). The classification accuracies of one generative (Gaussian distribution within the

NC and IC periods, resulting in a Bayes classifier) and one discriminant classifier (logistic regression) (Bishop, 2006) were compared.

The logit model was found to slightly outperform the Bayes classifier (average improvement of 6%). It was consequently integrated

into the HMM-based gating network for dynamic state detection.

Following (Valstar and Pantic, 2007), class priors Pðzt ¼ jÞ and

Bayes’ rule were combined to compute HMM emission probabilit

t

¼jjx Þ

. Dimensionality reduction was carried

ties Pðxt jzt ¼ jÞ / Pðz

Pðzt ¼jÞ

out before feeding neural features to the logit class decoder. 6fold cross-validation on the training data set was performed to

choose the optimal dimension of a latent subspace found by PLS

regression between xt and zt . ML estimates of the logit model

parameters were computed using the Iteratively Reweighted Least

Squares algorithm (Bishop, 2006).

3.3.2. SKF

A SKF (see Appendix B) was implemented for dynamic combination of K ¼ 2 Kalman filters, one specialized in NC periods and

the other in IC periods. The SKF state variable (i.e., response variable) ytSKF was chosen as the monkey’s wrist position and velocity,

as it was reported optimal for ECoG decoding (Pistohl et al., 2008).

Similarly to the MSLM state decoder, the dimension of the neural

features xt was reduced before application of the SKF. Dimensionality reduction permitted to limit the SKF computational cost. PLS

regression between xt and ytSKF was used to identify an informative

low-dimensional subspace. Subspace dimension was chosen by 6fold cross-validation on the training data set. OLS estimates of the

transition Ak and emission Ck matrices were computed on the

training data sets fXk ; Yk g. ML estimates of the variance matrices

Ck and Rk were found following (Aggarwal et al., 2013).

3.3.3. Thresholded Wiener filter

A Wiener filter was fit via PLS regression on the training data set

fX; Yg. The optimal number of PLS factors was estimated by 6-fold

cross-validation on the training data set. Trajectory thresholding

was used to integrate basic NC support into the decoder. Joint

^ti ; i ¼ 1; 2; 3 was performed using probit regresthresholding of y

n

o

^ z .

sion. It was identified on the training data set Y;

3.4. Performance indicators

Both state and trajectory decoding accuracy were assessed to

evaluate the hybrid decoders’ performance.

3.4.1. Performance indicators for state decoding

First, the positive or negative delay with which transitions were

predicted by decoders was estimated. A second indicator set characterized classification error (i.e., number of misclassified samples). The last indicator set focused on the number of false

activations/deactivations and on their durations.

353

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

jTN TPj

TN

: 36 18% and 14 10% respectively for the subdural and

FPþFN

epidural data sets), the classification error ERR ¼ TPþTNþFPþFN

(Mason et al., 2006) was additionally computed.

False activations and deactivations. The TPR, FPR and ERR are

sample-based indicators. They don’t take into account the dynamic

of misclassified samples. For example, several consecutive misclassified samples are liable to be less disturbing to users that the same

number of isolated misclassified samples. State decoding performance was thus considered in terms of false activations/deactiva

tions’ durations and occurrence number, where a ‘‘false activations” (respectively, a ‘‘false deactivation”) refers a block of consecutive NC samples misclassified as IC samples (respectively, a block

of IC samples mistaken for NC samples). A block of misclassified

samples was counted as one false activation or deactivation. Its

duration was computed as the number of samples divided by the

decision rate.

3.4.2. Performance indicators for continuous variable decoding

Accuracy of the estimates of continuous variables is typically

assessed via Pearson Correlation Coefficient (PCC) and the Normalized Root-Mean-Squared Error (NRMSE) (Spuler et al., 2015). The

PCC measures the amount of linear dependence between observed

^Þ

y

^ variables: PCCðy; y

^Þ ¼ corv ðy;

^

y and predicted y

r , where cov ðy; yÞ indiy

^

y

^ and ry refers to the standard

cates the covariance between y and y

1000

State transition temporal error (ms)

Transition detection The temporal resolution of transition detection was estimated by computing the average absolute value of the

temporal error (positive or negative delay) between true and estimated state transitions.

Confusion matrix-based performance indicators. The indicators

used to assess the performance of binary decoders (Mason et al.,

2006) mainly rely on the confusion matrix. The confusion matrix

gathers the number of NC states which are correctly (True Negatives, TN) or wrongly (False Positives, FP) labeled, and the number

of IC states which are correctly (True Positives, TP) or wrongly

(False Negatives, FN) labeled by the decoder. The True Positive Rate

TP

FP

TPR ¼ TPþFN

(sensitivity) and the False Positive Rate FPR ¼ FPþTN

are

widely used to assess the performance of asynchronous NC/IC

decoding. Both the TPR and FPR were applied to assess the state

decoding performance. Because the NC and IC classes were relatively well balanced in the considered data sets (average ratio

800

600

400

200

MSLM

WF

Subdural data set

ð:Þ0 denotes the transpose:

qﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ

PT

P

t 2

1 Tt¼1 yt .

kyk2 ¼

t¼1 ðy Þ and y ¼ T

where

^

y k2

NRMSE ¼ kky

k ,

y y

2

where

4. Results

The MSLM decoding performance was compared to SKF and WF

respective performances on the test data set of the 23 acquisition

sessions (2 and 2 monkeys, 7 and 16 sessions for the subdural

and epidural data sets, respectively). Significance of performance

differences was established using the Friedman test (Demšar,

2006) with the significance level a ¼ 0:05. Where appropriate, post

hoc pairwise decoder comparisons were performed using the Wilcoxon signed-rank test with Bonferroni correction.

4.1. State decoding performance

Fig. 2 shows the average temporal error with which a state transition was detected on the subdural and epidural data sets. Table 1

shows the respective performances of the MSLM, SKF and WF

decoders in terms of NC/IC samples classification accuracy. The

corresponding p-values are gathered in Table 2.

Epidural data set

Fig. 2. Temporal errors (ms) between observed and estimated state transitions for

the MSLM, WF and SKF decoders (average over 7 and 16 sessions for the subdural

and epidural data sets respectively). The state transition temporal error was

computed as the average absolute value of the positive or negative delay observed

between each true and estimated state transitions.

Table 1

State classification performance. Average False Positive Rate (FPR), True Positive Rate

(TPR) and Error (ERR) are displayed for the MSLM, WF and SKF decoders on subdural

and epidural ECoG data sets.

FPR (%)

TPR (%)

ERR (%)

Subdural (7 sessions)

MSLM

WF

SKF

4.0

6.0

20.2

96.5

94.5

84.6

3.9

6.0

20.5

Epidural (16 sessions)

MSLM

WF

SKF

9.7

17.6

63.0

84.5

83.3

65.3

11.9

16.3

50.5

Table 2

P-values for state classification performance. The significance of the differences

between the decoders’ respective performances was assessed using the Friedman test

with the significance level a ¼ 0:05. Post-hoc comparisons were performed using the

Wilcoxon signed-rank test with Bonferroni correction, i.e. aBonferroni ¼ 0:0167. Nonsignificant differences are indicated by ‘‘n.s.”.

2

deviation of y. The NRMSE measures the l -error between the vec0

^ ¼ ðy

^1 ; y

^2 ; . . . ; y

^T Þ0 ,

tors of observations y ¼ ðy1 ; y2 ; . . . ; yT Þ and y

SKF

Subdural (7 sessions)

Epidural (16 sessions)

MSLM/WF

MSLM/SKF

MSLM/WF

MSLM/SKF

FPR (%)

TPR (%)

ERR (%)

n.s.

0.0156

0.0023

<0.001

n.s.

0.0156

n.s.

0.0061

n.s.

0.0156

<0.001

<0.001

The state transition temporal resolution of the 3 decoders

results in misclassified samples located around true state transitions. In order to dissociate the decoder temporal resolution and

misclassification effects, samples located around transitions were

not considered when computing the confusion matrix. The number of discarded samples was taken equal to the maximal average

decoder transition delay plus its standard deviation, i.e. 6

(= 600 ms) and 10 (= 1 s) samples for the subdural and epidural

data sets respectively. It additionally permitted to take into

account possible manual state labeling inaccuracies in transition

surroundings. The MSLM state decoder outperformed the WF

state decoder for the FPR and ERR (average improvements of

25.0% and 36.1% respectively for the subdural data set, significant

improvements of 37.2% and 26.3% respectively for the epidural

data set). When compared to the SKF, the MSLM significantly

improved the FPR, TPR and ERR for both the subdural and epidural data sets. Table 3 shows false activations/deactivations in

354

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

Table 3

False activation/deactivation average frequency and duration for MSLM, WF and SKF

decoders on subdural and epidural ECoG data sets. False activations (respectively,

false deactivations) are blocks of consecutive NC samples misclassified as IC samples

(respectively, a block of IC samples mistaken for NC samples).

Subdural

Epidural

False activations

False deactivations

Number per

min

Duration

(ms)

Number per

min

Duration

(ms)

MSLM

WF

SKF

0.7

2.1

3.2

520

290

650

0.4

0.9

1.2

540

390

670

MSLM

WF

SKF

0.7

2.9

2.1

800

380

1490

0.9

2.0

1.3

870

410

650

terms of occurrence number and duration. The MSLM resulted in

significantly fewer false activations/deactivations than the WF

and SKF.

4.2. IC decoding performance

The relevance of the continuous IC kinematic models embedded

in the hybrid decoders was then assessed. As state labels are traditionally hidden during decoder application, the decoding accuracy

measured over IC states reflects both the quality of the IC continuous decoder and of the state decoder (e.g., a false negative results

in a high l2 -error). To decouple the impact of the state and continuous decoders, the MSLM, WF and SKF were applied on complete

test data sets fX; zg, i.e. the state sequence was not hidden during

application. Continuous performance indicators were then computed using exclusively the true IC samples. Fig. 3 shows the IC

continuous models decoding performance. Correlation and NRMSE

were not significantly different for the MSLM and WF. The MSLM

significantly outperformed the SKF in terms of NRMSE on the

epidural data set (average NRMSE improvement of 6.6% and

10.1% for the subdural and epidural data sets, respectively).

4.3. Overall decoding performance

Fig. 4 shows the global decoding performance of MSLM, SKF and

WF using continuous performance indicators. To discard perturbations induced by slight variations of monkeys’ tracked wrist position during NC periods, the target yt was taken to be the neutral

NC during NC periods. MSLM and WF performances were

value y

not significantly different. The MSLM significantly outperformed

the SKF for both correlation (average improvement of 19.6% and

16.2% on the subdural and epidural data sets, respectively), and

NRMSE (average improvement of 28.4% and 23% on the epidural

and subdural data sets, respectively). Examples of decoded trajectories are presented on Fig. 5. Fig. 6 illustrates the average influence of frequency, temporal and spatial features of the IC expert

and the gating network for Monkey B (6 sessions, 3 axes).

5. Discussion

In spite of proofs of concept in laboratory environments

(Collinger et al., 2013; Wodlinger et al., 2015), clinical applications

of movement-control BCIs remain rare (Mak and Wolpaw, 2009).

NC support is one of the crucial problems currently faced by the

BCI community (Graimann et al., 2010). NC support is particularly

important when users physically interact with the BCI effector, as

Epidural data set

3

2

2

PCC

PCC

Subdural data set

3

1

0

1

0

MSLM

WF

MSLM

SKF

*

4

4

WF

NRMSE

NRMSE

*

2

SKF

*

2

0

0

MSLM

WF

SKF

axis y

1

MSLM

axis y

2

WF

SKF

axis y

3

Fig. 3. IC decoding performance of the MSLM, WF and SKF decoders when the state sequence is known during application. The PCC and NRMSE between true and estimated IC

trajectories are averaged over 7 and 16 sessions, for the subdural and epidural data sets respectively. The vertical bar indicates the standard deviation associated to

each decoder. The significance of the differences between the decoders’ respective performances was assessed using the Friedman test with the significance level a ¼ 0:05.

Post-hoc comparisons were performed using the Wilcoxon signed-rank test with Bonferroni correction. Significant differences are indicated by a star (⁄).

355

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

Subdural data set

3

*

*

*

*

PCC

2

PCC

2

1

1

MSLM

WF

4

*

0

SKF

MSLM

WF

*

4

*

NRMSE

0

NRMSE

Epidural data set

3

2

SKF

*

2

0

0

MSLM

WF

SKF

MSLM

axis y

WF

axis y

1

SKF

axis y

2

3

Fig. 4. Overall decoding performance of the WF, MSLM and SKF decoders. The PCC and NRMSE between true and estimated trajectories are averaged over 7 and 16 sessions,

for the subdural and epidural data sets respectively. The vertical bar indicates the standard deviation associated to each decoder. The significance of the differences between

the decoders’ respective performances was assessed using the Friedman test with the significance level a ¼ 0:05. Post-hoc comparisons were performed using the Wilcoxon

signed-rank test with Bonferroni correction. Significant differences are indicated by a star (⁄).

position (mm)

y 1 -axis

y 2 -axis

position (mm)

MSLM

SKF

WF

300

250

200

150

100

0

y 3 -axis

position (mm)

-100

300

200

100

0

30

40

50

60

time (s)

70

30

40

50

60

70

30

time (s)

observed

40

50

60

70

time (s)

estimated

Fig. 5. Example of observed and estimated trajectories for the subdural data set (Monkey K). The projections of the monkey’s wrist trajectory onto the horizontal axes (y1 and

y2 ) and the vertical axis (y3 ) are indicated in solid black lines. Red trajectories represent the estimates yielded by the MSLM, SKF and WF decoders.

356

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

Fig. 6. Example of the average contribution of frequency, temporal and spatial modalities to the MSLM IC expert and to the MSLM state decoder (epidural data set, Monkey B,

average over 6 acquisition sessions). Contributions were assessed as the normalized summation of absolute values of models’ coefficients along each modality. Dashed blue

lines indicate the average weight value plus or minus the standard deviation.

is the case with neurally-driven ortheses. False activations are

likely to be particularly disturbing and stressful to users.

Dynamic data processing was chosen in the present article as a

tool liable to reduce BCI systems’ false activations. A hybrid logitHMM dynamic state decoder was integrated into a continuous trajectory scheme. In contrast with the Winner-takes-all strategy

used for model combination by most hybrid decoders, a probabilistic combination strategy was utilized. This statistically optimal

combination approach results in smoother transitions, which we

expect to be liable to facilitate effector control.

The training of the proposed MSLM decoder was here considered

in the case of complete training data sets. Parameter identification

was thus supervised with respect to both trajectory and state

sequences. Generalization to the case of training unsupervised with

respect to the state sequence is a prospective extension of the present study. It would permit to explore the use of multiple sub-states

for accurate modeling and classification of NC and IC periods, and

the combination of several IC experts to improve movement estimation. The respective relevance of Expectation-Maximization or

Variational Bayes training approaches will be investigated.

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

The proposed MSLM decoder was evaluated on epidural and

subdural ECoG data sets. It was compared to 2 hybrid decoders, a

SKF and a post-processed Wiener filter. Results were generally better when decoders were applied on subdural recordings rather

than on epidural recordings. The MSLM and the SKF presented similar IC trajectory decoding performances if the state sequence was

known during decoder application. On the other hand, in the practical case where state labels are hidden during decoder application,

the SKF’s overall correlation on subdural (respectively, epidural)

data sets was 19.6% (respectively, 16.2%) lower than the MSLM’s

one. The fact that the SKF was significantly outperformed by the

MSLM in terms of state detection (76% average increase of the classification error ERR) may have contributed to this degraded overall

decoding performance. The SKF state estimation is based on the

consistency of the NC- and IC-specific dynamic models with the

observed neural signals. It was here less accurate that MSLM state

detection, which directly exploits differences between the distributions of NC and IC neural features. The SKF tendency to drift during NC states was additionally observed on subdural recording (see

Fig. 5), and may have impacted the overall decoding performance.

These results suggest that the SKF-based NC support may be suboptimal for some applications, e.g. for neural signals with degraded

Signal-to-Noise Ratio or spatial resolution. A supplementary study

is thus required to investigate the performance of the SKF during

closed-loop experiments.

The MSLM IC expert and the WF performed equivalently when

they were applied to IC states (Fig. 3). Thus, excluding NC samples

when training the MSLM IC expert did not improve the model

between neural signals and kinematic parameters during IC periods. When compared to the WF, the MSLM state detection error

was lower by 36% and 26% for subdural and epidural ECoG data,

respectively. False activations and deactivations, i.e. blocks of consecutive false positives or false negatives, were longer (plus 380 ms

and 690 ms on average for the subdural and epidural data sets

respectively) but fewer (on average, minus 2:7 and 5:3 false activations/deactivations per minute for the subdural and epidural data

sets, respectively). We expect that a block of adjacent misclassified

samples is less disturbing to BCI users than a few isolated erroneous state estimates. Improvements came at the price of a slight

increase of the state detection delay (on average, plus 40 ms and

180 ms on the subdural and epidural data sets, respectively). It

was nevertheless essentially lower than the SKF delay (on average,

minus 160 ms and 230 ms on the subdural and epidural data sets,

respectively). This delay may explain why state detection improvement was not reflected in the MSLM and WF overall trajectory

decoding performances, which were not significantly different.

MSLM- and the WF-based state estimation mainly differ in two

aspects, which may both have contributed to MSLM significant

improvement of state detection in comparison to WF. First, the

MSLM makes use of a Markov hypothesis to model NC/IC state succession, while the WF relies on static state detection. Secondly,

state- and kinematic-related variables are independently modeled

by the MSLM decoder. By contrast, the WF infers state values from

the estimated trajectories. Several studies have reported that state

and kinematic parameters are encoded in different frequency

bands (Williams et al., 2013; Bundy et al., 2016). This observation

suggests the interest of building independent state and movement

decoders. The difference in neural feature contributions to the

state and movement decoders is illustrated by the modality

influence study presented in Fig. 6. Supplementary studies are

required to explore this topic.

The present offline study suggests that the proposed MSLM

decoder may be profitably used to introduce NC support into asynchronous kinematic decoders. Further investigations are necessary

to strengthen and possibly extend the present study. In particular,

the wrist trajectories used to assess decoders’ performance were

357

associated to a limited space exploration. The relevance of the

MSLM remains to be explored in the case of complex trajectories.

While piecewise linear static and dynamic decoders were considered in the present article, the use of other nonlinear decoders,

for example Generalized Linear Models (Eliseyev and Aksenova,

2014), Multilayer Perceptrons (Kim et al., 2006) or Support Vector

Machine Regression (Kim et al., 2006; Mehring et al., 2003), has

been proposed for the estimation of kinematic parameters from

neural signals. The respective performances of these nonlinear

decoders have yet to be compared for asynchronous ECoG

decoding.

Future work will consist in testing the proposed decoder during

neural control sessions in human subjects. ECoG-based BCI in

humans is the objective of CLINATEC ‘‘BCI and Tetraplegia” clinical

research protocol which Principal Investigator is Prof. A.-L. Benabid

(https://clinicaltrials.gov/show/NCT02550522).

The

wireless

64-channel ECoG implant WIMAGINEÒ has been specifically

designed for stable and long-term signal acquisition (Mestais

et al., 2015). CLINATEC BCI platform (Eliseyev et al., 2014) also

includes the 4-limb exoskeleton EMY (Morinière et al., 2015) and

the software environments required for real time processing of

the neural signals. As the forward algorithm used by the MSLM

for model gating is computationally efficient (Rabiner, 1989), we

expect the MSLM to be compatible with real-time BCI decoding.

Because of a change of context between open-loop and closedloop BCI settings, open-loop neural patterns may differ from

closed-loop patterns (Leuthardt et al., 2006; Jackson and Fetz,

2011; Jarosiewicz et al., 2013). The MSLM decoding performance

thus remains to be investigated during closed-loop experiments.

Acknowledgements

This work was supported in part by grants from the French

National Research Agency (ANR-Carnot Institute), Fondation

Motrice, Fondation Nanosciences, Fondation de l’Avenir, and Fondation Philanthropique Edmond J. Safra. The authors are grateful

to all members of the CEA-LETI-CLINATEC, and especially to A.

Eliseyev, G. Charvet, C. Mestais and Prof. A.-L. Benabid.

Appendix A. MSLM training

Let us consider a MSLM composed of K linear experts gated by a

HMM-based sequential decoder. Each linear expert i is

parametrized by Bi . Let A ¼ ðaij Þ be the transition matrix

aij ¼ Pðztþ1 ¼ jjzt ¼ iÞ; i; j ¼ 1; . . . ; K and p be the initial state

distribution pi ¼ Pðz1 ¼ zi Þ; i ¼ 1; . . . ; K. Finally, let bi gather the

parameters of the chosen emission distribution Pðxt jzt ¼ jÞ, e.g.

mean and variance in the case of Gaussian emissions.

Let fX; Y; zg be a complete training data set available for

Maximum-Likelihood estimation of the MSLM parameters

n

o

K

H ¼ fBi gKi¼1 ; p; A; fbi gi¼1 . Maximization of data complete loglikelihood Lc ðH; X; Y; zÞ ¼ ln PðX; Y; zjHÞ is equivalent to solving of

separate ML problems.

Lc ðH; X; Y; zÞ ¼ ln PðX; Y; zjHÞ;

ðA:1Þ

where

PðX; Y; zjHÞ ¼ PðX; Yjz; HÞPðzjHÞ

¼ Pðx1:T ; y1:T jz1:T ; HÞPðz1:T jHÞ:

ðA:2Þ

Following a typical HMM hypothesis, observations xt are

assumed to be temporally independent when conditioned on z1:T .

The same hypothesis is made for the dependent variable yt . Thus,

358

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

Pðx1:T ; y1:T jz1:T ; HÞ ¼

T

Y

P x1:T ; y1:T jz1:T ; H :

ðA:3Þ

t¼1

t

t

Because the emission distribution of x is only conditioned on z ,

Pðx1:T ; y1:T jz1:T ; HÞ ¼

T

Y

Pðxt ; yt jzt ; HÞ;

ðA:4Þ

t¼1

with Pðxt ; yt jzt ; HÞ ¼ Pðyt jxt ; zt ; HÞPðxt jzt ; HÞ.

Pðyt jxt ; zt ; HÞ corresponds to expert zt , and only depends on Bzt

and on noise distribution associated to the expert, i.e.

Pðyt jxt ; zt ; HÞ ¼ Pðyt jxt ; zt ; Bzt Þ. Pðxt jzt ; HÞ is the HMM emission

model corresponding to the state zt , and parametrized by bzt . Thus,

Pðxt jzt ; HÞ ¼ Pðxt jzt ; bzt Þ and

Pðx1:T ; y1:T jz1:T ; HÞ ¼

T

Y

Pðyt jxt ; zt ; Bzt ÞPðxt jzt ; bzt Þ:

ðA:5Þ

t¼1

The second term of (A.2) can be rewritten as

Pðz1:T jHÞ ¼ PðzT jz1:T 1 ; HÞPðz1:T 1 jHÞ

ðA:6Þ

Taking into account the first-order Markovian hypothesis and

the chosen parametrization,

Pðz1:T jHÞ ¼ PðzT jzT 1 ; AÞPðz1:T 1 jA; pÞ:

ðA:7Þ

By recurrence,

Pðz1:T jHÞ ¼ Pðz1 jpÞ

T 1

Y

Pðztþ1 jzt ; AÞ ¼ pz1

YT 1

t¼1

t¼1

azt ;ztþ1 :

ðA:8Þ

T 1

X

ln azt ;ztþ1 þ

þ

T

X

ln Pðx jz ; bzt Þ:

ðA:9Þ

t

t¼1

Let us introduce the following notation: dzt ;k ¼ 1 if zt ¼ k and

dzt ;k ¼ 0 otherwise, and xt ði; jÞ ¼ 1 if zt ¼ i and ztþ1 ¼ j, and

xt ði; jÞ ¼ 0 otherwise. Then (A.9) can be rewritten as

Lc ðH;X; Y;zÞ ¼

K

X

dz1 ;i lnðpi Þ þ

þ

T 1 X

K

X

K X

t¼1

i¼1

T X

K

X

i¼1

xt ði; jÞlnðaij Þ

j¼1

dzt ;i lnðPðyt jxt ; Bi ÞÞ þ

T X

K

X

t¼1 i¼1

dzt ;i ln Pðxt jbi Þ: ðA:10Þ

t¼1 i¼1

Switching the sum operators,

Lc ðH;X; Y;zÞ ¼

K

X

dz1 ;i lnðpi Þ þ

i¼1

þ

K X

T

X

K X

T 1

X

K X

i¼1

j¼1

xt ði; jÞlnðaij Þ

t¼1

dzt ;i lnðPðyt jxt ; Bi ÞÞ þ

i¼1 t¼1

K X

T

X

dzt ;i ln Pðxt jbi Þ: ðA:11Þ

i¼1 t¼1

Parameters of interest are thus decoupled in (A.11). Maximization of each term of (A.11) can be performed separately.

Appendix B. Switching Kalman filter

The Switching Kalman Filter is a hybrid decoder which probabilistically combines K Kalman Filters (Murphy, 1998). It permits

to infer a hidden trajectory yt from noisy measurements. The continuous response variable yt is composed by the trajectory coordinates and derivatives (velocity, acceleration, etc.). The SKF applies

to switching state-space models specified as follows:

ztþ1 ¼ Aswitch zt ;

y

tþ1

ðB:1Þ

¼ Azt y þ w ;

t

t

The initial probabilities Pðz1 ¼ kÞ; k ¼ 1 . . . K of the switching variable are gathered in p 2 RK .

In contrast with Kalman filtering, exact estimation of

^ t ¼ Eðyt jx1:t Þ is intractable because of the exponential number of

y

possible state sequences. The SKF therefore relies on approximate

solutions (Murphy, 1998).

Supervised Maximum-Likelihood training is possible when a

complete training data set fX; Y; zg is available. Maximization formula (Murphy, 1998) for EM-based SKF training are applied, taking

into account that the hidden switching state sequence is known.

ML estimates of the matrices Azt ¼k ; Czt ¼k ; Czt ¼k and Rzt ¼k are computed on the training data subsets fXk ; Yk g, where Xk and Yk gather

training samples observed at times t such that zt ¼ k. The switching transition matrix Aswitch is fit on the basis of transition frequen-

References

ln Pðyt jxt ; zt ; Bzt Þ

t¼1

t

The ‘‘transition” Eq. (B.2) describes the dynamic of the hidden

sequence yt 2 Rn , where Azt 2 Rn n is the transition matrix and

wt is a a Gaussian observation noise, Pðwt Þ Nð0; Czt Þ;

Czt 2 Rn n . The ‘‘emission” Eq. (B.3) characterizes the dependence

between the measurements xt and the hidden state value yt .

Czt 2 Rm n is the emission matrix, and v t is a Gaussian observation

noise, Pðv t Þ Nð0; Rzt Þ; Rzt 2 Rm m . Both the transition and emission equations are conditioned on the discrete switching variable

zt . Eq. (B.1), where Aswitch is the switching transition matrix,

describes the dynamic of zt which is assumed to be generated by

a first-order Markov process. The probability of the initial state

y1 is conditioned on z1 : y1 ¼ l0;z1 þ uz1 , with Pðuz1 Þ Nð0; P0;z1 Þ.

T

t¼1

T

X

ðB:3Þ

cies in the sequence ðzt Þt¼1 .

Finally, from A.2, A.5 and A.8,

Lc ðH; X; Y; zÞ ¼ ln pz1 þ

xt ¼ Czt yt þ v t ;

ðB:2Þ

Aggarwal, V., Mollazadeh, M., Davidson, A.G., Schieber, M.H., Thakor, N.V., 2013.

State-based decoding of hand and finger kinematics using neuronal ensemble

and LFP activity during dexterous reach-to-grasp movements. J. Neurophysiol.

109 (12), 3067–3081 <http://www.ncbi.nlm.nih.gov/pubmed/23536714>.

Ashmore, R.C., Endler, B.M., Smalianchuk, I., Degenhart, A.D., Hatsopoulos, N.G.,

Tyler-Kabara, E.C., Batista, A.P., Wang, W., 2012. Stable online control of an

electrocorticographic brain-computer interface using a static decoder. In:

Conference Proceedings: . . . Annual International Conference of the IEEE

Engineering in Medicine and Biology Society. IEEE Engineering in Medicine

and Biology Society. Annual Conference 2012, pp. 1740–1744. <http://

ieeexplore.ieee.org/xpls/abs{_}all.jsp?arnumber=6346285,

http://www.ncbi.

nlm.nih.gov/pubmed/23366246>.

Bashashati, A., Fatourechi, M., Ward, R.K., Birch, G.E., 2007a. A survey of signal

processing algorithms in brain-computer interfaces based on electrical brain

signals. J. Neural Eng. 4 (2), R32–R57 <http://www.ncbi.nlm.nih.gov/pubmed/

17409474>.

Bashashati, A., Ward, R.K., Birch, G.E., 2007b. Towards development of a 3-state selfpaced brain-computer interface. Comput. Intell. Neurosci. 2007, 84386 <http://

www.pubmedcentral.nih.gov/articlerender.fcgi?artid=2234253&tool=

pmcentrez&rendertype=abstract>.

Baxter, B.S., Decker, A., He, B., 2013. Noninvasive control of a robotic arm in multiple

dimensions using scalp electroencephalogram. In: International IEEE/EMBS

Conference on Neural Engineering, NER, pp. 45–47.

Bengio, Y., Frasconi, P., 1996. Input-output HMM’s for sequence processing. IEEE

Trans. Neural Netw. 7 (5), 1231–1249.

Bhattacharyya, S., Basu, D., Konar, A., Tibarewala, D.N., 2015. Interval type-2 fuzzy

logic based multiclass ANFIS algorithm for real-time EEG based movement

control of a robot arm. Robot. Auton. Syst. 68, 104–115. http://dx.doi.org/

10.1016/j.robot.2015.01.007.

Bishop, C.M., 2006. Pattern Recognition and Machine Learning, vol. 4. <http://www.

library.wisc.edu/selectedtocs/bg0137.pdf>.

Blokland, Y., Spyrou, L., Thijssen, D., Eijsvogels, T., Colier, W., Floor-Westerdijk, M.,

Vlek, R., Bruhn, J., Farquhar, J., 2014. Combined EEG-fNIRS decoding of motor

attempt and imagery for brain switch control: an offline study in patients with

tetraplegia. IEEE Trans. Neural Syst. Rehab. Eng. 22 (2), 222–229.

Bundy, D.T., Pahwa, M., Szrama, N., Leuthardt, E.C., 2016. Decoding threedimensional reaching movements using electrocorticographic signals in

humans. J. Neural Eng. 13 (2), 026021 <http://stacks.iop.org/1741-2552/13/i=

2/a=026021?key=crossref.7783bca9521da5800c8816f918b776b6>.

Cacciatore, T.W., Nowlan, S.J., 1994. Mixtures of controllers for jump linear and nonlinear plants. Adv. Neural Inform. Process. Syst. 6 (i), 719–726 <http://citeseer.

ist.psu.edu/cacciatore94mixtures.html>.

Chao, Z.C., Nagasaka, Y., Fujii, N., 2010. Long-term asynchronous decoding of arm

motion using electrocorticographic signals in monkeys. Front. Neuroeng. 3

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

(March),

3

<http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=

2856632&tool=pmcentrez&rendertype=abstract>.

Chin, C.M., Popovic, M.R., Thrasher, A., Cameron, T., Lozano, A., Chen, R., 2007.

Identification of arm movements using correlation of electrocorticographic

spectral components and kinematic recordings. J. Neural Eng. 4, 146–158.

Collinger, J.L., Wodlinger, B., Downey, J.E., Wang, W., Tyler-Kabara, E.C., Weber, D.J.,

McMorland, A.J.C., Velliste, M., Boninger, M.L., Schwartz, A.B., 2013. Highperformance neuroprosthetic control by an individual with tetraplegia. Lancet

381 (9866), 557–564 <http://www.pubmedcentral.nih.gov/articlerender.fcgi?

artid=3641862&tool=pmcentrez&rendertype=abstract>.

Demšar, J., 2006. Statistical comparisons of classifiers over multiple data sets. J.

Mach. Learn. Res. 7, 1–30 <http://dl.acm.org/citation.cfm?id=1248547.

1248548>.

Dietterich, T., 2009. Machine learning for sequential data: a review. Struct. Syntactic

Statist. Pattern Recogn. 2396, 227–246. http://dx.doi.org/10.1007/3-54070659-3_2 <http://www.springerlink.com/index/AV8L8HJL6YC2YA3M.pdf>.

Eberly, D., 2014. Derivative Approximation by Finite Differences Derivatives of

Univariate Functions, 1–7.

Eliseyev, A., Aksenova, T., 2014. Stable and artifact-resistant decoding of 3D hand

trajectories from ECoG signals using the generalized additive model. J. Neural

Eng. 11 (6), 066005 <http://www.ncbi.nlm.nih.gov/pubmed/25341256>.

Eliseyev, A., Mestais, C., Charvet, G., Sauter, F., Abroug, N., Arizumi, N., Cokgungor, S.,

Costecalde, T., Foerster, M., Korczowski, L., Morinière, B., Porcherot, J., Pradal, J.,

Ratel, D., Tarrin, N., Torres, N., Verney, A., Aksenova, T., Benabid, A.-l., 2014.

CLINATECÒ BCI platform based on the ECoG-recording implant WIMAGINEÒ and

the innovative signal-processing: preclinical results. In: 36th Annual

International Conference of the IEEE Engineering in Medicine and Biology

Society, pp. 1222–1225.

Eliseyev, A., Moro, C., Faber, J., Wyss, A., Torres, N., Mestais, C., Benabid, A.L.,

Aksenova, T., 2012. L1-penalized N-way PLS for subset of electrodes selection in

BCI experiments. J. Neural Eng. 9 (4), 045010 <http://www.ncbi.nlm.nih.gov/

pubmed/22832155>.

Flamary, R., Rakotomamonjy, A., 2012. Decoding finger movements from ECoG

signals using switching linear models. Front. Neurosci. 6 (March), 29 <http://

www.pubmedcentral.nih.gov/articlerender.fcgi?artid=3294271&tool=

pmcentrez&rendertype=abstract>.

Ghahramani, Z., 2001. An Introduction to Hidden Markov Models and Bayesian

Networks.

Graimann, B., Allison, B., Pfurtscheller, G., 2010. Brain-Computer Interfaces: A

Gentle Introduction, pp. 1–28. <http://link.springer.com/10.1007/978-3-64202091-9>.

Hasan, B.A.S., Gan, J.Q., 2011. Temporal modeling of EEG during self-paced hand

movement and its application in onset detection. J. Neural Eng. 8 (5), 056015.

Hochberg, L.R., Bacher, D., Jarosiewicz, B., Masse, N.Y., Simeral, J.D., Vogel, J.,

Haddadin, S., Liu, J., Cash, S.S., van der Smagt, P., Donoghue, J.P., 2012. Reach and

grasp by people with tetraplegia using a neurally controlled robotic arm. Nature

485 (7398), 372–375 <http://www.pubmedcentral.nih.gov/articlerender.fcgi?

artid=3640850&tool=pmcentrez&rendertype=abstract>.

Hortal, E., Planelles, D., Costa, A., Iáñez, E., Úbeda, A., Azorín, J.M., Fernández, E.,

2015. SVM-based brain-machine interface for controlling a robot arm through

four mental tasks. Neurocomputing 151 (P1), 116–121. http://dx.doi.org/

10.1016/j.neucom.2014.09.078.

Höskuldsson, A., 1988. PLS regression methods. J. Chemom. 2 (July), 211–228.

http://dx.doi.org/10.1002/cem.1180020306

<http://scholar.google.com/

scholar?hl=en&btnG=Search&q=intitle:PLS+regression+methods#0>.

Hotson, G., McMullen, D.P., Fifer, M.S., Johannes, M.S., Katyal, K.D., Para, M.P.,

Armiger, R., Anderson, W.S., Thakor, N.V., Wester, B.A., Crone, N.E., 2016.

Individual finger control of a modular prosthetic limb using high-density

electrocorticography in a human subject. J. Neural Eng. 13 (2), 026017. <http://

www.ncbi.nlm.nih.gov/pubmed/26863276,

http://stacks.iop.org/1741-2552/

13/i=2/a=026017?key=crossref.2d5ef8bc47143308a3a1ff00aabb53dd>.

Jackson, A., Fetz, E.E., 2011. Interfacing with the computational brain. IEEE Trans.

Neural Syst. Rehab. Eng. 19 (5), 534–541.

Jacobs, R.A., Jordan, M.I., Nowlan, S.J., Hinton, G.E., 1991. Adaptive mixtures of local

experts. Neural Comput. 3 (1), 79–87.

Jarosiewicz, B., Masse, N.Y., Bacher, D., Cash, S.S., Eskandar, E., Friehs, G., Donoghue,

J.P., Hochberg, L.R., 2013. Advantages of closed-loop calibration in intracortical

brain computer interfaces for people with tetraplegia. J. Neural Eng. 10 (4),

046012. <http://iopscience.iop.org/1741-2552/10/4/046012, http://iopscience.

iop.org/1741-2552/10/4/046012/pdf/1741-2552_10_4_046012.pdf>.

Kemere, C., Santhanam, G., Yu, B.M., Afshar, A., Ryu, S.I., Meng, T.H., Shenoy, K.V.,

2008. Detecting neural-state transitions using hidden Markov models for motor

cortical prostheses. J. Neurophysiol. 100 (4), 2441–2452.

Kim, K.H., Kim, S.S., Kim, S.J., 2006. Superiority of nonlinear mapping in decoding

multiple single-unit neuronal spike trains: a simulation study. J. Neurosci.

Methods 150 (2), 202–211.

Kim, S.-P., Sanchez, J.C., Erdogmus, D., Rao, Y.N., Wessberg, J., Principe, J.C.,

Nicolelis, M., 2003. Divide-and-conquer approach for brain machine

interfaces: nonlinear mixture of competitive linear models. Neural Netw.:

Off. J. Int. Neural Netw. Soc. 16 (5–6), 865–871 <http://www.ncbi.nlm.nih.gov/

pubmed/12850045>.

Kotecha, J.H., Djuric´, P.M., 2003. Gaussian particle filtering. IEEE Trans. Signal

Process. 51 (10), 2592–2601.

Leeb, R., Friedman, D., Müller-Putz, G.R., Scherer, R., Slater, M., Pfurtscheller, G.,

2007. Self-paced (asynchronous) BCI control of a wheelchair in virtual

environments: a case study with a tetraplegic. Comput. Intell. Neurosci. 2007,

359

79642 <http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=2272302&

tool=pmcentrez&rendertype=abstract>.

Leuthardt, E.C., Miller, K.J., Schalk, G., Rao, R.P.N., Ojemann, J.G., 2006.

Electrocorticography-based brain computer interface–the Seattle experience.

IEEE Trans. Neural Syst. Rehab. Eng.: IEEE Eng. Med. Biol. Soc. 14 (2), 194–198

<http://www.ncbi.nlm.nih.gov/pubmed/16792292>.

Li, Z., 2014. Decoding methods for neural prostheses: where have we reached?

Front. Syst. Neurosci. 8 (July), 1–6 <http://journal.frontiersin.org/article/10.

3389/fnsys.2014.00129/abstract>.

Li, Z., O’Doherty, J.E., Hanson, T.L., Lebedev, M.A., Henriquez, C.S., Nicolelis, M.A.L.,

2009. Unscented Kalman filter for brain-machine interfaces. PLoS ONE 4 (7),

e6243 <http://dx.plos.org/10.1371/journal.pone.0006243>.

Mak, J.N., Wolpaw, J.R., 2009. Clinical applications of brain computer interfaces:

current state and future prospects. EEE Rev. Biomed. Eng. 2, 187–199.

Mason, S., Kronegg, J., Huggins, J., Fatourechi, M., Schlögl, A., 2006. Evaluating the

Performance of Self-Paced Brain-Computer Interface Technology. Tech Report

(January), pp. 1–60. <http://ipl.ece.ubc.ca/bci_files/self_paced_tech_report2006-05-19.pdf>.

Mehring, C., Rickert, J., Vaadia, E., de Oliveira, S.C., Aertsen, A., Rotter, S., 2003.

Inference of hand movements from local field potentials in monkey motor

cortex. Nat. Neurosci. 6 (12), 1253–1254 <http://www.nature.com/doifinder/

10.1038/nn1158>.

Meila, M., Jordan, M.I., 1996. Learning fine motion by Markov mixtures of experts.

Adv. Neural Inform. Process. Syst. 8 (1567), 1003–1009. <http://papers.nips.

cc/paper/1063-learning-fine-motion-by-markov-mixtures-of-experts.pdf,

http://nfiles/1769/Meila?Jordan-1996LearningFineMotionbyMarkovMixturesofExperts.pdf, http://nfiles/1770/1063learning-fine-motion-by-markov-mixtures-of-experts.html>.

Mestais, C.S., Charvet, G., Sauter-Starace, F., Foerster, M., Ratel, D., Benabid, A.L.,

2015. WIMAGINE: wireless 64-channel ECoG recording implant for long term

clinical applications. IEEE Trans. Neural Syst. Rehab. Eng. 23 (1), 10–21.

Morinière, B., Verney, A., Abroug, N., Garrec, P., Perrot, Y., 2015. EMY: a dual arm

exoskeleton dedicated to the evaluation of brain machine interface in clinical

trials. In: IEEE International Conference on Intelligent Robots and Systems, pp.

5333–5338.

Murphy, K., 1998. Switching Kalman Filters. Tech. Rep. August. <http://citeseerx.ist.

psu.edu/viewdoc/download?doi=10.1.1.49.5703&rep=rep1&type=pdf>.

Pfurtscheller, G., Solis-Escalante, T., Ortner, R., Linortner, P., Muller-Putz, G.R., 2010.

Self-paced operation of an SSVEP-based orthosis with and without an imagerybased brain switch: a feasibility study towards a hybrid BCI. IEEE Trans. Neural

Syst. Rehab. Eng. 18 (4), 409–414.

Pistohl, T., Ball, T., Schulze-Bonhage, A., Aertsen, A., Mehring, C., 2008. Prediction of

arm movement trajectories from ECoG-recordings in humans. J. Neurosci.

Methods 167 (1), 105–114 <http://www.ncbi.nlm.nih.gov/pubmed/18022247>.

Rabiner, L., 1989. A tutorial on hidden Markov models and selected applications in

speech recognition. Proc. IEEE 77 (2), 257–286.

Rabiner, L., Juang, B., 1986. An introduction to hidden Markov models. IEEE ASSP

Mag. 3 (January). Appendix 3A.

Renals, S., Morgan, N., Bourlard, H., Cohen, M., Franco, H., 1994. Connectionist

probability estimation in HMM speech recognition. IEEE Trans. Speech Audio

Process. 2 (1), 161–174 <http://ieeexplore.ieee.org/xpls/abs_all.jsp?arnumber=

260359>.

Saa, J.F.D., de Pesters, A., Cetin, M., 2016. Asynchronous decoding of finger

movements from ECoG signals using long-range dependencies conditional

random fields. J. Neural Eng. 13 (3), 36017 <http://stacks.iop.org/1741-2552/13/

i=3/a=036017>.

Schalk,

G.,

Leuthardt,

E.C.,

2011.

Brain-computer

interfaces

using

electrocorticographic signals. IEEE Rev. Biomed. Eng. 4, 140–154.

Schalk, G., Miller, K.J., Anderson, N.R., Wilson, J.A., Smyth, M.D., Ojemann, J.G.,

Moran, D.W., Wolpaw, J.R., Leuthardt, E.C., 2008. Two-dimensional movement

control using electrocorticographic signals in humans. J. Neural Eng. 5 (1), 75–

84

<http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=2744037&

tool=pmcentrez&rendertype=abstract>.

Shimoda, K., Nagasaka, Y., Chao, Z.C., Fujii, N., 2012. Decoding continuous threedimensional hand trajectories from epidural electrocorticographic signals in

Japanese macaques. J. Neural Eng. 9 (3), 036015.

Solis-Escalante, T., Müller-Putz, G., Brunner, C., Kaiser, V., Pfurtscheller, G., 2010.

Analysis of sensorimotor rhythms for the implementation of a brain switch for

healthy subjects. Biomed. Signal Process. Control 5 (1), 15–20.

Spuler, M., Sarasola-Sanz, A., Birbaumer, N., Rosenstiel, W., Ramos-Murguialday, A.,

2015. Comparing metrics to evaluate performance of regression methods for

decoding of neural signals. In: Proceedings of the Annual International

Conference of the IEEE Engineering in Medicine and Biology Society, EMBS,

vol. 2015-November, pp. 1083–1086.

Srinivasan, L., Eden, U.T., Mitter, S.K., Brown, E.N., 2007. General-purpose filter

design for neural prosthetic devices. J. Neurophysiol. 98 (4), 2456–2475.

<http://jn.physiology.org/cgi/reprint/98/4/2456,

http://ovidsp.ovid.com/

ovidweb.cgi?T=JS&PAGE=reference&D=emed8&NEWS=N&AN=2007506292>.

Suway, S.B., Tien, R.N., Jeffries, S.M., Zohny, Z., Clanton, S.T., McMorland, A.J.C.,

Velliste, M., 2013. Resting state detection for gating movement of a neural

prosthesis. In: 2013 6th International IEEE/EMBS Conference on Neural

Engineering (NER), pp. 665–668 <http://ieeexplore.ieee.org/lpdocs/epic03/

wrapper.htm?arnumber=6696022>.

Valstar, M.F., Pantic, M., 2007. Combined support vector machines and hidden

Markov models for modeling facial action temporal dynamics. In: IEEE Int’l

Workshop on. HCI, LNCS, vol. 4796, pp. 118–127.

360

M.-C. Schaeffer, T. Aksenova / Journal of Physiology - Paris 110 (2016) 348–360

Velliste, M., Kennedy, S.D., Schwartz, A.B., Whitford, A.S., Sohn, J.-W., McMorland, A.

J.C., 2014. Motor cortical correlates of arm resting in the context of a reaching

task and implications for prosthetic control. J. Neurosci.: Off. J. Soc. Neurosci. 34

(17), 6011–6022 <http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=

3996218&tool=pmcentrez&rendertype=abstract>.

Velliste, M., Perel, S., Spalding, M.C., Whitford, A.S., Schwartz, A.B., 2008. Cortical

control of a prosthetic arm for self-feeding. Nature 453 (7198), 1098–1101

<http://www.ncbi.nlm.nih.gov/pubmed/18509337>.

Wang, P.T., Puttock, E.J., King, C.E., Schombs, A., Lin, J.J., Sazgar, M., Hsu, F.P.K., Shaw,

S.J., Millett, D.E., Liu, C.Y., Chui, L.A., Do, A.H., Nenadic, Z., 2013a. State and

trajectory decoding of upper extremity movements from electrocorticogram.

In: International IEEE/EMBS Conference on Neural Engineering, NER, 969–972.

Wang, W., Collinger, J.L., Degenhart, A.D., Tyler-Kabara, E.C., Schwartz, A.B., Moran,

D.W., Weber, D.J., Wodlinger, B., Vinjamuri, R.K., Ashmore, R.C., Kelly, J.W.,

Boninger, M.L., 2013. An electrocorticographic brain interface in an individual

with tetraplegia. PLOS ONE 8 (2), e55344 <http://www.pubmedcentral.nih.gov/

articlerender.fcgi?artid=3566209&tool=pmcentrez&rendertype=abstract>.

Wang, Z., Ji, Q., Miller, K.J., Schalk, G., 2011. Prior knowledge improves decoding of

finger flexion from electrocorticographic signals. Front. Neurosci. 5 (November),

127

<http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=3226159&

tool=pmcentrez&rendertype=abstract>.

Waterhouse, S.R., 1997. Classification and Regression using Mixtures of Experts (Ph.

D. Thesis).

Williams, J.J., Rouse, A.G., Thongpang, S., Williams, J.C., Moran, D.W., 2013.

Differentiating closed-loop cortical intention from rest: building an

asynchronous electrocorticographic BCI. J. Neural Eng. 10 (4), 046001 <http://

www.ncbi.nlm.nih.gov/pubmed/23715295>.

Wodlinger, B., Downey, J.E., Tyler-Kabara, E.C., Schwartz, a.B., Boninger, M.L.,

Collinger, J.L., 2015. Ten-dimensional anthropomorphic arm control in a human

brain-machine interface: difficulties, solutions, and limitations. J. Neural Eng.

12 (1), 016011 <http://www.ncbi.nlm.nih.gov/pubmed/25514320>.

Wood, F., Prabhat Donoghue, J., Black, M., 2005. Inferring attentional state and

kinematics from motor cortical firing rates. Conference Proceedings: Annual

International Conference of the IEEE Engineering in Medicine and Biology

Society. IEEE Engineering in Medicine and Biology Society. Conference, vol. 1,

pp. 149–152.

Wu, W., Black, M.J., Gao, Y., Bienenstock, E., Serruya, M.D., Donoghue, J.P., 2002.

Inferring hand motion from multi-cell recordings in motor cortex using a

Kalman Filter. In: SAB’02- Workshop on Motor Control in Humans and Robots:

On the Interplay of Real Brains and Artificial Devices, pp. 66–73.

Wu, W., Black, M.J., Mumford, D., Gao, Y., Bienenstock, E., Donoghue, J.P., 2004.

Modeling and decoding motor cortical activity using a switching Kalman filter.

IEEE Trans. Bio-med. Eng. 51 (6), 933–942 <http://www.ncbi.nlm.nih.gov/

pubmed/15188861>.

Yuksel, S.E., Wilson, J.N., Gader, P.D., 2012. Twenty years of mixture of experts. IEEE

Trans. Neural Netw. Learn. Syst. 23 (8), 1177–1193.