Application of Multivariate Empirical Mode Decomposition and Sample Entropy in EEG Signals via Artificial Neural Networks for Interpreting Depth of Anesthesia

EEG (Electroencephalography) signals can express the human awareness activities and consequently it can indicate the depth of anesthesia. On the other hand, Bispectral-index (BIS) is often used as an indicator to assess the depth of anesthesia. This study is aimed at using an advanced signal processing method to analyze EEG signals and compare them with existing BIS indexes from a commercial product (i.e., IntelliVue MP60 BIS module). Multivariate empirical mode decomposition (MEMD) algorithm is utilized to filter the EEG signals. A combination of two MEMD components (IMF2 + IMF3) is used to express the raw EEG. Then, sample entropy algorithm is used to calculate the complexity of the patients’ EEG signal. Furthermore, linear regression and artificial neural network (ANN) methods were used to model the sample entropy using BIS index as the gold standard. ANN can produce better target value than linear regression. The correlation coefficient is 0.790 ± 0.069 and MAE is 8.448 ± 1.887. In conclusion, the area under the receiver operating characteristic (ROC) curve (AUC) of sample entropy value using ANN and MEMD is 0.969 ± 0.028 while the AUC of sample entropy value without filter is 0.733 ± 0.123. It means the MEMD method can filter out noise of the brain waves, so that the sample entropy of EEG can be closely related to the depth of anesthesia. Therefore, the resulting index can be adopted as the reference for the physician, in order to reduce the risk of surgery.


Introduction
Accurate and non-invasive monitoring of depth of anesthesia (DOA) [1] is very desirable during surgery.To achieve this purpose, many of the techniques or devices have been examined or used as methods to indicate the DOA, such as heart rate, blood pressure, and electroencephalogram (EEG) signals [2,3].Among these methods, the analysis of EEG is very intuitive, because the main action of general anesthetic agents takes place in the brain.EEG normally measured non-invasively through scalp electrodes, reflects the spontaneous electrical activities of the human brain over a short period of time, according to the analysis of which the state of human can be determined.However, it is still very difficult to know the patient's level of consciousness at the current stage.So far, there are 95% spectral edge frequency (SEF95), median power frequency (MPF), spectral edge frequency, narcotrend index (NI), and bispectral index (BIS) [4][5][6][7].Except for bispectral analysis, these methods use linear computational algorithms.However, none of these methods have been proved to completely express the EEG message.EEGs exhibit significant complex behavior with strong nonlinear and dynamical properties [8].Linear methods cannot express or ignore the message in some cases.Therefore, nonlinear theory may be a better approach than traditional linear methods in characterizing the intrinsic nature of the EEG.
Originally, entropy was viewed as a thermodynamic property that is the measure of a system's thermal energy per unit temperature.While applying the concept, entropy can address the system randomness and predictability.About a decade ago, researchers try to explain the depth of anesthesia using entropy.If the entropy is small, it indicates that the patient status is in anesthesia.If the entropy is large, it expressed that the patient status is awake.Response entropy (RE) and state entropy (SE) [9,10] are based on spectral entropy [11].Spectral entropy uses FFT at first, however FFT is linear method so it might miss out some important imformation.ApEn (Approximate entropy) [12,13] for nonlinear physiological signals was developed in 1991.The sample entropy (SampEn) [14][15][16][17] is an improvement of ApEn with respect to computation and accuracy of signal regularity.Therefore, more information can be extracted form EEG using SampEn compared to ApEn.Sample entropy will report false results when there is noise in the EEG, therefore a filter is needed.Multivariate empirical mode decomposition (MEMD) is very powerful filter.It can filter out noise and does not undermine the original signal.Therefore, SampEn can be used in conjunction with MEMD as signal processing methods to interpret the EEG signals for monitoring depth of anesthesia to reduce the burden on physicians.
This paper is divided into five parts.First, the concept of depth of anesthesia is introduced.Then the SampEn and MEMD methodologies are described in the next section.In the third part, an experiment is designed in order to identify the best IMF combinations.Monitoring DOA by sample entropy is described in the fourth part.Finally, the results are discussed and the conclusions are drawn.

Materials
In this study, the EEG signals are collected from thirty patients, whose ages are ranged from 20 to 80 and under ear nose and throat (ENT) surgery with general anesthesia at the National Taiwan University Hospital (NTUH) of Taiwan.Subjects who had alcohol, smoking, medical illness or medication issues were excluded.The equipment in the operation room included a physiological monitor (IntelliVue MP60) and a portable computer.This equipment displays the patient's physiological signals, such as electrocardiographic (ECG), EEG, blood pressure (BP) and saturated percentage of oxygen (SpO2) in real time.The study is mainly aimed at single channel EEG signal analysis for interpretation of DOA using the BIS™ Quatro Sensor.It can only collect single channel EEG signal since only consciousness is measured during anesthesia.This study was also approved by Institutional Review Board and written informed consent was obtained for the permission of the patients.

Sample Entropy
We can imagine that the entropy is the index of the degree of confusion.When the signal is not changing, entropy is lower.When the signal has more confusion and no regularity, entropy is higher.Based on this concept, the entropy is very useful in that it can be applied to the analysis of EEGs.Sample entropy improves an approximate entropy so that the entropy of the results can be more sensitive.At the same time, we can get more information from the EEG.
First, given N data points from a time series x n x 1 , x 2 , … , x N to define SampEn, and set two parameters which are r and m.The former parameter (r) is a coefficient of tolerance, while the latter (m) is the dimension of the template vector (X ), and general formula of the k th X is x k , x k 1 , … , x k 1 m .Then, compute the total number of the distances that are less than R (R = r × SD (standard deviation of original data)).The distance function is as shown in Equation (1): Here, B r is used to express the i th total number, then calculate B r using Equation (2): Finally, we calculate B r that is similar to the computation of B r and calculate the sample entropy using Equation (3): In this study, various theoretical and clinical applications have been shown that 1 or 2, and 0.1~0.2SDare more suitable values to use.In this study, we tried parameters like 2 and 0.1SD.

Multivariate Empirical Mode Decomposition
EEG signals are very weak, so they are easily interfered by other signals, such as from electromyography (EMG), electrooculogram (EOG), or electrosurgical units (ESUs) [18].These signals are contaminating noise in the EEG, so they will cause erroneous results.Therefore, a method is required that can filter out noise and does not undermine the original signal.Huang et al. proposed empirical mode decomposition (EMD) [19,20], that can decompose the original signal into different intrinsic mode functions (IMFs), expressed as follows: (4) where is the original signal in time domain, is ith IMF, and is the residue.Hence, we can choose the different suitable IMFs combination and get rid of noise of IMFs to re-construct the signal.However, EMD has problems of mode mixing so they proposed the EEMD technique [21] to solve this problem.Unfortunately, EEMD is time-consuming and can add noise into the original signal.In 2010, Rehman and Mandic [22] propose an improved EMD method called multivariate empirical mode decomposition (MEMD).Later, in 2011 they proposed a noise-assisted MEMD method (N-A MEMD) [23], which is not only suitable for dealing with multichannel signals, but also solves the problem of mode mixing using white Gaussian noise added to different channels.Hence, the N-A MEMD method was used in this paper.In computation of N-A MEMD, the mean is calculated by means of the multivariate envelope curves, expressed as follows: where are the multivariate envelope curves of the whole set of direction vectors and is length of the vectors.Then, we compute the remainder r by r m .If the remainder fulfills the stoppage criterion, the remainder is a multivariate IMF.If not, the input will equal the remainder r and compute the remainder again.We repeat these steps until all of multivariate IMFs are found.Regarding the stoppage criterion, it is still the same as the original EMD proposed by Huang et al. [19] using decomposing signal until the signal becomes monotonic.

Artificial Neural Networks
The SampEn value range is from 0 to 3, but the surgeon does not know whether the patient is awake or unconscious when the SampEn value is 2, so an index is needed to indicate the patient's state.In order to model this patients' state, a commercial product (BIS machine) is used as the gold standard.Hence, the traditional method of linear regression and nonlinear method of ANN were compared.
ANN is a parallel computing model which is similar to the human nerve structure, so that it is also called the parallel distributed processing model or connectionist model [24].ANN uses repeated constant learning and error correction in order to achieve the best output, so that the whole system is like a brain that understands the new problems, analyzes, and finally sums up the best conclusion.In ANN there are three learning rules generally: supervised learning, unsupervised learning and reinforced learning.In supervised learning, a new weight is generated by the correspondence between the rules in the input and output values in the training.In this study, a back-propagation neural network (BPNN) is used for leaning to model SampEn and BIS values.

Analysis of Intrinsic Mode Functions
Ten IMFs obtained by N-A MEMD are shown in Figure 1.However, which IMFs are related to brainwaves and which IMFs are related to noise is still difficult to decide.In this study, the experiment is designed such that IMFs combinations which are needed for the next analysis are acquired by sample entropy and FFT.Based on the concept of entropy, the patient's brain activity is faster before the patient is injected with medicine.Therefore, the sample entropy value should be higher.During surgery, the patient's brain activity is slower due to the drug's effect.Hence, the sample entropy value should be lower.During recovery, the patient's brain activity gradually returns to normal, so the sample entropy value should raise gradually.EEG for thirty patients from NTUH was recorded and every patient's data was divided into three stages.Each stage only takes one minute.The first stage is before injecting drugs so the patients have no drug effects and their status should be awake.The second stage is during surgery when the EEG is the most gentle and undisturbed.At this time, because the role of anesthetic drugs, the state of the patient is now unconscious.The third stage is the recovery stage where one minute before the end of the collection of data is selected.At this time, the patient can move his body consciously.
Table 1 shows the frequencies of each IMF obtained by computing the frequency of each IMF using FFT.According to frequency ranges of the EEG signals from state and response entropies, which are from 0.8 Hz to 47 Hz and 0.8 Hz to 32 Hz respectively, IMF2, IMF3, IMF4, IMF5, and IMF6, were selected so there are a total of 31 possible combinations.The entropy of each stage and the combinations were calculated as shown in Table 2.If the entropy value of stage 2 is both less than entropy value of stage 1 and stage 3, the combination may be suitable for use.There are five combinations (i.e., IMF2, IMF2 + IMF3, IMF2 + IMF4, IMF2 + IMF3 + IMF4, and IMF2 + IMF3 + IMF6) which are qualified.Statistical approach was used to calculate the p value that entropy value of stage 1 and stage 2 and entropy value of stage 2 and stage 3.As shown in Table 3,   In addition, it was found that IMF2 frequency is in β wave (12.5~25Hz);IMF3 frequency is in α wave (7.5~12.5Hz);IMF4 frequency is in θ wave (3.5~7.5Hz);IMF5~9 frequency is in δ wave (1.5~3.5Hz).Therefore, it is similar to previous brain wave studies in that the combination of these (IMF2, IMF2 + IMF3 and IMF2 + IMF3 + IMF4) represent the patient awake to conscious state [25].To demonstrate the filtering effect of N-A MEMD, a sample example is given to illustrate the result.For example, a fifty-five year old aged patient who went to ENT surgery is given.Figure 5 shows that BIS and sample entropy are very close, except in between 60 to 70 min.During that time, the doctor used an electrosurgical unit (ESU).Hence, the BIS module judged that the signal was very noisy, so the value was −1.However, the N-A MEMD can filter the noise, and the sample entropy value can still express the patient's status.Therefore, the data is collected each thirty seconds before and after BIS lost the signal due to ESU.For this patient, BIS is only miss-classified once where it started at 3,756 s (i.e., 62.6 min) and lasted for 185 s as shown in Table 5 for patient 1.This phenomenon can be replaced by applying SampEn procedure, which was treated by N-A MEMD.Using this method, the value is still close to the 30 s before and after the BIS was equal to −1, and the operation of this patient was lasting for 102.08 min.Fourteen patients out of 30 patients who have several different times of missing BIS values were analysed as shown in Table 6.The most frequent record of BIS loss value is fifteen times, which is patient 14, for 750 s of the 229 min of the total operation time.

Receiver Operating Characteristic (ROC) Curve
Correlation coefficient and MAE can be used to compare different analysis methods, but there is a need for a quantitative value to explain the method to a degree of how good or bad it is.Therefore, in this study, the Receiver Operating Characteristic Curve [26,27] is used as an assessment method.Different cut-off points can be used to get sensitivities and specificities.Usually, a diagonal line is drawn on the chart as a benchmark, and the curve of the above on this line is better and below is not good.In this study, the area under the ROC curve (AUC) is used to assess these methods.Since the range of AUC values is from 0 to 1, the bigger values mean better results.When AUC is less than or equal to 0.5, the method has no discrimination.When AUC is greater than or equal to 0.7 and less than 0.8, the method has acceptable discrimination.When AUC is greater than or equal to 0.8 and less than 0.9, the method has excellent discrimination.When AUC is greater than or equal to 0.9, the method has outstanding discrimination.
In the study, the sample entropy from original EEG signal, after filtering the EEG signal, and via MEMD and ANN are compared with BIS as standard.In Table 7, "Original entropy" is the sample entropy from original EEG signal."Entropy via MEMD" is the sample entropy from after filtering EEG signal."ANN" is the sample entropy via MEMD and ANN.From this table, MEMD can improve accuracy of distinguishing the state of stupefaction or awake.Table 7 lists all the cases featuring mean AUC and SD of the methods of "ANN", "Entropy via MEMD" and " Original entropy" as 0.970 ± 0.028, 0.969 ± 0.028, and 0.733 ± 0.123, respectively.

Discussion and Conclusions
During surgical operations, inevitable noise such as EOG, EMG or ESU will interfere with EEG signals.Therefore, it is necessary to reduce the influence of noise using a filter first.For nonlinear signal processing, the N-A MEMD method can get the original signal without disrupting it.In this study, the original EEG signal was decomposed using N-A MEMD into many IMFs.Next, the main frequency of each IMF is calculated.Results show that frequencies of IMF2 to IMF6 are located in brainwave frequency range.Thirty-one IMF combinations were generated using the five significant IMFs.The combination (IMF2 + IMF3) is the best expression of the raw EEG.However, it is totally agreed that delta waves will occur during deep anesthesia.In this paper different combinations were compared as shown in Table 2. From this table, the entropy mean value of stage 2 is greater than the entropy mean value of stage 1 when considering IMF5 & IMF6 inside this analysis on the last six rows of Table 2.The results are not reasonable, so this combination is eliminated.This might be due to the type of the surgery [i.e., ear, nose, and throat (ENT) surgery] which does not require deep anesthesia, so the combination of IMF2 + IMF3 is considered to be the best.However, if the surgery lasts for a longer time, e.g., liver or heart transplant operations, the combination of IMF2 + IMF3 + IMF5 + IMF6 may be better.
The entropy value obtained using ANN having good results were achieved.In order to provide anesthesiologists valuable reference and evaluation, N-A MEMD methods should be used after further validations.Although EMD, EEMD, or CEEMD [28] are used to decompose this single channel signal, these three methods have some problems such as mode mixing, being time consuming, or too much noise being added to the original signals.Nowadays, there are many types of EMD and each one has merits and demerits in processing the signal.For example, eXtended-EMD (X-EMD) [29] can produce better results than MEMD and Turning Tangent EMD (2T-EMD) [30] in dealing with multivariate signal denoising.However, in 2011 Rehman et al. also proposed N-A EMD.This shows that the original MEMD is worse than X-EMD, but N-A MEMD improves the original MEMD and also solves the mode mixing problem.Therefore, the N-A EMD method applied in this study has proven to be very successful to deal with EEG filtering during surgery.However, in the near future it is hoped that other methods, such as Independent Component Analysis (ICA) [31] and Canonical Correlation Analysis (CCA) [32], can be compared to investigate the improvement in interpreting depth of anesthesia.
In Figure 5, it is different with BIS and entropy during 70 to 100 min.Because the BIS-index is poor in distinguishing the awake from unconscious patients [33], there is a need to cooperate with doctors in order to verify the method of this study.The collected signals are not only EEG and BIS, but also ECG, SPO2, BP, drug dosage, and operation of events and related information.According to the data and clinical experience, physicians can determine consciousness degree indexes more accurately than BIS.Then, the consciousness degree index can be used as the gold standard.Therefore, the aim of this paper is to establish the feasibility of this method at this stage, and the physician's judgment will be merged in near future.
there are three combinations' (IMF2, IMF2 + IMF3, and IMF2 + IMF3 + IMF4) of significant p values of stage 1 and stage 2 and entropy value of stage 2 and stage 3 are less 0.05.The selected combination (IMF2+IMF3) that has the least p value of stage 1 and stage 2 and entropy value of stage 2 and stage 3 are shown in Figure 2.

Figure 4 .
Figure 4.The training, validation, and testing of ANN results (red line is target and blue line is output).

Figure 5 .
Figure 5.The comparison of BIS and sample entropy.

Table 1 .
IMFs of frequency at each stage.

Table 2 .
Sample entropy of the combinations at each stage.

Table 5 .
The comparison of SampEn value and BIS value when BIS value of −1 in patient 1.

Table 6 .
The total time for BIS value is −1 during ESU disturbance.
Event & Time (s) Event no.Total event time (s) Operation time (min)

Table 7 .
AUC of all thirty patients.