A BCI based visual-haptic neurofeedback training improves cortical activations and classification performance during motor imagery

Nineteen right-handed adults (12 males and seven females, age range 21–28 years old, mean  ±  std.: 24.2  ±  2.27) participated in this experiment. All subjects are healthy without any history of neurological illness or limb movement disorder. Three subjects had no prior experience of MI-based BCIs. All of the subjects signed a consent form in advance. The procedure of the experiment was clearly explained to each subject before EEG data recording. This study was approved by the ethical committee of Tianjin University. See the supplementary data available online at (stacks.iop.org/JNE/16/066012/mmedia).

In this study, we are aiming to apply BCI based NFT to improve subjects' sensorimotor cortical activations and classification performance during MI. To this end, the experiment consisted of seven sessions containing four task stages (one screening session, two previous control sessions, two NFT sessions and two posterior control sessions) in total. Figure 1(a) illustrates the experimental paradigm procedure. Subjects could have an impermanent break between these sessions. For the first screening session, all the experimental scenarios would appear on the LCD screen. Subjects were asked to view and get familiar with the experimental paradigm. The screening session lasted about 5 min. The next two task sessions were the previous control sessions utilizing the classical Graz MI-BCI paradigm [9, 22]. And then, two visual- haptic NFT sessions were conducted lasting about 30 min. At last, there would be two posterior control sessions same as the previous control sessions. However, subjects were asked to apply the kinesthetic MI strategy trained in the NFT sessions.

Zoom In Zoom Out Reset image size

Figure 1(b) illustrates the experimental paradigm utilizing the classical Graz MI-BCI paradigm for previous and posterior control sessions, each consisting of 20 left hand MI (LH-MI) trials and 20 right hand MI (RH-MI) trials with randomized order. At the beginning of a trial, a green cross appeared in the middle of the screen and stayed for the preparation. After 1 s, a cue in the form of a green cross with left or right red arrow appeared for 1 s, indicating the subjects should imagine left or right kinesthetic grasping movement in the next 4 s with a green cross in the middle of the screen. Then the fixation cross disappeared and the screen became blank for a random time interval between 2 and 4 s. Subjects were asked to avoid electromyogram artifacts by restricting movements like eye blinking or swallowing during the task period.

Figure 1(c) illustrates the experimental paradigm for visual- haptic NFT sessions. Each session consisted of LH-MI training or RH-MI training with randomized order lasting about 15 min. At the beginning, a 'Relax' text appeared on the screen lasting for 30 s as the baseline of NFT. And then, at the beginning of each trial, a green cross appeared in the middle of the screen and stayed for the preparation. After 1 s, a cue in the form of a green cross with left or right red arrow appeared for 1 s, indicating the subjects should perform LH-MI or RH-MI training in the next 4 s. During this period, two virtual hands in the middle of the screen would get close to or away from the needles once per second according to real-time NFT parameters. Once the needle touched the hand (t₀), there would be an electrical stimulation (ES) as shown in figure 1(d) applying to left and right hand palms of the subject until the needle get away from the hand or this trial ended (t_s). Here, to induce better kinesthetic MI for all of the subjects, the ES protocol was utilized as a real-time proprioceptive feedback and haptic cue but not causing pain. At last, the screen became blank for a random time interval between 2 and 4 s. Subjects could have a break and self-regulation to wait for next trial.

Subjects were seated comfortably in a chair about 1 m away from a 23.5-inch LCD screen. EEG signals were acquired by a 64-channel SynAmps2 system (Neuroscan, Australia) with standard Ag/AgCl electrodes placed on the scalp according to the international 10–20 system. The reference electrode was placed on the nose and the ground electrode was placed on the forehead. The impedance for all electrodes was kept below 10 kΩ. The EEG signals were sampled at 1000 Hz. An online band-pass filter between 0.5 and 100 Hz and 50 Hz notch filter were enabled in the amplifier to filter the high-frequency noise, baseline drift and power line interference. For the preprocessing, the EEG raw data were downsampled at 200 Hz, then processed by the common average reference (CAR).

As aforesaid, the neurofeedback interface contained two virtual hands and two needles as shown in figure 2 (a RH-MI training trial as example). The angle range (0°, 15°, 30°, 45°) of virtual hand reflected a real-time feedback of the NFT parameter. If the training parameter was below the threshold (Goal_1, up to 15°), the virtual hands changed their angles to get away from the needles. The angles of the virtual hands increased whenever the feedback parameter stayed below the threshold for more than 2 s (Goal_2, up to 30°) and 3 s (Goal_3, up to 45°). Thus, the ultimate goal for the subjects was to keep away from the needles as far as possible. On the contrary, if the training parameter was above the threshold (the needles touched the hands, 0°), an ES (with biphasic current pulses of 300 μs duration, two self-adhesive ECG electrodes were placed at the obverse and reverse of each palm) would be applied to left and right hand palms, reminding the subject to have a self-regulation. The subjects were thus expected to apply better mental strategies of kinesthetic MI (imagining the haptic sensation, force or position of squeezing an object such as a ball, but not its scene in the mind) to achieve the goals, while no specific mental strategies were prescribed [14]. They were asked to utilize one strategy in a trial and could change it in the next trial if the current one was not successful to achieve Goal_1. We also collected subjects' mental strategies of kinesthetic MI after a NFT session and reminded them to use the same strategy in next sessions if it worked well. The threshold would be decreased by 10% in the next session if the subject had a relatively high performance. Theoretically, better mental strategies of kinesthetic MI after the NFT could achieve stronger sensorimotor cortical activations and higher BCI performance.

Zoom In Zoom Out Reset image size

For LH-MI or RH-MI training, the threshold of NFT parameter was here set to the mean power of lateralized relative ERD (lrERD) at alpha (8–13 Hz) and beta (14–29 Hz) rhythm bands respectively during the eyes-open relaxing baseline. The NFT parameter for the LH-MI or RH-MI training was the real-time lrERD power between C3 and C4 channels as shown in figure 3, which was calculated by the following equation:

Zoom In Zoom Out Reset image size

For RH-MI training parameter:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {\rm lrER}{{{\rm D}}_{{\rm C}3-{\rm C}4}}={ R}{{{ P}}_{{\rm ERD},{\rm C}3}}-{ R}{{{ P}}_{{\rm ERD},{\rm C}4}}.\nonumber \end{align} \tag{ 1 }$$

For LH-MI training parameter:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {\rm lrER}{{{\rm D}}_{{\rm C}4-{\rm C}3}}={ R}{{{ P}}_{{\rm ERD},{\rm C}4}}-{ R}{{{ P}}_{{\rm ERD},{\rm C}3}}\nonumber \end{align} \tag{ 2 }$$

where the relative ERD power (RP_ERD) at C3 or C4 channel using the power of n trials across specific time period(P_n) according to equation defined as follows [29]:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{P}_{relax}}=\frac{1}{{{T}_{relax}}}\sum\limits_{n\in {{T}_{relax}}}{{{P}_{n}}}\nonumber \end{align} \tag{ 3 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {{P}_{task}}=\frac{1}{{{T}_{task}}}\sum\limits_{n\in {{T}_{task}}}{{{P}_{n}}}\nonumber \end{align} \tag{ 4 }$$

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle R{{P}_{{\rm ERD}}}=\frac{{{P}_{task}}-{{P}_{relax}}}{{{P}_{relax}}}\times 100\nonumber \end{align} \tag{ 5 }$$

where P_relax and P_task are the average power spectra during the rest period (T_relax) and the task period (T_task). To this end, the averaged event-related power changes in a time-frequency domain could be visualized by the event-related spectral perturbation (ERSP), which provides detailed information for ERD/ERS patterns of different tasks. In our study, an ERSP of n trials was calculated according to equation defined as follows:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {\rm ERSP}(\,f,t)=\frac{1}{n}\sum\limits_{k=1}^{n}{\left( {{F}_{k}}{{\left(\,f,t \right)}^{2}} \right)}\nonumber \end{align} \tag{ 6 }$$

where F_k(f, t) indicates the spectral estimation at frequency f  and time t for the kth trial. The ERSP (dB) was computed through short-time Fourier transform (STFT) with a 256 points Hanning-tapered window from EEGLAB [30]. For the baseline-normalized ERSP, the mean power changes in a baseline period were subtracted from each spectral estimation. For the comparison between previous and posterior control sessions, the lrERD powers were averaged across all of the n trials (n  =  40) for  −2 to 5 s. For the NFT sessions, the lrERD powers were averaged across real-time EEG data (n  =  1) for every 1 s.

For further analysis of previous and posterior control sessions, the ERSP patterns from two key EEG electrodes (C3 and C4) were analyzed from  −2 to 5 s between 5 and 35 Hz, referring to the computing method of baseline-normalized ERSP. The mean power changes in a baseline period (2 s before MI) were subtracted from each spectral estimation. We also calculated the salient absolute alpha-ERDs (alpha_1: 8–10 Hz, alpha_2: 11–13 Hz) and beta-ERDs (beta_1: 15–20 Hz, beta_2: 22–28 Hz) respectively of all EEG channels with its topographical distributions according to the ERSP, which could indicate the cortical electrophysiological activations. Moreover, to investigate the quantitative ERD patterns for LH-MI and RH-MI training, we also calculated the lrERD powers of previous and posterior control sessions across all subjects.

In this work, a FBCSP algorithm was used for EEG features extraction [31]. FBCSP utilizes a bank of N_f band-pass filters to divide the EEG signals into components with different frequency bands, and then, constructs the projection matrix W_i (i  =  1, 2, ..., N_f) using the CSP algorithm [32] and extracts spatial features for each EEG component respectively. Support vector machine (SVM) was used for pattern classification. SVM is insensitive to overfitting and suitable for small training sets. A 2-class(LH-MI versus RH-MI) classification using SVM was realized by the LIBSVM software package for MATLAB (MathWorks Inc., Natick, MA, USA), a freely-available library of SVM tools [33].The multi-channel EEG data in the epoch from 0.5 to 3.5 s after MI onset were extracted for classification [34]. Then a leave-one-out strategy was used to the whole dataset obtained from previous and posterior control sessions respectively. The final classification accuracy based on the results obtained from all testing sets was given by:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle ac{{c}_{cv}}=\frac{1}{n}\sum\limits_{k=1}^{n}{acc(k)}\nonumber \end{align} \tag{ 7 }$$

where acc(k) represents the result of the kth fold accuracy and n  =  40 for the leave-one-out strategy.

To probe the NFT effects quantitatively for MI-BCI performance, we analyzed the 2-class MI-BCI performance of the previous and posterior control sessions. Firstly, the same frequency bands were selected to extract features using FBCSP. In this study, a total of six frequency bands (8–12, 12–16, ..., 28–32 Hz) covering the frequency range of 8–32 Hz were selected to calculate the classification performance. Besides, to probe whether there would be an improvement in MI performance after the NFT when separate ERD features were used intuitively, we also calculated the classification accuracies of the previous and posterior control sessions using the spectral powers of above selected six frequency bands. We chose 70% accuracy as the threshold for acceptable performance [35].

The PDC based brain networks were calculated for FC analysis. PDC could denote direct or cascade flows and connections in multivariate dynamic systems [36]. Consider that Y(k)  =  [y ₁(k),...,y _N(k)]^T, 1  ⩽  k  ⩽  n, is the N jointly stationary time series that can be adequately represented by a multivariate autoregressive (MVAR) model of order p  as:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle Y(n)=\sum\limits_{r=1}^{\,p}{{{A}_{r}}Y(n-r)+E(n)}\nonumber \end{align} \tag{ 8 }$$

where p  is the maximum number of lagged observations in the model, and A_r $\in$ R^N×N is the coefficient matrix at time lag r, with the frequency domain of the Fourier transform form:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle A(\,f)=I-\sum\limits_{r=1}^{\,p}{{{A}_{r}}{{e}^{-i2\pi fr}}}\nonumber \end{align} \tag{ 9 }$$

where I is an identity matrix. PDC values are defined as:

$$\begin{align} \newcommand{\e}{{\rm e}} \displaystyle {\rm PD}{{{\rm C}}_{j\to i}}(\,f)=\left| {{A}_{ij}}(\,f) \right|/\sqrt{\sum\limits_{k}{{{\left| {{A}_{kj}}(\,f) \right|}^{2}}}}.\nonumber \end{align} \tag{ 10 }$$

Therefore, PDC shows a ratio between the outflows from channel j  to channel i. A large PDC value indicates that there is a direct transmission between given channels in a given direction, whereas values close to 0 describe a lack of such a relationship. According to the above method, we could estimate the MVAR model and PDC [37]. In this study, we implemented PDC using the regression toolbox in MATLAB 2013b. EEG signals of 35 electrodes chosen from 64 electrodes were used to calculate the effective PDC connectivity networks. These 35 electrodes are overlying central and motor related cortical areas, involving FZ, F1–F6, FCZ, FC1–FC6, CZ, C1–C6, CPZ, CP1–CP6, PZ, P1–P6.

The ERD patterns and classification accuracies between previous and posterior control sessions were analyzed by paired-samples t-test across all subjects. With every two parameters of each subject from previous and posterior control sessions as two paired-samples, these tests were made with SPSS 22.0 (IBM SPSS Inc., Chicago, IL, USA).

Figure 4(a) shows the averaged time-frequency maps of C3 and C4 electrodes for previous and posterior control sessions. The maps present a long-lasting alpha-ERD (8–13 Hz) and beta-ERD (14–29 Hz) from MI task onset, especially clear for posterior control sessions. Also, the averaged ERSP waves of absolute alpha-ERD and beta-ERD are proposed, obtaining a consistent phenomenon as shown in figure 4(b). The contralateral dominance could be observed for LH-MI and RH-MI task, which have an obviously enhancement for the posterior control sessions.

Zoom In Zoom Out Reset image size

Furthermore, to investigate spatial distributions of the cortical activations for previous and posterior control sessions, the averaged topographical distribution maps of absolute ERDs at four salient frequency bands (alpha_1: 8–10 Hz, alpha_2: 11–13 Hz, beta_1: 15–20 Hz and beta_2: 22–28 Hz) are also shown across all subjects in figure 5. Additionally, the H-value (0 or 1) maps of paired-samples t-test (p   <  0.01) present significant differences at different frequency bands especially at alpha_2 and beta_2 ERDs covering sensorimotor cortical areas. These results also indicate the contralateral dominance could be observed more clearly for posterior control sessions.

Zoom In Zoom Out Reset image size

To quantitatively verify the effectiveness of the proposed NFT, the salient four ERD bands (alpha_1: 8–10 Hz, alpha_2: 11–13 Hz, beta_1: 15–20 Hz and beta_2: 22–28 Hz) were also selected respectively to compare their lrERD powers between previous and posterior control sessions. Figure 6 shows the averaged lrERD (according to C3 and C4 channels) powers of LH-MI and RH-MI task across all subjects. The paired-samples t-test yielded significant differences between previous and posterior control sessions at different frequency bands (pre versus post at beta_1 of LH-MI: t(18)  =  2.934, p   =  0.009, d  =  0.67; beta_2 of LH-MI: t(18)  =  6.131, p   =  0.000, d  =  1.41; alpha_1 of RH-MI: t(18)  =  5.787, p   =  0.000, d  =  1.33; beta_1 of RH-MI: t(18)  =  6.597, p   =  0.000, d  =  1.51; beta_2 of RH-MI: t(18)  =  6.552, p   =  0.000, d  =  1.50). The lrERD patterns indicate the relatively consistent trend with the above time-frequency phenomena. Therefore, the proposed visual-haptic NFT could indeed improve lateralized cortical activations significantly during left or right hand MI.

Zoom In Zoom Out Reset image size

Firstly, we compared the 2-class MI-BCI (LH-MI versus RH-MI) classification performance between previous and posterior control sessions to verify the effectiveness of the proposed transient NFT. As utilizing FBCSP (a total of six frequency bands: 8–12, 12–16, ..., 28–32 Hz) and SVM algorithms, figure 7 shows relatively high classification accuracies across all subjects and their mean accuracies. The performance in posterior control sessions achieves an 8.7% improvement relative to previous control sessions, reaching 84.6% in mean classification accuracy. The paired-samples t-test yields significant difference (pre versus post: t(18)  =  −5.051, p   =  0.000, d  =  −1.16). In addition, the MI-BCI performance has a relatively consistent improvement trend and reach above 70% after the proposed visual-haptic NFT for most of the subjects.

Zoom In Zoom Out Reset image size

Secondly, to further verify whether the proposed NFT could improve the ERD features based MI-BCI performance, the classification results were also calculated from separate ERD features (the same six frequency bands: 8–12, 12–16, ..., 28–32 Hz) with SVM algorithm and compared between previous and posterior control sessions. Figure 8 illustrates the classification performance of all subjects and mean accuracies for previous and posterior control sessions. The results show that the classification performance for most of the subjects could only be classified below 70% in previous control sessions, while half of the subjects above 70% in the posterior control sessions. Besides, the paired-samples t-test yielded significant difference of classification accuracy between previous and posterior control sessions (pre versus post: t(18)  =  −4.032, p   =  0.001, d  =  −0.97). The performance in posterior control sessions achieves a 9.2% higher accuracy relative to previous control sessions, reaching 71.9% in mean classification accuracy. These results could further demonstrate the feasibility of the proposed BCI based visual-haptic NFT to improve classification performance during MI.

Zoom In Zoom Out Reset image size

Since the lrERD patterns at alpha (8–13 Hz) and beta (14–29 Hz) rhythm bands were as the NFT parameters to improve the MI-BCI performance, correlation analysis was employed to verify whether there would be relationships or consistent trends between the lrERD and classification performance. After rejection of the outliers according to the criterion of Box-whisker plot in SPSS, Shapiro-Wilk test also indicated that the data were normally distributed overall previous and posterior control sessions. Thus, a 2-tailed Pearson correlation test was employed to examine their relationships. As shown in figure 9, the correlation analysis yielded significant negative correlations between lrERD patterns and classification accuracies: alpha (8–13 Hz) lrERDs of LH-MI and accuracies (r  =  −0.376, p   =  0.020), RH-MI (r  =  −0.453, p   =  0.005), beta (14–29 Hz) lrERDs of LH-MI and accuracies (r  =  −0.369, p   =  0.025), RH-MI (r  =  −0.378, p   =  0.019). These results indicate that the relationships between lrERD patterns and MI-BCI performance could be quantified through above method and it might be a promising approach to utilize MI induced lrERDs of electrophysiology to evaluate the effectiveness of motor training.

Zoom In Zoom Out Reset image size

In order to reveal the dynamic information flows of FC over activated brain regions, we integrated the PDC values at two frequency bands (alpha: 8–13 Hz and beta: 14–29 Hz) and time interval during MI tasks before and after NFT. All of the PDC values were normalized by subtracting the mean value and dividing the result by the standard deviation along different sessions (previous and posterior control sessions). To indicate the causal interaction between channels more clearly, the amount of effective connectivity was limited by setting a threshold which was 60% of the maximum PDC value of each subject. Besides, the statistical significance (unilateral t-test, p   <  0.01) of non-zero PDC values was assessed by means of a bootstrap approach using phase randomization according to the Theiler's method [37]. It was performed across all subjects for previous and posterior control sessions respectively. Figure 10 shows the effective connectivity networks of LH-MI and RH-MI tasks. The outflows are showed by blue arrows with their bases at 'source electrodes' and the tips pointing toward 'target electrodes'. There are significant increases of effective connectivity for posterior control sessions. This consistent phenomenon could be found at several electrodes covering sensorimotor regions. Additionally, neighboring electrodes also show more causal connects during MI after NFT. These results verify the effectiveness of the proposed BCI based visual-haptic NFT to enhance the FC patterns and indicate that the brain topology networks of motor training could be quantified through the PDC method.

Zoom In Zoom Out Reset image size