Reconstructing for joint angles on the shoulder and elbow from non-invasive electroencephalographic signals through electromyography

In this study, first the cortical activities over 2240 vertexes on the brain were estimated from 64 channels electroencephalography (EEG) signals using the Hierarchical Bayesian estimation while 5 subjects did continuous arm reaching movements. From the estimated cortical activities, a sparse linear regression method selected only useful features in reconstructing the electromyography (EMG) signals and estimated the EMG signals of 9 arm muscles. Then, a modular artificial neural network was used to estimate four joint angles from the estimated EMG signals of 9 muscles: one for movement control and the other for posture control. The estimated joint angles using this method have the correlation coefficient (CC) of 0.807 (±0.10) and the normalized root-mean-square error (nRMSE) of 0.176 (±0.29) with the actual joint angles.

Introduction

The field of Brain Machine Interface (BMI) has engaged in active research to help paralyzed patients regain some independence and to better integrate within societal activities. Brain-machine interface can be broadly divided into invasive and non-invasive modalities depending on how brain signals are harnessed. Invasive BMI, mainly targets motor related cortical areas. Non-human primates are often used to for example, harness the spikes and local field potentials from the primary motor cortex, known to interface with the spinal cord and containing signals useful to control arm movement. Such neural signal has been used to control external devices such as a robotic arm or a mouse cursor by reconstructing hand trajectories from the measured neural activities (Chapin et al., 1999; Wessberg et al., 2000; Serruya et al., 2002; Talylor et al., 2002; Carmena et al., 2003). More recently invasive BMI has been approved to be used in humans (Hochberg et al., 2006; Chadwick et al., 2011), including as well electrocorticography (ECoG) (Schalk et al., 2007; Sanchez et al., 2008).

In the case of non-invasive BMI, the state-of-the-art research uses motor imagery, a paradigm that classifies whether the subject performs left or right hand motor imagery using electroencephalography (EEG) signals (Ramoser et al., 2000; Wolpaw and McFarland, 2004). It was believed that the noisy EEG signal in non-invasive BMI would be insufficient to estimate three-dimensional hand movement (Lebedev and Nicolelis, 2006). However, recently Bradberry et al. (2010) succeeded in reconstructing the three-dimensional hand movements from the EEG signals while the subjects perform natural and self-initiated reaching actions.

In the present study, a new method using electromyography (EMG) signals is proposed that first reconstructs the EMG signals of the arm muscles from the source currents estimated from EEG signals, and then estimates joint angles on the shoulder and the elbow from the reconstructed EEG signals. By reconstructing the EMG signals of the arm muscles from EEG signals, it is possible to reconstruct not only kinematics-, but also dynamics-based information involving force generation. For example, impedance and joint torque can be obtained to build more realistic brain-machine interfaces, compatible with human motion execution. Furthermore, when reconstructing EMG signals from the EEG signals, by electrically stimulating the arm muscles of a paralyzed person, it is possible using a functional electrical stimulation (FES) system to engage the person in self-adaptive control of his/her arm.

In this study, source currents over 2240 vertexes were estimated from EEG signals of 64 channels through a hierarchical Bayesian method introducing a hierarchical prior (Sato et al., 2004). This method can effectively incorporate both structural and functional MRI data. In this hierarchical Bayesian method, the variance of the source current at each source location is considered an unknown parameter and estimated from the observed EEG data and prior information by using the Variational Bayesian (VB) method. The fMRI information was imposed as prior information on the variance distribution rather than the variance itself so that it gives a soft constraint on the variance. From the estimated source currents over 2240 vertexes, only 33 vertexes are selected, which is located in the left primary motor cortex contralateral to moving arm, to estimate the filtered EMG signals of 9 muscles by using a sparse linear regression method which can automatically select only useful features in estimating the filtered EMG signals. A modular artificial neural network was then used to reconstruct 4 joint angles on the shoulder and elbow from the estimated filtered EMG signals, which trains movement data and posture data in two different networks. This modular structure improves the accuracy of the estimation.

Materials and Methods

Experimental Task

Five healthy right-handed subjects (five men, Mage = 22.51, age range: 20–29 years) participated in the experiment. All five subjects do not have any experience of participating in the experiments of brain-machine interface study before. All participated subjects submitted a written form of consent before starting the experiment. The subjects performed a continuous arm-reaching task as shown in Figure 1A. The task consisted of pushing buttons in the following sequences: Hold-C-A-B, Hold-C-D-B, Hold-D-B-A, and Hold-D-C-A. Theses sequences are explained in greater detail below.

FIGURE 1

Figure 1. Experimental task. (A) A subject performs a continuous arm-reaching task while facing a touch screen displaying five lights and five buttons. (B) Sequential arm-reaching task (Hold-C-A-B sequence).

Here only the Hold-C-A-B sequence is explained (Figure 1B), since the others have similar patterns. First, a subject pushes the hold button for 1 s when the hold signals turns on. If the subject succeeds in pushing the hold button for 1 s, the C button turns on, and the hand of the subject has to move to the C button within 1 s. If the subject is successful in keeping the C button pressed for 1 s, the A button turns on. The hand of the subject is then supposed to push the A button within 1 s and keep it pressed for 1 ṡ. If the A button is successfully pressed for 1 s, the B button turns on. The hand of the subject then has to push the B button within 1 s and keep it pressed for 1 s.

When the subject succeeded in pressing all three buttons correctly, it was considered as success, and only successful trials were analyzed in this study. After pressing the three buttons, the subject takes a rest between 3 and 4 s, then, it goes to the next trial. The task of the next trial is decided randomly. By performing 10 trials for one of four tasks randomly, each subject conducted 40 trials within one set. A total of seven sets of experiment were conducted for one subject. The leave-one-out cross-validation method was used to analyze the measured data by using six sets for training data and one set for the test data.

fMRI Experiment

Figure 2 shows the fMRI task used to collect fMRI data as prior information to estimate cortical activity. One trial consisted of the execution task in which the participant moves the right index finger (e.g., up or down) every 1 s. This is followed by a resting period in which the participant takes a break for 15 s. Each participant conducted 24 trials of the fMRI task. The fMRI activity when participants take a rest (rest periods) was subtracted from the fMRI activity when participants moved their fingers (execution periods). All five participants conducted the fMRI task to get their individual prior information.

FIGURE 2

Figure 2. fMRI task. During 15 s, the participant moves the right index finger in the instructed direction every 1 s. The monitor then goes blank and the participant takes a break for 15 s while watching the monitor.

Estimation of Cortical Activities from EEG Signals

EEG signals were measured at 1 kHz sampling rate on 64 channels by using a biosemi system (Amsterdam, Netherlands). The measured EEG signals were taken baseline corrected (baseline data from −1 to 0 s) and band-pass filtered between 8 and 30 Hz using a fifth-order butterworth filter.

To estimate cortical activities from EEG signals, an inverse filter L (a matrix of dimensions 2240 × 64) in Equation 1 was used. By multiplying real-time EEG signals to the obtained inverse filter as in Equation 1, it is possible to quickly estimate the cortical activity.

$\begin{array}{cc}\begin{array}{l}L\left({\sum }_{\alpha }^{-1}\right)={\sum }_{\alpha }^{-1}\cdot G^{\prime }·{\left(G\cdot {\sum }_{\alpha }^{-1}\cdot G^{\prime }+{\beta }^{-1}{I}_M\right)}^{-1},\\ J(t)=L\left({\sum }_{\alpha }^{-1}\right)·E(t)\end{array}& (1)\end{array}$

Here, E(t) represents measured real-time EEG signals given by 64 × 256 Hz (sampling rate). J(t) denotes the estimated cortical activities over 2240 vertexes every second and is given by 2240 × 256 Hz entries. G (64 × 2240) is a lead field matrix which represents the impulse response of each source vector component at every measurement site (Baillet et al., 2001) and G′ denotes its transpose. The boundaries between brain, skull, and scalp were generated by using the Curry 5 software (Compumedics, USA). Here, the relative conductivities of the brain, skull, and scalp are 10.0125 and 1. I_M represents an identity matrix of M-by-M (M:number of sensors), β⁻¹ (64 × 64) corresponds to the inverse of the noise variance of the observed EEG signals. ∑⁻¹_α denotes the source covariance matrix, and is calculated as ∑⁻¹_α =diag(α⁻¹). Here, α⁻¹ (2240 × 2240) represents the source current variance which is considered unknown parameters in this study and estimated from the measured EEG data by applying a hierarchical prior on current variance.

Artifact dipoles were also incorporated in the estimation according to previous studies (Fujiwara et al., 2009; Morishige et al., 2009). Artifact dipoles were located at the center of the heart, the right shoulder, and wrist joints, the left and right eyeballs, and the carotid arteries, and estimated.

Estimation of current variance

In this study, the current variance α⁻¹ was estimated by the Automatic relevance determination (ARD) hierarchical prior (Neal, 1996).

$\begin{array}{cc}\begin{array}{l}P\text{}\left(J(t)|\alpha ,\beta \right)\propto \exp \text{}\left[-\frac{\beta }{2}{J}^{\prime }(t)·A·J(t)\right]\\ P\text{}\left({\alpha }_i\right)=\Gamma \text{}\left({\alpha }_i|{\alpha }_{0i},r_0\right),\\ P\text{}\left(\beta \right)=\frac{1}{\beta }\end{array}& (2)\end{array}$

where β is the inverse noise variance of the observed EEG signals, A = diag(α), and α is an I-by-1 vector whose component α_i is the inverse current variance corresponding to the i-th current dipole. Γ represents the Gamma distribution with mean α_0i and degree of freedom r₀. Intuitively, the hyper-parameter r₀ represents confidence of the hierarchical prior information. A prior current variance v_0i = α⁻¹_0i represents the prior information on current intensity. For large and small v_0i, estimated current J_i (t) tends to be large and small, respectively. These values were determined from the fMRI information:

$\begin{array}{cc}{v}_{0i}=v_{\text{base}}+\left(m_0-1\right)·v_{\text{base}}·{\left({\widehat{t}}_i\right)}^2,& (3)\end{array}$

where ${\widehat{t}}_i$ is a normalized T-value on the i-th vertex. Normalized T-values are computed by dividing the original T-values by the maximum of those T-values (thus ranging from 0 to 1).

v_base is a baseline of the current variance, which is estimated from the pre-movement interval (1.0–0.5 s before the movement initiation) of the EEG data by a Bayesian minimum norm estimation. A variance magnification parameter m₀, which is the other hyper-parameter, specifies the scaling between the current variances in the baseline and task periods. m₀ = 100 and r₀ = 10 were used.

Due to the hierarchical prior, the estimation problem becomes non-linear and cannot be solved analytically. Therefore, the VB method (Attias, 1999; Sato, 2001) is employed. In the VB method, J(t), α, and β are iterately updated until convergence.

Figure 3A depicts the fMRI activity while subject 1 conducts the Hold-C-A-B sequence task with the right arm. The left primary motor area is strongly activated. The fMRI information was used as the prior information to estimate cortical activities. Figure 3B shows the cortical activities of subject 1 estimated from the EEG signals for the Hold-C-A-B task. As expected, strong cortical activities are estimated in the left motor cortex. Meanwhile, several parts in the visual cortex are activated in the figure. The reason of the activation of the visual cortex is that while the subject performs the task of the experiment, he sees the target buttons emitting high intensity light.

FIGURE 3

Figure 3. (A) The fMRI activity while subject 1 performs the Hold-C-A-B sequence task. (B) Estimated cortical activities of subject 1 for the task of the Hold-C-A-B sequence (an average of 70 trials).

EMG Signal Processing

For all trials in this study, EEG, EMG signals, and the positions of the shoulder, the elbow, and the wrist of the subject were simultaneously measured. EMG signals were collected in the nine muscles involving four degrees of freedom (see Figure 4 and Table 1).

FIGURE 4

Figure 4. Four-degrees-of-freedom arm movement.

TABLE 1

Table 1. Muscles measured for EMG signals.

In order to measure the EMG signals, a silver/silver chloride surface electrode (NE-102, Nihon Kohden) was used. After differential amplification, each signal was sampled at 1 kHz with a 12-bit resolution. The signals were digitally rectified, averaged over 5 ms, and then filtered through a second-order low-pass filter with a cut-off frequency of ~3 Hz (Koike and Kawato, 1995).

$\begin{array}{cc}{f}_{\text{EMG}}(t)={\displaystyle \sum _{j=1}^n{h}_j}\text{EMG}(t-j+1),& (4)\end{array}$

$\begin{array}{cc}h(t)=6.44\times \left({\exp }^{-10.80t}-{\exp }^{-16.52t}\right)\text{},& (5)\end{array}$

The coefficient h_j in Equation 4 can be acquired by sampling h(t) in Equation 5 discretely. The resulting signal is very similar to the actual tension; consequently, it is called quasi-tension (Basmajian and DeLuca, 1985).

The method that uses a low-pass filter to estimate muscle tension shows good performance when the velocity of muscle contraction is slow. However, the method cannot estimate muscle tension precisely when the velocity of contraction is very high, and the method does not consider the non-linear characteristics of muscles, such as length and velocity. However, it is reasonable to assume that the output of the low-pass filter is similar to the actual tension (Mannard and Stein, 1973).

Kinematics

In order to measure the position of the shoulder, the elbow, and the wrist of the subjects, an infrared marker was attached on their arms and measured each position by using a 3D position measurement system (MacReflex, Qualisys). The sampling rate was 120 Hz. In order to calculate the joint angles of the four degrees of freedom in the shoulder and elbow from the positions measured, the inverse kinematics equations (Koike and Kawato, 1994) was used.

In Figure 5, if we set the transition matrix of θ₁ θ₂, …,θ₇ to A_x (θ₁), A_y (θ₂), …, A_z(θ₇) and the transition matrix of l₁ (the length of the upper arm), l₂ (the length of the fore arm), and l₃ (the length of the hand) to L_z (l₁), L_z (l₂), and L_z (l₃), we can represent the transition matrix of A_E, A_W, and A_H, which represents the relation from the elbow position E to the hand position H, like below,

$\begin{array}{cc}{A}_E=A_x\left({\theta }_1\right)A_y\left({\theta }_2\right)L_z(l_1)=\left[\begin{array}{cc}{C}_2& C_2{l}_{zl}\\ 0^T& 1\end{array}\right],& (6)\end{array}$

Here,

$\begin{array}{cc}{C}_2=C_x\left({\theta }_1\right)C_y\left({\theta }_2\right)=\left[\begin{array}{ccc}1& 0& 0\\ 0& c_1& -s_1\\ 0& s_1& c_1\end{array}\right]\times \left[\begin{array}{ccc}{c}_2& 0& s_2\\ 0& 1& 0\\ -s_2& 0& c_2\end{array}\right],& (7)\end{array}$

$\begin{array}{cc}{l}_{z1}={\left[\begin{array}{ccc}0& 0& -l_1\end{array}\right]}^T,& (8)\end{array}$

Equation 6 becomes

$\begin{array}{cc}{A}_E=\left[\begin{array}{cccc}{c}_2& 0& s_2& -s_2{l}_1\\ s_1{s}_2& c_1& -s_1{c}_2& s_1{c}_2{l}_1\\ -c_1{s}_2& s_1& c_1{c}_2& -c_1{c}_2{l}_1\\ 0& 0& 0& 1\end{array}\right],& (9)\end{array}$

The coordination of the elbow E (x_E, y_E, z_E) is represented in the 4th column that is,

$\begin{array}{cc}\left[\begin{array}{c}{x}_E\\ y_E\\ z_E\end{array}\right]=-l_1\text{}\text{}\left[\begin{array}{c}{s}_2\\ -s_1{c}_2\\ c_1{c}_2\end{array}\right],& (10)\end{array}$

Therefore,

$\begin{array}{cc}\left[\begin{array}{c}{x}_W\\ y_W\\ z_W\end{array}\right]=\left[\begin{array}{c}{x}_E\\ y_E\\ z_E\end{array}\right]+l_2\text{}\left[\begin{array}{c}-s_2{c}_4-c_2{s}_3{s}_4\\ s_1{c}_2{c}_4+(c_1{c}_3-s_1{s}_2{s}_3)s_4\\ -c_1{c}_2{c}_4+(s_1{c}_3-c_1{c}_2{s}_3)s_4\end{array}\right],& (11)\end{array}$

Finally, we can get the following Equations (12–15).

$\begin{array}{cc}\tan {\theta }_1=-\frac{y_E}{z_E},& (12)\end{array}$

$\begin{array}{cc}\sin {\theta }_2=-\frac{x_E}{l_1},& (13)\end{array}$

$\begin{array}{cc}\sin {\theta }_3=\frac{\left(x_E-x_W\right)\text{}/l_2-\sin {\theta }_2\cos {\theta }_4}{\cos {\theta }_2\sin {\theta }_4},& (14)\end{array}$

$\begin{array}{cc}\cos {\theta }_4=\frac{l_1^2+l_2^2-R^2}{2l_1{l}_2},& (15)\end{array}$

FIGURE 5

Figure 5. The structure of the joints of a human.

Estimation of EMG Signals from Estimated Cortical Activities

A sparse linear regression method (Toda et al., 2011) was used to estimate filtered EMG signals from the cortical activities estimated over 2240 vertexes.

$\begin{array}{cc}f{\text{EMG}}_i\text{}\left(t+\delta t\right)={\displaystyle \sum _{j=1}^{N_{\text{source}}}{w}_{ij}}\times J_j(t)+\text{bias},& (16)\end{array}$

Here, fEMG_i describes the i-th filtered EMG signal from the cortical activity on the j-th vertex (J_j). N_source denotes the number of vertexes used in estimating filtered EMG signals. In this study, since all subjects are right-handed, the cortical activities over 33 vertexes in the left primary motor cortex were used to estimate the filtered EMG signals. The weighting factor w_ij represents the strength influence from the cortical activity on the j-th vertex on muscle i-th muscle. δt is the delay between the cortical activity of the primary motor cortex and the EMG signals.

Modular Artificial Neural Network Model

In order to estimate joint angles from the filtered EMG signals, a modular artificial neural network (Jacobs et al., 1991) was used, as shown in Figure 6. Training the data of posture and movement in different networks will improve the accuracy of estimating joint angles compared to training the entire set of data in the same network, since the muscle tension is different in these two cases. Here, posture is defined as the state where the arm of the subject is in contact with a button on the screen, and movement is defined as the condition where the arm of the subject moves from one button to another. If training is done well, a gating network will select one of the two expert networks by its input signal. In this case, one of the two expert networks is used for posture control and the other is used for movement control. Since the gating network determines the output ratio for each expert network depending on its input signal, the sum of the outputs of the gating network should always be equal to 1.

FIGURE 6

Figure 6. Joint-angle estimation model with a modular architecture.

To achieve this, as shown in Equation 17, the output g_j of the gating network, which corresponds to the j-th expert network, is normalized by using the soft max activation function.

$\begin{array}{cc}{g}_j=\frac{e^{x_j}}{{\displaystyle \sum _{i=1}^N{e}^{x_i}}},& (17)\end{array}$

Here, x_i is the value determined by the input signal of the gating network and N is the total number of outputs of the gating network. The total output is calculated by multiplying the output of the gating network by the output of each expert network and summing the result, as given in Equation 18.

$\begin{array}{cc}\theta ={\displaystyle \sum _{i=1}^N{g}_i}{\widehat{\theta }}_i,& (18)\end{array}$

The gating network and each expert network are trained to maximize the likelihood function lnL (Equation 19) by the back propagation algorithm (Rumelhart et al., 1986).

$\begin{array}{cc}\ln L=\ln {\displaystyle \sum _{i=1}^N{g}_i}{e}^{\frac{-{\Vert \theta -{\widehat{\theta }}_i\Vert }^2}{2{\sigma }_i^2}},& (19)\end{array}$

The update of the weights of the gating network is calculated by a chain rule, as in Equation 20.

$\begin{array}{cc}\frac{\partial \ln L}{\partial x_i}={\displaystyle \sum _{i=1}^N\left(g(i|X,{\widehat{\theta }}_i)-g_i\right)}\text{},& (20)\end{array}$

Here, X is the input of the gating network, and the posteriori probability $g\left(i|X,{\widehat{\theta }}_i\right)$ is

$\begin{array}{cc}g\text{}\left(i|X,{\widehat{\theta }}_i\right)=\frac{g_i{e}^{\frac{-{\Vert \theta -{\widehat{\theta }}_i\Vert }^2}{2{\sigma }_i^2}}}{{\displaystyle \sum _{j=1}^N{g}_j}{e}^{\frac{-{\Vert \theta -{\widehat{\theta }}_j\Vert }^2}{2{\sigma }_j^2}}},& (21)\end{array}$

The update of the weights of each expert network is calculated by a chain rule as in Equation 22.

$\begin{array}{cc}\frac{\partial \ln L}{\partial {\widehat{\theta }}_i}={\displaystyle \sum _{i=1}^{N}g}\text{}\left(i|X,{\widehat{\theta }}_i\right)\frac{\theta -{\widehat{\theta }}_i}{{\sigma }_i^2},& (22)\end{array}$

Each network is trained by using the kick-out method (Ochiai and Usui, 1993).

The filtered EMG signals of the nine muscles were used as the input of each expert network model. The summed-squared velocity value of the four joint angles were used as the input of the gating network because when the value of the cortical activities in the primary motor cortex was directly used as the input of the gating network, the gating network could not distinguish between posture and movement. However, when using the summed-squared velocity of the four joint angles as the input, the gating network distinguished posture and movement correctly.

Analysis

The correlation coefficient (CC) was used to evaluate the similarity between actual and predicted signals. Accuracy was also evaluated using normalized root-mean-square error (nRMSE) between actual and predicted signals, defined as

$\begin{array}{cc}\text{nRMSE}=\sqrt{\frac{{\displaystyle \sum _{i=1}^n{\left(y_i^{\text{predicted}}-y_i^{\text{actual}}\right)}^2}}{n}}/\left(y_{\max }^{\text{actual}}-y_{\min }^{\text{actual}}\right)\text{},& (23)\end{array}$

where for each time i (i = 1, 2, …, n), y^predicted_i is the predicted signal and y^actual_i is the actual signal, and y^actual_max and y^actual_min are the maximum and minimum of actual signal, respectively.

Results

Estimation Result of Filtered EMG Signals from the Cortical Activities of the Primary Motor Cortex

The filtered EMG signals were estimated from the cortical activities in the primary motor cortex by using Equation 23. To determine the delay-time parameter, the intracortical microstimulation (ICMS) method (Heusler et al., 2000) was used and the delay time 17 ms was decided when the filtered EMG signals are estimated from the cortical activities of the primary motor cortex.

Of the 70 trials measured for each task (Hold-C-A-B, etc.), 60 trials were used for training data and 10 trials for the test data. The sparse linear regression method has an ability to automatically select only useful features in estimation among all extracted features. Therefore, this method is very strong against the overfitting problem. Figure 7A shows the weights of the selected features in estimating filtered EMG signals from the cortical activities in the left primary motor cortex while subject 3 performs the experimental task. In the case of Figure 7A, 20 vertexes are selected among 33 vertexes located in the left primary motor cortex to estimate filtered EMG signals.

FIGURE 7

Figure 7. Reconstruction of the filtered EMG signals from the cortical activities estimated on 33 vertexes in the left primary motor cortex. (A) The weights of the important features selected by the sparse linear regression to estimated filtered EMG signals from the cortical activities while subject 3 performs four tasks. (B) The filtered EMG signals estimated from the selected 20 features. Dotted lines (blue) represents the actual filtered EMG signals, and solid lines (red) show the reconstructed filtered EMG signals (normalized scale).

Figure 7B shows the filtered EMG signals of subject 1 estimated from the cortical activities over selected 20 vertexes in the left primary motor cortex. The estimated filtered EMG signals had a CC of 0.827 (±0.10) and nRMSE of 0.142 (±0.38) with the actual EMG signals. Table 2 shows the CC between the actual EMG signals and the reconstructed EMG signals of all of the 5 subjects participated in the experiment. The averaged CC and nRMSE of 5 subjects were 0.851 (±0.11) and 0.233 (±0.17).

TABLE 2

Table 2. The correlation coefficient (CC) and normalized root-mean-square error (nRMSE) between the actual EMG signals and the estimated EMG signals.

Estimation Result of Joint Angles from Filtered EMG Signals

After measuring 70 trials of the EMG signals and movement trajectories of the subject's arm, 60 trials were used as training data and one trial as test data. The number of training data samples was 1,080,720 (60 trials × 1 kHz × 4.503 s × 4 cases) and the number of test data samples was 180,120 (10 trials × 1 kHz × 4.503 s × 4 cases). In the case of the gating network, the network was trained by the summed-squared velocity value of the four joint angles. However, since this value cannot be used as test data, the velocity values from the filtered EMG signals were estimated. Figure 8 shows the four joint angles of subject 1 estimated from the cortical activities of the primary motor cortex. The CC and nRMSE between the estimated joint angles and the actual joint angles were about 0.817 (±0.10) and 0.212(±0.04). Table 3 depicts the CC between the actual joint angles of 5 subjects and the joint angles reconstructed by the modular artificial neural network model. The averages of the CC and nRMSE of the reconstructed joint angles were 0.807 (±0.10) and 0.176 (±0.29).

FIGURE 8

Figure 8. Estimated joint angles from the cortical activities of the primary motor cortex. Dashed lines represent the actual joint angles and solid lines show the estimated joint angles.

TABLE 3

Table 3. The correlation coefficient (CC) and normalized root-mean-square error (nRMSE) between the actual joint angles and the joint angles estimated by the modular artificial neural network model.

Discussion

In this study, the cortical activities on 2240 vertexes were estimated from the EEG signals of 64 channels using the hierarchical Bayesian method. Then, of the estimated cortical activities, only the cortical activities in the left primary motor cortex were used to reconstruct the EMG signals of nine muscles through the sparse linear regression method. When reconstructing EMG signals from the cortical activities, we could determine the delay time between the cortical activities and the EMG signals by searching the correlation of those two signals. However, the pattern of EMG signals has a simple waveform which has one or two peaks, and that of cortical activities is also similar. Thus, in this study, the ICMS method was used to decide the delay time. The delay time of 17 ms found from ICMS method was applied for all the subjects for estimating EMG signals from the cortical activities. In the future, we are going to study whether or not this delay time is effective for individuals in disease state. A modular artificial neural network model was used to estimate four joint angles on the elbow and the shoulder from the estimated EMG signals.

Why is it Important to Reconstruct EMG Signals from Brain Signals?

Morrow et al. (Morrow and Miller, 2003) succeeded in reconstructing the EMG signals of the distal forelimb muscles from the 50 M1 neurons of an non-human primate while performing a stereotyped precision grips task. Furthermore, Koike et al. (2006) estimated the EMG signals of seven arm muscles from the neural activities of 18 neurons of an non-human primate during an arm reaching task. Then, three joint angles (two at the shoulder and one at the elbow) were reconstructed from the estimated EMG signals. Similarly, most existing brain-machine interface studies reconstruct EMG signals from the neural activities of the primary motor cortex of non-human primates, by using invasive needle electrodes. In such cases, it is possible to obtain relatively clean brain signals.

When reconstructing EMG signals with non-invasive BMI technologies however, there are several difficulties because the skull, which is an insulator, is located between the brain and the sensors, thus introducing noise. Ganesh et al. (2008) succeeded in reconstructing the EMG signals of two antagonist muscles from fMRI signals measured in the primary motor cortex and pre-motor cortex. EEG signals have good time resolution, but its spatial resolution is poor. Consequently, it is difficult to estimate EMG signals with EEG signals. In this study, spatial resolution is improved by estimating cortical activities over 2240 vertexes from the EEG signals measured over 64 channels through the hierarchical Bayesian method. Among the features being abundant, the sparse linear regression method automatically selects only useful features in reconstructing EMG signals. The proposed method is very robust against the overfitting problem.

When EMG signals are reconstructed from the brain signals, there are several advantages: First, we can reconstruct not only position related information such as hand position but also force related information such as joint torque and stiffness from the estimated EMG signals (Koike and Kawato, 1993, 1994, 1995). For example, when we pick up an object, the brain stabilizes the posture of the arm by controlling muscle tensions. The stiffness is controlled by the co-contraction of the muscles. It is difficult to model this phenomenon by directly estimating the hand position because co-contraction causes different muscle patterns for the same posture. Similarly, when in addition to reconstructing the kinematics of hand motion, we obtain force information such as joint torque and stiffness from the brain signals, it is possible to control a robotic arm based on these information. In such cases it could also be possible to implement a brain-machine interface more compatible with features of the human arm.

Second, by using the estimated EMG signals as the command signals of the FES, we raise the possibility that a paralyzed person could in principle control his arm once we electrically stimulate his paralyzed muscles (Degnan et al., 2002; Uechi et al., 2004). Fagg et al. (2007), without modeling the characteristics of the musculoskeletal system, controlled arm movement by electrical stimulation of arm muscles through FES after reconstructing EMG signals from the neural activities in the primary motor cortex. Furthermore, Moritz et al. (2008), by facilitating the direct control of the stimulation of muscles from the neural activities of the primary motor cortex, made it possible for non-human primates to control bidirectional wrist torques from cortical cells. This research suggests that it may be possible to create more realistic neuro-prostheses. By modeling the musculoskeletal system, we may be able to extend non-invasive brain-machine interfaces to control anthropomorphic robotic devices.

Which Brain Part is Measured for Reconstructing EMG Signals?

The cortical region of choice to harness the control neural signal from seems to be important. In the case of studies using non-human primates, the neural activities of the primary motor cortex are mainly measured to reconstruct EMG signals (Nicolelis et al., 1998; Shoham et al., 2005; Wu et al., 2006). From the research result of Fried et al. (1991), the process of motor related information in the brain is that first the urge to move the arm occurs from the premotor cortex, then the occurred signals goes to the primary motor cortex via the supplementary motor area. The primary motor cortex is the final output part of motor related signals in the brain. The signal is transmitted to the arm muscles through the alpha motor neuron of the spinal cord, and finally it generates arm movement. Anatomically, since the primary motor cortex is linked to the muscles via one or more intermediate neurons, the neural activities of the primary motor cortex have high correlation with muscle activities. Figure 9 shows the CCs when reconstructing EMG signals from the cortical activities in several brain areas. In this study, it is found that when reconstructing EMG signals from the cortical activities estimated in the primary motor cortex, the highest CC is obtained.

FIGURE 9

Figure 9. The correlation coefficients when reconstructing EMG signals from the cortical activities estimated in different brain areas (M1, primary motor cortex; PMd, dorsal premotor cortex; PP, posterior parietal cortex; and All, using all brain areas).

Limitations and Future Work

There are some limitations with the use of modular neural networks for joint angle estimation. The estimated joint angles have a CC of 0.81 with the actual joint angles. The reason the modular artificial neural network model was used in estimating joint angles is because, in the case of isotonic movement, where force is outputted with a changing length of the muscle, the tension is different depending on the velocity that the muscle flexes or extends. In the case of muscle flexion, the tension decreases as the flex velocity increases. In the case of muscle extension, the tension increases as the extension velocity increases. The performance of estimating joint angles could be improved by training two networks with tension values, which change depending on the velocity, rather than training the data in the same network. One network was used for 0 velocity and the other for movement velocity. When joint angles are estimated from muscle tensions, the muscle tensions for posture have low values. In comparison, the muscle tensions for movement have significantly high values. If we trained these data in the same network, the network would determine that the error of posture data is much lower than that of movement data. Consequently, in the case of posture data, the estimated results are poor. In future work we will use different neural network structures for joint angle estimation. Recent advances in machine learning point at deep learning algorithms and neural networks (Salakhutdinov and Hinton, 2009, 2012) as a possibility for improving feature extraction to reconstruct the joint angles. We plan to explore these new avenues of research.

In this study, five normal subject's joint angles were estimated from EEG signals through EMG signals. In the case of individuals with spinal cord injuries where the pathway between primary motor cortex and muscles was disconnected, there was a necessity of identifying the relationship between EMG signals and joint angles of a normal subject. Then, the EEG signals of an individual with spinal cord injury is connected to EMG signals of the normal subject. In the future, we are going to study more about this topic with individuals with spinal cord injury. Furthermore, there is a possibility of using this proposed method in a study of post-stroke individual where primary motor cortex is not damaged.

Conflict of Interest Statement

The author declares that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

References

Attias, H. (1999). “Inferring parameters and structure of latent variables models by variation Bayes,” in Proceedings of 15th Conference on Uncertainty in Artificial Intelligence (Stockholm), 21–30.

Baillet, S., Riera, J. J., Marin, G., Mangin, J. F., Aubert, J., and Garnero, L. (2001). Evaluation of inverse methods and head models for EEG source localization using a human skull phantom. Phys. Med. Biol. 46, 77–96. doi: 10.1088/0031-9155/46/1/306

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Basmajian, J. V., and DeLuca, C. J. (1985). Muscles Alive. Williams and Wilkins.

Bradberry, T. J., Gentili, R. J., and Contreras-Vidal, J. L. (2010). Reconstructing three-dimensional hand movements from noninvasive electroencephalographic signals. J. Neurosci. 30, 3432–3437. doi: 10.1523/JNEUROSCI.6107-09.2010

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Carmena, J. M., Lebedev, M. A., Crist, R. E., O'Doherty, J. E., Santucci, D. M., Dimitrov, D. F., et al. (2003). Learning to control a brain-machine interface for reaching and grasping by primats. PLoS Biol. 1:e42. doi: 10.1371/journal.pbio.0000042

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Chadwick, E. K., Blana, D., Simeral, J. D., Lambrecht, J., Kim, S. P., Cornwell, A. S., et al. (2011) Continuous neural ensemble control of simulated arm reaching by a human with tetraplegia. J. Neural Eng. 8, 1741–1750. doi: 10.1088/1741-2560/8/3/034003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Chapin, J. K., Moxon, K. A., Markowitz, R. S., and Nicolelis, M. A. L. (1999). Real-time control of a robot arm using simultaneously recorded neurons in the motor cortex. Nat. Neurosci. 2, 664–670. doi: 10.1038/10223

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Degnan, G. G., Wind, T. C., Jones, E. V., and Edlich, R. F. (2002). Functional electrical stimulation in tetraplegic patients to restore hand function. J. Long Term Eff. Med. Implants 12, 175–188.

Pubmed Abstract | Pubmed Full Text

Fagg, A. H., Hatsopoulos, N. G., de Lafuente, V., Moxon, K. A., Nemati, S., Rebesco, J. M., et al. (2007). Biomimetic brain machine interfaces for the control of movement. J. Neurosci. 27, 11842–11846. doi: 10.1523/JNEUROSCI.3516-07.2007

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Fried, I., Katz, A., McCarthy, G., Sass, K., Williamson, P., Spencer, S., et al. (1991). Functional organization of human supplementary motor cortex studied by electrical stimulation. J. Neurosci. 11, 3656–3666.

Pubmed Abstract | Pubmed Full Text

Fujiwara, Y., Yamashita, O., Kawawaki, D., Doya, K., Kawato, M., Toyama, K., et al. (2009). A hierarchical Bayesian method to resolve an inverse problem of MEG contaminated with eye movement artifacts. Neuroimage 45, 393–409. doi: 10.1016/j.neuroimage.2008.12.012

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Ganesh, G., Burdet, E., Haruno, M., and Kawato, M. (2008). Sparse linear regression for reconstructing muscle activity from human cortical fMRI. Neuroimage 42, 1463–1472. doi: 10.1016/j.neuroimage.2008.06.018

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Heusler, P., Cebulla, B., Boehmer, G., and Dinse, H. R. (2000). A repetitive intracortical microstimulation pattern induces long-lasting synaptic depression in brain slices of the rat primary somatosensory cortex. Exp. Brain Res. 135, 300–310. doi: 10.1007/s002210000530

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Hochberg, L. R., Serruya, M. D., Friehs, G. M., Mukand, J. A., Saleh, M., Caplan, A. H., et al. (2006). Neuronal ensemble control of prosthetic devices by a human with tetraplegia. Nature 442, 164–171. doi: 10.1038/nature04970

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Jacobs, R. A., Jordan, M. I., Nowlan, S. J., and Hinton, G. E. (1991). Adaptive mixtures of local experts. Neural Comput. 3, 79–87. doi: 10.1162/neco.1991.3.1.79

CrossRef Full Text

Koike, Y., and Kawato, M. (1993). Virtual trajectories predicted from surface EMG signals. Soc. Neurosci. 12:543.

Koike, Y., and Kawato, M. (1994). Estimation of arm posture in 3D-space from surface EMG signals using a neural network model. IEICE Trans. Fundam. E77-D, 368–375.

Koike, Y., and Kawato, M. (1995). Estimation of dynamic joint torques and trajectory formation from surface electromyography signals using a neural network model. Biol. Cybern. 73, 291–300. doi: 10.1007/BF00199465

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Koike, Y., Hirose, H., Sakurai, Y., and Iijima, T. (2006). Prediction of arm trajectory from a small number of neuron activities in the primary motor cortex. Neurosci. Res. 55, 146–153. doi: 10.1016/j.neures.2006.02.012

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Lebedev, M. A., and Nicolelis, M. A. (2006) Brain-machine interfaces: past, present and future. Trends Neurosci. 29, 536–546. doi: 10.1016/j.tins.2006.07.004

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Mannard, A., and Stein, R. (1973). Determination of the frequency response of isometric soleus muscle in the cat using random nerve stimulation. J. Physiol. 229, 275–296.

Pubmed Abstract | Pubmed Full Text

Morishige, K., Kawawaki, D., Yoshioka, T., Sato, M. A., and Kawato, M. (2009). “Artifact removal using simultaneous current estimation of noise and cortical sources,” in Lecture Notes in Computer Science: ICONIP 2008 (Berlin: Springer-Verlag), 335–342.

Moritz, C. T., Perlutter, S. I., and Fetz, E. E. (2008). Direct control of paralyzed muscles by cortical neurons. Nature 456, 639–642. doi: 10.1038/nature07418

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Morrow, M. M., and Miller, L. E. (2003) Prediction of muscle activity by populations of sequentially recorded primary motor cortex neurons. J. Neurophysiol. 89, 2279–2288. doi: 10.1152/jn.00632.2002

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Neal, R. M. (1996). Bayesian Learning for Neural Networks. New York, NY: Springer. doi: 10.1007/978-1-4612-0745-0

CrossRef Full Text

Nicolelis, M. A. L., Ghazanfar, A. A., Stambaugh, C. R., Oliveira, L. M. O., Laubach, M., Chapin, J. K., et al. (1998). Simultaneous encoding of tactile information by three primary cortical area. Nat. Neurosci. 1, 621–630. doi: 10.1038/2855

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Ochiai, K., and Usui, S. (1993). Improved kick out learning algorithm with delta-bar-delta-bar rule. IEEE Int. Conf. Neural Netw. 1, 269–274. doi: 10.1109/ICNN.1993.298568

CrossRef Full Text

Ramoser, H., Muller-Gerking, J., and Pfurtscheller, G. (2000). Optimal spatial filtering of single trial EEG during imagined hand movement. IEEE Trans. Rehabil. Eng. 8, 441–446. doi: 10.1109/86.895946

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Rumelhart, D., Hinton, G., and Williams, R. (1986) Learning representations by back-propagation errors. Nature 323, 533–536. doi: 10.1038/323533a0

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Salakhutdinov, R., and Hinton, G. (2009). Deep Boltzmann Machines. NIPS Neural Information Processing Systems. Vancouver, BC: MIT Press.

Salakhutdinov, R., and Hinton, G. (2012). An efficient learning procedure for deep Boltzmann machines. Neural Comput. 24, 1967–2006. doi: 10.1162/NECO_a_00311

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Sanchez, J. C., Gunduz, A., Carney, P. R., and Principe, J. C. (2008) Extraction and localization of mesoscopic motor control signals for human ECoG neuroprosthetics. J. Neurosci. Methods 167, 63–81. doi: 10.1016/j.jneumeth.2007.04.019

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Sato, M. A. (2001). Online model selection based on the variational Bayes. Neural Comput. 13, 1649–1681. doi: 10.1162/089976601750265045

CrossRef Full Text

Sato, M. A., Yoshioka, T., Kajihara, S., Toyama, K., Goda, N., Doya, K., et al. (2004). Hierarchical Bayesian estimation for MEG inverse problem. Neuroimage 23, 806–826. doi: 10.1016/j.neuroimage.2004.06.037

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Schalk, G., Kubanek, J., Miller, K. J., Anderson, N. R., Leuthardt, E. C., Ojemann, J. G., et al. (2007) Decoding two-dimensional movement trajectories using electrocorticographic signals in humans. J. Neural Eng. 4, 264–275. doi: 10.1088/1741-2560/4/3/012

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Serruya, M. D., Hatsopoulos, N. G., Paninski, L., Fellows, M. R., and Donoghue, J. P. (2002). Instant neural control of a movement signal. Nature 416, 141–142. doi: 10.1038/416141a

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Shoham, S., Paninski, L. M., Fellowss, M. R., Hatsopoulos, N. G., Donoghue, J. P., and Normann, R. A. (2005). Statistical encoding model for a primary motor cortical brain-machine interface. IEEE Trans. Biomed. Eng. 52, 1312–1322. doi: 10.1109/TBME.2005.847542

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Talylor, D. M., Tillery, S. I., and Schwartz, A. B. (2002). Direct cortical control of 3D neuroprostheic devices. Science 296, 1829–1832. doi: 10.1126/science.1070291

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Toda, A., Imamizu, H., Kawato, M., and Sato, M. A. (2011). Reconstruction of two-dimensional movement trajectories from selected magnetoencephalography cortical currents by combined sparse Bayesian methods. Neuroimage 54, 892–905. doi: 10.1016/j.neuroimage.2010.09.057

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Uechi, M., Naito, Y., Sato, M., and Koike, Y. (2004). Stiffness teaching and motion assist system using functional electrical stimulation and electromyogram signals. J. Rob. Mechatron. 16, 446–455.

Wessberg, J., Stambaugh, C. R., Kralik, J. D., Beck, P. D., and Laubach, M. (2000). Real time prediction of hand trajectory by ensembles of cortical neurons in primates. Nature 408, 361–365. doi: 10.1038/35042582

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Wolpaw, J. R., and McFarland, D. J. (2004). Control of a two-dimensional movement signal by a noninvasive brain-computer interface in humans. Proc. Natl. Acad. Sci. U.S.A. 101, 17849–17854. doi: 10.1073/pnas.0403504101

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text

Wu, W., Gao, Y., Bienenstock, E., Donoghue, J. P., and Black, M. J. (2006). Bayesian population decoding of motor cortical activity using a Kalman filter. Neural Comput. 18, 80–118. doi: 10.1162/089976606774841585

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text