MULTICHANNEL SPECTRAL PATTERN SEPARATION - AN EEG PROCESSING APPLICATION Tomasz M. Rutkowski1,∗,† , Andrzej Cichocki1 , Toshihisa Tanaka2,1† , Danilo P. Mandic3,† , Jianting Cao4,1,† , Anca L. Ralescu5,∗ 1
Brain Science Institute RIKEN, Saitama, Japan
[email protected] 2 Tokyo University of Agriculture and Technology, Tokyo, Japan 3 Imperial College London, London, UK 4 Saitama Institute of Technology, Saitama, Japan 5 University of Cincinnati, OH, USA ABSTRACT A problem of information separation in multichannel recordings is important in engineering applications such as brain computer/machine interfaces (BCI/BMI). Whereas this problem is not entirely new, engineering approaches connecting the mental states of humans and the observed electroencephalography (EEG) recordings are still in their infancy, mostly due to problems with electrophysiological denoising. The electrophysiological signals captured in form of the EEG carry brain activity in form of the neurophysiological components which are usually embedded in much higher power electrical muscle activity components (electromyography - EMG; electrooculography - EOG; etc.). In this paper we present an approach to remove muscular interference caused by eye-movements from EEG recorded during auditory experiments in an eight channel recording setting. This is achieved by analyzing the correlation of the oscillatory modes within a multichannel signal in the Hilbert domain. Simulations in a real world auditory BCI setting support the analysis. Index Terms— biomedical signal processing, electroencephalography, electrooculography, electromyography, spectral analysis 1. INTRODUCTION The problem of separating brain electrical activity, captured in form of the EEG, from other electric signals in human body (EOG, EMG, ECG, etc.) and environmental interference (electrical devices) is a serious obstacle to many neuroscience experiments (see Figure 1 with EEG recordings contaminated by EOG). Solutions to those problems are a prerequisite to the analysis of the information processing mechanism of the brain [1, 2, 3]. This work focuses on the removal of the muscular artifact, which carries significant power, from useful EEG data. EEG signals recorded on the scalp levels are usually highly contaminated by noise due to a very low level of neurophysiological signals and to high power of general electrophyciological signals, such as EMG or EOG, as well as the presence of different devices in the recording environment which cause electromagnetic interference. To tackle these problems, the current study proposes to use ∗ The
author was supported in part by AOARD grant. author was supported in part by JSPS & Royal Society under the Japan-UK Research Cooperative Program 2007 & 2008. † The
978-1-4244-2354-5/09/$25.00 ©2009 IEEE
373
empirical mode decomposition (EMD) [4], a new technique to decompose EEG signals. Once the components from each recorded EEG channel (see Section 2. for recording and channels setup details) are obtained we aim to compare them in the spectral domain in order to separate common activities propagated across the human scalp from neurophysiological ones. The aim of this paper is therefore twofold: (i) to introduce a novel computational framework based on Empirical Mode Decomposition (EMD), convenient for simultaneous data conditioning and information separation within the EEG data, and (ii) to provide an insight into possible clustering of different electrical activities originating in different sources of human body, in particular, to separate eye movements from neurophysiological signals. 2. METHODS To this end, experiments were conducted in Advanced Brain Signal Processing Laboratory of RIKEN Brain Science Institute, Japan. Subjects participated in a spatial audio source direction localization experiment in a surround sound 7.1 speakers system environment. The EEG electrodes were connected to the head channels Fp1, Fp2, C3, C4, F3, F4, T7, and T8, as in extended 10/20 EEG recording systems [5]. Such experimental paradigm causes subjects often to move their eyes causing ocular interference in recoded EEG. For this reason as reference additionally two ocular channels vEOG, hEOG capturing vertical and horizontal eye movements, were recorded. 2.1. Empirical Mode Decomposition (EMD) EMD utilizes empirical knowledge of oscillations intrinsic to a signal in order to represent them as a superposition of components, called as intrinsic mode functions (IMF), with well defined instantaneous frequencies. To obtain an IMF from a single channel EEG it is necessary to remove local riding waves (abrupt changes in time frequency representation) and asymmetries, which are estimated from local envelopes of minima and maxima of the waveform. The technique of finding IMFs corresponds thus to the separation of band limited semi-orthogonal components from recorded EEG. It also corresponds to eliminating riding-waves from the signal, which ensures that the instantaneous frequency will have no fluctuations caused by an asymmetric wave form. In each cycle, the IMF is defined by zero crossings and involves only one mode of oscilla-
ICASSP 2009
Recorded EEG
Fp1
100 0
Fp2
−100
C3
0
100
0
0
40 20 0 −20 −40
0
1500
50 0 −50
F3 F4
0
1500
−50 40 20 0 −20 −40 40 20 0 −20 −40
0
1500
0
0
−50
−50
200
200
0
0
0
1500
1500
50
50
0
0
0
0
1500
100 0 −100
−50
0
1500
50
0
0 1500
0
1500
100 0 −100
1500 [samples]
(ii) at any point, the mean value of the envelope defined by the local maxima and the envelope defined by the local minima is zero.
1500
Since IMF represents an oscillatory mode within a signal; its periods, which are defined by zero crossings, correspond to the only one mode of oscillation. Both the amplitude and frequency of this oscillation may vary over time, in other words, the oscillation is not necessarily stationary nor narrow-band. The process of extracting an IMF from a signal x(t) is called “the sifting process” [4] and consists of the following steps:
0
1500
0
1500
(i) determine the local maxima and minima of x(t); (ii) generate the upper and lower signal envelope by connecting those local maxima and minima respectively by an interpolation method (e.g., linear, spline, piece-wise spline [4, 6]) (in this paper the linear method was chosen);
−50
1500
(iii) determine the local mean m1 (t), by averaging the upper and lower signal envelope; 0
(iv) subtract the local mean from the data: h1 (t) = x(t) − m1 (t).
1500
0
1500
0
1500
Ideally, h1 (t) satisfies the criterion of an IMF, however, typically this procedure needs to be repeated until the first IMF is extracted. In order to obtain the second IMF we applied the sifting process to the residual ε1 (t) = x(t)−IMF1 (t), obtained by subtracting the first IMF from x(t); the third IMF is in turn extracted from the residue ε2 (t) and so on. The decomposition is complete when two consecutive sifting results are similar; the empirical mode decomposition of the signal x(t) may be written as:
0 0
1500
−50 50
100 0 −100 0
0
(i) the number of extrema and the number of zero crossings are either equal or differ at most by one;
50
100 0 −100
100 0 −100
−50
100
−100 1500 0
0
1500
0
50
−50
0
100 0 −100
−200 1500 0
0
1500
50
0 1500
0
20 0 −20
−50
−200
T7
1500
50
0
T8
40 20 0 −20 −40
1500
50
vEOG
1500
50 0 −50 0
hEOG
0
−100 1500 0
0
tion, thus not allowing complex riding waves. Notice that an IMF is not limited to be a narrow-band signal, as it would be in the classic Fourier or wavelets decompositions. In fact, an IMF can be both amplitude and frequency modulated simultaneously, as well as non-stationary or non-linear. The EMD decomposes a signal in hand into a number of IMFs [4] (oscillatory modes); which satisfy the following two conditions:
Denoised EEG
50
0 −100 1500 0
100
−100
C4
EOG interferences 100
x(t) =
0 0
1500 [samples]
−50
n X
IMFk (t) + εn (t),
(1)
k=1
0
1500 [samples]
Fig. 1. Illustration of the muscular interference problems for eight EEG channels: Fp1, Fp2, C3, C4, F3, F4, T7, and T8. The reference recordings capturing vertical (vEOG) and horizontal (hEOG) eye movements recorded simultaneously and depicted in two bottom panels. The results of the proposed approach (see Section 2) of EEG and muscular artifacts separation are presented in the middle and right columns. The middle column presents “clean” muscle activity, which is an expected result for EOG and an interference in case of EEG. The right column presents only the electrophysiological activity. The spectral analysis of the above signals is presented in Figure 2.
374
where n is the number of extracted IMFs, and the final residue εn (t) is either the mean trend or a constant. Note that the IMFs are not guaranteed to be mutually orthogonal, but often are close to orthogonal; it is also noteworthy that IMFs are adaptive, that is, two independent realizations of a signal with the same statistics may have a different number of IMFs. 2.2. Hilbert-Spectral Clustering of EMD Components In order to compare all IMFs extracted from the analyzed channels (two EOG and eight EEG in this paper) we propose to cast them separately to Hilbert spectra domain in order to capture the detailed content (intrinsic frequency tracks) carried by all of them. The amplitude and phase ridge traces of all IMFs (note that adaptive nature of EMD may result in different numbers of IMFs in each channel) are combined together and correlated. From the IMFs the corresponding time–frequency representations can be produced by applying the
0
Fp1
10
Spectra of original EEG
−2
20 40 60 80 100
0
10
−2
10 20 40 60 80 100
0
10
−2
10
Fp2 C3 C4
20 40 60 80 100
0
20 40 60 80 100
0
20 40 60 80 100
−2
10
20 40 60 80 100
0
10
−2
10
Using the above procedure in a single channel mode, the EEG signals from chosen electrodes could be decomposed separately, thus forming subsets of IMFs, from which low frequency drifts and high frequency spikes can be removed. The most “interesting” part of EEG is usually in the middle of frequency range. To analyze multichannel EEG signal sets recorded synchronously in a single experiment we propose to decompose all channels separately preventing possible information leakage among the channels.
0
10
−2
10
0
20 40 60 80 100
0
0
10
−2
10
10 F3
0
ωj denotes an instantaneous frequency. The Hilbert transform allows us to observe the variable amplitude and the instantaneous frequency in a form of very sharp and localized functions of frequency and time (in contrast to Fourier expansion, for example, where frequencies and amplitudes are fixed for their bases). Such an approach is very suitable for the analysis of non-stationary EEG and modeling of common/sychronized activities within certain channels.
0
0
−2
10
0
20 40 60 80 100
0
10
−2
10
0
20 40 60 80 100
0
20 40 60 80 100
0
20 40 60 80 100
0
20 40 60 80 100
For this end, use correlations coefficients of Hilbert amplitude and frequency traces as “a distance measure” in order to capture spectral similarity across the IMFs. Once the cross–correlation analysis is performed for all Hilbert transformed IMFs from all analyzed channels, a hierarchical cluster analysis using a set of dissimilarities for the n objects to be clustered is performed [7] (using “R” package [8]) for amplitude and frequency ridges separately. Initially, each vector representing amplitude or frequency ridges values is assigned to its own cluster and then the algorithm proceeds iteratively, at each stage joining the two most similar clusters. Such procedure continues until there is just a single cluster. At each stage distances between clusters are recomputed by the Lance–Williams dissimilarity update formula with a single linkage clustering method. This method is closely related to the minimal spanning tree concept and it adopts a “friends of friends” strategy for clustering [7] (see results in Figures 2 and 3).
0
10 F4
20 40 60 80 100
−2
10
0
−2
10
0
20 40 60 80 100
0
10
0
−2
10
10 T7
0 0
10
−2
10
0
20 40 60 80 100
0
0
10
−2
10
10 T8
k=1
10 0
−2
10
0 vEOG
Separated spectra
0
10
−2
10
10
20 40 60 80 100
0
10
0
−2
10
10
A result of such procedure in the frequency domain is presented in Figure 3 where the sets of clusters are visualized. The first set is for distances below 0.4 and it defines a relatively compact cluster of similar components across the EEG channels. Those components are classified as similar and originating from very strong EOG interference.
−2
10
0 hEOG
Hilbert transform to each component [4], time-frequency representation: „ Z « n X IMFk (t) exp i ωk (t)dt , (2) R(t) =
20 40 60 80 100
0
10
−2
10
0
20 40 60 80 100
0
20 40 60 80 100 frequency [Hz]
0
10
−2
10
0
20 40 60 80 100 frequency [Hz]
Fig. 2. Comparison of the original contaminated spectra of recorded EEG (the left column) and separated EOG with clean EEG (the right column) within each channel with utilization of the method described in Section 2. The red/dashed–line spectra represent interferences which were removed from the signals while black/solid–lines represent the remaining signals of interest. It is interesting to note that the proposed approach with EMD analysis allows us to keep underlying EEG activity in the frequency bands dominated by EOG. This would be not possible to achieve with classical filtering methods - the whole band would be classified as interference.
375
The reconstruction (summation within each EEG channel separately) for those IMFs is presented in the middle panel of Figure 1 and it clearly confirms the hypothesis of EOG related activity. The remaining set of clusters with distances above 0.4 and with rising distances (not as compact cluster as for previous EEG case) represent IMFs carrying neurophysiological signals only and the resulting reconstruction is presented in the right panel of Figure 1. Spectral analysis of the proposed separation is visualized in Figure 2 where it is worth to note that the proposed method spares low frequency spectral content of neurphysiological signals under very strong EOG interference (see the range 0Hz to about 10Hz which is not completely “zeroed”, what could happen with classical filtering approaches). A comparison with contemporary EEG denosing techniques (due to space constrains not discussed here) as independent component analysis (ICA) shown also strength of the proposed result due to nonlinear, nonstationary and convolutive nature of signal mixing environment (brain, muscle and skin tissues).
3. CONCLUSIONS
(1−correlation) as a distance measure of IMF frequency ridges 0.0
0.2
C3_imf1 T7_imf1
0.4
0.6
A framework to separate interfering muscle activity (EOG in this paper) from EEG has been presented. This has been achieved by proposing a novel EEG decomposition technique, which allows a flexible sub-band signal decomposition while preserving the nonlinear and non-stationary features of the signals which is very crucial for brain activity analysis. The so obtained components from each EEG channel processed separately have been further transformed to the Hilbert domain and compared within amplitude and phase domains using the clustering technique in order to identify those similar (spectrally correlated) across channels. The resulting reconstruction has allowed us to separate common non-EEG related interferences from underlying brain activity in the data–driven signal processing approach without information leakage between channels. The proposed approach was tested in several EEG recording sessions in a multiple subjects confirming the presented here results. This is a step forward in EEG signal processing applications which could be useful primarily for creating user friendly brain– machine–interfaces that would be less susceptible to common interferences resulting from human body activity.
C4_imf1 vEOG_imf2 Fp2_imf2
F3_imf1 F3_imf2 Fp1_imf2 hEOG_imf1 F3_imf4 C3_imf2 T8_imf2 C3_imf3 C3_imf4 F4_imf3 T7_imf3 F4_imf1 T8_imf4 vEOG_imf4 T7_imf2 T8_imf1 C4_imf3 F3_imf3 T7_imf4 F4_imf2 hEOG_imf2 hEOG_imf3 vEOG_imf3 T8_imf3 C4_imf2 Fp2_imf7 Fp2_imf4 Fp2_imf3 Fp1_imf3 T8_imf5 C4_imf4 F4_imf5 C4_imf5 hEOG_imf4 vEOG_imf6 Fp1_imf4 hEOG_imf5 F3_imf8 vEOG_imf5 Fp1_imf5 T7_imf8 hEOG_imf7 T8_imf6 Fp1_imf6 F4_imf7 hEOG_imf6 Fp1_imf9 C4_imf7 C4_imf8 C3_imf7 T7_imf7 C3_imf9 Fp2_imf9 F3_imf7 F3_imf6 T7_imf9 F4_imf8 vEOG_imf7 T8_imf7 C3_imf6 C3_imf10 C3_imf5 hEOG_imf8 Fp2_imf8 C3_imf8 C4_imf6 Fp2_imf6 F4_imf6 Fp2_imf5 Fp1_imf7 T7_imf6 F3_imf5 T7_imf5 Fp1_imf8 F4_imf4 vEOG_imf1 Fp1_imf1 Fp2_imf1
4. REFERENCES
Fig. 3. Clustering of Hilbert spectral frequency ridges based on IMFs, revealing similarity among several components from different channels. Note the groups of clusters resulting from crosscorrelation analysis of frequency ridges. The components having spectral distances to the similar components in other electrodes lower than 0.4 were considered as common EOG activity. The resulting reconstruction of neurophysiological signals and “clean” EOG is presented in right column of Figure 1.
376
[1] C. A. Joyce, I. F. Gorodnitsky, and M. Kutas, “Automatic removal of eye movement and blink artifacts from eeg data using blind component separation.,” Psychophysiology, vol. 41, no. 2, pp. 313–325, 2004 Mar. [2] T. M. Rutkowski, A. Cichocki, A. L. Ralescu, and D. P. Mandic, “Emotional states estimation from multichannel EEG maps,” in Advances in Cognitive Neurodynamics ICCN 2007 Proceedings of the International Conference on Cognitive Neurodynamics, Rubin Wang, Fanji Gu, and Enhua Shen, Eds. 2008, Neuroscience, pp. 695–698, Springer Berlin & Heidelberg. [3] D. Looney, L. Li, T. M. Rutkowski, D. P. Mandic, and A. Cichocki, “Ocular artifacts removal from EEG using EMD,” in Advances in Cognitive Neurodynamics ICCN 2007 Proceedings of the International Conference on Cognitive Neurodynamics, Rubin Wang, Fanji Gu, and Enhua Shen, Eds. 2008, pp. 831– 834, Springer Berlin & Heidelberg. [4] N. Huang, Z. Shen, S. Long, M. Wu, H. Shih, Q. Zheng, N. C. Yen, C. Tung, and H. Liu, “The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis,” Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, vol. 454, no. 1971, pp. 903–995, 1998. [5] E. Niedermeyer and F. L. Da Silva, Eds., Electroencephalography: Basic Principles, Clinical Applications, and Related Fields, Lippincott Williams & Wilkins, fifth edition, 2004. [6] T. M. Rutkowski, A. Cichocki, and D. P. Mandic, Signal Processing Techniques for Knowledge Extraction and Information Fusion, chapter Information Fusion for Perceptual Feedback: A Brain Activity Sonification Approach, pp. 261—274, Signals and Communication. Springer, April 2008. [7] F. Murtagh, “Multidimensional clustering algorithms,” COMPSTAT Lectures, vol. 4, 1985. [8] R Development Core Team, R: A Language and Environment for Statistical Computing, R Foundation for Statistical Computing, Vienna, Austria, 2008, http://www.R-project.org.
