Copyright © 2008 Kianoush Nazarpour et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

A novel blind signal extraction (BSE) scheme for the removal of eye-blink artifact from electroencephalogram (EEG) signals is proposed. In this method, in order to remove the artifact, the source extraction algorithm is provided with an estimation of the column of the mixing matrix corresponding to the point source eye-blink artifact. The eye-blink source is first extracted and then cleaned, artifact-removed EEGs are subsequently reconstructed by a deflation method. The a priori knowledge, namely, the vector, corresponding to the spatial distribution of the eye-blink factor, is identified by fitting a space-time-frequency (STF) model to the EEG measurements using the parallel factor (PARAFAC) analysis method. Hence, we call the BSE approach semiblind signal extraction (SBSE). This approach introduces the possibility of incorporating PARAFAC within the blind source extraction framework for single trial EEG processing applications and the respected formulations. Moreover, aiming at extracting the eye-blink artifact, it exploits the spatial as well as temporal prior information during the extraction procedure. Experiments on synthetic data and real EEG measurements confirm that the proposed algorithm effectively identifies and removes the eye-blink artifact from raw EEG measurements.

I. Introduction

The electroencephalogram (EEG) signal is the superposition of brain activities recorded as changes in electrical potentials at multiple locations over the scalp. The electrooculogram (EOG) signal is the major and most common artifact in EEG analysis generated by eye movements and/or blinks [1]. Suppressing eye-blink over a sustained recording course is particularly difficult due to its amplitude which is of the order of ten times larger than average cortical signals. Due to the magnitude of the blinking artifacts and the high resistance of the skull and scalp tissues, EOG may contaminate the majority of the electrode signals, even those recorded over occipital areas. In recent years, it has become very desirable to effectively remove the eye-blink artifacts without distorting the underlying brain activity. In this regard, reliable and fast, either iterative or batch, algorithms for eye-blink artifact removal are of great interest for diverse applications such as brain computer interfacing (BCI) and ambulatory EEG settings. Various methods for eye-blink artifact removal from EEGs have been documented that are mainly based on independent component analysis (ICA) [1, Chapter 2], linear regression [2], and references therein. Approaches, such as trial rejection, eye fixation, EOG subtraction, principal component analysis (PCA) [3], blind source separation (BSS) [4–6], and robust beamforming [7] have been also documented as having varying success. A hybrid BSS-SVM method for removing eye-blink artifacts along with a temporally constrained BSS algorithm have been recently developed in [5, 6]. Moreover, methods based on [8] adaptive and spatial filters [9] have also been presented in the literature for eye-blink removal. It has been shown that the regression- and BSS-based methods are most reliable [1, 2, 5–7, 10], despite no quantitative comparison for any reference dataset being available.

Statistically nonstationary EEG signals yield temporal and spatial information about active areas within the brain and have been effectively exploited for localizing the EEG sources and the removal of various artifacts from EEG measurements. For instance, in [11] PCA is utilized to decompose the signals into uncorrelated components where the first component, the component with highest variance, represents eye-blink artifact. However, the use of PCA introduces nonuniqueness due to an arbitrary choice of rotation axes. Although this nonuniqueness may be resolved by introducing reasonable constraints, recently, ICA has been applied to eliminate this problem by imposing the statistical independence constraint which is stronger than the orthogonality condition exploited by PCA [12]. However, the eye-blink component should be identified manually or in an automatic correction framework [5] if one uses ICA. In these conventional methods, usually prior concepts such as orthogonality, orthonormality, nonnegativity, and in some cases even sparsity have been considered during the separation process. However, such mathematical constraints usually do not reflect specific physiological phenomena. In essence, there are two different approaches for incorporating prior information within the semiblind EEG source separation (extraction); firstly, the Bayesian method [13] which introduces a probabilistic modeling framework by specifying distributions of the model parameters with respect to prior information. Often the probabilistic approach is too complicated, analytically and practically, to be implementable specifically in high-density EEG processing; slow convergence drawback should also be highlighted. The second more feasible approach proposes expansion of conventional gradient-based minimization of particular cost functions by including rational physiological constraints. Theoretically, widely accepted temporally or spatially constrained BSS (CBSS) [5, 14–16] algorithms are the outcome of above-mentioned methodology. However, CBSS methods still suffer from extensive computational requirements (unlike blind source extraction methods, i.e., [17]) of source separation and severe uncertainties regarding the accuracy of the priors.

Simple and straightforward priors, such as the spectral knowledge of ongoing EEGs or spatial topographies of some source sensor projections, can be realistically meaningful in semiblind EEG processing. In this regard, an interesting work on topographic-time-frequency decomposition is proposed in [18] in which, however, two mathematical conditions on time-frequency signatures, namely, minimum norm and maximal smoothness, are imposed. It has been shown that these conditions may provide a unique model for EEG measurements. Consolidating [18], recently in [19] the space-time-frequency (STF) model of a multichannel EEG has been introduced by using parallel factor analysis (PARAFAC) [20]. More recently, we have utilized the STF model for the first time in single trial EEG processing for brain computer interfacing, where spatial signature of selected component is employed as a feature vector for classification purpose [1, 21].

In this paper, a novel physiologically inspired semiblind signal extraction technique for removing the eye-blink artifacts from single trial multichannel EEGs is presented. Our SBSE method is based on that introduced in [17], while by investigating the STF signatures of extracted factor(s) by PARAFAC, the eye-blink factor is automatically selected and its spatial distribution is exploited in the separation procedure as a prior knowledge. The main advantages of our method are as follows:

(1) in the BSS- and CBSS-based methods [4, 6, 15, 16, 22–24], identification of the correct number of sources is an important issue and requires high computational costs. However, the simplicity of our method is due to using the spatial a prior information to guarantee that the first extracted source is the one of interest, that is, the eye-blink source. Therefore, there is no need to extract other sources which significantly reduces the computational requirements. EEGs are then reconstructed in a batch deflation procedure; (2) unlike methods presented in [4, 5], there is no need to compute objective criteria for distinguishing between eye-blink and spurious peaks in the ongoing EEGs; (3) unlike the regression-based methods [25], the proposed method does not need any reference EOG channel recordings (typically three channels); (4) there is no need to separate the dataset into training and testing subsets as in [6]. As long as, by using any primitive method we identify an eye-blink event, the presented method can be utilized to remove the artifact from EEGs.

This paper is organized as follows. In Section 2, we present the SBSE method and compare its performance to that of an existing spatially constrained BSS algorithm presented in [16]. Afterwards, we briefly review the fundamentals of the PARAFAC method in Section 2.2 and suggest our effective procedure to identify the spatial signature of the eye-blink relevant factor. The results are subsequently reported in Section 3, followed by concluding remarks in Section 4.

2. Algorithm Development

Eye-blink contaminated EEG measurements at time are assumed as zero-mean real mutually uncorrelated sources , where denotes the vector transpose, mixed by an real full column rank matrix where generally is the th column of and specifically is the column of corresponding to the eye-blink source . The vector of time mixture samples is given as (1) where is the additive white Gaussian zero-mean noise. We assume that the noise is spatially uncorrelated with the sensor data and temporally uncorrelated. Since the sources are presumed to be uncorrelated, the time lagged autocorrelation matrix can be calculated as (2) for , where is the index of the maximum time lag, that is, and denotes the statistical expectation operator. In (2), is the time lagged autocorrelation value of .

2.1. Semiblind Eye-Blink Signal Extraction

The vector in (1), that is, recorded EEGs, is a linear combination of the columns of the mixing matrix, that is, the s, weighted by the associated source and contaminated by sensor noise . Therefore, the most straightforward way to extract the th source, the eye-blink artifact , is to project onto the space in orthogonal to, denoted by , all of the columns of except , that is, . Hence, by defining a vector and , and adopting the notation of an oblique projector [17, 26], we may write (3) where is an estimate of one source, say , and denotes the space in orthogonal to , that is, . In (3), represents the oblique projection of onto the space . Then, can be extracted using the spatial filter as (4) in which the scalar has been omitted and has been dropped from both sides of (3). In second-order statistics-based BSE [17], both and are unknown and in order to extract one source the following cost function has been proposed: (5) where , is a column vector and denotes the squared Euclidean norm. We employ multiple time lags instead of a single time lag which minimizes the chance, in practice, of the time-lagged autocorrelation matrices employed having duplicate eigenvalues and, hence, leading to failure in the extraction process [5]. The cost function utilized in (5) exploits the fact that for BSE, should be collinear [27] with incorporating the coefficients which provides with the proper scaling. The trivial answer for (5) is . This solution has been avoided by imposing the condition . Successful minimization of (5) leads to the identification of , which extracts the source of interest (SoI) in (4).

The main advantage of using (5) for BSE over other conventional BSE methods which incorporate higher order statistics [12] is that it is indeed computationally simple and efficient for extraction of nonstationary sources. However, fundamentally in BSE, in the course of extraction, it is not possible to tune the algorithm to extract the SoI as the first extracted source in order to significantly decrease the processing time which is essential in real-time applications. Therefore, some prior knowledge should be incorporated into the separation process to extract only the SoI. To this end, we consider an auxiliary cost function (6) where is a column vector and is prior spatial information of the eye-blink source, that is, the estimation of , provided by PARAFAC (Section 2.2).

By minimizing coupled with (5) in a Lagrangian framework, that is, , we effectively extract the SoI as the first extracted source. Moreover, as it will be shown in Section 3, including has significant incremental effect in the minimization and results in faster convergence of . In mathematical terms the novel cost function is (7) where is the Lagrange multiplier. In (7), the values are free parameters to scale during an iterative solution to (7) and .

Essentially, there are two approaches in using the spatial priors which vary the degree of freedom of the optimizing process, that is, (7). In the optimizing procedures, we can either strictly minimize the difference between and iteratively as much as possible regardless of the probable errors while estimating or on the other hand, by employing a milder approach and allowing in the optimization process to deviate from the prior vector by an -norm-bounded threshold. In mathematical terms, in soft constraining, we consider as the estimation error where ; is a known positive constant. For the majority of spatially constrained BSS applications, that is, [16, 22] and references therein, the latter conservative approach is preferable to strict ones, even if is accurately estimated. However, to the authors' belief, for eye-blink artifact removal from EEGs hard constraining the extraction algorithm is sufficient since sparsely occurring eye-blink is the dominant source superimposed on EEGs. Therefore, the estimation of is trustworthy. We, in this paper, have explored the former approach and assumed that the estimation of by the PARAFAC-based STF model is accurate enough. We have also experimentally found that although the introduction of in (7) does not have any rotational effect on , it does result in better minimization of . The interested reader is referred to [7] in which we have realized a conservative method for the eye-blink artifact removal from the EEGs.

The solution to (7) is found by alternatively adjusting its parameters, that is, an alternating least squares (ALS) method. We iteratively update each of the four unknown vectors till convergence. Firstly, we fix and and update . Taking the gradient of with respect to leads to an optimal analytical solution for as (8) where denotes replacing by . Thereafter, we fix , , and and update . As in [17], utilizing the property that , the gradient of with respect to becomes (9)

The update rule for is as (10) where . Then, fixing , , and , we adjust while ensuring . Consider (11) and is adjustable by (12)

For updating , the rest of the variables are fixed, that is, , and and we proceed by minimizing (7) with respect to , that is, (13) is updated as (14) We retain as a vector instead of a scalar to present a consistent formulation. Finally, in order to solve (7) for the Lagrange multiplier, that is, , we define vector as a vector whose elements are all zero except for the th component which is one, that is, . Considering that in (12), can be easily updated by putting after each iteration. Therefore, we assign a new value for as (15)

The performance of the proposed semiblind signal extraction procedure has been evaluated through a comparison with the spatially constrained blind signal separation (SCBSS) algorithm proposed in [16, 22] for a set of synthetic mixtures of analytic sources.

Four signal sources, namely, two sinusoids of frequencies of 10 Hz and 12 Hz representing brain rhythmic waves, a spiky source standing for eye-blink artifact and a white Gaussian distributed signal as the background brain activity have been synthetically mixed. The mixing matrix (generated randomly from a standardized normal distribution) used in this paper is (16)

The source waveforms and the mixtures are presented in Figures 1(a) and 1(b). The source signals have been selected as such in order to cover the range of sub-Gaussianity to super-Gaussianity. The original mixtures have been plotted in Figure 1(b) in solid blue lines, where and are highly affected by the spiky source, . Here, the objective is to visually compare our proposed method with that of [16] in which a spatially constrained blind source separation (SCBSS) method based on FastICA [12] has been suggested for eye-blink artifact removal. In Figure 1(b), the outcome of our semiblind signal extraction method has been plotted in red solid lines which has effectively removed the signal from the mixtures. It is also worth considering the clean artifact free parts of the mixtures which have been reconstructed perfectly. Moreover, the outputs of the established method of [16] in artifact removal from EEGs have been shown in solid green lines. Evidently, the outcome of our method does overlap that of [16]. The correlation coefficient of two discrete random variables and over a fixed interval is mathematically defined as: (17) where is the number of time samples. Figure 1(c), demonstrates averaged values between segments of cleaned mixtures (after removing ) and original mixtures by using proposed method and that of [16, 22]. CC values of about unity show that SBSE method provide similar results as to SCBSS.

Figure 1: Simplified scalp EEG measurements; brain source signals in (a) and mixed recordings (b). (a) shows four synthetic sources, namely, and which represent brain rhythmic activities, for background white noise, and the eye-blink artifact source. (b) illustrates the mixed signals in solid blue lines, that is, , where and are highly contaminated by the eye-blink source, . The artifact removed mixtures have been also plotted by using our proposed method, plotted in solid red, and that of [16] in solid green lines. Evidently, our proposed method presents reasonably similar performance to that of the semiblind separation method in [16]. In (c), the averaged CC values between the segments of cleaned mixtures (after removing ) and the original mixtures by using SBSE method and SCBSS algorithm in [16] have been depicted. CC values of about unity again justify that the SBSE method provides similar results as to SCBSS.

In these simulations, we have presumed that spatial distribution (signature) of the source of interest, , is known in advance. This assumption helps to validate our SBSE method comparing to [16, 22] regardless of how accurate various existing methods perform in estimating the aforementioned vector.

Moreover, through simulation studies we have found consistent faster convergence of our optimization scheme, as reported in Section 3, as compared to that in [17] which highlights that incorporating auxiliary cost function into extraction process significantly upgrades the performance. Next, we establish how PARAFAC is utilized to provide the required a prior information.

2.2. PARAFAC

PARAFAC is a widely accepted tool in extracting disjoint multidimensional phenomena with application to food science, communications, and biomedicine [7, 10, 19–21, 28, 29, 30, 31]. In this paper, by exploiting PARAFAC, we decompose the eye-blink contaminated EEG measurements in order to extract the factor relevant to the eye-blink artifact for use within the SBSE. The resulting spatial signature of the eye-blink-related factor, that is, is exploited to formulate (7). The spatial signatures of this factor is directly related to the level of eye-blink contamination for each electrode and is thereby comparable to the column of the mixing matrix that propagates the point source eye-blink artifact into the EEG channels. Physiologically, this assumption is rational since eye-blink is attenuated while propagating from frontal to central and occipital areas of the brain.

In our approach, the multichannel EEG data are transformed into time-frequency domain. This gives the two-way EEG recording, that is, the matrix of space(channel)-time, an extra dimension and yields a three-way array of space-time frequency. In other words, for EEG channels, we compute the energy of the time-frequency transform for time instants and frequency bins. By stacking these matrices (of size ) and adopting the Matlab matrix notation, we set up the three-way array and introduce it to PARAFAC.

Conventional methods, for instance, PCA or ICA, analyze such data by unfolding some dimensions into others, reducing the multiway array into matrices. However, the aforementioned unfolding procedures make the interpretation of the results ambiguous since they remove specific information endorsed by those dimensions. Consequently, rather than unfolding these multiway arrays into matrices, we exploit PARAFAC to explore the space-time-frequency (STF) model of EEG recordings. The key idea behind this research is in considering EEGs as superposition of neural electro-potentials. EEGs may be represented by using the linear models which are defined in three domains, that is, space, time, and frequency, in order to simultaneously investigate their spatial, temporal, and spectral dynamics [1, 7, 10, 19, 21, 30]. Here, we have assumed that each distinct local EEG activity (on the scalp) is uncorrelated with the activities of the neighboring areas of the brain. EEGs can be modeled as sum of the distinct components where each distinct component is formulated as the product of its basis in space, time, and frequency domains. The interested reader is referred to [28, 29, 32] for further mathematical details of the PARAFAC model, the uniqueness conditions, and its robust iterative fitting which are out of the scope of this paper.

Complex wavelet transform

To setup a three-way array, in the present study, a continuous wavelet transform is utilized to provide a time-varying representation of the energy of the signals over all channels. The complex Morlet wavelets , with , and , are used here in which the tradeoff ratio is 7, to create a wavelet family. This wavelet configuration is known to be optimized in EEG processing [19]. The time-varying energy of a signal at a specific frequency band is the squared norm of the convolution of a complex wavelet of the signal , that is, , where stands for the convolution product and the modulus operator is denoted by .

In mathematical terms, the factor analysis is expressed as where is the factor loading, the factor score, the error, and the number of factors. Similarly, the PARAFAC for the three-way array is presented by unfolding one modality to another as (18) where is the factor score corresponding to the second modality and is the Khatri-Rao product and denotes the Kronecker product [33]. Equivalently, the th matrix corresponding to the th slice of the second modality of the 3-way array is expressed as (19) where is a diagonal matrix having the th row of along the diagonal. ALS is the most common way to estimate the PARAFAC model. In order to decompose the multiway array to parallel factors the cost function (normally the squared error) is minimized as in [20] (20) Here, is the three-way array of wavelet energy of multichannel EEG recordings and , and denote the spatial, temporal, and spectral signatures of , respectively. In this paper, the trilinear alternating least squares (TALSs) method [34] is used to compute the parameters of the STF model. We in [7], inspired by [30], have proposed a novel computationally simple method for STF modeling of EEG signals in which in order to reduce the complexity present in the estimation of the STF model using the three-way PARAFAC, the time domain is subdivided into a number of segments and a four-way array is then set to estimate the space-time-frequency-time/segment (STF-TS) model of the data using the four-way PARAFAC. Subsequently, the STF-TS model is shown to approximate closely the classic STF model, with significantly lower computational cost.

In summary, our method consists of the following stages. Given an artifact contaminated EEG data, we (1)bandpass filter the EEGs between 1 Hz and 40 Hz,(2)set up the three-way array, that is, , as stated in Section 2.2,(3)execute PARAFAC and select the eye-blink artifact relevant factors as will be fully described in Section 3,(4)exploit the spatial signature of the eye-blink artifact factor in SBSE cost function (7),(5)reconstruct the artifact removed EEGs in a deflation framework. See the appendix.

3. Results

We applied the SBSE algorithm to real EEG measurements. The database was provided by the School of Psychology, Cardiff University, UK, and contains a wide range of eye-blinks and, therefore, gives a proper evaluation of our method. The scalp EEG was obtained using 28 Silver/Silver-Chloride electrodes placed at locations defined by the 10–20 system [1]. EEGs have been recorded to provide a reference dataset specifically for the purpose of evaluating different artifact removal methods from one healthy subject and contains numerous eye-blinks, eye movements, and motion artifacts. The data were sampled at 200 Hz, and bandpass filtered with cut-off frequencies of 1 Hz and 40 Hz. In order to reduce the computational costs of the PARAFAC modeling, we selected 16 channels out of the above-mentioned 28 channels as illustrated in Figure 2.

Figure 2: The result of the proposed eye-blink artifact removal method for a sample of real EEG signals recorded from the selected 16 electrodes. In (a), the EOG is evident just after the time 2 seconds and more prominent on the frontal electrodes, that is, FP1 and FP2. However, in (b), the same segment of EEG after being corrected for eye-blink artifact using the proposed algorithm is illustrated. Note the small spike-type signals, indicated by arrows, right after the first eye-blink are precisely retained after eye-blink artifact removal.

Each EEG segment was transformed into the time-frequency domain by means of the complex wavelet transform where the frequency band from 2 Hz to 25 Hz with resolution of 0.1 Hz has been considered. This three-way array is then introduced to PARAFAC where the number of factors is selected as one or two, as highlighted in the following experiments, identified by using the method of core consistency diagnostic (CORCONDIA) [35]. Automatically, PARAFAC identifies the most significant factors with CORCONDIA values greater than a set threshold, that is, 85% [35], within each recording. Two sample results are demonstrated here in order to elaborate the potential of our method.

3.1. Experiment 1

Figure 2(a) shows EEG measurements which are contaminated with two eye-blinks at approximate times of two and half and five seconds. The effects of the eye-blinks are evident mostly in the frontal electrodes, namely, FP1, FP2, F3, F4, F7, and F8. However, central, C3 and C4 and occipital O1 electrodes are also partly affected. Implementation of PARAFAC on this measurement results in the STF model, the spectral, temporal, and spatial signatures which are depicted in Figures 3(a) to 3(c). Although there are two eye-blinks, CORCONDIA suggests the number of factors to be one as in Figure 3(d). This value is rational since both of the eye-blinks originate from a certain vicinity (frontal lobe of the brain) and occupy the same frequency band and there is no significant brain background activity. By using spatial distribution of the extracted factor as a prior information, eye-blink artifacts are effectively removed. In order to minimize (7) initial values of the vectors , , , and independently drawn from standardized normal distributions , is initialized to 5 and is set to the spatial signature of the extracted factor. Figure 4 compares the average value of over 50 independent runs. Two scenarios have been devised by varying the number of time lags, that is, . Note that in [17], . Evidently, in both scenarios, performance of proposed SBSE method is superior to that of the method in [17]. After approximately 10 iterations, the extracting vector is identified. Furthermore, by incorporating the prior knowledge, it is guaranteed that extracts the eye-blink source. The effect of the eye-blink is then removed from the multichannel EEG using the batch deflation algorithm in [36]. The impressive issue on the resolution of the proposed algorithm is that it does not affect the very low amplitude spike-type signals right after first eye-blink, indicated by arrows, during extraction process, Figure 2.

Figure 3: The extracted factor by using PARAFAC; (a) and (b) illustrate, respectively, the spectral and temporal signatures of the extracted factors and (c) represents the spatial distribution of the extracted factor which has been considered as the a prior knowledge during extraction procedure, (d) shows that the number of factors suggested by CORCONDIA to be one since the bars corresponding to and are less than the threshold, that is, 85%.

Figure 4: The averaged (over 50 independent runs) convergence characteristics, , of the SBSE and BSE of [17] are depicted for two values of , that is, 10 in (a) and 25 in (b). In both subplots the solid and dashed curves correspond, respectively, to the proposed SBSE and BSE of [17].

3.2. Experiment 2

Performance of the method with same initial values for another set of EEGs from the database is demonstrated in Figure 5 where in left subplot, the truncated 4 seconds of EEG recordings before and after eye-blink removal processing are plotted. Figure 5(b) illustrates averaged correlation coefficients between artifact removed channel signals and original contaminated ones with their corresponding standard deviations over 25 independent runs. As expected, CC values corresponding to the signals recorded from frontal electrodes are relatively low showing these signals are significantly altered; artifact removed. However, values corresponding to other channel signals, that is, parietal, central, temporal, and occipital, are almost unity demonstrating that our algorithm does not affect clean EEG measurements.

Figure 5: The results of the proposed eye-blink artifact removal method for a set of real EEG signals recorded from 16 electrodes; (a) shows the eye-blink contaminated EEGs in red and the artifact corrected EEGs in blue where the eye-blink artifact is evident just before time 2 seconds and more prominent on the frontal electrodes, that is, FP1 and FP2. However, in (b), averaged CC values between the artifact corrected channel signals and the original contaminated EEGs with their corresponding standard deviations over 25 independent runs are plotted. CC values corresponding to the frontal channel signals are relatively lower than the values corresponding to other channel signals which are almost unity, (b) illuminates how our algorithm reconstructs the artifact-freed EEGs faithfully without affecting clean signals coming from nonfrontal areas.

The STF model of this recording is introduced by PARAFAC. In contrast to previous experiments, CORCONDIA suggests since PARAFAC identified a significant brain background activity during occurrence of eye-blink. Figures 6(a) to 6(d) illustrate the estimated signatures of 16-channel EEG signal contaminated by eye-blink. The first component (factor 1) of the STF model demonstrates the eye-blink-relevant factor. (1) It mainly occurs in the frequency band of around 5 Hz while the other factor exists in the entire band and represents the ongoing activity of the brain or perhaps a broadband white noise-like component, Figure 6(a). (2) The temporal signature of the first factor definitely shows a transient phenomenon such as eye-blink while that of Factor 2 consistently exists in the course of EEG segment, Figure 6(b). (3) Unlike in Figure 6(d), in Figure 6(c), the spatial distribution of the extracted factor is confined to the frontal area, which clearly demonstrates the effect of eye-blink. The other factor shows the background activity of the brain as it spreads all over the scalp.

Figure 6: The extracted factors by using PARAFAC; (a) and (b) illustrate, respectively, the spectral and temporal signatures of the extracted factors; (c) and (d) present the spatial distribution of the factors, respectively. Evidently, factor 1 demonstrates the eye-blink phenomenon as it occurs in frequency band of around 5 Hz (a), it is indeed transient in the time domain (b) and it is confined to the frontal area.

Hence, we employ spatial distribution of the first extracted factor in the SBSE.

3.3. Performance Evaluations

In order to provide a quantitative measure of performance for the proposed artifact removal method, the CC values of the extracted eye-blink artifact source and the original and the artifact removed EEGs are computed, see Figure 7.

Figure 7: The averaged CC values (and their corresponding standard deviations) between the extracted eye-blink and the restored EEGs before and after artifact removal of different channels in (a) and (b), respectively. The experiments have been performed for 20 different eye-blink contaminated EEG recordings. Note that the scales are different by 10³.

The values reported in Figure 7 have been computed as follows. For each of the 20 different artifact contaminated EEGs, we executed our proposed algorithm. The aforementioned CCs for each run were then computed between the extracted eye-blink and the EEGs before and after the artifact removal. These values have subsequently been averaged and shown in Figure 7. Furthermore, their corresponding standard deviations have also been reported. As expected, the CC values have been significantly decreased by using the proposed method. Simulations for 20 EEG measurements demonstrate that the proposed method can efficiently identify and remove the eye-blink artifact from the raw EEG measurements.

As a second criterion for measuring the performance of the overall system, we selected a segment of EEG, called and the reconstructed EEG which does not contain any artifact, and measured the waveform similarity by (21) When the value of is zero, the original and reconstructed waveforms are identical. From the 20 sets of EEGs, the average waveform similarity was as low as  dB (standard deviation  dB). These results suggest that the observations have been faithfully reconstructed.

3.4. Robustness

As indicated earlier, in soft constrained blind source extraction (separation [16]) schemes, even if the estimation of is slightly biased, the optimization algorithm takes that into account and accommodates it during the extraction of the source of interest. However, as indicated in Section 2.1, in this paper a hard approach has been taken where the algorithm strictly minimizes the cost function, in (7) regardless of the probable errors or biases while estimating .

Interestingly, the scenario is not actually as restricted as it seems; that is, even if there is a small deviation in the from the actual which sounds quite rational, SBSE is able to accommodate that without any need for further formulations as in [16]. The truth lies in the alternating least squares approach in updating , that is, (11) where SBSE tries to estimate the best set of and simultaneously both ideally orthogonal to in order to minimize the cost function (7). Therefore, even if is utilized instead of the , as the result of STF modeling and PARAFAC in the cost function (7), the optimization process results in converging to the originally estimated , that is, . In the sequel the results of a series of experiments with different are presented in order to consolidate the proposed SBSE method for EB artifact removal. Let us start with an experiment where instead of is introduced to SBSE where is computed as (22) where is a vector of 16 elements ideally drawn from a zero-mean and unit-variance normal distribution, that is, . Using (22), the norm of is highly likely to be less than 0.6. Therefore, if , it is probable that SBSE compensates for the deviation of from and extracts the EB artifact. For instance in Figures 8 and 9, an example has been provided where . In Figure 8(a), obtained by PARAFAC is depicted which should be used in (7). Figure 8(b) shows the perturbed by which has been replaced in (7) instead of and introduced to SBSE. Finally, in Figure 8(c), the resulting after the alternative least squares optimization has been illustrated. Evidently, Figure 8(c) is quite similar to Figure 8(a).

Figure 8: In (a), is depicted, (b) shows the deviated by which has been put in (7) instead of , (c) illustrates the resulting after ALS optimization procedure.

Figure 9: The result of the artifact removal from EEGs depicted in Figure 5(a). EEG traces plotted in red color are the original artifact contaminated signals. EEGs in blue color are the resulting artifact removed signals using . Traces in black are the resulting artifact restored EEGs by using instead of .

The result of the artifact removal is depicted in Figure 9. EEG traces in red are the original artifact contaminated recordings. Traces in blue are the resulting artifact removal using the original estimate of , that is, , by PARAFAC. EEG plots in black, which entirely overlap with the blue ones, are the resulting artifact restored EEGs by using the artificially perturbed , that is, put in (7).

Thereafter, instead of is introduced to SBSE. The vector is computed in the same way as by keeping the coefficient as 0.1 in (22), norm . Since Figure 10(b), is significantly different in steering direction from Figure 10(a), SBSE may not compensate for the deviation . In Figure 10(a), resulted by PARAFAC is depicted which should have been put in (7). Figure 10(b) shows the perturbed by which has been replaced in (7) instead of and introduced to SBSE. Finally, in Figure 10(c), the resulting after the alternative least squares optimization has been illustrated. The vector plotted in Figure 10(c) does not converge to the vector plotted in Figure 10(a).

Figure 10: In (a), is depicted, (b) shows the deviated by which has been put in (7) instead of , (c) illustrates the resulting after ALS optimization procedure.

The result of the artifact removal is depicted in Figure 11. Again as Figure 9, the EEG traces in red are the original artifact contaminated recordings. Traces in blue are the resulting artifact removal using the original estimate on , that is, , by PARAFAC. However, EEG plots in black show an absolute failure in artifact removal procedure by .

Figure 11: The result of the artifact removal from EEGs depicted in Figure 5(a). EEG traces plotted in red color are the original artifact contaminated signals. EEGs in blue color are the resulting artifact removed signals using . Traces in black are the resulting of the unsuccessful artifact removal procedure by using instead of .

It can be concluded that in order that the SBSE presents a robust performance even if is perturbed by a norm bounded small deviation, its direction should not be changed. That is, if the bias is fairly distributed over the elements of , since a normalized version is used in the formulations, based on our experience, it is highly unlikely that SBSE does not compensate for it.

4. Concluding Remarks

It is generally accepted that the eye-blink artifact can be removed from EEGs by using the BSS- and regression based methods for multichannel EEGs data with/without the reference EOG electrodes. However, nowadays this challenging topic is tended to be solved by a semiblind method rather than in a totally blind signal processing framework [5, 7, 10, 15, 16, 22]. Notwithstanding these recently published semiblind approaches, we propose an analytic and rational method to acquire the prior information, that is, the spatial signature of the eye-blink signal, from the EEG measurements. Therefore, we do not follow the conventional heuristic approaches such as that of [15] where an approximation of the temporal structure of the eye-blink source signal is included in ICA. Furthermore, to the best of our knowledge, there has not been any method specifically based on semiblind signal extraction for eye-blink artifact removal from EEGs. The presented method is computationally simpler than the spatially constrained blind source separation method of [16, 22] since there is no need to estimate all the columns of the mixing matrix in (1).

The vector of spatial distribution of the eye-blink factor has been identified using PARAFAC. For the first time in this work, we have utilized the vector of spatial signature of the eye-blink factor resulted by the STF modeling of EEGs as the estimation of the column vector of the mixing matrix that introduces the eye-blink source to the EEGs. This assumption is rational since the eye-blink can be considered as a strong point source which is merely attenuated while propagating from frontal area to the central and occipital parts of the brain. This spatial distribution of the eye-blink factor then has been incorporated to our SBSE algorithm. The EEGs are processed using the time-lagged second-order SBSE algorithm and the artifact is autonomously extracted; then, the EEGs are reconstructed in a deflation framework. Based on our experiments, the proposed SBSE algorithm consistently removes the eye-blink artifacts from the EEG signals.

Appendix

The Deflation Method

In order to achieve EB-free EEG recordings, , after the extraction of the EB source using (4), we apply the deflation procedure which eliminates the previously extracted signal, , from the recording mixtures, that is, : (A.1) where, as in [36, Section 5.2.5], can be estimated either adaptively or simply after minimization of the mean square cost function with respect to : (A.2) This results in the following efficient batch one-step formula to estimate : (A.3) where is achieved by (8 ). In fact, is an estimation of , the th column of , neglecting arbitrary scaling and permutations of columns ambiguities.

Acknowledgments

This work is supported in part by The Leverhulme Trust, UK, and Cardiff University, UK. The authors would like to acknowledge Dr. Edward Wilding at the School Psychology, Cardiff University, UK, for the provision of the dataset.

References

1. S. Sanei and J. A. Chambers, EEG Signal Processing, John Wiley & Sons, New York, NY, USA, 2007.
2. A. Schlögl, C. Keinrath, D. Zimmermann, R. Scherer, R. Leeb, and G. Pfurtscheller, “A fully automated correction method of EOG artifacts in EEG recordings,” Clinical Neurophysiology, vol. 118, no. 1, pp. 98–104, 2007.
3. M. Fatourechi, A. Bashashati, R. K. Ward, and G. E. Birch, “EMG and EOG artifacts in brain computer interface systems: a survey,” Clinical Neurophysiology, vol. 118, no. 3, pp. 480–494, 2007.
4. K. H. Ting, P. C. W. Fung, C. Q. Chang, and F. H. Y. Chan, “Automatic correction of artifact from single-trial event-related potentials by blind source separation using second order statistics only,” Medical Engineering and Physics, vol. 28, no. 8, pp. 780–794, 2006.
5. L. Shoker, S. Sanei, and J. A. Chambers, “Artifact removal from electroencephalograms using a hybrid BSS-SVM algorithm,” IEEE Signal Processing Letters, vol. 12, no. 10, pp. 721–724, 2005.
6. L. Shoker, S. Sanei, W. Wang, and J. A. Chambers, “Removal of eye blinking artifact from the electro-encephalogram, incorporating a new constrained blind source separartion algorithm,” Medical and Biological Engineering and Computing, vol. 43, no. 2, pp. 290–295, 2005.
7. K. Nazarpour, Y. Wongsawat, S. Sanei, J. A. Chambers, and S. Oraintara, “Removal of the eye-blink artifacts from EEGs via STF-TS modeling and robust minimum variance beamforming,” in Proceedings of the 29th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC '07), pp. 6211–6214, Lyon, France, August 2007.
8. S. Puthusserypady and T. Ratnarajah, “H∞ adaptive filters for eye blink artifact minimization from electroencephalogram,” IEEE Signal Processing Letters, vol. 12, no. 12, pp. 816–819, 2005.
9. T. D. Lagerlund, F. W. Sharbrough, and N. E. Busacker, “Spatial filtering of multichannel electroencephalographic recordings through principal component analysis by singular value decomposition,” Journal of Clinical Neurophysiology, vol. 14, no. 1, pp. 73–82, 1997.
10. K. Nazarpour, S. Sanei, and J. A. Chambers, “A novel semi-blind signal extraction approach incorporating PARAFAC for the removal of the removal of eye-blink artifact from EEGs,” in Proceedings of the 15th International Conference on Digital Signal Processing (DSP '07), pp. 127–130, Cardiff, UK, July 2007.
11. A. C. K. Soong and Z. J. Koles, “Principal-component localization of the sources of the background EEG,” IEEE Transactions on Biomedical Engineering, vol. 42, no. 1, pp. 59–67, 1995.
12. A. Hyvärinen, J. Karhunen, and E. Oja, Independent Component Analysis, John Wiley & Sons, New York, NY, USA, 2001.
13. C. Phillips, J. Mattout, M. D. Rugg, P. Maquet, and K. J. Friston, “An empirical Bayesian solution to the source reconstruction problem in EEG,” NeuroImage, vol. 24, no. 4, pp. 997–1011, 2005.
14. M. A. Latif, S. Sanei, J. A. Chambers, and L. Shoker, “Localization of abnormal EEG sources using blind source separation partially constrained by the locations of known sources,” IEEE Signal Processing Letters, vol. 13, no. 3, pp. 117–120, 2006.
15. C. J. James and O. J. Gibson, “Temporally constrained ICA: an application to artifact rejection in electromagnetic brain signal analysis,” IEEE Transactions on Biomedical Engineering, vol. 50, no. 9, pp. 1108–1116, 2003.
16. C. W. Hesse and C. J. James, “On semi-blind source separation using spatial constraints with applications in EEG analysis,” IEEE Transactions on Biomedical Engineering, vol. 53, no. 12, pp. 2525–2534, 2006.
17. X.-L. Li and X.-D. Zhang, “Sequential blind extraction adopting second-order statistics,” IEEE Signal Processing Letters, vol. 14, no. 1, pp. 58–61, 2007.
18. T. Koenig, F. Marti-Lopez, and P. A. Valdés-Sosa, “Topographic time-frequency decomposition of the EEG,” NeuroImage, vol. 14, no. 2, pp. 383–390, 2001.
19. F. Miwakeichi, E. Martínez-Montes, P. A. Valdés-Sosa, N. Nishiyama, H. Mizuhara, and Y. Yamaguchia, “Decomposing EEG data into space-time-frequency components using parallel factor analysis,” NeuroImage, vol. 22, no. 3, pp. 1035–1045, 2004.
20. R. Bro, “PARAFAC. Tutorial and applications,” Chemometrics and Intelligent Laboratory Systems, vol. 38, no. 2, pp. 149–171, 1997.
21. K. Nazarpour, S. Sanei, L. Shoker, and J. A. Chambers, “Parallel space-time-frequency decomposition of EEG signals for brain computer interfacing,” in Proceedings of the 14th European Signal Processing Conference (EUSIPCO '06), Florence, Italy, September 2006.
22. C. W. Hesse and C. J. James, “The fastICA algorithm with spatial constraints,” IEEE Signal Processing Letters, vol. 12, no. 11, pp. 792–795, 2005.
23. W. Lu and J. C. Rajapakse, “Approach and applications of constrained ICA,” IEEE Transactions on Neural Networks, vol. 16, no. 1, pp. 203–212, 2005.
24. N. Ille, Artifact correction in continuous recording of the electroand magnetoencephalogram by spatial filtering, Ph.D. thesis, University of Mannheim, Mannheim, Germany, 2001.
25. L. C. Parra, C. D. Spence, A. D. Gerson, and P. Sajda, “Recipes for the linear analysis of EEG,” NeuroImage, vol. 28, no. 2, pp. 326–341, 2005.
26. R. T. Behrens and L. L. Scharf, “Signal processing applications of oblique projection operators,” IEEE Transactions on Signal Processing, vol. 42, no. 6, pp. 1413–1424, 1994.
27. H. S. M. Coxeter and S. L. Greitzer, Collinearity and Concurrence, The Mathematical Association of America, Washington, DC, USA, 1967.
28. N. D. Sidiropoulos, G. B. Giannakis, and R. Bro, “Blind PARAFAC receivers for DS-CDMA systems,” IEEE Transactions on Signal Processing, vol. 48, no. 3, pp. 810–823, 2000.
29. N. D. Sidiropoulos, R. Bro, and G. B. Giannakis, “Parallel factor analysis in sensor array processing,” IEEE Transactions on Signal Processing, vol. 48, no. 8, pp. 2377–2388, 2000.
30. Y. Wongsawat, S. Oraintara, and K. R. Rao, “Reduced complexity spacetime-frequency model for multi-channel EEG and its applications,” in Proceedings of the IEEE International Symposium on Circuits and Systems (ISCAS '07), pp. 1305–1308, New Orleans, La, USA, May 2007.
31. T. Acar, Y. Yu, and A. P. Petropulu, “Blind MIMO system estimation based on PARAFAC decomposition of higher order output tensors,” IEEE Transactions on Signal Processing, vol. 54, no. 11, pp. 4156–4168, 2006.
32. S. A. Vorobyov, Y. Rong, N. D. Sidiropoulos, and A. B. Gershman, “Robust iterative fitting of multilinear models,” IEEE Transactions on Signal Processing, vol. 53, no. 8, pp. 2678–2689, 2005.
33. C. R. Rao and S. Mitra, Generalized Inverse of Matrices and Its Applications, John Wiley & Sons, New York, NY, USA, 1971.
34. Y. Rong, S. A. Vorobyov, A. B. Gershman, and N. D. Sidiropoulos, “Blind spatial signature estimation via time-varying user power loading and parallel factor analysis,” IEEE Transactions on Signal Processing, vol. 53, no. 5, pp. 1697–1710, 2005.
35. R. Bro and H. A. L. Kiers, “A new efficient method for determining the number of components in PARAFAC models,” Journal of Chemometrics, vol. 17, no. 5, pp. 274–286, 2003.
36. A. Cichocki and S. Amari, Adaptive Blind Signal and Image Processing: Learning Algorithms and Applications, John Wiley & Sons, New York, NY, USA, 2005.