Keywords

1 Introduction

An electrooculogram (EOG) is the recording of electrical activity originating from blinking as the movement of the eyeball itself, which is obtained by using electrodes placed on the skin around eyes. The an electrical potential difference is created between the cornea and the ocular fundus and the time course of the potential subsequently monitored with typical values ranging from 250 to 1000 \(\upmu \)V [1, 2]. Both blinking and movement of the eye is a source of noise when neuronal activity is recorded via EEG. This is largely due to the fact that the EEG potential is less than one tenth of the EOG, typically ranging from 5 to 100 \(\upmu \)V. The estimation of pure EOG and EEG signals is essential for parallel verification of changes in visual attention and event related neuronal activities. Both signals are generated from the same biological systems and interfere with each other as a result. To separate these two types of signals, identification of the individual characteristics of each signal in the time-frequency domain is key. Once these characteristics are identified, the unique shape of EOG and EEG signals under various conditions can be determined.

Signal feature extraction methods using linear analysis in time-frequency domains include Fast Fourier transform, wavelet transform, eigenvectors and so on [3]. For signal source separation, blind source separation (BSS), independent component analysis (ICA),the principal component technique (PCA) and many other methodologies have been used to extract independent signal features. However every technique has advantages and disadvantages and we have yet to reach the stage in which real time analysis can be employed as a single method. For example, PCA is known as a sophisticated technique to suppress artifacts. It specifies principal components (PC) to reconstruct the overall data structure and to remove components with small amplitudes and irregular changes. However, there are cases involving signals that have very low amplitudes as compared to those from other sources that are, nevertheless, necessary for reconstruction of the signals. It is difficult to specify that any or more PCs represent the artifacts as the identification of PCs requires prior knowledge of the artifacts [4, 5]. Therefore research trends have shifted towards decomposition by ICA with higher order statics to determine independence in signals. Since ICA is based on the measurement of signal independence, the noise of the input is amplified by ICA, inhibiting the detection of true EOG components by the spread of Gaussian noise over the components, the method is not recognized as a complete solution to discreminate between artifact components and signal components [6–9].

Morphological Component Analysis, however is not restricted by the concept of independence, which was recently developed and has attracted attention in sparse signal processing as a means to decompose signals and images [10]. Due to the sparseness in bio-signal representation provided by various signal sources, the signal are decomposed with a redundant basis or a mixed overcomplete dictionary [11]. An overcomplete dictionary is described as a collection of different types of mathematical basis functions representing evoked potentials and their associated background noise respectively. Yong et al. [14] demonstrated the removal of artifacts from EEG signals using MCA with a Dirac \(\delta \) function basis (DIRAC), discrete wavelet transform with a Daubechies filter (DWT) and a discrete cosine transform (DCT). The results were subsequently compared with other traditional methods of blind source separation such as AMUSE and EFICA.

In this paper, we focused on a decomposition method for EOG signals using Morphological Component Analysis with DIRAC, UDWT (undecimated discrete wavelet transform) and DCT basis, which was inspired by the approach of Yong et al. [14]. EOG and EEG are similar biological signals and their interference with each other induces the appearance of similar components in both signals. In principle, artifacts are derived from measurement equipment, cables used for signal transmission and amplification, power-supply instabilities, and so on. When biological signals are contaminated, it is difficult to define the signal/noise ratio or to identify artifacts. The proposed method in this has the advantages of facilitating the separation of different types of signals, including the measurement of artifacts from the desired signal(EOG), thereby making analysis of EOG possible.

This paper in organized into the following sections: the concept of MCA is explained in Sect. 2, the EOG measurement is explained in Sect. 3, experimental results are given in Sect. 4 and finally the conclusion and discussion are described in Sect. 5.

2 Morphological Component Analysis Concepts

Biomedical signal decomposition has been performed using various methods such as PCA, Wavelets, ICA and so on [3]. Blind Source Separation (BSS) like ICA relies on the independence of signals, this concept can be extend to a sparse representation in linear analysis to treat bio-medical signals and time course of natural phenomena. BSS algorithms estimate the vector m under the condition that the underlying sources is \(X \in \mathbb {R}^{m*N}\) and the observed signal is taken from the k channel or sensor, \(S \in \mathbb {R}^{k*N}\). The vector length N is the number of samples taken during the entire recording time T, which is given by \(T=N \cdot f_s\) where \(f_s\) is the sampling rate. Assuming that the relevant signals are recorded from channels as independent sources, the mixing of each independent source can be classified by

$$\begin{aligned} Signal (amplitude) =electrical coupling \cdot electrical potential \end{aligned}$$
(1)

and can then be written as:-

$$\begin{aligned} S=BX+W \end{aligned}$$
(2)

where B is the \(k \times m\) coefficient mixing matrix the sources X and W is assumed to be noise, which is generated from power supply, electronic devices for amplification and so on. B and X are the unknown values to be estimated by this method. Various algorithms are organized in the context of different assumptions. The MCA is recently established as a general form to use various combinations of such dictionaries, depending on target applications [10].

We assume that the target EOG signal can be represented by a linear superposition of several signal sources with a time constant [12]. According to the sparsity and morphological diversity of the EOG signal, i.e. the assumption that each of the m sources \(\{ S_1, \cdots , S_m \}\) is sparse in an overcomplete dictionary D, the raw vector \(S_i\) is modeled by a linear combination of p morphological components:

$$\begin{aligned} S_i=\sum _{k=1}^{p}S_{i,k}=\sum _{k=1}^{p} \beta ^i_{k}\phi _{k} \end{aligned}$$
(3)

where \(S_{i}\) is ith source, \(S_{i,k}\) is the time series of the kth morphological component and the \(\beta ^i_{k}\) are the coefficients corresponding to dictionary \(\phi _{k}\).

The dictionary D is a collection of mathematical basis functions representing parameterized waveforms \(\phi _{m}\) with parameter number m. The parameter m determines a set of indexing frequencies, such as time-frequency dictionaries [13, 14]. \(\phi _{m}\) plays a role in discriminating between different signal components. For instance, a morphological component might be categorized as sparse in one particular dictionary but not in another. To find the sparsest of all signals within an augmented dictionary containing all values of \(\phi _{m}\) , we use the following equation:

$$\begin{aligned} \{\beta _{1}^{opt},\ldots ,\beta _{m}^{opt} \} = \arg \min _{\{\beta _{1},\ldots ,\beta _{m} \}} \sum _{i=1}^{m}\parallel \beta _{i} \parallel _{1}+ \lambda \parallel S - \sum _{i=1}^{m} \beta _{i}\phi _{i} \parallel _{2}^{2}. \end{aligned}$$
(4)

To solve this equation(4), we used a numerical solver called the Block-Coordinate Relaxation Method [10] with an appropriate value of \(\lambda \).

We hypothesized that the EOG signal could be decomposed into three categories: irregular spike activities, slow EOG changes and EEG related activities with multi-frequency components. We used the DIRAC, UDWT and DCT dictionaries respectively.

3 Experiment and Measurements

3.1 Experiments Procedure

In our experiment, electrodes were placed on the skin around eyes and designated as channels \(\{ V_z, V^R_u, V^R_d, V^L_u, V^L_d, H^R,H^L \}\). The suffixes u and d represents the upper and lower sides of the eye respectively and suffixes R and, L represents the right and left eyes respectively. The electrode for channel Vz is placed between the eyebrows. All potentials were recorded with a resolution of 0.1 \(\upmu V\), a sampling interval of 1 microsecond, i.e. and a sampling rate of 1 / 1000 s. The subjects were instructed to see a fixed point in the screen and spontaneous eye blinks were monitored.

3.2 EOG Signal Decomposition Using MCA

The MCA concept used for modeling of the EOG signal. The decomposition of EOG signal S as

$$\begin{aligned} S=\sum _{i=1}^{m} \beta _{i}\phi _{i} \end{aligned}$$
(5)

Here the EOG signal can be classified by linear combination of morphological components. It is described by

$$\begin{aligned} S= \sum _{i=1}^{m_1}\beta ^{DIRAC}_i \cdot \phi _i^{DIRAC} +\sum _{i=1}^{m_2} \beta _i^{UDWT} \cdot \phi ^{UDWT}_i +\sum _{i=1}^{m_3}\beta ^{DCT}_i \cdot \phi ^{DCT}_i \end{aligned}$$
(6)

where \(\beta _{DIRAC}\), \(\beta _{UDWT}\) and \(\beta _{DCT}\) are the coefficient vectors corresponding to the complete dictionaries of \(\phi _{d}\), \(\phi _{UDWT}\) and \(\phi _{DCT}\), denoting respectively the Dirac basis, Undecimated discrete wavelet transform and discrete cosine transform. \(m=m_1+m_2+m_3\).

Fig. 1.
figure 1

An example of the MCA decomposition. The original signal is obtained from \(V^L_u\), 6.5 s (6500 samples), including a eye blink. In Block-Coordinate Relaxation Method, we used 3 as \(\lambda \) to be the noise removal condition according to the result of Starck et al. [10]. The top, middle and bottom panels represent the original signal, coefficients of three dictionaries (DIRAC, UDWT and DCT) and the reconstructed signal from the three components. UDWT clearly pursuits reproduces the envelop of the target signal by showing the a large of amplitude change around 2.5 s to be an envelop and there is no after the removal of spike-like activities. The EOG potential is known to change slowly rather than EEGs [1, 2]. This result indicates the combination of DIRAC, UDWT and DCT dictionaries provides an effective method to extract pure EOG signals from artifacts with multi timescales. cc and iter represent respectively correlation coefficient and the number of iterations.

Full size image
Fig. 2.
figure 2

The time evolution of the correlation coefficient (cc) between the original signal and the reconstructed signal (Fig. 1; bottom) during iterations of Block-Coordinate Relaxation Method [10]. Target data is the same as shown in Fig. 1. The decomposition starts to extracts the signal envelope by using UDWT with a cc value higher than 0.7 and without in the absence of any DIRAC and DCT coefficients in DIRAC and DCT, having cc value larger than 0.7. The development of the cc value is attributed develops due to the amplification of the DCT coefficients in DCT fitteding to fast high-frequency cyclic oscillations as well as, and then DCT pursuits remaining spikes. After the value of cc reached a maximumthe peak point of cc, the DIRAC components represent become an over-fitting curve to pursuit the original signal. Once it exceeds the maximum cc, the time evolution goes down drastically.

Full size image
Fig. 3.
figure 3

The resultant frequency distribution of for each part, which is obtained from the reconstructed signals in the individual DIRAC, UDWT and DCT dictionaries of DIRAC, UDWT and DCT and itstheir respective coefficients (Fig. 1; middle) In the power spectrum, UDWT corresponds to the slowest component, DIRAC spreads over all the frequencies with bumps and DCT corresponds to lower frequencies. The DCT frequency range is from the 1- to 9 Hz and, which is covers the ranges of the alpha, delta and theta rhythms of EEG signals.

Full size image

4 Results

In the experiment, EOG signals were obtained to identify the high-amplitude baseline changes attributed to blinking. In this analysis, we constructed a dictionary that included information from the DIRAC, UDWT and DCT dictionaries. This dictionaries were selected according the characteristics of signal. The categorization of morphological components depends on the particular dictionary used, as many dictionaries do not agree on whether a component should be labeled sparse. Consequently these components represent only one type of signal feature. By applying the proposed MCA method with the DIRAC, UDWT and DCT basis, we were able to successfully decompose the EOG signal into it’s separate features. Figure 1 shows that the original signal was successfully decomposed into three parts. This signal was then reconstructed with a high correlation coefficient(0.9841) following 450 iteration of the Block-Coordinate Relaxation Method [10], the numerical method used to apply MCA. The DIRAC components represent spike-like activities, UDWT represent components that change shape slowly, which are believe to correspond to pure EOG signals and the DCT components represent EEG signals mixed with cyclic background noise. As shown in Fig. 2, the correlation coefficient between the original and reconstructed signal monotonically increases, reaching a maximum value and then drastically decreases, thereby upsetting the balance between the three parts.

The power spectrum of the individual parts, shown in Fig. 3, indicates that the UDWT portion is closely correlated with the slowest components and decreases smoothly from the beginning of the measurement. DIRAC portion starts to be a flat line, which develops a bumping distribution over frequencies. A hill-like distribution in the lower (1–9 Hz) frequency ranges appears in the DCT portion of the spectrum, representing alpha, delta and theta rhythms. However noise generated by electronic devices, power lines and other type of noise cannot easily be singled out in spectrum. It is possible that baseline changes arising from power supply instability of the power supply are embedded in the UDWT components causing an overshoot for measurement devices relying on an AC amplifier. This would explain spiking activity as being a mixture of biological artifacts and device-driven artifacts. Higher frequency EEG signals such as gamma, which ranges between 25 to 100 Hz, do not appear in DCT components but may however, be accounted for in the DIRAC results. Specific tendencies will be clarified in further analysis, whcih will include detail analyses of EEG channels. The correlation coefficient between the original and reconstructed signals was measured at 0.989 as averaged over channels \(\{ V_z, V^R_u, V^R_d, V^L_u, V^L_d, H^R\) and \(H^L \}\). This demonstrates that decomposition via the MCA method ia an effective tool to analyze various features of biological signals and further, to remove artifacts from the target signal.

5 Conclusion

Using MCA, we were able to successfully discriminate between signals and background noise in signal decomposition. The signals seem to be EEG signals as well as other artifacts, as determined by differences in their morphological characteristics. In this paper, we analyzed biological signals involving various types of affects, which are influenced by the performance of tasks as well the conditions under which it is performed. The DIRAC dictionary decomposed the signal into spike-like activities, while the UDWT dictionary highlighted slower movements and finally, the DCT dictionary determined that background signals arise from EEG or pure tone signals. By choosing an appropriate value of \(\lambda \) and maximizing the iteration period, it was determined that the iterative algorithm works well for the reconstruction of the signal coefficients. This result attests to the respective effectiveness of artifact removal raw signals and further that the main component of EOG changes smoothly overtime.