Abstract
Brain-computer interface (BCI) systems translate the human neurophysiological activities into commands through EEG analysis. Improving the BCI performances leads to faster and easier use and less fatigue. In this study, we proposed a new prepossessing approach to increase the robustness of a steady-state visual evoked potential (SSVEP) based BCI. Inspiring from the known properties of the SSVEP frequency components, the goal was to enhance the signal quality by making it more convenient to be interpreted by the decision-making step. We first investigated the potential to detect the deteriorating periods based on the physiological properties of the SSVEP. The proposed system localizes the intervals which can obscure the SSVEP frequencies by a new algorithm founded on the processing and the analysis of the instantaneous phase. The piecewise linear regression allows a sampler comprehension of the phase signal. Then, these intervals are filtered by the moving average filter to enhance the SSVEP quality. Finally, the decision making is made by the canonical correlation analysis (CCA) algorithm. The results of experiments, using real EEG signals from five subjects, show that the proposed approach significantly increases the performances in terms of accuracy and information transfer rate by about 7.3% and 3.85 bits/min, respectively, in case of 2 s segment length. On the other hand, the spatial filtering methods of the literature weaken the system performances.
1. Introduction
The brain-computer interfaces (BCIs) are one of the last means to regain the touch with the surrounding when no neuromuscular pathway is available, such as in the case of people suffering from locked-in syndrome. They offer a new neurophysiological link between a person (generally severely handicapped) and a computer or any other machine and provide new options for communication and control [1, 2]. BCI has a variety of paradigms including P300, motor imagery, and steady-state visual evoked potential (SSVEP).
The SSVEP-based BCI is one of the most successful interfaces. In fact, it provides a high information transfer rate (ITR) compared to other BCIs with no need for training and less required electrodes [3, 4]. The SSVEP signal is generated in the occipital area of the brain when a subject gazes at a visual target flickering at a constant frequency beyond 1 Hz [5]. It is composed of near sinusoidal components at the stimulation frequency and its multiples also called harmonics. Thus, it has a closely constant spectral distribution centered on the fundamental frequency.
For an easy, natural, and convenient use of the BCI, we need to provide a system with the highest possible accuracy and ITR. Hence, the question arises of how to improve the performances of an SSVEP-based BCI. The answer to the question is quite straightforward. We need to ameliorate the quality of the SSVEP signal. In fact, the quality of the EEG is the key to the communication speed of BCI. An accurate target identification using data with lower signal-to-noise ratio (SNR) requires longer signal epochs leading to the decrease of the ITR. Thus, we need efficient and effective algorithms for preprocessing.
Several approaches and methods have been used in the preprocessing step of the SSVEP-based BCI to enhance the EEG quality in order to improve the decision making by ensuring that the system is not controlled by artifactual and noncerebral signals. One of the first adopted approaches is the time domain filtering used to eliminate frequency bands that may include noise and artifacts. For instance, in [6], a band-pass infinite impulse response (IIR) filter from 22 Hz to 48 Hz is used to retain the pertinent parts of the EEG signal. Moreover, a filter bank combining a set of band-pass filters was used to select the frequency bands of interest in order to estimate their spectral powers in the following processing steps [7]. Despite the simplicity and ease of implementation of the time domain filtering approach, its success requires that relevant and irrelevant signals be located in different frequency bands [8, 9]. In the case of overlapping, this filtering will not be able to reduce the noise level without attenuating the SSVEP frequencies. The spatial filtering approach has also gained a lot of attention. More specifically, the common averaging rereferencing (CAR) method has been used to reduce the levels of noise and artifacts injected by the reference electrode [7]. It is based on the fact that each event-related potential’s (ERP) component is positive over some areas of the head and negative over some other areas. Thus, if the head model is assumed to be spherical, the sum of all ERP’s components is zero (neutral) over the entirety of the head. Consequently, averaging all the EEG channels leads to the creation of a zero potential reference [10]. More popular methods are based on the blind source separation (BSS) paradigm. For example, the AMUSE method allows to estimate the independent components (ICs) with automatic ranking. The harmful components related to the noise and the artifacts will be removed, and then, a projection back to the scalp level creates the artifact-free signal [11]. In fact, the first and the last components will be removed. The first component containing the slowest brain activity is supposed to be related to electrooculogram (EOG) artifacts and other slow artifacts while the last component containing the fastest activity is supposed to be related to muscle artifacts and other high-frequency noise [12].The independent component analysis (ICA) also allows to decompose the EEG signal into different ICs. To differentiate between the relevant and the irrelevant components, Wang et al. [13] calculated the normalized amplitude at the stimulus frequency. The ICs with the highest normalized amplitude are considered to be SSVEP activities. The correlation analysis was also used to separate the ICs. For example, the ICs which have the highest correlation with the signal from AF3 electrode are considered as related to the EOG artifacts. On the other hand, the ICs which have the highest correlation with the signal from Oz electrode are considered as related to the SSVEP activities [14]. Multivariate data analysis (MVA) algorithms were also used in the SSVEP preprocessing step. For example, the minimum energy combination (MEC) is a spatial filtering method dedicated only to enhancing the SSVEP quality. It uses the principal component analysis (PCA) and the orthogonal projection to accentuate the SSVEP frequencies and attenuate the remaining frequencies [15]. Moreover, the partial least square (PLS) regression is another MVA method used to maximize the covariance between the EEG and SSVEP model signals [16]. The common spatial pattern (CSP) algorithm was also used to find a linear combination of the EEG data to increase the power of the stimulation frequency [17, 18]. It creates a divergence between the distributions of dataset from two different classes. The first class, called target subset, is obtained by filtering the epochs with band-pass filter centered on the stimulation frequency. The second class, called nontarget subset, is obtained by filtering the epochs with band-pass filter centered on each nontarget SSVEP frequency.
All these presented preprocessing methods use the time, frequency, and statistical properties of the noise and artifacts. Another promoting paradigm is the deployment of the physiological property. In fact, in case an artifact is identified by its shape, its removal becomes easier with a higher probability to conserve the useful components. Within this context, Zhang et al. [19] proposed an approach to eye blink artifacts removal. They start by modeling a generic skeleton of the artifacts with straight lines based on the ballistic properties of the eye blink generation, and then signals from the channel Fp1 are examined in order to fit the skeleton to each individual blink. Next, the detected artifacts are removed from all channels. Reducing the levels of noise and artifacts is essential for better performances; however, preserving the original signal shape with minimal distortion is fundamental in several applications. To this end, several works have been interested in the moving average filter (MA) [20, 21]. Besides its simplicity, its ability to maintain the characteristics of the curve allows to keep the significant information which leads to a better interpretation.
All the presented preprocessing methods filter the signal regardless of its quality level which can deteriorate it when no enhancement is needed. Also, most of them require a large number of channels. In this paper, we would like to propose a new method that segments the signal into two groups according to the quality level and then filters only the periods that require an enhancement. The proposed method should be able to work even with a small number of channels (one channel).The idea behind our method is to improve the performances by benefiting from the sought signal properties. Our proposed paradigm is motivated by the fact that the stimulation frequencies which reflect the shape of the sought signal are preknown. Thus, it is possible to use the physiological information about the SSVEP signal to control the preprocessing stage. We proposed a new preprocessing method based on the piecewise linear regression of the instantaneous phase to detect the deteriorating periods. Moreover, the moving average filter was integrated for the improvement of these periods. For decision making, we use the canonical correlation analysis (CCA) to detect the target stimulus. The proposed algorithm was evaluated using an open database for SSVEP-based BCI. This database was collected with a 4-target BCI from 5 subjects. For comparison, we use all the spatial filtering methods previously presented. Results show that the ability of our method to choose when to apply the filtering allows it to reach better performances.
The rest of this paper is organized as follows. Section 2 presents the used database. The proposed method and its blocks are detailed in Section 3. The experimental results and the discussion are given in Section 4, and then the work is concluded in Section 5.
2. SSVEP Database
In the present study, the system performance analysis is assessed using a freely available online database [22]. Five subjects, aged 23, 25, 31, 42, and 46, respectively, participated in the experiment. EEG signals were acquired from 16 channels (O2, AF3, AF4, P4, P3, F4, Fz, F3, FCz, Pz, C4, C3, CPz, Cz, Oz, and O1 according to the 10/20 international system) by means of the g. USBAmp EEG signal acquisition system at a sampling rate of 256 Hz. Each subject was asked to sit 1 m far from 4 green LEDs of a 1 cm diameter and to gaze once at each LED for 30 s. Each LED flickers with a constant frequency different from the other. The four frequencies are 5, 6, 7, and 8 Hz, respectively.
3. Proposed Method for SSVEP Quality Enhancement
The SSVEP signal is composed of different frequency components proportional to the stimulation (fundamental frequency) and its harmonics. Each component represents one scale of oscillation. The noise, the artifacts, and the background activities inject some additive deteriorating quantities to the time series samples. These alterations in the amplitude over a time period can cause some changes in the frequency components and lead to the appearance of new ones with non-SSVEP sources which can distort the final decision. The objective of our proposed method is to enhance the SSVEP signal by selecting the periods responsible for the decline of the results (e.g., decisions of the system) and then by applying a filtering on them to restore the maximum amount of the SSVEP components. The filtering will not be applied on the advantageous periods. Therefore, a higher probability to conserve the useful information will be ensured.
Figure 1 illustrates the proposed SSVEP enhancement method and its integration in the whole BCI system.
After acquiring the required amount of the EEG data from the scalp, the multichannel EEG signal is preprocessed by the new prepossessing method. The enhanced signals are then sent to the CCA method to find the frequency at which the subject gazes. The proposed enhancement method is composed of two blocks. The first block allows to detect the segments which can obscure the SSVEP frequencies by creating non-SSVEP components. The second block is the filtering used to reduce the level of these components.
Figure 2 details the different steps of the proposed SSVEP enhancement method.
The multichannel EEG signal is treated channel by channel. Each single channel is preprocessed separately by the proposed method, and then the enhanced multichannel EEG signal is compared with the reference signals by the CCA method.
The noise, the artifacts, and the background activities are totally or mostly random, and thus, they can weaken the SSVEP components or create new components only for short periods of time. Consequently, we can characterize the non-SSVEP component as an abrupt change in frequency which can last for a short time period. Conversely, an SSVEP component has a constant frequency lasting for a longer period. Obviously, more abrupt changes should be found if the levels of noise and artifacts are higher. Based on these assumptions, we can distinguish the useful periods from the harmful ones by the sharpness of the frequency changes and the lengths of the intervals during which the change is maintained. To simplify the estimation of these two criteria, we proposed to analyze the signal in the instantaneous phase domain. In fact, the frequency components are transformed into sawtooth signals in the phase representation. This technique allows to represent the problematic linearly.
Several methods have been proposed in the literature to estimate the instantaneous phase of a time signal. In this work, we used the Hilbert transform which is one of the most robust and most useful methods. Figure 3 shows an example of the instantaneous phase found using the Hilbert transform applied to an EEG signal acquired when the subject was gazing at a visual stimulation flickering at 5 Hz.
The abrupt changes in the frequency lead to sudden changes in the phase signal. Thus, the abrupt changes in the slopes of the instantaneous phase, as shown in Figure 3, are caused by the noise and artifacts injected in the EEG signal. The rapidity of these changes is due to the random nature of the noise and artifacts, which means that the occurrence of these changes is proportional to the SNR of the SSVEP.
To analyze the slopes of the instantaneous phase, we have used the piecewise linear regression (PLR) [23, 24]. In fact, using the PLR, the signal is split into several portions that have linear trends and a straight line is fitted to each interval. The difference between the line and the actual signal in each period should be below an accepted tolerance which depends on a threshold value that must be set according to the application. This threshold has an impact on the breakpoint positions and on the lengths of lines. The idea behind our method is to automatically specify the threshold value in such a way the length of each straight line indicates whether this part of the signal is relevant (e.g., reflects the SSVEP signal) or irrelevant (e.g., obscures the SSVEP signal). In order to guarantee the finest separation between the periods and since the noise and artifacts are abrupt and quick, the proposed algorithm sets the threshold value in a way that the irrelevant intervals are indicated by straight lines composed of only two samples. The adjacent irrelevant intervals form one irrelevant segment that will be filtered separately. For this task, we choose to use the MA filter. In fact, the MA filter is applicable even with a small number of samples. Moreover, it is optimal for reducing random noise while retaining a sharp step response. Thus, it can improve the trends of a noisy signal to make more significant correlation coefficients with the CCA method.
The following section explains in more detail the main steps of the proposed method.
3.1. Threshold Value Determination Algorithm
Figure 4 presents the different steps of the threshold value determination algorithm.
The threshold value determination algorithm is composed of three steps beginning by estimating the SNR of the SSVEP. In fact, as mentioned above, the level of noise and artifacts directly influences the occurrence of the abrupt changes and so the occurrence of breakpoints. Breakpoints are unknown points where two successive lines joined. They represent positions where slopes of functions with linear trends start changing more sharply than a predefined threshold. Our goal is to deduce the appropriate amount of breakpoints based on the SNR. The SNR estimation technique used in our proposed method is the one presented in the work of Kus et al. [25]. This technique measures the activity level at a given stimulation frequency with respect to the level of activity in adjacent frequencies according to the following equation:where n = 6, , and is obtained by evaluating the spectral power at frequency by means of periodogram with Tukey window with parameter equal to 0.1.
The power of the fundamental frequencies is the highest in the SSVEP signal followed by the first harmonics then the next harmonics. Thus, it is easier to analyze the fundamental frequencies than the harmonics. For this end, we chose to express the SNR of the SSVEP signal as the average of SNR values at fundamental frequencies:where is the set of used stimulation frequencies and is the number of stimuli.
After estimating the SNR of the SSVEP, this quantity is used by the second step to specify the adequate percentage of breakpoints (BP) according to the following equation:where and are two system parameters estimated using the genetic algorithm (GA).
The GA is one of the most performant optimization methods. It takes its advantages from the Darwinian evolution theory and genetics. The goal of the GA, in the present work, was to specify the optimal values of A and B that maximize the accuracy of the whole BCI system (calculated after the decision made by the CCA block). The GA starts with a random initial population of individuals and then analyzes in each generation the founded fitness values (accuracy values), and thereby, it evolves to the optimal solution through an iterative computation process [26, 27]. As a consequence, a system parameters optimization phase is required before the normal use of the BCI system. Once, the offline preparatory phase determines the values of A and B, the BCI system became ready to start the normal real-time use phase. It is important to notice that the value of differs from one segment to the other one which makes the system adaptive. During the system evaluation experiments, each time we picked trials from one subject as the test set and combined trials from all remaining subjects as the set needed to execute the GA.
In the next step (specify the appropriate threshold value), the proposed method explores all possible percentages of breakpoints and picks the closest one to the value fixed earlier by the GA and then specifies the threshold that will be used for the remaining analysis.
Algorithm 1 details the different instructions followed to get the appropriate threshold value. The highest breakpoints percentage (100%) is found when the number of breakpoints is equal to the half number of the EEG samples (when the threshold is minimal, 0.000001 in our case). This algorithm starts by verifying if the threshold is equal to the minimal or the maximal value; otherwise, it splits the threshold range into 20 intervals and then analyzes them one by one. This technique reduces the computational time.
|
3.2. Piecewise Linear Regression
The piecewise linear approximation aims to segment a time series of samples into straight line with the highest possible quality of approximation. One possible way to form this problem is to consider a threshold value specified by the user and representing the tolerated amount of error, where the goal is to produce the new representation in such a way the maximum error for any segment does not exceed the threshold. As is smaller than , this representation can be used to facilitate the analysis and the interpretation of the data. Several techniques have been proposed to create the breakpoints between the lines. The sliding windows algorithm is the most popular one due to its simplicity, its intuitiveness, and its compatibility with online applications. Using the sliding windows, each segment length increases from one sample to the next one while the error boundary was not reached. Once the error exceeds the threshold, a new breakpoint is created and the process will be repeated with the next segment.
The objective of the simple linear regression is to analyze the relation between two quantitative variables: one explanatory also called independent variable and one response also called dependent variable (conventionally, the and coordinates in a Cartesian coordinate system). The sought relation is expressed by a straight line that should fit as accurately as possible the points in the Cartesian coordinate to form a prediction of the new dependent variable as a function of the independent variables [23]. This linear function is given bywhere and are the parameters (coefficients) of the model indicating, respectively, the intercept term with the y-axis and the slope of the regression line. These parameters are found by minimizing the sum of squared deviations between the dependent variables and the corresponding values from the regression line (also called the sum of squares error).
The resolution of this optimization problem is given by
The solution of this problem is as follows:where and are the average of and , respectively, is the variance of , and is the covariance of and .
In piecewise linear regression, the simple linear regression is used to create the straight lines and the sum of squares error is used as threshold to create the breakpoints [24].
3.3. Moving Average Filter
The MA filter is a variety of finite impulse response (FIR) filter which smooths the signal to remove the random variations. It is useful to improve the time domain signal information. The current output sample is calculated by averaging the previous samples where is the length or the order of the filter. A temporal moving window is applied to the signal to calculate the successive local averages. The output of the MA filter with order is expressed as follows:where is the input original signal.
3.4. Standard Canonical Correlation Analysis-Based Method
After the signal preprocessing step, a decision-making step is needed to provide the system output command. One of the most efficient methods to recognize the frequency of the target at which the subject gazes is the CCA-based method. This method has been used in many works due to its simplicity, its reliability, and its high accuracy. Firstly introduced by Lin et al. [6, 28], this method integrates the CCA to quantify the similarity between the multichannel EEG and models of optimal SSVEP signals. Therefore, correlation coefficients are estimated between the multichannel EEG signals and reference signals ( is the number of stimuli). Each reference is formed from several time series corresponding to one stimulation frequency and its harmonics.
As shown in the previous equation, consists of a set of sine and cosine signals where and are two constants representing the i-th stimulation frequency and the number of harmonics taken into account, respectively, and is the count up timer. Hence, where is the number of time points.
The stimulation frequency that produces the maximum coefficient value is selected and identified as the target frequency.
4. Results
SSVEP-based BCIs have gained a research interest due to their advantages: portability and ease of use, higher SNR as well as faster ITR. To respect these pros, we chose to limit the number of channels to 6 and use only the nearest channels to the occipital area (O2, Oz, O1, P3, Pz, and P4). Also, the length of the EEG segments is reduced to 2 s. Each 30 s of multichannel EEG signal is split into nonoverlapped 2 s segments and then each segment is analyzed apart.
4.1. Detection of Irrelevant Intervals
As presented in the previous section, the first block of our proposed method allows to detect the irrelevant intervals which deteriorate the system performances. Figure 5 shows an example of the results of this block.
(a)
(b)
(c)
Before this experiment, the system has already executed the GA in a preparation phase in order to determine the best values of the parameters A and B.
Figure 5(a) shows the results of the PLR applied on an instantaneous phase signal estimated from an original EEG segment. Circles on the curve indicate the beginnings of straight lines, while crosses indicate the endings. These beginnings and endings can be considered as features extracted from the phase signal as they create a new representation easier to be interpreted. In fact, just calculating the difference between the two points allows to know if the interval they bound is relevant or irrelevant.
Figure 5(b) illustrates, for the same data, the results of the irrelevant period recognition step. Results show that the short lengths of the PLR lines allow to distinguish between the relevant and the irrelevant intervals without latency regardless of their lengths. The type of the interval depends only on the instants of the changes in the slope.
Figure 5(c) illustrates the instantaneous phase estimated from the preprocessed signal found with our proposed algorithm. From the trend of the curve, it is clear that the phase signal is more stable and the slopes are steadier with fewer changes. This result proves that there are fewer changes in the frequency components of the SSVEP signal.
4.2. System Performance Evaluation
Two metrics have been proposed in the literature to evaluate the BCI systems which are the accuracy and the ITR. The accuracy rate is calculated by dividing the number of correct decisions by the number of trials. The ITR reflects the speed of the system and takes into account the values of accuracy and execution time. Wolpaw et al. [29] proposed a way to express the ITR as given in the following equation:where is given in bits/symbol, is the number of targets, and is the classification accuracy. A more common ITR derived from the first one is given as follows [30]:where is given in bits/min and CTI is the command transfer interval representing the average time needed to convert a brain feature activity into a command.
The proposed method was compared with the spatial filtering methods presented in the Introduction. Experiments were conducted using a PC computer with 8 Go RAM and an Intel Core i7, 4.6 GHz processor. The software codes were implemented using MATLAB R2016a with the EEGLAB plugin for ICA methods. All methods of EEG signal preprocessing with the description of the obtained signals are presented in Table 1. These signals were used in the CCA step for decision making.
Methods that belong to the BSS family use all the 16 channels to create the ICs, and then, only the 6 previously mentioned channels are used in the CCA step.
Table 2 shows the system accuracy rate found for the five subjects using all preprocessing methods. Segments of 2 s length were used in this experiment.
Results show that the preprocessing methods of the literature ameliorate the SSVEP quality in some cases. However, the majority of the results were declined. In fact, in case of CAR methods, errors can be due to a limited electrode density, incomplete electrodes cover, as only the top of the head can be covered by electrodes, and the inaccurate modeling of the head by a simple sphere. Moreover, the BSS methods are not able to distinguish the EEG segments that are contaminated by artifacts from the artifact-free segments. Thus, removing ICs from every decomposition result can lead to the removal of important ICs. Also, an IC can include a mixture of components related to artifacts and other components related to SSVEP. Thus, the rejection of entire ICs in a systematic way can lead to the loss of important information included in these ICs. CSP and PLS methods allow to increase the average accuracy; however, their performances are not stable for all subjects. Both of them largely decline the accuracy using EEG of subject 2.
It is obvious from the obtained results that the proposed method performs better than the standard CCA and the spatial filtering ones. Indeed, it increases the average system accuracy rate by about 7.3% and the average ITR by about 3.85 bits/min. This proves the effectiveness of the proposed method in the enhancement of the SSVEP frequency component and the reduction of the level of noise and artifacts. Moreover, the accuracy rises for all subjects without exception even when the precision is low, as in case of the subject 3 where the accuracy is below 65%; hence, we can see that our method remains stable even facing a high level of noise and artifact.
Table 3 shows a comparison of the average time needed for the processing of segments of 2 seconds length.
Results show that our proposed method performs better than the ICA-based methods. In fact, the estimation of the mixing matrix and its inverse is very time-consuming. However, CSP, PLS, and CAR methods are far better than our method in terms of computation time. The average calculation time was about 0.28 s for the 2 s segment length. This latency between the end of the segment acquisition and the decision making can weaken the natural use of the BCI. Nevertheless, to improve the computational costs, a more performant computer with more powerful resources can be used. Also, the codes can be optimized with software solutions.
Table 4 shows the system performances for different segment lengths. Since the user cannot use the BCI during the calculation time, this period was added to the segment length to calculate the ITR.
Results show that the proposed method improves the decision of the CCA method even after reducing the length of the segments, which proves its stability and robustness. Indeed, as the system can select irrelevant intervals of reduced length and apply filtering on them, it remains functional with precision regardless of the length of the segments. Segments with 2 s length allow to reach the best accuracy and ITR values making the system closer to an online use outside the laboratory. However, to ensure the real-time use, the average calculation time needed to make a decision must be improved.
5. Conclusion
SSVEP-based BCI is one of the most promoting systems to bring BCIs outside the laboratory for everyday use. The high accuracy and ITR make it a good choice for control and communication applications. In the present study, we proposed a novel preprocessing method to enhance the SSVEP quality. The proposed approach is a real-time signal processing chain composed of the Hilbert transform, the piecewise linear regression, the Tukey periodogram estimator, and the moving average filter with the use of the genetic algorithm for parameter optimization. Results show that the proposed system boosts the CCA-standard system and reaches an average accuracy rate of 84.29% and an average ITR of 29.57 bits/min for a segment length of 2 s.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was carried out under the MOBIDOC scheme, funded by the European Union (EU) through the EMORI program and managed by the National Agency for the Promotion of Scientific Research (ANPR).