Materials and Methods - 適用於中風復健之即時低複雜度腦機介面開發

Materials

The datasets include our collected datasets from Taoyuan Chang Gung Hospital with 4 normal subjects, 2 stroke subjects, and BCI Competition IV dataset 1_b, 1_f, and 1_g.[4][18]

Our Dataset Participants

The normal 4 subjects were 2 males, 2 females, with age ranged from 23 to 33. All participants labeled as TP, WR, WH, YC are right-handed, and are the research group members of our own laboratory. All of them are the first time users of BCI.

The first stroke subject is Mr. Cheng, aged 66, with his left brain injured. The second stroke subject is Mr. Lee, aged 46, with his right brain injured. These two patients have been rehabilitating for a while, and are in the late phase of rehabilitation treatment. The subjects are not severely damaged; both can perform easy instructions by their own will, but not fluently. All these subjects are the first time users of BCI.

Table 3-1: Summary of subjects

Experimental Paradigm

Figure 3-1 shows the timing scheme of the experiment. The 0~4 seconds is the rest period, with blank screen. Our rest period is longer, so that it has the equal length of movement period, for the convenience of analysis tool. Then at the 4s, the screen shows a fixation cross to hint the subjects to concentrate for the coming cue. Finally, the 5~9s is the movement period, the screen will show an arrow pointing left or right with equal chance.

For each subject, we collect 3 to 4 runs, with 80 trials (normal subjects), or 60 trials (stroke subjects) each run. The movement tasks in our test are listed as below.

(1) FH_fast: fast finger tapping, about 4 times per second.

(2) FH_comf: finger tapping with the comfortable speed, about 1~2 times per second.

(3) GR: grasp, we assume that grasp involves more muscle movement, and larger movement than finger tapping, it might result in stronger EEG waveform.

(4) MI: motion imagery, imagine raising your hand.

(5) MIGR: motion imagery, imagine hand grasping.

(6) TMP_fast: thumb press, this is for the stroke subject, because he had difficulty Figure 3-1: Experiment Paradigm

tapping his finger, we modify the movement from finger to thumb.

(7) PF: upper limb pushing forward, this is also for the stroke patient, because there is a push forward action in the rehabilitation movement, also the doctor said that pushing movement have stronger EEG waveform, we assume it is easier for the movement detection. The summary of subject’s task are listed in Table 3-2.

The task for each subject is slightly different though the whole collection series, we keep modifying the task movement to get better result, and the collection of one subject is finished within a day.

Table 3-2: Summary of subjects' task

Experimental Setup

The EEG was collected by the Mitsar EEG-201 amplifier (http://www.mitsar-medical.com/eeg-machine/eeg-amplifier-compare/) with 19 channels showed in Figure 2-3. The sampling rate is at 250Hz, and the recorded EEG waveform is passed to the compatible software, called WinEEG. The participants were asked to sit on a chair toward the screen and to fix their body in a room with weak light. The instruction were presented by visual cue from the computer screen with the software, called PSYTASK.

Also, the subjects were instructed to minimize eye blinking and to avoid physical movement throughout the EEG recording progress.

TP WR WH YC mrc mrl

BCI Competition IV Dataset Participants

The BCI competition IV dataset 1_b, 1_f, and 1_g were provided by Berlin BCI group [18], the reference link about the detail information is at:

http://www.bbci.de/competition/iv/. The subjects were all healthy.

Experimental Paradigm

The BCI competition IV datasets, there are two kinds set of data. One is for the calibration (training), and one is for evaluation (testing). For calibration datasets, in the first two runs, arrows pointing left, right, or down were presented as visual cues on a computer screen, with each trial 8 seconds. First, the screen will show 2s fixation cross to inform the subjects to be attention. Then, cues were displayed for a period of 4s during which the subject was instructed to perform the cued motor imagery task. Finally, 2s of blank screen shown in the center of the screen. Then, for the evaluation datasets, there are 4 runs, the motor imagery tasks were cued by soft acoustic stimuli (words left, right, and foot) for periods of varying length between 1.5 and 8 seconds. The end of the

Figure 3-2: BCI Competition IV dataset 1 paradigm [18]

motor imagery period was indicated by the word stop. Intermitting periods had also a varying duration of 1.5 to 8s. Like shown in Figure 3-2.[18]

Experimental Setup

The BCI Competition IV datasets are recorded by 64 channels, 1000Hz sampling rate, and 2 classes (+idle state).

Methods

After the EEG collection by the WinEEG software, we transfer the EEG signal to MATLAB platform though the format called EDF plus. The reason that we analyze the waveform on MATLAB platform is that we are familiar with the toolbox called

EEGLAB (http://sccn.ucsd.edu/eeglab/)[6] and BCILAB

(http://sccn.ucsd.edu/wiki/BCILAB)[7] these two toolbox are based on MATLAB, and we have more transparency and plasticity to modify the code we want, and subtly perform the instruction we gave.

From EDF format to MATLAB, first we cut the data out, then use the EEGLAB to import new data, then we add the labels of each channel, and use the default channel locations of EEGLAB. Last, we load the event types (left / right) and latencies through .txt file exported from WinEEG. After all these steps, we save it as .set file of EEGLAB for the future processing by BCILAB toolbox.

Then we introduce the BCILAB toolbox to do the process. BCILAB is a powerful toolbox based on EEGLAB, and it is also compatible with MATLAB. It comprise all the key components for BCI analyzing, like Figure 3-3. All the components are constructed in module form, also with GUI, and each module has different parameters.

We can choose different modules combine together, and define the parameters to fulfill the BCI processing that we want. BCILAB also combines some machine learning module function to classify the prediction, and finally, it has many ways (mean square

error, kullback-leibler divergence, negative log-likelihood …) to evaluate the performance.

Preprocessing

The purpose of preprocessing is to make the raw EEG data more clear, and reject the artifacts (eye blink, background noise …). The most common way includes band-pass filtering, independent component analysis (ICA) [1][8][9], principle component analysis (PCA)[8], or the simplest way is to manual removal of fuzzy EEG raw data, which needs some experience, simple but lack of validity.

However, we don’t do any preprocessing in our methods, and there are several reasons:

1. For the band pass filtering, in section 3.2.2, we will mention our approach: FBCSP (filter bank common spatial pattern), it already separates the raw EEG data to different frequency bands, so there is no need to do the preprocessing of band pass filtering.

2. For ICA and PCA, the algorithm transform the raw EEG data into the components (with the same number of channel numbers), and we have to manually reject the “noise”

component, or manually select just the C3/C4 component we want, and in some case it Figure 3-3: BCI flow diagram

is hard and has a lot uncertainty to pick up the right component we want; in addition, although there are some algorithms that can select the component automatically, but they require too much computation resources, which is not suitable for the purpose of this thesis. Also, the EEG data varies with time; hence, the ICA varies with time greatly.

It takes quite amount of time to keep calculating the ICA with the updating income of data while in the online application. Due to the above reasons, we don’t do the ICA/PCA in the preprocessing part, and we still get quite satisfactory results.

Feature Extraction

For the training part, I choose the first 40% of the data to become the training data, and the rest 60% are for the testing data. The target marker for training comes from the event (0 represents the rest, 1 represents the movement) that labels on the data. The rest event is at 0sec, and the movement event is at the 5sec as shown in the Figure 3-1, so we extract the epoch from 0.5sec to 3.5sec for each marker event.

Then we choose a feature extraction approach to transform the EEG data to the most suitable feature for further classification by machine learning. There are many approaches in BCILAB, the most common approaches are Log-Bandpower[15], Common Spatial Pattern (CSP)[3][10][11], and Filter-Bank Common Spatial Pattern (FBCSP)[5][12]-[14].

After choosing the approach, we combine the classification to train the model, then we can use the model to test the testing data, and get the accuracy. The following Table 3-3 shows the result. All the three approaches are with 19 channels, LDA classifier, 0.5~3.5 sec epoch extracted. For Log Power and CSP approaches, the algorithm filtered the 7~30Hz, and FBCSP filtered 0.5~3; 4~7; 8~12; 13~30; 31~42Hz by default parameters.

Table 3-3: Accuracy analysis respect to different approaches Approach

mcr_TMP_fast 65.28% 79.17% 79.17%

mrl_FH 88.89% 94.44% 93.06%

mrl_GR 75.00% 88.89% 84.72%

mrl_MIGR 93.06% 98.61% 98.61%

mrl_PF 87.50% 91.67% 94.44%

Average 76.70% 82.88% 86.95%

As we can see from the Table 3-3, FBCSP outperforms the rest algorithms greatly;

however, the FBCSP is the 5 more times the computation and calculation time of CSP algorithm because of the 5 times more frequency band to compute. Therefore, in order to cost down the computation time and complexity to reach our purpose and goal of immediate online classification with almost 80% accuracy rate, we try out to modified 1. Channels selection 2. Frequency bands selection parameter to fit our requirement.

From Figure 3-4, we can see that CSP and FBCSP both are good approaches for stroke patients, but FBCSP performs better for normal patients, this is because of that normal subjects have more clear frequency bands than stroke subjects, we will mention this furthermore in 3.2.4.

Figure 3-4: Feature extraction accuracy analysis

Channel Selection

(a) (b)

In Figure 3-5, provided by [16], (a) is the normal subject’s right hand grasping fMRI plot, while (b) is the stroke patient subject’s right hand grasping fMRI plot. (c) is the normal subject’s elbow flexion fMRI plot; on the other hand, (d) is the stroke patient subject’s elbow flexion fMRI plot. We can clearly see that for normal subjects,

Figure 3-5: Normal vs. Patient fMRI plot for right hand grasping & elbow flexion

the involved area in brain concentrate around the Cz, C3 (right hand movement, triggered by contra-lateral brain). However, for stroke subjects, the response area spread all over the brain, but still includes C3, Cz, C4, and even Fz. Therefore, according to Figure 3-5, we try to cost down the channel computation in our FBCSP algorithm, but still maintaining the accuracy.

Table 3-4 is the accuracy analysis respect to different channel number, with FBCSP algorithm, training epoch extracted from each event’s 0.5~3.5 sec. 9 channels refer to central 9 channels (F3, Fz, F4, C3, Cz, C4, P3, Pz, P4), and the 6 channels refer to front-central 6 channels (F3, Fz, F4, C3, Cz, C4). The interesting result shows that while we choose only C3, Cz, C4, the accuracy still maintains over 80%, but with great channel number reduction. Also proves the results in Figure 3-5 that the activation area involved while hand motion or MI (motion imagery) are mostly in C3, Cz, C4.

From Figure 3-6, we also tried C3, Fz, C4 these 3 channels to analyze, and we found that comparing to C3, Cz, C4, the normal subjects (TP, WR, WH, YC) showed great progress in the accuracy, but the dataset mrc_GR and mcr_TMP_fast drop dramatically. Due to the above reason, and we will show in the next chapter that normal subject with Fz channel in training will reduce the response time of online classification.

Hence, we finally decided to use 4 channels (Fz, C3, Cz, C4) instead of 3 channels (C3, Cz, C4 or C3, Fz, C4) to train the data for the convincing results.

Table 3-4: Accuracy analysis respect to different channel numbers Channel

Figure 3-6: Channel selection accuracy analysis

Frequency Band Selection

(a)

(b)

(c)

(d)

Figure 3-7: Event Related Spectral Power of different subjects

Figure 3-7 shows the Event Related Spectral Power (ERSP) of different subjects, the x-axis of plot is time from 0~10,000(ms), the y-axis is the frequency from 0~50(Hz), and the meter aside is the power of that channel. Plot (a), (b) are the representative subjects of normal subjects with different tasks shown in row, and with different hand movement’s C3/C4 channel shown in column. We clearly see that, at time 4~10 sec, there are strong Event Related De-synchronization (ERD, with power decrease) in 8~12Hz (μband) marked by the red circle, which strongly proves that while motion or even motion imagery (MI) there are strong ERD in C3/C4’sμband.

Plot (c), (d) show the two stroke patient subjects’ ERSP. The energy band scatters at all frequencies, not clear as normal subjects’, but still, we can see that there are strongerμband ERD at the un-damaged brain area no matter what hand moves. In plot (c), subject mrc’s damaged area is at left hemisphere of brain, so it shows C4 (right hemisphere of brain) has stronger ERD; likewise, plot (d), subject mrl’s damaged area is at right hemisphere of brain, so it shows C3 (left hemisphere of brain) has stronger ERD.

The default frequency bands of BCILAB’s FBCSP algorithm are 0.5~3Hz; 4~7 Hz; 8~12 Hz; 13~30 Hz; 31~42Hz, these 5 bands, including the μ rhythm(8~12Hz), and β(13~30Hz), which are known for the most related bands while motion &MI.

According Figure 3-7 above, we select different bands to analyze, to see whether we can have a reduction of these frequency bands or not. We have several combination of bands listed below:

(A): 0.5~3;4~7;8~12;13~30;31~42Hz, default frequency bands of FBCSP (B): 0.5~3;4~7;8~12;13~30Hz, exclude the 31~42Hz

(C): 4~7;8~12;13~30;31~42Hz, exclude the 0.5~3Hz

(D): 4~7;8~12;13~30Hz , exclude both lowest & highest frequency bands (E): 8~12;13~30;31~42Hz, exclude the lowest 2 bands

(F): 0.5~3;8~12;13~30Hz, lowest band &μ, β bands (G): 8~12;13~30Hz , μ, β bands

(H):4~8;8~12;12~16… 36~40Hz, four frequencies bands each, from 4~40Hz

The result of testing these combination of bands are listed in Table 3-5 and Figure 3-8, with FBCSP algorithm, 0.5~3.5 seconds epoch extracted, 4 channels (Fz, C3, Cz, C4) trained for the model, and LDA for the classifier. The result shows that the combination of (D) had best result in normal subjects, and quite satisfying result in stroke subjects. The result was not surprising, because we can see from Figure 3-7 that exclude the 8~12Hz, there are still some ERD in the 4~7 Hz. Thus, finally we decide to choose 4~7Hz, 8~12Hz, 13~30Hz to be our training frequency bands.

Figure 3-8: Frequency band selection accuracy analysis

Table 3-5: Accuracy analysis respect to different frequency band

Frequency

Summary

After all the analysis above, we find few parameters listed below.

1. Algorithm: FBCSP

2. Channels selection: Fz, C3, Cz, C4 (4 channls) 3. Frequency band selection: 4~7Hz, 8~12Hz, 13~30Hz

These parameters are to fit our application: fast, less computation, but still convincing accuracy results shown in Figure 3-9(the complexity is compared to 19 channels 5 frequency bands FBCSP), for online instant react BCI system. Furthermore, we hope to develop this BCI system into a portable, wearing device by the means of ASIC IC design. First, the epoch is extracted respect to each event at 0.5~3.5 seconds using the 4 channels (Fz, C3, Cz, C4) EEG data, and the prior 40% data are used for training, the rest 60% are for the testing. Then we select the FBCSP algorithm, with frequency bands cut at 4~7; 8~12; 13~30 Hz to be our feature extraction approach.

Finally, the features are fed into the LDA classifier to get the training model. The BCILAB uses the module called bci_train to train the data by the approach we selected.

After the training, the model compares the event 0.5~3.5 seconds’ data that labeled by the event marker to classifies this epoch’s type (rest or movement), gives it the percentage to rest event or motion event, and uses 5-fold cross validation and mean square error (MSE) to evaluate the performance. Then, we can use the rest 60% data for the testing, same as evaluating the training data, getting the results of mean square error. The results of training and testing are listed below in Table 3-6.

We see that there are still great variation among the subjects, some tasks have better performance for some subjects, while others don’t, and each kind of task gets different performances too, shown in Figure 3-10 (the task accuracy are the average of all subjects who have done that kind of task). Nonetheless, our total average result had

over 80% of accuracy both in normal subjects and stroke subjects. From [17], we can know that accuracy from 43~58% are accuracy by chance, which means our results are not at chance, and convincing.

Figure 3-10: Tasks summary analysis Figure 3-9: Performance analysis of final result

Table 3-6: Final result of datasets training & testing performance training testing

Online BCI Implementation

Motivation

Our research topic is a project that cooperate with Chang Gung Hospital. We wish to develop a BCI device to help stroke patients’ rehabilitation by the means of the online instant BCI response triggering the robot arm that helps the stroke patient to move their hand, enhancing the neural network communication from brain to limb. Hence, the online BCI classification is the crucial part to fulfill our goal. Different from the previous chapter that does the offline analysis, in this chapter, we want to develop an online analysis method with quick response, high accuracy, based on the parameters and algorithm of the pervious chapter.

Methods

The approach is the same as that in Chapter 3, we choose the following parameters:

1. Algorithm: FBCSP

2. Channels selection: Fz, C3, Cz, C4 (4 channls) 3. Frequency band selection: 4~7Hz, 8~12Hz, 13~30Hz

Then we can train the model by the module, bci_train in BCILAB as described in the Chapter 3. After that we take the rest 60% data and the training model to another BCILAB module called onl_simulate. This module has the parameter “sampling rate”.

It will use this sampling rate to update the simulate result with the same length of training model. We choose 20Hz to be our sampling rate in our approach. The training phase paradigm is shown in Figure 3-1, and the online simulated paradigm is shown in Figure 4-1.

After the onl_simulate module process, we get a vector of cognitive state at various specified time points, called predictions, p(i), which is the simulated result with respect to time. p(i) is also the percentage of movement intention: 1 refers to movement, and 0 refers to rest, like the prediction line in the top picture of Figure 4-3. However, the prediction results are a continuous, non-smooth line with many glitches, which will result in abrupt movement for the robot arm. To avoid this, we will transform it to a discrete (just 0 and 1) smooth line like the modified classification result showed in the bottom picture of Figure 4-3, to better control the robot arm.

Figure 4-2: Flow diagram of online BCI Figure 4-1: Online simulate paradigm

The whole algorithm for smooth movement is listed below, also the flow diagram is shown in Figure 4-2. First, we set a variable, called “threshold(i)”, which is the combination of constant value, 0.4, and the root mean square of previous 10 points (in our 20Hz sampling rate case, 10 points refers to 0.5 seconds data). As we set the threshold(i) to be the combination of previous data and a constant, the threshold(i) can be updated with time, which gains more flexibility and accommodation to conquer the variation of data though time.

threshold(i) = 0.4 × 80% + √^∑^𝑖^𝑖−9₁₀^𝑝(𝑖)²× 20% (4-1)

Then we define a counter, “C(i)”, with the range from 0~18, that counts how many times the p(i) is larger than threshold(i), or smaller than threshold(i). We accumulate the number that p(i) is larger than threshold(i), in such way, the glitches in p(i) can be eliminated, and therefore, we can get a smoother waveform like the middle plot in Figure 4-3. This flow is shown as below.

在文檔中適用於中風復健之即時低複雜度腦機介面開發 (頁 21-77)