Changes in EEG frequency characteristics during sevoflurane general anesthesia: feature extraction by variational mode decomposition

Mode decomposition is a method for extracting the characteristic intrinsic mode function (IMF) from various multidimensional time-series signals. Variational mode decomposition (VMD) searches for IMFs by optimizing the bandwidth to a narrow band with the H1\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathcal{H}}^{1}$$\end{document} norm while preserving the online estimated central frequency. In this study, we applied VMD to the analysis of electroencephalogram (EEG) recorded during general anesthesia. Using a bispectral index monitor, EEGs were recorded from 10 adult surgical patients (the median age: 47.0, and the percentile range: 27.0–59.3 years) who were anesthetized with sevoflurane. We created an application named EEG Mode Decompositor, which decomposes the recorded EEG into IMFs and displays the Hilbert spectrogram. Over the 30-min recovery from general anesthesia, the median (25–75 percentile range) bispectral index increased from 47.1 (42.2–50.4) to 97.4 (96.5–97.6), and the central frequencies of IMF-1 showed a significant change from 0.4 (0.2–0.5) Hz to 0.2 (0.1–0.3) Hz. IMF-2, IMF-3, IMF-4, IMF-5, and IMF-6 increased significantly from 1.4 (1.2–1.6) Hz to 7.5 (1.5–9.3) Hz, 6.7 (4.1–7.6) Hz to 19.4 (6.9–20.0) Hz, 10.9 (8.8–11.4) Hz to 26.4 (24.2–27.2) Hz, 13.4 (11.3–16.6) Hz to 35.6 (34.9–36.1) Hz, and 12.4 (9.7–18.1) Hz to 43.2 (42.9–43.4) Hz, respectively. The characteristic frequency component changes in specific IMFs during emergence from general anesthesia were visually captured by IMFs derived using VMD. EEG analysis by VMD is useful for extracting distinct changes during general anesthesia.


Introduction
The theoretical background behind using electroencephalogram (EEG) recording to evaluate the effects of general anesthesia in anesthetized patients involves observation of the changes in the electrical activity of cerebral cortical neurons that follow consciousness changes occurring in response to the dose of an anesthetic agent [1].As a measure of biological information reflecting the effects of general anesthesia, EEG is a surface recording of the summed electrophysiological activity of the cerebral cortex, which varies over time according to surgical stimuli and general anesthetic doses that affect the level of consciousness.More specifically, the EEG during wakefulness consists mainly of α, β, and γ waves with frequencies higher than 8 Hz.Under general anesthesia using anesthetics such as propofol and inhalational anesthetics, α wave components of approximately 10-14 Hz and θ and δ wave components with higher amplitude appear.These EEG changes have been neuroscientifically explained by anesthetics increasing the sensitivity of GABA A receptors, promoting inhibitory neurotransmission, suppressing the formation of action potentials in the cerebral cortex and thalamus and hyperpolarizing them, and then generating a burst firing rhythm of thalamocortical projection neurons [2][3][4].The processed EEG monitor analyzes these changes in a clinical setting and replaces them with EEG-index values.Anesthesiologists follow changes in EEG values over time, considering factors such as surgical invasion and anesthetic agent.Commercially available processed EEG monitors, such as the bispectral index (BIS) monitor (Medtronic, Minneapolis, MN, USA) and SEDline (Masimo, Irvine, CA, USA), provide processed EEG-derived parameters such as BIS and PSi (patient status index), and density spectral arrays that display two-dimensional color power spectra of frequency components over time.A fundamental challenge in measuring the effects of general anesthesia on consciousness is developing technology that objectively captures the constantly changing EEG components in general anesthesia and replaces them with more specific parameters.
Commercially available processed EEG monitors generally display the results of Fourier-transform-based spectral analysis and a spectrogram showing spectral changes in time, thereby providing information on EEG changes during anesthesia that anesthesiologists can quickly and visually grasp.However, the Fourier transform is not the only method of analyzing the frequency component of signals, and other analysis methods that capture temporal changes in signals in an easy-to-understand manner have been reported.In particular, mode decomposition has attracted attention in various fields because it allows capture of nonlinear biological signals as multidimensional time-series data and the extraction of essential information on feature structures hidden therein [5].One such technique, mode decomposition combined with the Hilbert transform, has recently been used to analyze various nonstationary signals to obtain instantaneous frequencies and amplitudes.In mode decomposition, waveforms are decomposed into several characteristic modes called intrinsic mode functions (IMFs), and if all the IMFs are summed, the original waveform can be reproduced.When applied to EEG analysis, mode-decomposed waveforms are EEG components with specific narrow frequency band characteristics.Features can be extracted from EEG frequency changes by adopting a Hilbert transform that obtains instantaneous frequencies and amplitudes.Among the available mode decomposition methods, empirical mode decomposition (EMD) is a component of the Hilbert-Huang transform proposed by Huang et al., and has recently been reported as a useful method for analyzing EEG during general anesthesia [6][7][8][9][10][11].
EMD has been reported to be an ad-hoc algorithm that lacks mathematical theory [12], and there have therefore been various proposals to replace EMD with mode decomposition methods using mathematically more robust approaches.One such proposal is variational mode decomposition (VMD) [12].In VMD, the central frequency is iteratively estimated as the centroid of the mode's power spectrum to adaptively determine the relevant band motivated by narrowband characteristics.VMD finds a modal ensemble that reconstructs a given input signal by narrowband optimization of the modal bandwidth with the H 1 norm while each IMF preserves the online-estimated central frequency.This study aimed to develop an adaptive algorithm for the time-frequency component analysis of EEGs of patients under general anesthesia through mode decomposition using the VMD method, and to explore the possibility of developing EEG parameters related to the effects of general anesthesia on consciousness.To analyze an EEG during general anesthesia, we applied the VMD method and Hilbert transform in the form of a modified Hilbert-Huang transform, and investigated changes in the time-frequency amplitude characteristics of an EEG recorded through the time-series progression of general anesthesia.
2 Theoretical background, methods, and materials

Anesthesia management and data acquisition
All experiment protocols involving humans were conducted in accordance with the principles of the Declaration of Helsinki.
The current study was approved (No.ERB-C-1074) by the institutional review board for human experiments at the Kyoto Prefectural University of Medicine (IRB of KPUM), and for this non-interventional and noninvasive retrospective observational study, the requirement for informed patient consent was waived by the IRB of KPUM; patients were provided with an opt-out option, of which they were notified in the preoperative anesthesia clinic.EEG data were collected from 10 patients whose demographics and clinical characteristics are summarized in Table 1.In our facility, a BIS monitor is routinely used for adult and pediatric patients undergoing general anesthesia surgery.The anesthesiologists in charge of management did not receive notice of the study, but planned the anesthesia methods for scheduled surgeries under our facility's standard care protocol without any feedback regarding the online analysis of processed EEG signals.According to our facility's standard protocol, patients were not premedicated before anesthesia induction.Anesthesia was induced through small doses of fentanyl (1 µg kg −1 per dose) with continuous intravenous infusion of remifentanil (0.125-0.25 µg kg −1 min −1 ) followed by a one-shot intravenous injection of propofol (2 mg kg −1 ) and rocuronium (0.8-1.0 mg kg −1 ).Tracheal intubation was then conducted, and the anesthesia was maintained through inhalation of sevoflurane (2-3%), additional maintenance doses of fentanyl (1 µg kg −1 per dose), and rocuronium (0.2 mg kg −1 at intervals of 20-30 min), and continuous intravenous infusion of remifentanil (0.125-0.25 µg kg −1 min −1 ).At the emergence phase from general anesthesia, sugammadex (4 mg kg −1 ) was administered before extubation.To acquire EEG data, we used the software application EEG Analyzer [13,14], through which raw EEG signals were recorded as text files on a personal computer (Surface Pro 4, Microsoft Co., Redmond, WA, USA) via the RS-232 interface of a VISTA A-3000 BIS monitor (VISTA Application revision 3.22, Platform revision 2.03, Medtronic, Minneapolis, MN, USA) with a BIS Quatro sensor mounted on the frontal regions according to the manufacturer's recommendations.Note that the BIS VISTA records EEG at a sampling rate of 16,384 samples/sec with a default bandpass filter set of 0.25-100 Hz, and raw EEG waveforms (16-bit raw EEG signals at a sampling frequency of 128 Hz, eight packets/second) with processed EEG parameters, such as EMGlow, BIS value, burst suppression ratio (SR), electromyography (EMG) parameter EMGlow, and signal quality index (at one packet/second) were acquired using the binary output mode of the BIS monitor through its serial output.We analyzed the changes in correlations over time between the raw EEG waves and parameters such as the observed BIS, SEF 95 , total power, and EMG low .

Algorithm for VMD and the Hilbert transform
VMD aims to divide a real-valued input signal f into k narrowband sub-signals u k (IMFs) characterized by a discrete central frequency ω k with a particular sparsity property.VMD involves three steps, according to the original description published by Dragomiretskiy et al. [12].
1. Compute the analytic signal associated with each mode according to formula (1) using the Hilbert transform, where * stands for convolution.This results in each mode having a purely positive spectrum: where u k (t) is the k th IMF and δ(t) stands for a pulse signal.2. Demodulate the analytic signal to the baseband by multiplying with an exponential function adjusted to the respective estimated central frequency ωk: 3. Estimate the bandwidth using the H 1 Gaussian smooth- ness of the demodulated signal; i.e., the squared L 2 norm of the slope of the demodulated analytic signal.For each mode u k , compute the associated analytic signal using the Hilbert transform to obtain the one-sided frequency spectrum.The constrained variational evaluation of a given signal x(t) is described as where {u k }≔ {u 1 ,……,u k }, and {ω k }≔ {ω 1 ,……,ω k } are shorthand notations for the set of all modes and their central frequencies, respectively.The sum ∑ k ∶= ∑ K k=1 of all mode functions u k is equal to the specified time series signal f(t).
The solution to this constrained variational problem introduces both a quadratic penalty term and a Lagrangian multiplier λ in the mathematical optimization method to convert the problem to that of an extended Lagrangian function L without constraints: As a result, a mathematical analysis can be adopted using an optimal extended Lagrangian function that "Finds the one that minimizes f(x) among x ∈ ℝ n that satisfies the constraint g(x) = 0": (2) We solve the unconstrained problem using the quadratic penalty term α and the Lagrangian multiplier (double rise) λ in Eq. ( 3): (5) min The solution to this minimization problem is found as a saddle point of the extended Lagrangian L in a series of iterative suboptimizations called the alternate direction method of multipliers (ADMM).This process adopts VMD in the first step of computing u 1 k , 1 k , and 1 .These are further updated to a degree (n) in the next step (u n+1 k , n+1 k , and n+1 ).Fur- ther computations can be reduced to algorithms that optimize directly in the frequency domain.The analyzed signal is decomposed into modes ( u k ) using the central frequency ( k ) in Eq. ( 7) of the ADMM algorithm for VMD, as reported by Dragomiretskiy et al. [12].Supplementary Video 1 shows the VMD suboptimization of IMFs by ADMM using Python (ver.3.8) code (EEG data from case #6).
Algorithm: optimization concept for VMD [12]: Update ̂for all ≥ 0: Update ̂: until convergence: For an EEG signal x(t), VMD decomposes the signal into a series of IMFs, Cn (n = 1, 2,..., N), where N is the number of IMFs.After VMD, a signal x(t) is written as Applying the Hilbert transform to the IMF components, it follows that in which where (t) and a(t) are the instantaneous frequency and amplitude of the IMF in obtaining a time-frequency distribution for signal x(t) and the Hilbert amplitude spectrum,

H(x, t).
In summary, in the VMD algorithm, each mode is updated in the frequency domain and the central frequency is recomputed at each iteration.This places the new ω k at the centroid of the power spectrum of the corresponding mode.This average carrier frequency is the frequency of the least-squares linear regression of the instantaneous phase observed in the mode.
The VMD analysis software EEG Mode Decompositor (ver 2.5, Mac OS X [Apple silicon] and MS-Windows 64 bit version, https:// anesth-kpum.org/ blog_ ts/?p= 4018) (Fig. 1) was developed for EMD [11] and updated by porting the original Python VMD algorithm code (VMD_Python; Shenmusmart, https:// github.com/ shenm usmart/ VMD_ python, [15]) based on the report published by Dragomiretskiy et al. [12] to Processing [16] with its integrated development environment (IDE) (Supplementary Program Code 1).The Apache Commons Math Library (The Apache Software Foundation, https:// commo ns.apache.org) for the Fourier transform, and Digital Signal Processing in Java (JDSP, https:// jdsp.dev/ index.html) [17] for the Hilbert transform, were incorporated into the Processing code.For drawing graphs, we used our SwGraph library that runs in Processing.The EEG Mode Decompositor displays the original EEG wave (8 s), IMFs, FFT spectrums, FFT spectrograms (color density spectral arrays [DSAs]) calculated using a Blackman window, the Hilbert spectrum of the IMFs, and Hilbert spectrograms (11) of the IMFs (Fig. 1).Supplementary Program Code 3 shows the Python code used for VMD analysis.Although the relationship between the isolated IMF and the narrow band's center frequency depends on the power distribution of the overall EEG frequency components, we wanted to observe significant changes in specific narrow-frequency bands in the 0-64 Hz range, such as δ, θ, α, β (low-β, high-β), and γ waves, so we decided to separate the EEG into six IMFs as a default mode.However, we also tried VMD decomposition into two to five IMFs in the analyses performed in this study.

Data processing and statistics
Microsoft Excel for Mac (ver.16.16.5,Microsoft Corp., Redmond, WA, United States) was used to create graphs.Changes in various EEG parameters between two-time points, 30 min before and at the time of emergence, were statistically compared by conducting Kruskal-Wallis tests with SPSS Statistics (28.0.1.0,IBM Corp., Armonk, NY, USA), and p-values < 0.05 were considered significant.The RMSE and R 2 were calculated as estimators for the regression analyses.

Results
In the first part of this section, we outline the EEG analysis by the VMD method in the three phases of general anesthesia induction, maintenance, and awakening, using recorded EEG from case #6.In the second part, we discuss burst suppression and VMD in deep anesthesia.In the third part, we focus on the VMD analysis of EEG during the emergence phase in general anesthesia using recorded EEG from all 10 cases.

Case analysis of the VMD application
Using the results of case #6 (Table 1, 27-year-old woman, ovarian tumor/laparoscopic oophorectomy, sevoflurane anesthesia), we present a typical example of VMD adaptation decomposing the EEG into six IMFs (as a default mode of EEG Mode Decompositor) from induction to emergence (for 3 h 37 min) during sevoflurane anesthesia (Supplementary Video 2) [19].Figure 2 presents the time-course changes in the BIS, SEF 95 , and Hilbert spectrum of each IMF at 2-min intervals for 10 min around the induction of anesthesia (− 3 min to + 7 min), 10 min after induction (+ 57 min to + 1 h 7 min), and 10 min before awakening (+ 3 h 27 min to + 3 h 37 min).At the time of induction of anesthesia, the BIS changed from 95 to 30 in 1 min, and SEF 95 decreased from 26 to 7 Hz (owing to the intravenous administration of fentanyl and remifentanil).Along with the hypnotic effect, the peak value of the Hilbert spectrum at IMF-4, IMF-5, and IMF-6 shifted from the region below 10 Hz to 10-14 Hz.During the maintenance of general anesthesia, the BIS was approximately 40, SEF 95 was approximately 16 Hz, and the α wave rhythm of approximately 11-14 Hz was detected in the IMF-4 to IMF-6 region.During the emergence phase, the BIS increased from 40 to 98 for 10 min, SEF 95 gradually increased from 15 to 27 Hz, and the peak value at IMF-5 and IMF-6 increased from approximately 12 Hz to over 15 Hz.In IMF-5 and IMF-6, the peak power decreased, the power distribution appeared in the high-frequency band above 20 Hz, and there was dispersion over a wide area.
Next, for the EEG of case #6, 64-s epochs of the EEG (recorded at 128 Hz) were analyzed in three stages: the maintenance period (at + 3 h 26 min), the transition period from maintenance to awakening (at + 3 h 35 min), and the emergence period (at + 3 h 36 min).The EEG for the first 4 s of the epoch is plotted in Fig. 3A.The power spectrum of the epoch was analyzed by adopting the multitaper method and a short-time Fourier transform (Fig. 3B).The spectrogram of this epoch is shown in Fig. 3C.Supplementary Program Code 2 shows the Python code used for spectral analysis.
During the maintenance period of anesthesia, spindle waves (which are recognized as α-wave rhythms) with a peak at 11 Hz, and δ waves below 4 Hz, were observed.In the transition period, the α wave rhythm fluctuated between 8 and 14 Hz then gradually disappeared.The power of 20-30-Hz β waves increased.With awakening, the α wave rhythm disappeared completely, the power decreased at all frequencies, and waves appeared at 30 Hz and higher frequencies.
The VMD of the EEG epoch is shown with its six IMFs in Fig. 4A.These IMFs are almost orthogonal components with no overlap in frequency, as shown in Fig. 4B.The corresponding Hilbert spectrum and the marginal spectrum are shown in Fig. 4B and C, respectively.We applied the VMD algorithm to decompose the related 4-s EEG into six IMFs, applied the Hilbert transform to each IMF, obtained the During the maintenance period of anesthesia, IMF-1 to IMF-6 were in the region of 18 Hz or lower.There were δ waves below 4 Hz in IMF-1 to IMF-3, α wave rhythms in IMF-4.However, in the transition period, IMF-6 largely shifted to the 30-Hz band, and with awakening, the peak of IMF-6 rose to approximately 34 Hz.Except for IMF-1, all IMFs moved to higher frequency bands.

Burst suppression and the VMD application
Next, in our data set, we observed the relationship between burst suppression (which can be observed in EEG because of the strong anesthetic effect) and mode decomposition.Note that it was reported that suppression ratio values > 40% are linearly correlated with BIS values from 30 to 0 [20].Of the 10 cases in this analysis, a burst suppression event, as judged by the BIS monitor, was recorded in four patients.In two of these cases the events occurred when the BIS value exceeded 70 during the emergence phase from general anesthesia, and they were considered to be erroneous evaluations.In the remaining two events, SR was 8 for one minute at a BIS value of approximately 50 in case #3, and SR was 8-17 for approximately 3 min immediately after propofol-bolus infusion at the time of induction of anesthesia in case #4.The burst suppression event and VMD analysis in case #4 are summarized in Supplementary Fig. 1.At the induction of anesthesia, immediately after the propofol administration, when the BIS value rapidly decreased from 97 to 23, burst suppression was recorded with a suppression ratio of 8%-18% for approximately 3 min.In this state, IMF-4 captured the cycling up and down in amplitude of the burst suppression wave for several minutes.However, because the VMD was analyzed using the EEG signal from an 8-s epoch, the IMF-4-derived Hilbert spectrum calculated from a spectral distribution of one EEG epoch appeared unaffected by the short-time-duration burst suppression event.

The analysis of VMD application during the emergence phase from sevoflurane-inhalational anesthesia in 10 cases
Next, for each of the 10 patients, we analyzed the change in each parameter's average value and standard deviation over the 30 min before emergence.The median BIS increased from 47.1 (25th-75th percentile range, 42.2-50.4)to 97.4 (96.5-97.5)(p < 0.001), SEF 95 increased from 14.5 (13.2-15.9)Hz to 22.6 (20.3-24.1)Hz (p = 0.001), and the total power decreased from 64.2 (61.5-65.8)dB to 58.9 (56.9-62.2) dB (p = 0.179).EMG low showed a significant increase from 27.6 (26.5-27.9)dB to 52.The current algorithm leaves it up to the user to decide how many IMFs to separate in VMD.In VMD, online optimization is performed around frequency bands with significant powers, so reducing the number of discrete IMFs may not significantly increase the bandwidth of IMFs.To understand how the number of IMFs selected for the decomposition influences each IMF's spectrogram, we here tried VMD decomposing into two to five IMFs from induction to emergence (for 3 h 37 min) for case #6 (Supplementary Video 3-6, played at 230 times the actual speed).Mode Decompositor screen).During the maintenance phase of anesthesia, the Hilbert spectrogram of IMF-3 derived from the VMD decomposition into three IMFs successfully extracted the α wave-peak power shown in the DSA.As a result, the intensity of the extracted α wave power in the Hilbert spectrograms was increased because of the reduction in the number of IMFs to be resolved.Supplementary Fig. 2 shows the Hilbert spectrogram of each IMF derived from VMD decomposition into six IMFs and DSAs for 30 min before emergence in all ten cases.IMF-1-IMF-2, Fig. 6 Time course changes of the central frequency and the total power of of IMFs for 30 min before emergence from general anesthesia.The total power was displayed on a log scale (in decibels, dB; obtained from power (P, µV 2 ) by dB = 10 log 10 P).Data are medians (dark blue lines), 25th-75th percentile ranges (light blue area), and max-min ranges (grey area).*p < 0.05 against the value at − 30 min by Kruskal-Wallis test IMF-3, IMF-4, IMF-5, and IMF-6 well-extracted the δ, θ, α, low β, and high β range, respectively.The increases in the γ frequency band power were covered by IMF-4 to IMF-6 towards the emergence from anesthesia.
We performed Gaussian process regression analysis to examine whether the fluctuations in the central frequency band of each IMF correlated with the changes in the BIS value for the 10 min up until emergence in the 10 cases, as shown in Supplementary Fig. 3. Blue dots indicate measurements, light-blue bands indicate the mean plus a one standard deviation (68%) credible interval, and orange lines indicate the posterior prediction distribution of the mean (for 30 virtual cases).The RMSE, mean absolute error (MAE), root mean squared logarithmic error (RMSLE), and R 2 of the posterior predictive distribution for the mean were calculated.The results showed IMF-4 to have the strongest correlation with the BIS, with an RMSE of 12.65, MAE of 10.30, RMSLE of 0.212, and R 2 of 0.507.This demonstrates that the response recorded in the BIS during the emergence phase did not simply reflect the changes in IMF characteristics.

Discussion
Spectral analysis is often applied to measurements to understand the frequency characteristics of target signals.Fourier analysis based on trigonometric functions is commonly used in spectral analysis.However, because Fourier analysis assumes stationarity and periodicity of the data, it is not possible to accurately understand the frequency characteristics from the analysis results when dealing with highly nonstationary data.Although time-frequency analysis such as short-time Fourier analysis and wavelet analysis is effective for understanding the frequency characteristics of such highly nonstationary data, these methods also have drawbacks associated with the uncertainty principle in time-frequency analysis.An analysis window function such as a Hann (Hanning), Hamming, Blackman, or multitaper function is used in short-time Fourier analysis.When the width of the finite interval of the window function is reduced, the characteristics of the signal can be examined within a short period, but the resolution is reduced.Conversely, increasing the width of the finite interval of the window function increases the frequency resolution but decreases the time resolution.In other words, there is a trade-off relationship (in analogy to Heisenberg's uncertainty principle) between the frequency and time resolution, and the frequency and time cannot be specified at the same time.
Mode decomposition of signals is technically similar to the procedure used in radio and television broadcasting.In radio and television broadcasting, the video and audio electric signals sent from multiple broadcasting stations are captured through an antenna or a network cable and are separated into their specific frequency bands using decoders such as radios and televisions.Then, the audience selects the specific channel to listen to or watch.An EEG is recorded on the scalp as a sum of physiological electrical activities derived from a group of nerve cells in multiple anatomical regions of the cerebrum.Decomposing EEG waves into modes is similar to decoding television or radio waves into specific transmission bands (Fig. 8).However, there is a difference between broadcast radio waves and EEG waves.In EEG waves, unlike radio broadcasting, the broadcasting station is unknown (it may be an anatomical part of the brain), as is the content of the broadcasts (electrical activity of brain nerve cell groups) and under what situations the broadcasts occur.EEG signals are decomposed into modes according to a particular mathematical rule without a clear indication of the broadcasting station or the program schedule.
Here, where we consider the effects of general anesthesia with anesthetics such as sevoflurane, desflurane, and propofol, the emergence of sleep spindles in the 10-14-Hz band (located in the α frequency band that appears with general anesthesia) is a major characteristic.Synchronous membrane potential amplitude in thalamic neurons and cerebral cortical pyramidal cells is considered to be the mechanism of sleep spindle generation.Another feature that general anesthesia enhances is the appearance of slow waves belonging to the δ wave region.In contrast, the awakening process from general anesthesia is characterized by the disappearance of sleep spindles and the appearance of fast waves in the β and γ regions.It is currently of interest whether these features occurring during general anesthesia can be effectively extracted from EEG by mode decomposition.
The EMD mode decomposition method was proposed by Huang et al. [6] in 1998.Originally, the Hilbert-Huang transform method adopted EMD and obtained the instantaneous frequency and amplitude using the Hilbert transform for the extracted IMF.The analysis of EEG under general anesthesia using the EMD method has been reported [8][9][10], including our most recent publication about the application of the Hilbert-Huang transform using EMD [11].Although the EMD method can perform mode decomposition with a simple algorithm, the decomposition signal depends on detecting extreme points, interpolating extreme points to the carrier envelope, and the stopping criterion imposed by repeated operations.It was suggested that there is room for theoretical development and improvement regarding the robustness of the decomposition method owing to the heavy dependence on, and fragility of, the mathematical theory.For the above reasons, we applied a modified Hilbert-Huang transform combining VMD with the Hilbert transform to create a frequency analysis technique different from EEG analysis during general anesthesia.
The VMD mode decomposition method was proposed by Dragomiretskiy et al. in 2014 [12].In VMD, the bandwidth of an IMF is reduced to a narrow band with the H 1 norm while maintaining the function's central frequency, which is estimated online by solving a variational problem (a problem of finding a procedure that minimizes a process) to find the IMF to be optimized and reconstructed.Because the problem is solved mathematically as a constrained optimization problem, the solution to the formula has two main aspects: the objective function and the constraints.The goal of VMD is to minimize the bandwidth (objective function) at the central frequency of each mode, subject to the sum of all modes equaling the original signal.Past studies applying mode decomposition to EEG under general anesthesia have been limited to the EMD method, including our recent publication [11]; no studies have reported on mode decomposition applying the VMD method to EEG analysis during general anesthesia.The actual mathematical correspondence of the VMD method needs to be by tuning an arithmetic loop using the ADMM method, which requires more complex arithmetic processing than the EMD method.We ported the Python programming code to the Java application language Processing in the present example.Within Processing's IDE, we performed continuous EEG analysis for each 8-s epoch in the VMD.In this study, the EEG during general anesthesia was decomposed into six IMF modes, and the Hilbert transform was performed on each IMF to obtain the instantaneous frequency and amplitude.
In particular, VMD extracted the sleep spindle component of the α wave, which characterizes the effect of general anesthesia and captures the changes in hypnosis and arousal due to general anesthesia.Compared with the EMD method, the VMD calculation algorithm, which is based on mathematical optimization theory, is complex and takes a little longer to run on a computer.However, we confirmed that the algorithm for time-frequency analysis could be constructed for practical use without any execution problems using the current level of personal computing power.In the EMD method, the extracted IMF proceeds from the high-frequency band picked up by IMF-1 to the low-frequency band in order.In other words, the IMF-1 of EMD becomes the narrowband part with the highest center frequency.In the EMD during the maintenance period of general anesthesia, sleep spindles are likely to be extracted from IMF-1 to IMF-2.The distribution ranges (> 15 Hz) of the center frequencies in the Hilbert spectrum of IMF-1 and IMF-2 were greater with EMD than with VMD.
The optimal number of separated IMFs in VMD requires further investigation.In the VMD method decomposing into six IMFs, the frequencies are distributed from IMF-1 in the low-frequency band to IMF-6 in the high-frequency band.The sleep spindles during general anesthesia are extracted into mainly IMF-4.The distribution ranges of the center frequencies in the Hilbert spectrum were less than 10 Hz in IMF-3 and IMF-4, which are narrower than (> 15 Hz) those in EMD.During emergence from general anesthesia, the central frequencies of IMF-2 to IMF-6 increase and the total powers of IMF-4 to IMF-6 decrease because of the effects of β and γ waves that appear simultaneously.
During general anesthesia and deep sedation procedures, anesthesiologists should consider changes in the level of consciousness to avoid the possible disadvantages of too light or too deep anesthesia/sedation when administering anesthetics and sedatives.It would be helpful to quickly capture EEG changes and evaluate them to provide information on the action of the anesthetic on the consciousness level.In typical analysis used to capture changes in an EEG related to the so-called depth of anesthesia, many parameters depend on trade secrets, although the recent reverse engineering of the BIS monitor clarified the algorithm for the calculation of BIS values as openibis [21].According to past reports [22], the fluctuation of the BIS value does not correlate with the pharmacological efficacy of the general anesthetic.Calling the BIS value calculated by the BIS monitor a 'depth of anesthesia index' is controversial, and we should be careful when considering the concept of the depth of anesthesia and its association with pharmacologically correlative effects on various anesthetic levels [23].All the same, during general anesthesia with sevoflurane, the resolving of the EEG into IMFs with VMD was excellent for extracting the specific EEG characteristics enhanced by general anesthesia.Compared with displaying the entire frequency image in the DSA, the VMD analysis may provide anesthesiologists with better monitoring of the state of general anesthesia.We believe that new EEG analysis methods, including the VMD method applied here, will lead to better strategies for measuring the anesthetic action on various anesthetic levels, levels that anesthesiologists wish to target for safe patient management during invasive procedures.

Research involving human participants
The data used as an example in the data analysis of this review were approved and the requirement for written informed consent was waived by the IRB of KPUM.
Consent to participate and for publication For this non-interventional and noninvasive retrospective observational study, the requirement for informed patient consent was waived by the IRB of KPUM; patients were provided with an opt-out option, of which they were notified in the preoperative anesthesia clinic.

Fig. 1
Fig.1The application software mode decompositor (v.2_5N) used in the analysis of the mode decomposition of an electroencephalogram (EEG).The EEG, intrinsic mode functions (IMFs), fast Fourier

Fig. 7
Fig. 7 The DSA and the Hilbert spectrogram (HSG).The DSA (left graphs), the HSG of IMF-2 derived from VMD intotwo IMFs, the HSG of IMF-3 derived from VMD into three IMFs, and the HSG of

Fig. 8
Fig. 8 Mode decomposition of signals is technically similar to decomposition of radio and television broadcasts.In radio and television broadcasting, the video and audio signals sent from multiple broadcasting stations are separated into their original frequency bands using the decoders of receivers such as radios and televisions.A specific channel is selected to watch a broadcast program.An EEG recorded on the scalp records the combined physiological electrical activities of groups of nerve cells in multiple anatomical regions of the cerebrum.Decomposing EEG waves into modes is similar to decoding television or radio waves into specific transmission bands.However, there is a difference between broadcast radio waves and EEG waves.In EEG waves, unlike radio broadcasting, it is unknown which broadcasting station (anatomical part of the brain) broadcasts what kind of content (electrical activity of brain nerve cell groups) and when.In other words, unlike radio broadcasting, radio signals are decomposed into modes according to a particular mathematical rule without a clear indication of the broadcasting station and the program schedule

Table 1
Patient characteristics and anesthesia management parameters Data are shown as a median [25-75 percentile], or b mean ± sd.ht: height, bw: body weight, TUR-Bt: transurethral resection of bladder tumor, TKA: total knee arthroplasty