Real-Time Depth of Anaesthesia Assessment Based on Hybrid Statistical Features of EEG

Yi Huang

^1,*

Peng Wen

¹,

Bo Song

¹ and

Yan Li

School of Engineering, University of Southern Queensland, Toowoomba 4350, Australia

School of Mathematics, Physics and Computing, University of Southern Queensland, Toowoomba 4350, Australia

Author to whom correspondence should be addressed.

Sensors2022,22(16), 6099;https://doi.org/10.3390/s22166099

Submission received: 7 July 2022 /Revised: 8 August 2022 /Accepted: 12 August 2022 /Published: 15 August 2022

(This article belongs to the SectionBiomedical Sensors)

Downloadkeyboard_arrow_down

Browse Figures

Review Reports Versions Notes

Abstract

This paper proposed a new depth of anaesthesia (DoA) index for the real-time assessment of DoA using electroencephalography (EEG). In the proposed new DoA index, a wavelet transform threshold was applied to denoise raw EEG signals, and five features were extracted to construct classification models. Then, the Gaussian process regression model was employed for real-time assessment of anaesthesia states. The proposed real-time DoA index was implemented using a sliding window technique and validated using clinical EEG data recorded with the most popular commercial DoA product Bispectral Index monitor (BIS). The results are evaluated using the correlation coefficients and Bland–Altman methods. The outcomes show that the highest and the average correlation coefficients are 0.840 and 0.814, respectively, in the testing dataset. Meanwhile, the scatter plot of Bland–Altman shows that the agreement between BIS and the proposed index is 94.91%. In contrast, the proposed index is free from the electromyography (EMG) effect and surpasses the BIS performance when the signal quality indicator (SQI) is lower than 15, as the proposed index can display high correlation and reliable assessment results compared with clinic observations.

Keywords:

EEG;depth of anaesthesia;real-time;machine learning

1. Introduction

Monitoring the depth of anaesthesia (DoA) during surgery is challenging [1,2]. When anaesthetic agents are applied, the response of the central nervous system could be reflected in the electroencephalography (EEG) [3]. The EEG signals exhibit low-voltage high-frequency and high-voltage low-frequency characteristics in the awake state and deep anaesthesia state, respectively [4]. Thus, the neural activities in the central nervous system could be revealed by EEG. During surgery, the main problem is prescribing precision doses of appropriate anesthetic agents for each patient. An excessively high dosage of anaesthetic agents may lead to coma and postoperative complications, and on the other hand, if the dosage of anaesthetic agents is insufficient, the patient may suffer from pain and waking up [5]. Thus, the real-time DoA assessment can help both patients and anaesthetists.

In the past two decades, the Bispectral Index (BIS) has been widely used for assessing and monitoring the DoA. BIS was developed by Aspect Medical Systems in 1992 [6,7], and is one of the most popular commercial products in the DoA monitoring markets. The BIS index was developed from several complex parameters in the time domain and frequency domain, which are integrated as a one dimension index ranging from 0 (unconscious state) to 100 (awake state) [8,9,10]. Besides, there are some other electrophysiological monitors, such as the Narcotrend index [11], M-Entropy [12], and patient state index [13], for tracing the DoA. However, these DoA monitors still have some limitations, such as inconsistent states between clinical observation and device monitoring index, insensitivity at switch points between consciousness and unconsciousness, and not being accurate across patients [14,15]. Moreover, these real-time monitor systems in clinical use have poor performance when the signal quality index (SQI) is low.

A wide range of features, including frequency domain and time domain features, have been proposed for monitoring DoA over the years. Nguyen-Ky et al. [14] proposed to measure the DoA using a wavelet transform method. In this method, the EEG signal was decomposed into different levels to extract the desired features, and an eigenvector of wavelet coefficients was calculated to develop a new index. Wavelet-weighted median frequency and wavelet coefficient energy entropy methods were used to develop new indexes with high rates of agreement with BIS [16]. Several studies commonly use permutation and sample entropy to create a correlated index with BIS [17,18]. Sarkela et al. [19] used spectral characteristics of the EEG burst suppression as features to propose an automatic method for burst suppression detection and segmentation. Lashkari and Boostani [20] introduced an improved instantaneous frequency (IF) by applying a Kalman filter to develop a noticeable index correlated with the BIS index. In addition, several different types of regression models and classifiers were used to monitor the trending of anaesthetic states, such as an artificial neural network [21,22], an adaptive neuro-fuzzy system [18], a genetic algorithm with a support vector machine (GA-SVM) [23], random forest [24], and convolutional neural network [25]. Neural networks have the advantage of being more flexible in classifying EEG signals. However, these classifiers require a large amount of data to train a robust model, most of the datasets in the research field of anaesthesia are not open to the public, and it is even more challenging to find the public DoA datasets with observation comments and labels from the attending anaesthetists [25].

This study used the threshold value for the wavelet denoising method to pre-process the raw EEG signal data. After feature extraction and selection, five features were evaluated as the input of the Gaussian process regression model to classify different states of anaesthesia. All of the experiments in this study were carried out using MATLAB 2018b and a workstation with an Inter^® Core™ i9-10900K CPU @ 3.7 GHz processor.

The first section of the paper provided a brief introduction to the work.Section 2 described the details of the datasets and data acquisition. The pre-processing, features extraction, states of anaesthesia classification, and real-time application are also introduced in this section.Section 3 presents the experiments and results obtained. The discussion on the EMG effect and signal quality were evaluated inSection 4.Section 5 concluded the paper.

2. Materials and Methods

The proposed threshold wavelet denoising method started with applying a sliding window technique. The sliding window size was 10 s, having a 9 s overlap with the previous window. Then, the statistical features of the denoised EEG signal were extracted. The regression model was trained by the machine learning method to classify the depth of anaesthesia. Finally, the new index,

N_{D o A}

, was proposed to monitor the DoA in real-time.Figure 1 describes the workflow of the proposed methodology in this study.

2.1. Subjects and EEG Recording

The EEG and BIS datasets are directly from our collaborating hospitals, including Toowoomba Hospital, Brisbane Prince Charles Hospital, Queensland, Australia, and our industry partners, such as Shenzhen Delica Medical Equipment Co., Ltd., Shenzhen, China. The EEG and BIS data were obtained using a BIS

{VISTA}^{TM}

monitoring system, version 3.22. All intravenous dosing and intraoperative events were recorded on data log files by the attending anaesthetist. The data log files also contain the BIS index, raw EEG data, signal quality indicator (SQI), impedance, electromyography (EMG), and monitor log of errors. Among the two channels of raw EEG, this study used the one with more consistent EEG quality. In addition, the sampling frequency of the recording data was 128 Hz with a 16-bit signed integer. In this study, a total of 73 subjects with corresponding data were used and the summary of patient demographics is shown inTable 1.

2.2. Signal Denoising and Pre-Processing

It is difficult to accurately discriminate the anaesthetic response signal from the raw EEG signal because the noises, such as eye movements, muscle activities, and other artefacts, may corrupt the recorded signals, especially in the awake state. Therefore, all raw EEG signals were pre-processed to reduce the noise. Firstly, the mean and standard-deviation-based threshold methods are used to remove outliers [26]. Then, a wavelet threshold method based on entropy is used to remove the low-amplitude noise and spike noise interference.

Two threshold methods are used to conduct the denoising process, where

T_{h}

is the threshold value. Donoho [27,28] proposed a threshold wavelet method to remove the noise using a threshold value:

T_{h} = σ \sqrt{2 \log N}

(1)

σ = M A D / 0.6745

(2)

where

σ

is the noise standard deviation,N is the size of the wavelet coefficient arrays, andMAD is the media of absolute detail coefficient of each level in the wavelet transform.

The function of soft-threshold

S_{j}^{'}

has a sign function related to the detail coefficients

C_{j}

S_{j}^{'} = {\begin{matrix} s i g n (C_{j}) (| C_{j} | - T_{h}), i f | C_{j} | > T_{h} \\ 0, i f | C_{j} | \leq T_{h} \end{matrix}

(3)

The function of hard-threshold

H_{j}^{'}

has a value equal to 0 if the detail coefficient is smaller than the thresholdT. Otherwise, it is equal to

C_{j}

H_{j}^{'} = {\begin{matrix} C_{j}, i f | C_{j} | > T_{h} \\ 0, i f | C_{j} | \leq T_{h} \end{matrix}

(4)

An improved method was proposed using Stein’s unbiased risk estimate [28]. An adaptive wavelet energy entropy method was used to filter noise in [29]. In this study, the threshed

T_{h}_{n e w}

used a permutation entropy method based on Stein’s unbiased risk estimate [28].

H_{p e} (m) = - \sum_{j = 1}^{m!} p_{j} \ln p_{j}

(5)

wherem is the dimension of time series signals and

p_{j}

is the possibility of each permutation.

The total energy of permutation entropy is

P E_{t o t a l} = - \sum_{k = 1}^{K} \sum_{j = 1}^{m!} p_{j} \ln p_{j} = \sum_{k = 1}^{K} H_{p e} (m)

(6)

The relative permutation entropy energy is

r_{p e} = \frac{H_{p e} (m)}{P E_{t o t a l}}

(7)

The new adaptive threshold

T_{h}_{n e w}

was proposed as

T_{h}_{n e w} = (\log r_{p e} \sqrt{\frac{\log (n \log (n))}{2}} - a) \times b

(8)

wheren is the window length of an EEG segment anda andb are two constants (9 and 6, respectively) that are calculated offline empirically. In this study, a sliding window technique with a step of one second (s) and a fixed length (L = 10 s) was used to read and process the raw EEG signals. Thus, there are 9 s signals that are overlapped between two sequential sliding windows.

2.3. Features’ Extraction

Five features were selected based on the pre-processed data to evaluate their differences in different anaesthetic states. The features are the sample entropy (SE), fuzzy entropy (FE), permutation entropy (PE), Hurst exponent (HE), and power spectral density (PSD), and their calculation methods are provided below.

Entropy is the parameter to measure the degree of irregularity of a time series signal. In this study, SE, FE, and PE are calculated to produce a new depth of anaesthesia index [17,30].

S E = - \ln \frac{A^{m} (r)}{B^{m} (r)}

(9)

where

A^{m}

is the number of template pairs having

d [X_{m + 1} (i), X_{m + 1} (j)] < r

and

B^{m}

is the number of template pairs having

d [X_{m} (i), X_{m} (j)] < r

. Tolerancer was chosen as 0.2 ×std (std is the standard deviation of the time series signals), and the dimension m was chosen as

F E (X, m, n, r, N) = \ln \emptyset^{m} (n, r) - \ln \emptyset^{m + 1} (n, r)

(10)

whereN is the length of the time series. Empirically, tolerancer was chosen as 0.15 ×std,n was selected as 2, and dimensionm was chosen as 2 to produce the best performance.

P E = \frac{\sum_{j = 1}^{J} P_{j} \ln P_{j}}{\ln (J)}

(11)

where embedding dimensionm and time delay

τ

are 4 and 1, respectively.

The Hurst exponent, a numerical method, was used to predict the trend in EEG signals. The time series signalsY with the length ofN were divided into components

y_{m} = {y_{1, n}, y_{2, n}, \dots, y_{m, n},}

and each component has the same length

m ϵ {N, \frac{N}{2}, \frac{N}{4}, \dots}

. The following procedure was used to calculate the range response in the Hurst method.

Calculating the mean value of each component:

m e a n = \frac{1}{m} \sum_{i = 1}^{m} y_{i, n}

(12)

2.: Creating a mean-adjusted series:

X_{i, n} = y_{i, n} - m e a n, f o r i = 1, 2, \dots, m

(13)

3.: Calculating the cumulative deviate series:

Z_{i, n} = \sum_{t = 1}^{i} X_{i, n}

(14)

4.: Computing the rangesR:

R_{n} = \max {Z_{1, n}, Z_{2, n}, \dots, Z_{m, n}} - \min {Z_{1, n}, Z_{2, n}, \dots, Z_{m, n}}

(15)

5.: Computing the standard deviation $S_{n}$ :
$S_{n} = \sqrt{\frac{1}{m}} \sum_{i = 1}^{m} {(y_{i, n} - m e a n)}^{2}$
(16)

6.: Calculating the rescaled range:

{(R / S)}_{n} = \frac{1}{(\frac{N}{m})} \sum_{n = 1}^{[\frac{N}{m}]} (R_{n} / S_{n})

(17)

7.: Calculating the rescaled range:

R_{r} = \min_{n = 1, 2, \dots} {m e a n (R_{n})}

(18)

Multiple signal classification (MUSIC), as the method of an eigenvector-based frequency estimator, is used. The feature extraction was from the power spectral density (PSD) of wavelet coefficients, and the method is described below based on [14,31]:

1.: $A_{j}$ and $D_{j}$ were used as input signals for the eigenvector method:

S D : {\begin{matrix} D_{j} \to P (D_{j}) = E D_{j} \\ A_{j} \to P (A_{j}) = E A_{j} \end{matrix}

(19)

where

j

was set as six-level wavelet decomposition, the order of 16 Daubechies wavelet filter was used [32], and the principal eigenvector was 6.

2.: Means of $E A_{j}$ and $E D_{j}$ are as follows:

M e a n : {\begin{matrix} D_{j} \to M (E D_{j}) = m e a n (E D_{j}) \\ A_{j} \to M (A D_{j}) = m e a n (A D_{j}) \end{matrix}

(20)

3.: Standard deviation (STD) of $E A_{j}$ and $E D_{j}$ are as follows:

S T D : {\begin{matrix} D_{j} \to S (E D_{j}) = s t d (E D_{j}) \\ A_{j} \to S (A D_{j}) = s t d (A D_{j}) \end{matrix}

(21)

4.: Derived from steps 2 and 3:

D : {\begin{matrix} M_{j} \to \frac{1}{2} {\log (m e a n (M (E D_{j}))) + \\ l o g (m e a n (M (A D_{j})))} \\ S_{j} \to \frac{1}{2} {\log (m e a n (S (E D_{j}))) + \\ l o g (m e a n (S (A D_{j})))} \end{matrix}

(22)

5.: Deriving the feature:

F e a = \frac{k_{1} \times M_{j} + k_{2} \times S_{j}}{k_{3}}

(23)

The constant parameters

k_{1}

k_{2}

, and

k_{3}

were chosen as 28, 90, and 3, respectively.

These five features were extracted individually after pre-processing the raw EEG signals. Then, they were used as the input data to train a robust model to classify different states of anaesthesia.

2.4. Regression Models and Evaluation Measures

In this study, a Gaussian process regression model, as a nonparametric Bayesian approach, was applied. This kind of regression model has its advantage on small datasets and provides reliable uncertain predictions. We also utilized linear regression, SVM, and regression tree models for comparison and evaluation. The dataset with 73 subjects was divided into the training and testing sets. A total of 60 subjects were randomly selected for training, and the aforementioned five features of each subject were calculated for training and evaluating regression models to find the best performance models based on different combinations of features.

The training models were evaluated using R-square (

R^{2})

, root-mean-square-error (RMSE), mean-square-error (MSE), and mean absolute error (MAE).

R^{2} = 1 - \frac{\sum_{i} {(y_{i} - f_{i})}^{2}}{\sum_{i} {(y_{i} - \bar{y})}^{2}}

(24)

where

y_{i}

is the data set value (BIS value),

\bar{y}

is the mean of the BIS values, and

f_{i}

is the calculated value from regression results. The range of the

R^{2}

values are between 0 and 1. The higher value of

R^{2}

indicates a higher correlation between the BIS value and the extracted feature values, and vice versa.

M S E = \frac{1}{n} \sum_{i = 1}^{n} {(Y_{i} - \bar{Y_{i}})}^{2}

(25)

M A E = \frac{\sum_{i = 1}^{n} | Y_{i} - \bar{Y_{i}} |}{n}

(26)

wheren is the number of the data points,

Y_{i}

is a set of the observed values, and

\bar{Y_{i}}

is a set of the regression results.

After the best performance model was selected and evaluated using the testing datasets, the Bland–Altman method was used to assess the degree of agreement between the proposed index,

N_{D o A}

, and the BIS index. In addition, the Pearson correlation coefficient was also applied to assess the correlation of

N_{D o A}

with the BIS index and clinical observations.

r = \frac{\sum^{} (x_{i} - \bar{x}) (y_{i} - \bar{y})}{\sqrt{\sum^{} {(x_{i} - \bar{x})}^{2} {(y_{i} - \bar{y})}^{2}}}

(27)

2.5. Real-Time DoA Monitor under Sliding Window Framework

In order to have a real-time response, the sliding window technique was applied and the window size was chosen as 10 s, and the overlap between two adjacent windows was 9 s. In addition, the tuned index,

D o A_{t u n e d}

, was calculated using Equation (27), from the mean value between previous and present data.

D o A_{t u n e d} = 0.8 \times I n d e x_{p r e} + 0.2 \times I n d e x_{c u r}

(28)

where

I n d e x_{p r e}

is the mean value of the DoA index from the last four seconds,

I n d e x_{c u r}

is the current DoA index, and

D o A_{t u n e d}

presents the real-time index second by second. Using this sliding window method, the DoA index could be updated every second. Consequently, the time delay for monitoring the DoA could be reduced. Meanwhile, the starting time of processing the recording raw EEG data was from the fifth second, and there was a four-second delay.

3. Results

3.1. Raw EEG Signal Pre-Processing and Features’ Selection

Figure 2 presents the sample denoising results with the proposed threshold

T_{h n e w}

using the permutation entropy method. The raw EEG signal (from patient ID: L112161431) has low amplitude noise (between 2000 s and 2010 s) and spike noise (between 2980 s and 2990 s), as shown inFigure 2a–c, respectively. Compared with the denoised EEG signal using the threshold

T_{h}

, as shown inFigure 2d,e, the denoised EEG signal using the proposed threshold

T_{h n e w}

, as shown inFigure 2f,g, has much less noise. After using both denoising threshold methods, the EEG signal has a similar amplitude. However, the proposed threshold method in this study removed almost all the noise, as shown inFigure 2d–g.

There are EEG signals from two channels (CH1 and CH2) in the datasets. In this study, CH2 data were selected because their correlation between the features and anaesthetic states was higher than that of CH1. The five features of fuzzy entropy, permutation entropy, sample entropy, power spectral, and Hurst range response were extracted from the denoised CH2 EEG signal. All of the features were measured using coefficients of determination (

R^{2}

), which reflects how close a feature is to different anaesthetic states. The range of

R^{2}

values is between 0 and 1. A higher value of

R^{2}

indicates a higher correlation between the BIS and the extracted feature, and vice versa. As an example of the patient ID: L01040838, the values of

R^{2}

are 0.7497, 0.6974, 0.7034, and 0.6214 for fuzzy entropy, permutation entropy, sample entropy, and power spectral, respectively. The feature of Hurst range response has a relatively lower correlation with the BIS compared with the other four features. However, the Hurst range response has an obvious state change point from consciousness to unconsciousness at 762 s, as shown inFigure 3.

3.2. Regression Model Training

A total of 60 subjects were randomly selected as the training set from 73 subjects. The total EEG recording time for the 60 subjects was 385,116 s, with the corresponding BIS values for each second. The remaining 13 subjects were used to evaluate the training regression models. Fivefold cross-validation, as a resampling procedure, was used to protect against overfitting by partitioning the data sets.

The MATLAB regression learner app [32] was used for training with four different models, including linear regression, trees, SVM, and Gaussian process models.Table 2 lists the best performance for each type of training model. For example, fine Gaussian SVM has the greatest value of

R^{2}

and smallest value of RMSE among linear SVM, quadratic SVM, cubic SVM, medium SVM, and coarse SVM. The results inTable 2 show that the squared exponential Gaussian process regression model (SEGP-RM) has a better overall performance. Thus, the SEGP-RM was selected for training the robust model.Figure 4 andFigure 5 show the fivefold cross validation results of the SEGP-RM. As a result,

R^{2}

is 0.76, while RMSE, MSE, and MAE are 6.18, 38.22, and 5.15, respectively. The response plot is used to explore the correlation between the predicted values and true values, and the correlation coefficient is 0.89. The predicted plot versus actual plot inFigure 4 shows how well the regression model makes predictions for different response values. A perfect regression model has a predicted response equal to the true response, so all of the points (blue points in the below figure) lie on a diagonal line (red line in the below figure). The vertical distance from the line to any point is the prediction error for that point.

3.3. Performance of the Training Models and Results of the Nex Index

To further validate the performance of the trained models, the remaining 13 subjects were used for evaluation. The metric, Pearson correlation coefficient, which is defined in Equation (26), is used to evaluate the four trained regression models listed inTable 2.

Figure 6 indicates the Pearson correlation coefficients of the 13 testing subjects. Patient number 14 represents the average value of all testing data results. The correlation coefficient of the SEGP-RM is higher than the other three regression models, robust linear, fine trees, and fine Gaussian SVM. The highest correlation coefficient is 0.940 for patient No. 8 (subject ID: L01131248), and the lowest one is 0.580 for patient No. 3 (subject ID: L01060841). Meanwhile, the average value of correlation coefficients from these 13 testing subjects is 0.812.

The difference between

N_{D o A}

and BIS is defined as

d = N_{D o A} - B I S

. The mean value of the difference between

N_{D o A}

and BIS is

d i f f = m e a n (d)

and the standard deviation of the difference is

s d = s t d (d)

. The Bland–Altman method suggests that 95% of the

d

should lie between

d i f f + 2 s d

and

d i f f - 2 s d

(two black lines inFigure 7a) if the data

d

have a normal distribution. InFigure 7b, the pattern matches well with the distribution of

d

. The Bland–Altman plot of the data from 13 testing subjects, calculates the mean difference

d i f f = 2.1348

as the ‘bias’ (red dash line). For the 95% agreement limitation, the upper limitation is

d i f f + 2 s d = 21.46

and the lower one is

d i f f - 2 s d = - 17.19

. The agreement rate between

N_{D o A}

and BIS is 94.91%, as shown in the scatter plot inFigure 7a.

4. Results and Discussion

4.1. Length of Sliding Window

During surgery and clinical operations, the whole anaesthetic process is controlled by anaesthetists. The surgery’s incisions started when the patient’s anaesthetic states were considered at the middle anaesthesia stage (BIS index between 40 and 60). The real-time technique proposed in this paper only has four seconds’ delay at the beginning of recording data (awake stage), which may not affect clinical data processing and anaesthetic state tracing in clinical applications. The sliding window size affects the real-time results of the anaesthesia monitoring and efficiency. If the sliding window size is too big, the results may cause a delay in the clinical operations, and too small of a window size may affect the feature extractions during the transitions of two anaesthesia states. In this study, a fixed-length sliding window of 10 s was selected.

4.2. The Effect of Signal Quality

The signal quality indicator is an index recorded with BIS values in real time, which measures the quality of the raw EEG signals. Artifacts, such as eye movement, muscle activity, head motion, and the electrical knife operations, affect the impedance calculation and result in poor signal quality. We found that the BIS values were not displayed on the monitor screen when SQI was lower than 15. Meanwhile, we reviewed the raw data log files, and the corresponding BIS values are −3276.8 when the signal quality is poor (SQI lower than 15). However, our proposed real-time monitor index,

N_{D o A}

, can still provide the corresponding anaesthetic states in the case of poor signal quality.

Figure 8 demonstrates a sample case with poor signal quality (subject ID: L01121339). The BIS values are −3276.8 during the time from 369 to 371 s, from 1708 to 1710 s, and from 1790 to 2017 s. However, the proposed index,

N_{D o A}

, computed and provided the DoA values clearly during the period.Figure 9 presents another case with poor signal quality (subject ID: L01210841). The subject was a 17-year-old male patient weighing 46 kg. The intraoperative EEG monitoring in surgery started at 08:41:41 and finished at 13:37:27. The detail of anaesthesia medications for this patient is reported inTable 3. In this case, the SQI index was below 55 in the awake stage (from 1 to 610 s). Furthermore, the invalid BIS index (−3276.8) could not display on the screen monitor as poor signal quality (SQI < 15) during these periods from 89 s to 132 s and 389 s to 412 s, while the proposed

N_{D o A}

was able to estimate the DoA. Compared with the BIS index, the proposed method has a better correlation and agreement with clinical observations during periods of poor signal quality. Moreover, two additional cases with the same phenomenon could be observed inFigure 10 andFigure 11.

4.3. The Effect of EMG Signals

EMG signals, caused by muscle activities, were recorded as a spark from forehead muscles areas. These signals can have a negative effect on increasing the value of the BIS index, and the neuromuscular blocking agents can reduce or eliminate this side effect [33]. InFigure 12 there are two spikes with the BIS index higher than 70 in a general anaesthesia state from 2560 s to 2580 s and from 2870 s to 2890 s. The index of EMG signals also reflected this fact. Meanwhile, the SQIs are between 90 and 100 during these two periods, which indicates reliable and accurate signal quality with the BIS index. Therefore, the EMG activities during these periods might result in two spikes for the BIS indexing. However, as observed from the waveform of the

N_{D o A}

, the trend of the proposed index does not change in the general anaesthesia state.

5. Conclusions

This paper studies the real-time monitoring of anaesthesia states using hybrid statistical features and machine learning methods. Firstly, a wavelet threshold method was used to remove the low-amplitude and spike noise interference. Then, hybrid statistical features were extracted from the denoised EEG signals. Next, the regression models were trained using these extracted features by applying the machine learning methods, linear regression, SVM, regression trees, and Gaussian process regression. The squared exponential Gaussian process regression model outperformed the other models. This yields a better result in accuracy and a higher degree of agreement rates. In addition, the proposed index,

N_{D o A}

, has a higher correlation with anaesthetists’ observations than the BIS index during periods of poor signal quality.

The results were evaluated by comparing the correlation and agreement with the BIS index (current gold standard) and the clinical observations recorded by the attending anaesthetists. The results of

N_{D o A}

are highly correlated with the BIS index. The highest correlation coefficient,

R^{2},

is 0.940, and the average is 0.812. Furthermore, the agreement measured by the Bland–Altman method presented a higher degree of agreement between the proposed

N_{D o A}

and the BIS index. The scatter plot shows an agreement of 94.91%.

Author Contributions

Y.H. presented the project idea and completed the modelling, experiments, and the writing of this manuscript, while P.W., B.S., and Y.L., being supervisors, contributed to the design of the study, the completion of the project, and the editing of this manuscript. All authors read and approved the final. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

The study was conducted in accordance with the Declaration of Helsinki, and approved by the Institutional Ethics Committee of University of Southern Queensland (No: H09REA029, 2009/016).

Informed Consent Statement

The medical condition and procedure of data collection were explained to all subjects/patients. Consent, including all ethical issues, was obtained from all patients. The study was approved by the University of Southern Queensland Human Research Ethics Committee (No: H09REA029) and the Toowoomba and Darling Downs Health Service District Human Research Ethics Committee (No: TDDHSD HREC 2009/016).

Data Availability Statement

The data and materials used in this study are available at the University of Southern Queensland under the research data management policy.

Conflicts of Interest

The authors declare no conflict of interest.

References

Bard, J.W. The BIS monitor: A review and technology assessment.AANA J.2001,69, 477–483. [Google Scholar]
Constant, I.; Sabourdin, N. Monitoring depth of anesthesia: From consciousness to nociception. A window on subcortical brain activity.Pediatr. Anesthesia2014,25, 73–82. [Google Scholar] [CrossRef] [PubMed]
Jospin, M.; Caminal, P.; Jensen, E.W.; Litvan, H.; Vallverdu, M.; Struys, M.M.R.F.; Vereecke, H.E.M.; Kaplan, D.T. Detrended Fluctuation Analysis of EEG as a Measure of Depth of Anesthesia.IEEE Trans. Biomed. Eng.2007,54, 840–846. [Google Scholar] [CrossRef] [PubMed]
Musizza, B.; Ribarič, S. Monitoring the Depth of Anaesthesia.Sensors2010,10, 10896–10935. [Google Scholar] [CrossRef]
Pandit, J.J.; Cook, T.M.; Jonker, W.R.; O’Sullivan, E. A national survey of anaesthetists (NAP5 Baseline) to estimate an annual incidence of accidental awareness during general anaesthesia in the UK.Anaesthesia2013,68, 343–353. [Google Scholar] [CrossRef] [PubMed]
Moca, V.V.; Scheller, B.; Mureşan, R.C.; Daunderer, M.; Pipa, G. EEG under anesthesia—Feature extraction with TESPAR.Comput. Methods Programs Biomed.2009,95, 191–202. [Google Scholar] [CrossRef]
Diykh, M.; Miften, F.S.; Abdulla, S.; Saleh, K.; Green, J.H. Robust approach to depth of anaesthesia assessment based on hybrid transform and statistical features.IET Sci. Meas. Technol.2020,14, 128–136. [Google Scholar] [CrossRef]
Ahmadi, B.; Negahbani, E.; Amirfattahi, R.; Zaghari, B.; Mansouri, M. Extraction of BIS index sub-parameters in different anesthetic and sedative levels. In Proceedings of the 2008 9th International Conference on Signal Processing, Beijing, China, 26–29 October 2008; pp. 2665–2668. [Google Scholar] [CrossRef]
Nguyen-Ky, T.; Wen, P.; Li, Y. An Improved Detrended Moving-Average Method for Monitoring the Depth of Anesthesia.IEEE Trans. Biomed. Eng.2010,57, 2369–2378. [Google Scholar] [CrossRef]
Voss, L.; Sleigh, J. Monitoring consciousness: The current status of EEG-based depth of anaesthesia monitors.Best Pract. Res. Clin. Anaesthesiol.2007,21, 313–325. [Google Scholar] [CrossRef]
Kreuer, S.; Bruhn, J.; Larsen, R.; Bialas, P.; Wilhelm, W. Comparability of Narcotrend™ index and bispectral index during propofol anaesthesia.Br. J. Anaesth.2004,93, 235–240. [Google Scholar] [CrossRef]
Bruhn, J.; Myles, P.S.; Sneyd, R.; Struys, M.M.R.F. Depth of anaesthesia monitoring: What’s available, what’s validated and what’s next?Br. J. Anaesth.2006,97, 85–94. [Google Scholar] [CrossRef] [PubMed]
Drover, D.; Ortega, H.R. Patient state index.Best Pract. Res. Clin. Anaesthesiol.2006,20, 121–128. [Google Scholar] [CrossRef]
Nguyen-Ky, T.T.; Wen, P.; Li, Y.; Gray, R. Measuring and Reflecting Depth of Anesthesia Using Wavelet and Power Spectral Density.IEEE Trans. Inf. Technol. Biomed.2011,15, 630–639. [Google Scholar] [CrossRef]
Diykh, M.; Li, Y.; Wen, P.; Li, T. Complex networks approach for depth of anesthesia assessment.Measurement2018,119, 178–189. [Google Scholar] [CrossRef]
Afrasiabi, S.; Boostani, R.; Koochaki, S.; Zand, F. Presenting an effective EEG-based index to monitor the depth of anesthesia. In Proceedings of the 16th CSI International Symposium on Artificial Intelligence and Signal Processing (AISP 2012), Shiraz, Iran, 2–3 May 2012; pp. 557–562. [Google Scholar] [CrossRef]
Liu, Q.; Ma, L.; Fan, S.-Z.; Abbod, M.F.; Shieh, J.-S. Sample entropy analysis for the estimating depth of anaesthesia through human EEG signal at different levels of unconsciousness during surgeries.PeerJ2018,6, e4817. [Google Scholar] [CrossRef] [PubMed]
Shalbaf, A.; Saffar, M.; Sleigh, J.W.; Shalbaf, R. Monitoring the Depth of Anesthesia Using a New Adaptive Neurofuzzy System.IEEE J. Biomed. Heal. Inform.2017,22, 671–677. [Google Scholar] [CrossRef]
Särkelä, M.; Mustola, S.; Seppänen, T.; Koskinen, M.; Lepola, P.; Suominen, K.; Juvonen, T.; Tolvanen-Laakso, H.; Jäntti, V. Automatic Analysis and Monitoring of Burst Suppression in Anesthesia.Int. J. Clin. Monit. Comput.2002,17, 125–134. [Google Scholar] [CrossRef]
Lashkari, A.; Boostani, R. A Kalman-based instantaneous frequency estimation for anesthetic depth measurement. In Proceedings of the 2017 22nd International Conference on Digital Signal Processing (DSP), London, UK, 23–25 August 2017; pp. 1–4. [Google Scholar] [CrossRef]
Gu, Y.; Liang, Z.; Hagihira, S. Use of Multiple EEG Features and Artificial Neural Network to Monitor the Depth of Anesthesia.Sensors2019,19, 2499. [Google Scholar] [CrossRef]
Shalbaf, R.; Behnam, H.; Sleigh, J.W.; Steyn-Ross, A.; Voss, L.J. Monitoring the depth of anesthesia using entropy features and an artificial neural network.J. Neurosci. Methods2013,218, 17–24. [Google Scholar] [CrossRef]
Liang, Z.; Huang, C.; Li, Y.; Hight, D.F.; Voss, L.J.; Sleigh, J.W.; Li, X.; Bai, Y. Emergence EEG pattern classification in sevoflurane anesthesia.Physiol. Meas.2018,39, 045006. [Google Scholar] [CrossRef]
Peker, M.; Arslan, A.; Sen, B.; Celebi, F.V.; But, A. A novel hybrid method for determining the depth of anesthesia level: Combining ReliefF feature selection and random forest algorithm (ReliefF+RF). In Proceedings of the 2015 International Symposium on Innovations in Intelligent SysTems and Applications (INISTA), Madrid, Spain, 2–4 September 2015; pp. 1–8. [Google Scholar] [CrossRef]
Liu, Q.; Cai, J.; Fan, S.-Z.; Abbod, M.F.; Shieh, J.-S.; Kung, Y.; Lin, L. Spectrum Analysis of EEG Signals Using CNN to Model Patient’s Consciousness Level Based on Anesthesiologists’ Experience.IEEE Access2019,7, 53731–53742. [Google Scholar] [CrossRef]
Li, X.; Li, D.; Liang, Z.; Voss, L.J.; Sleigh, J.W. Analysis of depth of anesthesia with Hilbert–Huang spectral entropy.Clin. Neurophysiol.2008,119, 2465–2475. [Google Scholar] [CrossRef] [PubMed]
Donoho, D.L. De-noising by soft-thresholding.IEEE Trans. Inf. Theory1995,41, 613–627. [Google Scholar] [CrossRef]
Donoho, D.L. Nonlinear wavelet methods for recovery of signals, densities, and spectra from indirect and noisy data.Proc. Symp. Appl. Math.1993,47, 173–205. [Google Scholar] [CrossRef]
Sang, Y.-F.; Wang, D.; Wu, J.-C. Entropy-Based Method of Choosing the Decomposition Level in Wavelet Threshold De-noising.Entropy2010,12, 1499–1513. [Google Scholar] [CrossRef]
Liang, Z.; Wang, Y.; Sun, X.; Li, D.; Voss, L.J.; Sleigh, J.W.; Hagihira, S.; Li, X. EEG entropy measures in anesthesia.Front. Comput. Neurosci.2015,9, 16. [Google Scholar] [CrossRef]
Eriksson, A.; Stoica, P.; Soderstrom, T. Markov-based eigenanalysis method for frequency estimation.IEEE Trans. Signal Process.1994,42, 586–594. [Google Scholar] [CrossRef]
MATLAB.9.9.0.1524771 (R2020b); The MathWorks Inc.: Natick, MA, USA, 2020. Available online:https://www.mathworks.com/ (accessed on 15 May 2022).
Whitham, E.M.; Lewis, T.; Pope, K.J.; Fitzgibbon, S.P.; Clark, C.R.; Loveless, S.; DeLosAngeles, D.; Wallace, A.K.; Broberg, M.; Willoughby, J.O. Thinking activates EMG in scalp electrical recordings.Clin. Neurophysiol.2008,119, 1166–1175. [Google Scholar] [CrossRef]

Figure 1. The diagram of the new DoA index development.

Figure 2. Results of denoising raw EEG data, including the low amplitude (time is between the 2000 s and 2010 s) and spike noise (time is between 2980 s and 2990 s). (a) Raw EEG data, patient ID: L12161431. (b) Raw EEG data with low amplitude noise. (c) Raw EEG data having spike noise. (d) EEG signal after denoising low amplitude with threshold

T_{h}

. (e) EEG signal after denoising spike with threshold

T_{h}

. (f) EEG signal after denoising low amplitude with proposed threshold

T_{h n e w}

. (g) EEG signal after denoising spike with proposed threshold

T_{h n e w}

T_{h}

. (e) EEG signal after denoising spike with threshold

T_{h}

. (f) EEG signal after denoising low amplitude with proposed threshold

T_{h n e w}

. (g) EEG signal after denoising spike with proposed threshold

T_{h n e w}

Figure 3. Relationship between the Hurst range response value and the BIS value.

Figure 4. Predicted response from the SEGP-RM vs. actual values (BIS).

Figure 5. Response plot of the SEGR-RM for the 60 training subjects.

Figure 6. Pearson correlation coefficients of the testing results and BIS values for 13 subjects (patient 14 represents the average results of the 13 testing subjects).

Figure 7. (a) Bland-Altman plot between

N_{D o A}

index and BIS. (b) Distribution plot between new DoA index and BIS.

Figure 7. (a) Bland-Altman plot between

N_{D o A}

index and BIS. (b) Distribution plot between new DoA index and BIS.

Figure 8. The proposed index,

N_{D o A}

, can monitor the anaesthetic states change, but BIS cannot when SQI is lower than 15 (subject ID: L01121339).

Figure 8. The proposed index,

N_{D o A}

, can monitor the anaesthetic states change, but BIS cannot when SQI is lower than 15 (subject ID: L01121339).

Figure 9. (a) Comparison between

N_{D o A}

and BIS in the case of poor signal quality (subject ID: L01210841). (b) Zoom out during the period between 0 and 650 s. (c) Zoom out during the period between 16,940 and 17,748 s.

Figure 9. (a) Comparison between

N_{D o A}

and BIS in the case of poor signal quality (subject ID: L01210841). (b) Zoom out during the period between 0 and 650 s. (c) Zoom out during the period between 16,940 and 17,748 s.

Figure 10. Comparison between

N_{D o A}

and BIS in the case of poor signal quality (subject ID: L01131131).

Figure 10. Comparison between

N_{D o A}

and BIS in the case of poor signal quality (subject ID: L01131131).

Figure 11. Comparison between

N_{D o A}

and BIS in the case of poor signal quality (subject ID: L01131400).

Figure 11. Comparison between

N_{D o A}

and BIS in the case of poor signal quality (subject ID: L01131400).

Figure 12. (a) EMG causes the BIS index to increase at two spikes, and the trend of

N_{D o A}

does not change during these periods. Subject ID: L01041002. (b) Zoom out during the period between 2500 and 3000 s.

Figure 12. (a) EMG causes the BIS index to increase at two spikes, and the trend of

N_{D o A}

does not change during these periods. Subject ID: L01041002. (b) Zoom out during the period between 2500 and 3000 s.

Table 1. Summary of patient demographics.

Age (year)	Weight (kg)	Height (cm)	Gender (M/F)
2–83	11–101	90–194	39/34

Table 2. Training model regression of the anaesthetic states.

Training Model	$R^{2}$	RMSE	MSE	MAE
Robust Linear	0.43	11.21	125.58	8.81
Fine Tress	0.6	9.41	88.53	6.25
Fine Gaussian SVM	0.73	6.61	43.67	5.12
Squared Exponential Gaussian Process Regression	0.79	5.66	32.04	4.84

Table 3. Details of anaesthesia medication.

Subject	$E v e n t$	Timestamp	Time (s)
L01210841	E1. Add Dexmedetomidine Hydrochloride (30 ug)	08:25:00	-
	E2. Intraoperative EEG monitoring	08:41:41	0
	E3. Add Dezocine (6 mg)	08:55:00	799
	E4. Add Tropisetron Hydrochloride (2 mg)	09:05:00	1399
	E5. Add Ketorolac (30 mg)	09:16:00	2059
	E6. Surgery End	13:37:27	17,746

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Huang, Y.; Wen, P.; Song, B.; Li, Y. Real-Time Depth of Anaesthesia Assessment Based on Hybrid Statistical Features of EEG.Sensors2022,22, 6099. https://doi.org/10.3390/s22166099

AMA Style

Huang Y, Wen P, Song B, Li Y. Real-Time Depth of Anaesthesia Assessment Based on Hybrid Statistical Features of EEG.Sensors. 2022; 22(16):6099. https://doi.org/10.3390/s22166099

Chicago/Turabian Style

Huang, Yi, Peng Wen, Bo Song, and Yan Li. 2022. "Real-Time Depth of Anaesthesia Assessment Based on Hybrid Statistical Features of EEG"Sensors 22, no. 16: 6099. https://doi.org/10.3390/s22166099

APA Style

Huang, Y., Wen, P., Song, B., & Li, Y. (2022). Real-Time Depth of Anaesthesia Assessment Based on Hybrid Statistical Features of EEG.Sensors,22(16), 6099. https://doi.org/10.3390/s22166099

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further detailshere.

Article Metrics

Article Access Statistics

For more information on the journal statistics, clickhere.

Multiple requests from the same IP address are counted as one view.

Movatterモバイル変換

Article Menu

Real-Time Depth of Anaesthesia Assessment Based on Hybrid Statistical Features of EEG

Abstract

1. Introduction

2. Materials and Methods

2.1. Subjects and EEG Recording

2.2. Signal Denoising and Pre-Processing

2.3. Features’ Extraction

2.4. Regression Models and Evaluation Measures

2.5. Real-Time DoA Monitor under Sliding Window Framework

3. Results

3.1. Raw EEG Signal Pre-Processing and Features’ Selection

3.2. Regression Model Training

3.3. Performance of the Training Models and Results of the Nex Index

4. Results and Discussion

4.1. Length of Sliding Window

4.2. The Effect of Signal Quality

4.3. The Effect of EMG Signals

5. Conclusions

Author Contributions

Funding

Institutional Review Board Statement

Informed Consent Statement

Data Availability Statement

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI