A Wavelet-based Energetic Approach for the Analysis of Electroencephalogram

: Electroencephalography (EEG) is the recording of electrical activity along the scalp produced by the firing of neurons within the brain. The main application of EEG is in the case of epilepsy, as epileptic activity can create clear abnormalities on a standard EEG study. EEG signals, like many biomedical signals, are highly non-stationary by their nature. Wavelet analysis has found a prominent position in the investigation of biomedical signals for its ability to analyze such signals, in particular EEG signals. Wavelet transform is capable of separating the signal energy among different frequency bands (i


Introduction
pilepsy is a common brain disorder associated with abnormal neuronal activity.About 1% of the world population suffers from epilepsy and 30% of epileptic patients are not cured by medicine (NINDS 2001).Its major manifestation is the epileptic seizure which may involve a discrete part of the brain (partial) or the whole cerebral mass (generalized).EEG is the recording of electrical activity along the scalp produced by the firing of neurons within the brain.The main application of EEG is in the case of epilepsy, as epileptic activity can create clear abnormalities on a standard EEG study.Careful analysis of EEG records can provide valuable insight and improved understanding of the mechanisms causing epileptic disorders.
The Fourier transform (FT) has been the traditional method applied to the analysis of time series signals for decades.Fourier coefficients are determined by the entire signal support and frequencies are not localized in time, since the infinite basis functions are used in FT.Consequently, Fourier analysis provides only globally time-averaged information, whereas it lacks any local behavior within the signal.Hence, it is suitable for extracting frequency information from stationary signals.
EEG signals like many biomedical signals are non-stationary, and Fourier transform does not give an effective analysis for such signals.For non-stationary signals one method to partly overcome this difficulty is the usage of short-time Fourier transform, STFT (the windowed Fourier transform) in which the signal is multiplied by a sliding window of limited extent, considering the signal as quasi-stationary for such a short period.In essence, STFT extracts several frames of the signal to be analyzed with a window that moves with time.If the time window is sufficiently narrow, each frame extracted can be viewed as stationary so that Fourier transform can be used.With the window moving along the time axis, the relation between the variance of frequency and time can be identified.However the compromise between the temporal and frequency resolution, established by the window size, is the same for all frequencies.
The wavelet transformation is well-suited to representing various aspects of EEG signals such as trends, discontinuities, and repeated patterns where other signal processing approaches fail.Wavelet is an effective time/frequency analysis tool for analyzing transient signals.In the wavelet transform (WT) case, WT employs a windowing technique with variable-size windows.Wavelets are mathematical functions that cut up data into different frequency components, and then study each component with a resolution matched to its scale.They have advantages over traditional Fourier methods in analyzing physical situations where the signal contains discontinuities and sharp spikes.The fundamental idea behind wavelets is to analyze according to scale.
Wavelets are functions that satisfy certain mathematical requirements and are used in representing data or other functions.Approximation using superposition of functions has existed since the early 1800's; Joseph Fourier discovered that he could superpose sines and cosines to represent other functions.However, in wavelet analysis, the scale used to look at data has a special role.Wavelet algorithms process data at different scales or resolutions.If we look at a signal with a large window, we can notice gross features, and if we look at a signal with a small window, we can catch small features.This makes wavelets interesting and useful.
Wavelets have been traced all the way back to Alfred Haar in 1910; however, the starting point of their modern history coincides with two publications in the 1980s by Mallat (1989) and Daubechies (1990).Mallat (1989) identified the important concept of multi-resolution analysis which is the corner stone of modern wavelet theory, while Daubechies (1990) constructed the first orthogonal wavelet bases that were compactly supported.
Wavelets have become a tool of choice for scientists, leading to efficient solutions in time and space frequency analysis problems, as well as a number of other applications.Nowadays wavelet methods of analysis and representation have a significant impact on the science of medical imaging and the diagnosis of disease and screening protocols.Wavelet applications in medicine include, but are not limited to, capillary pressure, coronary artery disease, auditory nerve models, blood flow velocity, ECG timing, distortions and noise detection of abnormalities, heart rate variability, cardiac arrhythmias, ECG data compression, evoke potentials, epileptic seizures and epileptogenic foci, classification of EEG, pathological sounds, ultrasounds, and vibrations (blood flow heart and lung sounds), medical imaging (ultrasonic, magnetic resonance, optical images, computed tomography and others) electromyography signals and sleep apnea (Aldroubi and Unser, 1996).

E
It may be emphasized that wavelet transform has been extensively applied to EEG with different purposes: analysis and characterization of epileptic activity (such as spikes, slow waves, polyspikes, sharp waves, etc.), in the perspective of obtaining clues on the processes underlying the onset of an epileptic attack (Attelis et al., 1997;Bhandari et al., 2007;Gigola et al., 2004;Rosso et al., 2006), development of algorithms for the prediction and on-line/off-line automatic detection of epileptic seizures (Attelis et al., 1997, Gigola et al., 2004, Latka et al., 2003) in order to assist clinitians in monitoring hospitalized patients and reviewing EEG recordings as well as to improve the quality of life of epileptic subjects.
In Inouye et al. (1990) a change of power spectrum in alpha frequency before the spike and wave complexes was reported.Fisher et al. (1992) reported that frequencies above 40 Hz are poorly visualized on conventional EEG scalp recordings, and that they have recorded frequency components up to 150 Hz in digitally recorded EEGs.High frequency increases were largely localized to the region of the seizure focus.In Inouye et al. (1994) structural changes were observed just before spike occurance during seizure.In Magosso et al. (2009) wavelet methods were applied to EEG data obtained from epileptic patients suffering from drug resistance temporal lobe seizures acquired at Bellaria Hospital (Bologna).In this particular paper, Magosso et al. (2009) analyzed the energy distribution of the EEG to determine if it is altered among the different scales of the wavelet representation and exhibits distinct patterns of energy redistribution.In the present work, we exploit the same approach for our data obtained at Kocaeli University's Medical Hospital in Turkey to characterize the epileptic attack in quantitative terms and to obtain indications concerning the genesis of the seizure propagates among adjacent regions of the brain.We have obtained similar results for digitally recorded video EEG recordings for epileptic patients as Magosso et al. (2009).We have shown that energy distribution of the EEG has altered among the different scales of wavelet representation at seizure onset and during the seizure.
The paper is organized as follows.The first section gives an introduction.In the second section, wavelets and energy computation from wavelet coefficients will be introduced, while in the third section analysis of epileptic EEG will be discussed.In the last two sections, data, methodology, application to EEG and results are presented.

Wavelet analysis and energy computation
The wavelet is a quickly vanishing oscillation function localized both in frequency and time.In both continuous and discrete forms of wavelet analysis, the signal is decomposed into scaled and translated versions   where a and b are the scale and translation parameters, respectively, with , ab and 0.
2 fL  (the space of finite energy signals), then the continuous wavelet transform (CWT) of the signal is defined as where ,  means the inner product and the superscript * means complex conjugate.The factor 1 a is used to normalize the energy so that it stays at the same level for different values of a and b.The wavelet function   , ab t  becomes narrower when a is increased, and displaced in time when b is varied.Therefore, a is called the scaling parameter which captures the local frequency content and b is called the translation parameter which localizes the wavelet basis function at time tb  and its neighborhood.CWT at every possible scale a and translation b provides a redundant representation of the signal; hence CWT requires a heavy burden of computations compared to the discrete wave transform (DWT).DWT provides sufficient information both for analysis and synthesis of the original signal, with a significant reduction in the computation time.The filters of different cutoff frequencies can be used to analyze the signal at different scales.The signal is passed through a series of high pass filters to analyze the high frequencies, and it is passed through a series of low pass filters to analyze the low frequencies (Addison, 2002).
The DWT is obtained by discretizing the parameters a and b.We may choose 2, j a The DWT can be written as


 form an orthonormal basis for 2 () L  .In that case, the original signal can be reconstructed from the resulting wavelet coefficients accurately and efficiently without any loss of information (Strang, 1996); the number of wavelet coefficients will be the same as the number of data points in the original signal.In other words, the DWT discards all redundant information in CWT by employing a set of orthogonal basis functions.Furthermore, the DWT may be interpreted in terms of a multi-resolution analysis, where a hierarchy of approximation and details of the signal is constructed in nested subspaces of 2 () L  .Multi-resolution scheme provides an effective way of implementing DWT (Mallat, 1989).In this scheme, the square integrable space 2 () L  is decomposed into a direct sum of the subspaces , where  indicates direct sum, then the subspaces j V are a multi-resolution approximation of the square integrable space 2 L and can be obtained from dilation and translation of a single scale function , ( ). jk t  Thus, the subspaces j W are the orthogonal complementary of the subspaces j V : for all .
j  Given a signal   , ft its multi-resolution decomposition formula at level H is defined as is the mother wavelet, while  is a companion function, called the scaling function: It expresses the synthesis of the original signal from wavelet coefficients.In each level j, the series in ( 6) has the property of complete oscillation (Chui, 1992), which makes the decomposition useful in applications to time localization of events.
By using down-sampling operations along with low-pass and high-pass filtering, signal decomposition as in ( 5) can be efficiently implemented (fast pyramid algorithm).It is important that the two filters are related to each other, and they are known as a quadrature mirror filter.The high-pass filters are associated with the wavelet functions , ( ), jk t  while the low-pass filters are associated with the scaling functions , ( ). jk t


At each decomposition step this pair of filters decomposes the signal into low-frequency components (approximation coefficients), and high-frequency components (details coefficients).Filtering is applied first to the original signal, and then recursively, to the approximation series only.At every iteration, the output of each filter is down-sampled by a factor of 2 (decimation) halving the data each time.The down-sampling is done to give speed to the algorithm, reducing the computation at each iteration geometrically (after J iteration the number of samples being manipulated shrinks by 2 M ).
Let the signal be of finite length, say N and 2. M N  Let s T denote the sampling time.Theoretically, wavelet decomposition may involve all scales from negative to positive infinity as seen in ( 6).In practice, wavelet decomposition is applied to discrete signals.Since there is not enough resolution between two consecutive samples to construct finer details at scales 0, j  the sample values of the signal are generally taken as the signal approximation coefficients at scale 0 j  and the analysis is limited only to positive scales (Siddiqi, 2004).The approximated and detailed signals at scale j will have only 2 Mj  samples each, because of the down sampling operation.Coefficients at each scale j are placed at instant   , 2 0, , 2 1 .
the range of scales that can be investigated is 1, jM  since the decomposition can proceed only until the individual details contain a single coefficient.
If the decomposition is done for all resolution levels M, the wavelet expansion will be 21 /2 , 10 where k starts from 0 since we assume, without loss of generality, that the signal starts from where 0 f is the frequency of the mother wavelet.When the family is an orthonormal basis in 2 () L  , the concept of energy is linked with the usual notions derived from the Fourier theory, and the sum of the square of the coefficients of the series is the energy of the function, i.e. the energy series associated with coefficient series is given by Energy coefficients as computed by ( 9) have a different localization and density over the temporal axis, depending on the scale: at scale j, the coefficients are placed at instants 2 , 0, , 2 1.

j M j s k T k
  Thus, to study and to compare the temporal evolution of energy at different scales, it is necessary to recover for the halved time resolution at each scale due to the down-sampling operation.In the following, two methods which allow uniformly time distributed atoms of energy to be obtained across all scales will be intoduced.In the last section these methods will be applied to the analysis of the EEG data.

Grouping and spreading-out energy coefficients
In this method the original signal will be divided into non-overlapping temporal windows of fixed length . Here 2 J is the number of signal samples falling within the window, and an atom of energy within each window at every resolution level j will be computed.Two different cases have to be considered:

Case a)
Jj  : In this case, at each scale j ( 1, , ), jM  the assigned window contains 2 Jj  energy coefficients.
The atom of energy within the window n ( 0, 1, , 1) W nN  can be computed by grouping all the energy coefficients falling within the window: The centre point of the corresponding time window will be set according to (12), so the first window will be centered on 0, and half of the coefficients contained in it are null.Thus for each scale, a series of where the function   round x rounds x to the nearest integer.
Therefore the total energy within a time window n can be obtained by and the relative energy associated with the resolution j in the time window n can be obtained by , , , ˆjn

Averaging energy coefficients
In this method atoms of moving average energy will be computed.Series of moving average energy can be computed according to the following equation ,, where the functions   fix x rounds x to the nearest integer towards zero.
The total means energy within a time window n   and the relative mean energy is defined as , , , ˆjn

Analysis of epileptic EEG
EEG is a record of electrical potential generated by cerebral cortex nerve cells (Latka et al., 2003).Careful analysis of EEG records can provide valuable insight and improved understanding of the mechanisms causing epileptic disorders.Wavelet transformation, which is known as mathematical prism or microscope, is being exploited for the analysis and proper understanding of the EEG records (Attelis et al., 1997;Magosso et al., 2009;Adeli et al., 2003;Furati et al., 2006;Gigola et al., 2004;Mallat, 1989).Siddiqi et al. (2009) have analyzed EEG's of epileptic and normal persons through energy distribution among different approximate and detailed levels.Any redistribution may indicate change in the characteristics of the EEG signal which in turn may represent specific events in the course of seizure.It may be remarked that if the information of a signal in general, or an EEG signal in particular, carrying specific information consists of different components, each individual component makes a different contribution to the total signal.For instance, each person's voice is of different wavelet frequency scale characteristics; that is, the individual frequencies make a different contribution to the total energy of the voice.Therefore the signal can be distinguished by the characteristics of its energy spectrum.
To set up the correlation between target signal and the wavelet energy spectrum we must build up a model database (Bhandari et al., 2007;Furati et al., 2006;Gencay et al., 2002;Iske and Randen, 2006;Manchanda et al., 2007;Percival and Walden, 2000;Rivera, 2003;Siddiqi, 2004;Siddiqi et al., 2007;Yue et al., 2005;Yue and Tao, 2006) combining a series of specific models that are distinguished from one another.In the model database each model should have a group model signal.The energy spectrums of these model signals are obtained by applying wavelet transform.It may be emphasized that wavelet transform has been extensively applied to EEG with different purposes: analysis and characterization of epileptic activity (such as spikes, slow waves, polyspikes, sharp waves, etc.), in the perspective of obtaining clues on the processing underlying the onset of an epileptic attack (Gigola et al., 2004;Attellis et al., 1997;Rosso et al., 2006).
Development of algorithms for the prediction and on-line/off-line automatic detection of epileptic seizures (Attellis et al., 1997;Gigola et al., 2004;Latka et al., 2003) in order to assist clinitians in monitoring hospitalized patients and reviewing EEG recordings as well as to improve the quality of life of epileptic subjects.
In Magosso et al. (2009) wavelet methods have been applied to EEG data obtained from epileptic patients suffering from drug resistance temporal lobe seizures acquired at Bellaria Hospital (Bologna).The objective of this study is to analyse whether the energy distribution of the EEG is altered among the different scales of the wavelet representation and exhibit distinct patterns of energy redistribution.This information may be exploited to characterize the epileptic attack in quantitative terms and to obtain indications concerning the genesis of the seizure propagates among adjacent regions of the brain.

Data and methodology
In this study we collected data of 30 epileptic patients suffering from temporal lobe seizures at Medical School's Hospital of Kocaeli University.Data were recorded by the neurology laboratory of the hospital by using the International 10-20 recording system (see Figure 1) using a digital EEG recording device.The scalp EEG recording was used since most hospitals have limited or no usage of intracerebral EEG recordings.Using scalp EEG recordings we get the same major modifications of frequency and of energy distribution as with intracerebral EEG recording.The signals were sampled at 200 Hz and stored on a 32-64 channel computerized video-EEG system.Each patient was hospitalized for EEG video monitoring for several hours, and seizures in their EEG were detected visually.The data were cut into small pieces containing a few minutes before and after the seizures.The multi-resolution wavelet analysis was applied to all channels using the Daubechies order 4 wavelet (Db4).Each EEG signal was decomposed into seven resolution levels, in order to consider all frequency bands which are commonly considered in the analysis of EEG signals.

Application and concluding remarks
Figure 2 shows the 7-level decomposition of T6-O2 channel from one patient, with seizure activity starting at 115 t  s.The signal(s), the seven levels of details (D1-D7) and the residual approximation (A7) are shown in the Figure .The approximation and detail records are reconstructed from the scaling coefficients and wavelet coefficients, respectively.The original signal is the superimposition of details D1-D7 and approximation A7.
Figure 3 shows the power spectral density (PSD) of the details and approximation estimated with the Welch's method.The frequency components captured by details move from the high-frequencies towards the low-frequencies as scale increases from 1 to 7 (35-100 Hz at scale 1 vs. 0.75-1.5 Hz at scale 7), whereas approximation A7 contains all the residual lower frequency information (<0.75 Hz) of the signal for our data.From Figure 3 the spectra clearly indicate the frequency content captured by each detail and the approximation.(where s T is the sample time) and the energy of all details were computed within each window.In particular, all energy atoms at the same resolution level, contained within the same window, were summed up to have energy at that particular level (see ( 14)).Our results show that EEG recordings of epileptic patients exhibit a rearrangement of relative energy among the different frequency bands just before and during seizure.Our analysis confirms this result specifically for patients with temporal lobe epilepsy.However this relocation of the energy is not the same at all channels.

0
t  .Scalograms are the graphical representation of the square of the wavelet coefficients for the different scales.They are an isometric view of the sequence of the wavelet coefficients versus wavelength.A scalogram clearly shows more details, identifies the exact location/time and detects low-frequency cyclicity of the signal.The scalogram surface highlights the location (depth) and scale (wavelength) of dominant energetic features within the signal.The combination of the various vectors of coefficients at different scales (wavelengths) forms the scalogram.The depths (location/time) with the largest (strongest) coefficients indicate the position where the particular wavelength change is taking place.The scalogram provides a good space-frequency representation of the signal.The scale 'a' is related to the frequency by 0 f a f

Figure 3 .
Figure 3. Power spectrum density (PSD) of the details and approximation.An important aspect in the analysis of EEGs during epilepsy is the energy redistribution among the different details; this redistribution may indicate changes in the characteristics of EEG signals which, in turn, may represent specific events in the course of seizure.In order to characterize the temporal evolution of the energy redistribution of EEG signals, the signal was divided into moving average windows each with a duration 2 J s LT (where s T is the sample time) and the energy of all details were computed within each window.In

Figure 4 .Figure 5 .
Figure 4. Total energy of the signal(s) and the relative energies contained in each detail D4-D7.

Figure 6 .
Figure 6.2-D scalogram of T6-O2 channel of the signal(s).In Figure 5 snapshots of the relative energy distribution in all EEG channels for a patient with temporal lobe epilepsy just before seizure ( 65 t  s), at seizure beginning ( 105 t  s), and during several seconds of seizure with reference to detail 4 can be seen.Different colors indicate different relative energy; dark blue for colored pictures (light grey for black/white pictures) show values close to zero, while dark red (light grey) shows approximately 0.6 or higher values.Figures 6 and 7 show, respectively, 2-D and 3-D scalograms of T6-O2 channel of the signal(s) and confirm the seizure.The depths (location/time) with the largest (strongest) coefficients indicate the position where the particular wavelength change is taking place.
Extending decomposition over all resolution levels, the complete wavelet expansion can be obtained as are scaled and translated versions of the original scaling function  . t 