Single channel approach for filtering electroencephalographic signals strongly contaminated with facial electromyography

Abstract

Eliminating facial electromyographic (EMG) signal from the electroencephalogram (EEG) is crucial for the accuracy of applications such as brain computer interfaces (BCIs) and brain functionality measurement. Facial electromyography typically corrupts the electroencephalogram. Although it is possible to find in the literature a number of multi-channel approaches for filtering corrupted EEG, studies employing single-channel approaches are scarce. In this context, this study proposed a single-channel method for attenuating facial EMG noise from contaminated EEG. The architecture of the method allows for the evaluation and incorporation of multiple decomposition and adaptive filtering techniques. The decomposition method was responsible for generating EEG or EMG reference signals for the adaptive filtering stage. In this study, the decomposition techniques CiSSA, EMD, EEMD, EMD-PCA, SSA, and Wavelet were evaluated. The adaptive filtering methods RLS, Wiener, LMS, and NLMS were investigated. A time and frequency domain set of features were estimated from experimental signals to evaluate the performance of the single channel method. This set of characteristics permitted the characterization of the contamination of distinct facial muscles, namely Masseter, Frontalis, Zygomatic, Orbicularis Oris, and Orbicularis Oculi. Data were collected from ten healthy subjects executing an experimental protocol that introduced the necessary variability to evaluate the filtering performance. The largest level of contamination was produced by the Masseter muscle, as determined by statistical analysis of the set of features and visualization of topological maps. Regarding the decomposition method, the SSA method allowed for the generation of more suitable reference signals, whereas the RLS and NLMS methods were more suitable when the reference signal was derived from the EEG. In addition, the LMS and RLS methods were more appropriate when the reference signal was the EMG. This study has a number of practical implications, including the use of filtering techniques to reduce EEG contamination caused by the activation of facial muscles required by distinct types of studies. All the developed code, including examples, is available to facilitate a more accurate reproduction and improvement of the results of this study.

Keywords: EEG; EMG; adaptive filtering; facial electromyography; signal decomposition.

Figure 1

Participants were instructed to make a variety of facial expressions by activating muscles whose electrical activity corrupts the electroencephalogram. The facial expressions were performed with both open and closed eyes. In the neutral condition (A) there was no muscular contraction, whereas in the other conditions the following muscles were activated: Frontalis (B), Masseter (C), Orbicularis Oculi (D), Zygomatic (E), and Orbicularis Oris (F).

Figure 2

Each muscle was contracted 15 times in one of three timing patterns: long (3 s), medium (1 s), and short (0.5 s). Each pattern of contraction was repeated five times randomly. A 2-s neutral period followed each contraction. This protocol was executed with open and closed eyes.

Figure 3

Filtering EEG signals corrupted by EMG via a series of steps. Linear and non-linear trends are eliminated, as well as outliers. The signal is then decomposed using one of the methods outlined, and the resulting components are thresholded. The thresholding process requires the identification of noise periods in the signal, which are provided by a binary signal generated by an EMG burst detector (Figure 4). Once the components have been thresholded, the filtered signal is reconstructed, producing a reference signal, i.e., EEG or EMG reference signal, which can be used as a reference for one of the indicated adaptive filters. Various characteristics are estimated to evaluate the filtering process at distinct stages. Note that when the method EMD-PCA is used it is not necessary to execute the soft-thresholding stage.

Figure 4

Sequence of required steps for detecting EMG bursts. First, the input signal is preprocessed by removing linear and nonlinear trends, and then the resulting signal is decomposed using EMD. The estimated components are soft-thresholded using a priori knowledge of the signal's noise level. Signal filtration is achieved by summing the thresholded components. The EMG envelope is determined by estimating the signal's energy, and bursts are detected using a threshold. As a result of this step, a binary signal is generated in which low levels indicate noise and high levels indicate EMG activity.

Figure 5

Set of time (GL, GH, GXin, and GXout) and frequency (p^r, p^f, $f_{med}^r$, and $f_{med}^f$) domain features used to evaluate the performance of distinct methods for adaptive filtering EEG corrupted by facial EMG. The features compare the signals before and after adaptive filtering and consider the regions in which there is the presence and absence of muscle activity.

Figure 6

Overview of the analysis required for characterization of contamination caused by different facial muscles, as well as a comparison of the performance of decomposition methods based on the features GH and GL. The analysis takes into account data grouping by participants, muscles, EEG sensors, and filtering methods. For each group, a similarity measure based on the normalized Euclidian distance can be estimated between a pair of vectors representing GH and GL estimates for varying a parameter used in the soft-thresholding of the signal components. The similarity metrics are used to create spatial brain maps that depict the contamination of EMG levels at various areas. Statistical analyses are carried out for GH, GL, and similarity measures.

Figure 7

Examples of typical GL and GH feature vectors obtained for high (Fp2–F8) and low (O1-O2) levels of EMG contamination on EEG signals. By decomposing the raw EEG signals with EMD, the reference signals were obtained. The contamination level of the EEG signal can be captured by the distance between the feature vectors, in the sense that an increased distance is related to a lower signal to noise ratio. The examples demonstrate typical collected signals for the open and closed eyes scenarios.

Figure 8

Typical EMG and EEG signals collected during the experimental trials. The EMG signals from the left and right Frontalis are shown. These signals were filtered to remove linear and non-linear trends. The EMG bursts were detected and then the binary signals oscillating from two levels were generated. Simultaneously collected EEG signals are shown for Fp2–F8 (high contamination) and O1-O2 (low contamination) locations. The binary signals are placed over the EEG signals for the indication of the periods in which there was EMG contamination.

Figure 9

Using the international 10-20 system, the topological maps illustrate how the studied muscles contaminate distinct brain regions. Lighter colors represent the most contaminated locations, whereas darker colors denote the least contaminated areas. The presented results are for subjects from 1 to 5. The colors represent the similarity measure between the GL and GH features.

Figure 10

Topological maps for subjects from 6 to 10.

Figure 11

Typical GL and GH feature vectors estimated using different decomposition techniques for Subject 1. Each plot consists of six vector pairs, one pair for each method. The outcomes are presented for individual EEG sensors and muscles.

Figure 12

(A) Box plot of GL for distinct decomposition algorithms, regardless of subject and sensor location. The smaller the value of GL, the more appropriate the filtering method. The dashed lines represent the best result obtained with the SSA approach. Statistically significant differences between methods are represented by labels. All possible combinations of two were evaluated. For example, the EEMD method is represented by label “b” and was statistically different from CiSSA (label “a”), EMD (label “c”), EMD-PCA (label “d”), and SSA (label “e”). (B) Box plot of GH for distinct decomposition algorithms, regardless of subject and sensor location. The larger the value of GH the more suitable is the method for filtering. (C) Box plot of the mean normalized Euclidean distance between GL and GH for each muscle, independent of subjects and EEG sensors. The larger the value of this metric, the more contamination is caused by the muscle. (D) Box plot of the mean normalized Euclidean distance between GL and GH for each muscle and subject, independent of the EEG sensor.

Figure 13

Evaluation of distinct adaptive filtering methods based on the time-domain features. The assessment is independent of the decomposition method and specific to the type of reference signal. (A–D) show results referent to GH and GL. (E,F) Present the results related to GXin and GXout.

Figure 14

A typical EEG signal corrupted by facial EMG. The EEG signal is from the Fp2-F8 region because it is most affected by facial electromyography. EMG bursts can be seen on the detrended EEG signal. The EEG signal was used as a reference signal, estimated from EMD, and then filtered using the RLS filter. The residue, which is the difference between the detrended and filtered EEG data, clearly shows the EMG activity that was eliminated from the signal. The inset plots at the top indicate the selection of two EEG regions contaminated by EMG. The filtered signal is shown over the contaminated signal in red. For each region, the detrended EEG, reference signal, filtered EEG, and residue are shown.

Figure 15

Median frequency and its power of the electroencephalogram (EEG) together with its components for the whole signal (ENTIRE) and the two experimental conditions (OPEN EYES and CLOSED EYES). A contrast between the raw non-filtered signal with the filtered signal is presented. The asterisks show the pair of variables in which the variable associated to the filtered signal was significantly reduced in comparison to the non-filtered signal, i.e., raw signal.

Figure 16

Execution time of distinct decomposition (A) and adaptive filtering (B) methods as function of the number of samples. The box plots show the central trend and dispersion of execution times (C,D).