BISPECTRAL ANALYSIS OF SCALP ELECTROENCEPHALOGRAMS: QUADRATIC PHASECOUPLING PHENOMENON IN DETECTING BRAIN TUMOR

Since the power spectral analysis of a non-Gaussian process generated by a nonlinear mechanism e.g., E EG does not provide much information on the underlying nonlinear dynamics due to the lack of phase information, the higher-order statistics such as th e bispectra are used to better understand the under lying nonlinear dynamics e.g., the quadratic phase coupli ng henomena. The quadratic phase couplings have be en observed in the EEG by the researchers over a decad e for many diagnostic applications such as epilepsy , sleep, mental states. This study discusses the use of bispectral analysis of the EEG recorded from the posterior region of the head of the brain tumor pat ient in quantifying the quadratic phase couplings t o indicate the presence of the tumor. The Bicoherence I d x (BCI) or simply the Bicoherence (BIC) has be en used for the purpose. Self-couplings (around 27-52% ) in the [8-13] Hz (alpha) band and phase couplings (around 23-42%) in the [1-8] Hz (delta-theta) band have been observed for the normal subjects while on ly self-couplings (around <6.5% and around 40-53%) hav e been seen in both bands for the brain tumor patients. Significant lowering of coupling strength s (from 38.15% (±12.76%) to 3.51% (±3.28%)) in the alpha band and mild increase of them (from 32.76% ( ±18.73%) to 45.49% (±17.49%)) in the delta band have been observed for the brain tumor patients. Th e Power Ratio Index (PRI) based on the power spectr um is only statistically inferior (p>0.05) to the BIC in discriminating the brain tumor case from the nor mal one.


INTRODUCTION
Diagnosis and subsequent (early) treatment are either missed or delayed in 69% of the brain tumor cases due to the fact that the most of the brain tumor symptoms are highly misleading and around 26% of these cases suffer a delay of more than a year before proper diagnosis (MFBTRI, 2013). Once the brain tumor symptoms are found, the advanced neuroimaging techniques such as MRI and CT or biopsy are not immediately suggested due to the following facts: They are either costly or invasive or do involve risks like hazardous radiation, especially in case of children, pregnant women and patients with implant devices (Black, 2010). Since an early treatment increases the survival rate, a better method that does not involve much cost, risks or complexity is required to reduce the delay in the diagnosis of the brain tumors (Black, 2010). One such option is the use of the scalp Electroencephalograms (EEGs) (Fattal-Valevski et al., 2012). Qualitative investigations on the use of scalp EEG for the diagnosis of brain tumors are numerous (Accolla et al., 2011;Bagchi et al., 1961;Hartman and Lesser, 2012;Walter, 1936;Watemberg et al., 2002). However the quantitative AJAS Science Publications works such as automated classification or localization using the quantitative (scalp) EEG (qEEG) features are very few (Chetty and Venayagamoorthy, 2002;Habl et al., 2000;Karameh and Dahleh, 2000;Murugesan and Sukanesh, 2009;Nagata et al., 1985;Silipo et al., 1999). Some of these literatures Silipo et al., 1999) are interesting and they suggest the use of nonlinear analysis of the scalp EEG for the purpose.
The power spectrum does not provide much information about the non-Gaussian processes and the processes generated by the nonlinear mechanisms due to the lack of phase information (Nikias and Raghuveer, 1987). The conditions of non-Gaussianness and nonlinearity are often met in many practical applications, e.g., the biological signals and the biological systems generating them. Hence the analysis of such processes requires the use of the higher-order spectra such the bispectra as they provide more information such as the deviation from Gaussianity and linearity, phase couplings. In this study the application of bispectrum, which unveils the quadratic phase coupling phenomena arising from the second-order nonlinearity, has been considered.
The bispectrum has been extensively used to analyze the EEG time series to discriminate various pathological conditions, metal tasks, physical (sleepy, anesthetic) states (Hosseini et al., 2010;Saikia and Hazarika, 2011;Swarnkar et al., 2010;Venkatakrishnan et al., 2011;Yufune et al., 2011). It has been shown in these papers that the quadratic phase couplings in EEG do change in accordance with the physical and physiological conditions of the brain. As the self-emitting oscillations of the complicated neural network of the brain (Buzsaki and Draguhn, 2004) can be considered non-Gaussian (Nikias and Raghuveer, 1987), the application bispectrum to EEG is quite intuitive. As the presence of a brain tumor affects the brain both physically and physiologically, the expectation of changes in the (nonlinear) dynamics of the brain and subsequently the use of bispectrum for quantifying them are intuitively understandable.
In this study, the bispectral analysis of the scalp EEG recorded from the posterior (parietal and occipital) regions of the heads of both the brain tumor patients and normal subjects has been presented. The Bicoherence Index (BCI) or simply Bicoherence (BIC) is used to quantify the quadratic phase couplings in the EEG records. The results obtained have been tested for their statistical significance in characterizing the presence of the brain tumor. A power spectral feature, namely the Power Ratio Index (PRI) has also been analyzed for comparison.

Acquisition of Data
Nineteen-channel Common Average Reference (CAR) montage EEG records, namely Fz-REF, Cz-REF,  Pz-REF, FP1-REF, FP2-REF, F3-REF, F4-REF, C3-REF, C4-REF, P3-REF, P4-REF, O1-REF, O2-REF,  F7-REF, F8-REF, T3-REF, T4-REF, T5-REF and T6-REF in the standard 10-20 electrode system, with eyes closed, from 67 brain tumor patients and 42 normal subjects are obtained at a sampling rate of 256 Hz for about 15-20 min. The age is approximately uniformly distributed ranging from 14 to 53 years with a mean of 31.36 years and a median of 30 years for the brain tumor patients and ranging from 16 to 53 years with a mean of 33.34 years and a median of 34 years for the normal subjects. The majority of the brain tumor cases fall in the age group between 20 and 40 years.

Preprocessing
A 50 Hz FIR notch filter and an IIR EMG filter are used to eliminate the 50 Hz power line interference and the EMG artifacts, respectively. Certain artifacts such as eye movements, eye rolling and essential tremors are still present in the EEG records. Generally the removal of these artifacts by the traditional method of filtering is not possible since they fall in the useful bandwidth of cerebral EEG.
The method of Independent Component Analysis (ICA) is the best solution for this situation (Mammone et al., 2012;. A recently proposed ICA technique known as the Modified Wavelet ICA (MwICA)  was used to remove these artifacts.
From this set of "clean" EEG records, ten min (10×60×256 = 1,53,600 data points) of only the seven posterior channels, namely, Pz, P3, P4, O1, O2, T5 and T6, ('REF' has been dropped for convenience), of each EEG record are retained for the proposed analysis. All these 7-channel, 10 min EEG records are then bandpassfiltered to 1-40 Hz as this is the band of interest generally for most of the clinical applications (including this proposed analysis) (Picot et al., 2009).

Basics of Bispectrum and Bicoherence
The concept of the bispectrum and its properties are extensively presented in the literature (Nikias and Raghuveer, 1987).
The bispectrum, B x (ω 1 , ω 2 ) of a discrete, real, stationary (univariate) random process, X (t) = {x(t)}, t = 0,±1,±2,... with zero mean is defined as the two-AJAS Fourier transform of X(t) and * represents the complex conjugate. The properties that are worth mentioned here are: (i) the bispectrum is generally complex and has the magnitude and phase, (ii) the complex plane of the bispectrum has 12 symmetric regions, (iii) the bispectrum of a (linear) Gaussian process is everywhere zero, (iv) the bispectrum of a linear non-Gaussian process is everywhere constant and (v) the bispectrum of a nonlinear process is varying and peaking due to the phase couplings arising from the nonlinearity. The property (ii) implies that there is a "principal" region only in which the bispectrum needs to be computed. This principal region is the triangular region, 2 1 0, ω ≥ ω and 1 2 ω + ω ≤ π . A quadratically nonlinear system is the one with second-order nonlinearity e.g., a Volterra-Wiener system of order 2. The phenomenon, where the phases have the same relations as the frequencies, is called the quadratic phase coupling. Generally, in case of a frequency triplet in the output of a system, there are three possible cases: (i) when they are not harmonically related and are with random phases (i.e., no phase-coupling), the bispectrum is identically zero everywhere, (ii) when they are harmonically related with same phase relationship as frequency (i.e., phase-coupling), the bispectrum peaks due to the quadratic phase coupling and (iii) when they are harmonically related with random phases, a situation known as the frequency-coupling, the bispectrum may still peak under certain frequency and amplitude combinations. In all the three cases, the power spectrum invariably shows three peaks at these three frequencies. Due to this ambiguity, the estimated bispectral peak, especially from a finite data set, needs to be statistically tested for its significance in quantifying the presence of phase coupling. This ambiguity is overcome in the normalized bispectrum i.e., bicoherence (Elgar and Guza, 1988).

Estimation of Bispectrum and Bicoherence
The bispectral estimation methods are broadly classified into two: The conventional or non-parametric methods and parametric methods. Two conventional or non-parametric methods, namely, the direct and indirect methods are generally found in the literatures. The parametric methods provide higher frequency resolution than the conventional methods and however suffer from the difficulty in determining the model order parameter. The conventional methods are popular for their ease of implementation. In this study, the direct class of the conventional methods is used. The direct method, which is based on the Fast Fourier Transform (FFT), is extensively described in (Kim and Powers, 1979;Nikias and Raghuveer, 1987). The bicoherence is the normalized version of the bispectrum (Kim and Powers, 1979). Given a discrete-time sequence, x(n), n = 0,1, 2,..., N -1, where N is the cardinality of the sequence, the sequence is segmented into K short sequences, x i (n) = x[n + iM-(i +1)(D-1)] , n = 0,1, 2,...M -1 and i = 0,1, 2,..., K-1 of M samples each with an overlap of D samples and each segment, x i (n) is centered The bicoherence is then given by Equation 1 and 2: Where: ( ) The X i (k), i = 0,1, 2,...,K-1 in Equation 3 are the discrete Fourier transforms of x i (n). Each segment, x i (n) can be windowed, if desired, in order to reduce the spectral leakage. The possible window functions include Hamming, Hann and Blackman. For a given N, the overlapping not only reduces the variance in the estimated bispectrum, as it increases the number of segments (i.e., K) for the ensemble averaging, but also provides a chance for a larger segment length (i.e., M) thereby improving the frequency resolution. A segment size of M samples gives a frequency resolution of ∆f = f s /M where f S is the sampling rate of the time series in Hz and an ensemble average over K segments yields a variance proportional to N/K 2 . Hence, given a set of f S and N the values of M and D must be properly selected so that the desired frequency resolution is achieved and at the same time the variance in the estimated bispectrum is as low as possible.
The estimation of the bispectrum requires the EEG signal to be stationary (Siu et al., 2008). Since the EEG signals are non-stationary, they are often analyzed in segments in order to ensure the stationarity using a criterion of weak stationarity, which requires that the statistical parameters up to certain order remain (practically) constant over the entire period of the segment (Blanco et al., 1995). The most popularly used weak stationarity is the second-order stationarity which requires the second-order statistics, mean and standard deviation (std), to remain constant at some prescribed significant level (e.g., 5%). The test for detecting the stationary segment length was carried out as follows. Each channel record was split into overlapping segments of length 1s (256 samples) called the bins. The overlap was 0.5s (128 samples). The means (bin means) and stds (bin stds) of each of these bins were computed. Then the means and stds of the bin means and bin stds were computed and the constancy of the bin means and bin stds was tested by checking whether their values lied within the 95% confidence interval of their means i.e., within mean±2std for various lengths of consecutive bins. The percentage of number of consecutive bins that met this criterion was calculated. All the consecutive bins of lengths, 2s, 3s and 4s remained (weakly) stationary. For the bin length of 5s and above, the number of consecutive bins that failed the stationarity test was above 6%. The bin length of 4s is chosen for this analysis as this ensures both the quasi-stationarity and the sufficiency of the data points to estimate the bispectrum (Miller et al., 2004). Thus there are totally 10,050 (67 cases×10 min×60 sec/4 sec) epochs of seven posterior channels of brain tumor EEG and 6,300 (42 cases×10 min ×60 sec/4 sec) epochs of seven posterior channels of normal EEG for the analysis.
From each 4-sec (1024 data points i.e., N = 1024) epoch of each posterior channel EEG record, the bicoherence is estimated with K = 30 overlapping segments of size M = 256 (1 sec) and overlap D = 230 (0.8984 sec i.e., 90% overlap). The largest (peak) Bicoherence (BIC) values in the [1-8] Hz (delta-theta) band and the [8-13] Hz (alpha) band from each epoch of each posterior channel EEG record are then selected and averaged over the all the epochs (10,050 from brain tumor case and 6,300 from normal case) to obtain the final channelwise value for the proposed analysis.

Test for Gaussianity and Linearity
The Hinich test is the statistical hypothetical testing against the null hypothesis that the time series under consideration is a Gaussian process generated by a linear mechanism. The Hinich test is extensively formulated in (Hinich, 1982). In this test, two statistics, namely the S statistic and the R statistic derived from the routine bicoherence and the non-central Chi-square distribution, are used as the measure of deviation from the Gaussianity and the measure of degree of linearity, respectively. The test makes use of the facts that the bicoherence is zero for a Gaussian process and constant for a linear process. The hypothetical test for the linearity is considered only if the null hypothesis on the Gaussianity is rejected. The HOSA toolbox (freely available at http://www.mathworks.com/matlabcentral/fileexchange/ 3013) contains a MATLAB routine called glstat to calculate these statistics.
The S and R statistics required for the Hinich test are estimated using the bicoherence for each 4-sec epoch. It should be noted that though the Hinich test is performed on the epoch basis i.e., 10,050 epochs of brain tumor EEG and 6,300 epochs of normal EEG are analyzed individually, the results are used, on the average, to statistically determine the Gaussianity and linearity of the entire channel EEG record. Nagata et al. (1985) used a parameter, called the power ratio index (PRI), to correlate the epicenter of the power dysfunction to the locality of the tumors in 15 patients with malignant brain tumors. He defined the PRI as the ratio of absolute power in the [1-8] Hz (deltatheta) band to that in the [8-30] Hz (alpha-beta) band. This parameter was also calculated from all the 4-second EEG epochs for the sake of comparison.

Statistical Analysis of BIC and PRI Values
The objective of the statistical analysis of the BIC and PRI values obtained is to assess the ability of the AJAS

Science Publications
proposed parameter in characterizing the presence of the brain tumor. Since the size of the data set is large enough (>30) in both the brain tumor and normal cases for the assumption of the normality of the estimated BIC and PRI values according to the central limit theorem (Montgomery and Runger, 2010), the two-tailed z-test is used to statistically test the channelwise and overall mean BIC values in discriminating the brain tumor patient from the normal subject at 5% significance level and to estimate the 95% confidence intervals of these mean differences. The confidence interval is here rather than the p-values following the suggestions from (Gardner and Altman, 2000).

Implementation
The entire analysis was carried out using the MATLAB @ (The MathWorks Inc., USA) software package on a personal computer.

RESULTS
The results of the test for determining the stationary segments using the second-order weak stationarity criterion are shown in Fig. 1a and b from two exemplary EEG records, one being the brain tumor and the other, normal. The dotted line shows the 95% confidence interval of the respective second-order statistic (mean or std) for two sets of 7 consecutive bins i.e., two 4-sec segments. The markers (asterisks) indicate the values of the respective statistic.
All the epochs of the seven posterior channel EEG records from both the brain tumor patients and normal subjects fail the Hinich test for the Gaussianity. Around 94% of the epochs fail the linearity test in the brain tumor case while around 81% of the epochs fail the linearity test in the normal case.
The contour plots of the bispectra of seven posterior-channel scalp EEG records from a normal subject and a brain tumor patient are shown respectively in Fig. 2 and 3 where f1 and f2 are the bifrequency values. Strong self-couplings (Raghuveer, 1988) are seen in the [8-13] Hz (alpha) band in the normal case while strong self-(some phase-) couplings are seen in the [1-8] Hz (delta-theta) band in the brain tumor cases.
The channel-wise distributions of all the (largest) bicoherence (BIC) values (circles) and their mean values (triangles) in the [1-8] Hz (delta-theta) band and in the [8-13] Hz (alpha) band across all the normal subjects and across all the brain tumor patients are shown respectively in Fig. 4 and 5 In contrast, only self-couplings are seen in both the bands for the brain tumor case.
The 95% confidence intervals of the channel-wise mean (largest) BIC values for the normal and brain tumor cases in the [1-8] Hz (delta-theta) band and the [8-13] Hz (alpha) band are shown, respectively in Fig. 6  and 7 where a value, 0.0 equals 0% coupling strength and 1.0 equals 100% coupling strength. The coupling strength in the [1-8] Hz (delta-theta band) varies from about 40% to 53% for the brain tumor patients while it is in the range around 23-42% for the normal subjects Meanwhile, the coupling strength in the [8-13] Hz (alpha) band is very low (<6.5%) for the brain tumor case and very high (about 27% to 52%) for the normal case, especially in and around the occipital region.
The statistics of the overall mean (largest) BIC values in the [1-8] Hz (delta-theta) band and the [8-13] Hz (alpha) band are shown in Table 1. The mean difference in the [8-13] Hz (alpha) band is statistically more significant (p<0.0001) than that in the [1-8] Hz (delta-theta) band (p<0.001). The 95% confidence interval on the mean difference of the largest BIC values in the [1-8] Hz (delta-theta) band is quite large, almost 11% more than the observed difference. However, this quantity in the [8-13] Hz (alpha) band is only approximately onequarter (more precisely 23%) of the observed difference. Table 2 shows the statistical analysis of the Power Ratio Index (PRI) using the two-tailed z-test. The 95% confidence interval on the mean difference between the PRI for the brain tumor case and that for the normal case is statistically insignificant (p<0.05).

DISCUSSION
The normal cases are selected with an age distribution similar to the age distribution of the brain tumor cases since the EEG characteristics are generally dependent on the age (Nunez and Srinivasan, 2006). This avoids any ambiguity on the results obtained due to the age-dependency of the scalp EEG.
The results of the Hinich test for the Gaussianity and linearity suggest that an EEG record with eyes closed can most often be considered as a non-Gaussian process generated by a (possibly quadratically) nonlinear mechanism at a 95% confidence level (or 5% significance level). This complies with the result presented in (Pradhan et al., 2012) and concludes that the attempt to quantify the phase couplings in these records is meaningful. The EEG time series remains non-Gaussian all the time and nonlinear most of the time. It may also be noted that the nonlinear structure in the EEG record is more dominant in the brain tumor case than in the normal case. Such dominance of the nonlinear character of the scalp EEG has been already observed under certain other pathological conditions as well (Diambra et al., 2001). This might be due to the fact that the presence of a tumor adds inhomogeneity (spatial non-uniformity) to the locality of the tumor (O'Connor and Robinson, 2005).
The reason for the selection of the posterior channels is the fact that the alpha activity is dominant in the frequency range 8-13 Hz during wakefulness in the posterior region of the head, particularly over the occipital region with eyes closed (Nunez and Srinivasan, 2006). The strong self-couplings in the [8-13] Hz band for a normal case indicate the presence of strong alpha rhythms. The presence of strong alpha rhythms in the frequency range 8-13 Hz is generally an indication of the healthy brains in 95% of the cases. Such strong selfcouplings in the alpha rhythms for normal subjects with eyes closed have already observed by Barnett et al. (1971). These self-couplings are almost lost in the brain tumor cases. The noteworthy association of the alpha rhythms with a variety of physiological, cognitive and behavioral states of the brain is once again evident from these results.
Since the phase coupling is considered as a sort of (nonlinear) phase synchrony (Kim and Powers, 1979;Siu et al., 2008), the strong phase (self-) coupling in the alpha band suggests a strong phase synchrony in the alpha band for the normal case. Such a strong phase synchrony is seen in the delta-theta band for both the normal and brain tumor cases.
For the brain tumor patients the [1-4] Hz (delta) band exhibits slightly stronger phase synchrony while the phase synchrony in the [8-13] Hz (alpha) band seems to be significantly lost. Decreased coupling strengths in the high-frequency band and increased ones in the low-frequency band are often indicative of low brain activity (e.g., levels of sedation, depth of sleep) (Barnett et al., 1971;Roustan et al., 2005;Venkatakrishnan et al., 2011;Yufune et al., 2011). Thus the presence of tumor reduces the brain activity.
The absence of phase coupling, which is found in the [1-8] Hz (delta-theta) band for the normal case, under the presence of a brain tumor is an indication of loss of quadratically nonlinear interaction in this band. These findings support the nonlinear interdependency between the brain regions and its contribution to the alpha rhythms (Breakspear and Terry, 2002).
The physiological model proposed for the corticothalamic dynamics in assessing the characteristics of the EEG associated with the thalamic tumors in (O'Connor and Robinson, 2005) and subsequent discussions therein are the evidences for the effect of the supratentorial tumors on the corticothalamic inputs and corticocortical networks which are responsible for the generation of the alpha rhythms (Schaul, 1998). The finding presented in this study not only supports this model but might also be helpful in extending this model to other supratentorial tumors on the assumption of non-Gaussianity for the model response, the EEG.
The final note is that the use of the bispectral index in the assessment of the level of consciousness during the administration of the anesthetic agent for the brain tumor patients has become questionable and is under serious study (e.g., http://clinicaltrials.gov/show/NCT01060631).
The anticipated correlation between the quadratic phase coupling phenomenon and the presence of the brain tumor is possibly the outcome of the presented analysis.
The Power Ratio Index (PRI) is not as significant in discriminating the brain tumor EEG from the normal one as in assessing the tumor locality (Nagata et al., 1985). The statistical results on the Power Ratio Index (PRI) based on the power spectrum once again prove that the bispectrum is far more informative than the power spectrum.
The strengths of the proposed analysis are the robustness of the bispectrum against the observational noises such as (white or colored) Gaussian noises (Nikias and Raghuveer, 1987) and the fact that the discriminating factor is the coupling strength in the alpha band which is less or not at all corrupted by the noises whereas the delta-theta and beta-gamma bands are always sensitive, respectively to lowfrequency artifacts such as movements, EKG and high-frequency artifacts such as EMG (Nunez and Srinivasan, 2006). The major setback of the method is its computational expenses. However, with latest advancements in the processor technology, the computational speed is not an issue.

CONCLUSION
From the Results and Discussion, it can easily be concluded that the bispectral analysis of the posterior channel EEG records can be useful for the diagnosis of (supratentorial) brain tumors. The absence of phase couplings in the [1-8] Hz (delta-theta) band and subsequently self-couplings in the 8-13 Hz (alpha) band at the posterior, especially occipital region of the head is indicative of a structural lesion in the supratentorial region of the brain.

ACKNOWLEDGEMENT
The researchers thank Dr. J. Mohanasundaram MD, Dean, Madras Medical College (MMC), Chennai and Dr. A. Sundaram MD, Vice-Principal, MMC, Chennai for having approved the collection of the data required for this research study from MMC and Dr. V. Sundar, Prof., Dept. of Neurosurgery, MMC, Chennai for having rendered his help in getting the approval from the Ethical Committee of MMC. The authors also thank Dr. S. Arunan, Asst. Prof., Dept. of Neurology and Dr. K. Thirumaran, Asst. Prof., Dept. of Neurosurgery for having rendered their help in selecting the cases (brain tumor patients and normal subjects) and recording their scalp EEGs. The first author is indebted to Dr. N. Kumaravel, Prof and Head, Dept. of ECE, Anna University, Chennai who has been a consistent encouragement to all his research activities.