M O N I T O R I N G T H E A N E S T H E T I C - I N D U C E D UNCONSCIOUSNESS (HYPNOSIS) U S I N G W A V E L E T A N A L Y S I S OF T H E E L E C T R O E N C E P H A L O G R A M by Tatjana Zikov Bachelor of Science, University of Belgrade, Belgrade, Yugoslavia, 1995 A thesis submitted in partial fulfillment of the requirements for the degree of Master of Appl ied Science in T H E F A C U L T Y O F G R A D U A T E S T U D I E S Department of Electrical and Computer Engineering We accept this thesis as conforming to the required standard T H E U N I V E R S I T Y OF BRITISH C O L U M B I A October 2002 © T a t j a n a Zikov, 2002 In presenting t h i s thesis i n p a r t i a l f u l f i l l m e n t of the requirements f o r an advanced degree at the University of B r i t i s h Columbia, I agree that the Library s h a l l make i t f r e e l y a v a i l a b l e for reference and study. I further agree that permission for extensive copying of t h i s thesis f o r s c h o l a r l y purposes may be granted by the head of my department or by his or her representatives. It i s understood that copying or p u b l i c a t i o n of this thesis for f i n a n c i a l gain s h a l l not be allowed without my written permission. Department of The U n i v e r s i t y of B r i t i s h Columbia Vancouver, Canada October n 1002. ^ 1-Abstract This thesis presents a novel method for the discrimination between various levels of anesthetic-induced unconsciousness (i.e. hypnosis) in the context of intra-operative anesthesia. This is achieved by means of the wavelet analysis of the patient's Electroencephalogram (EEG). While monitors of hypnosis are already available, they rely on complex signal processing algorithms, and in some instances, on multivariate statistical analysis of a large quantity of data. To date, the BIS® monitor (Aspect Medical Systems, Inc., MA) provides the most accepted solution. It displays an index calculated based on bispectral analysis, and has been the subject of many clinical trials over the past 5 years. As such, it is widely known and accepted by the medical community, particularly in the United States. Wavelet analysis can be viewed as a generalization of Fourier analysis that introduces time localization in addition to frequency decomposition of a signal. This major advantage of wavelets makes them particularly suitable for the analysis of non-stationary signals such as the E E G . In this thesis, we present a simple wavelet-based technique that calculates, in real time, an index of hypnosis (so-called WAV Index) based on a patient's electroencephalogram. Wavelet analysis significantly reduces the computational complexity to perform this task. Further, neither a large number of patients nor an extensive amount of clinical E E G data is needed for the index derivation. In addition, a technique for the de-noising of the electroencephalogram to remove ocular artifacts is proposed. These artifacts, which result from ocular activity (eye blinks and movements), strongly perturb the E E G signal by superimposing large spikes and step-like patterns. Hypnotic states such as the awake state and light sedation are typically corrupted by these waveforms. The solution proposed successfully removes these artifacts from the E E G , thus providing a clean signal for the wavelet analysis. This technique further works in the case of head movement artifacts. Finally, clinical validation has been carried out to assess the performance of the proposed index. While being particularly well-correlated to the bispectral index, the proposed index has a faster response to fast transients and it is more consistent in terms of patient intra-variability. ii Table of Contents Abstract ii List of Tables vi List of Figures vii List of Acronyms x Acknowledgements • xiii Dedication xiv 1 Introduction 1 1.1 Anesthesia versus Depth of Anesthesia 1 1.1.1 Early Concepts 1 1.1.2 Restatement of Concepts 3 1.1.3 A Hypothesis of Anesthesia 5 1.2 Monitoring the State of Anesthesia . . 6 1.2.1 Prevention of Awareness 6 1.2.2 Can We Monitor Unconsciousness? 7 1.2.3 An Ideal Monitor 8 1.3 Intraoperative Assessment of Consciousness 8 1.3.1 Clinical Signs 9 1.3.2 Other Approaches 10 1.4 The E E G as a Measure of the Anesthesia Adequacy 10 1.4.1 Quantitative Analysis of the Spontaneous E E G 11 1.4.2 Evoked Potentials 17 1.5 Bispectral Analysis of the EEG for Monitoring the State of Anesthesia 18 1.5.1 Bispectral Analysis 18 1.5.2 Bispectral Index Scale (BIS) Monitor 19 1.6 Summary and Project Statement 24 2 Wavelets 25 2.1 Wavelet Transform: Theory and Implementation 26 2.1.1 Continuous Wavelet Transform 27 iii TABLE OF CONTENTS iv 2.1.2 Discrete Wavelet Transform 28 2.1.3 Redundant (Over-complete) and Stationary Wavelet Transform 39 2.1.4 Wavelet Packets 40 2.2 Wavelets in Biomedical Applications 40 2.3 Motivation for Wavelet Analysis of the EEG 43 2.4 Summary 44 3 WAV Index: Methodology and Implementation 47 3.1 Experimental Data Sets 47 3.1.1 Raw EEG 48 3.1.2 Processed Data 49 3.2 Feature Extraction Methodology 49 3.2.1 Training Data Sets 50 3.2.2 Methodology: A Mathematical Approach 51 3.2.3 Best Wavelet Selection 53 3.3 Implementation 54 3.3.1 Functional Block Diagram 57 3.3.2 Analysis Algorithm 57 3.4 Preliminary Results 65 3.5 Summary 67 4 Pre-Conditioning of the E E G Signal 69 4.1 Introduction 69 4.2 Ocular Artifacts 73 4.3 Wavelet-Based De-Noising 75 4.3.1 Wavelet Thresholding 76 4.3.2 Wavelet Estimators and the Choice of Wavelet Transform 77 4.3.3 Methodology 77 4.4 Results 78 4.5 Application to the WAV Index 82 4.6 Summary 83 5 Clinical Validation 85 5.1 Arthroscopy/Menisectomy (A/M) Observational Study 85 5.1.1 Protocol 85 TABLE OF CONTENTS v 5.1.2 Results 86 5.1.3 Survey 90 5.2 Electro-Convulsive Therapy (ECT) Study 90 5.2.1 Protocol 94 5.2.2 Results 96 5.3 Summary 96 6 Conclusion 98 6.1 Significance of the Work 98 6.2 Contributions 99 6.3 Future Work 102 Bibliography 104 A Bases and Frames 115 A . l Bases 115 A . 2 Frames 116 A.3 Unconditional Bases 117 B Effect of Amplitude Normalization 118 List of Tables 5.1 Arthroscopy/Menisectomy clinical study: patient's data and titration. Doses are expressed in fig (sufentanil, fentanyl, remifentanil) or mg (propofol) 92 5.2 WAV Index vs. bispectral index: survey of the 22 arthroscopy/menisectomy cases. A positive time shift indicates that the WAV Index leads the BIS. The values inside parentheses indicate the slope of the index at the observed point 93 vi List of Figures 1.1 Suppression of reflex responses to noxious stimuli by general anesthetics (adapted from [1]) 4 1.2 The state of anesthesia: relationship between unconsciousness and analgesia (adapted from [2]). . . . 5 1.3 Various EEG patterns (from clinical data) corresponding to increasing doses of anesthetics: (a) awake; (b) general anesthesia; (c) slow delta waves; and (d) burst suppression 12 1.4 Power spectrum distribution for awake (a) and anesthetized (b) patients 15 1.5 Bicoherence values of signals y(k) and z(k) from Eqs. 1.2 and 1.3 (adapted from [3]) 20 1.6 The Bispectral Index Scale (adapted from [4] and the web site of Aspect Medical Systems, Inc.). . . 21 2.1 The short-time Fourier and wavelet transforms (adapted from [5]): (a) different modulates of an elementary STFT window; (b) tiling of the time-frequency plane by the STFT; (c) different scales of the analyzing wavelet; and (d) tiling of the time-frequency plane by the WT. . 26 2.2 The CWT magnitude of a chirp signal . . . . . . . . . . 28 2.3 Haar MRA: (a) scaling function; (b) wavelet; and (c) two different resolutions (adapted from [6]). . . 32 2.4 Scaling function and wavelet vector spaces 32 2.5 Perfect reconstruction filter bank for the,fast wavelet algorithm (adapted from [7]) 35 2.6 Daubechies 4 filters: (a) low and high pass filter coefficients and their corresponding scaling and wavelet functions; (b) power spectra of scaling (solid) and wavelet (dashed) functions (adapted from [6]) 37 2.7 Three-level DWT (L — 3): (a) two-band filter bank structure; (b) analysis tree; and (c) filter bank frequency responses 38 2.8 Three-level wavelet packet decomposition: (a) analysis tree (approximations and details); (b) filter bank frequency responses 41 2.9 The three-harmonic signal sout with a drifting phase: (a) Simulink synthesis; (b) phase change over time 44 2.10 The power spectrum of the signal sout over time calculated for 4-second epochs 45 2.11 The absolute value of wavelet coefficients of the signal sout over time 45 3.1 Surgical procedure data: (a) raw electroencephalogram; (b) bispectral index 48 vii LIST OF FIGURES viii 3.2 Median Frequency (MEF) and Spectral Edge Frequency (SEF) indexes based on the EEG signal from Fig. 3.1.a 50 3.3 Training data sets: (a) awake; and (b) anesthetized 51 3.4 Best wavelet selection for Daubechies filters: (a) discrimination parameter D; (b) mean of the vari-ability va; (c) mean of the variability v^ ; and (d) variance of the variabilities v a and vw for the band Bdt 55 3.5 Daubechies 14 decomposition - PDFs of the detail coefficients averaged over the awake and anes-thetized training data sets: (a) Bdt; (b) Bd2; (c) Bdz; and (d) Bd4 56 3.6 Functional block diagram of the system for estimating the hypnotic state based on wavelet analysis. . 58 3.7 Pre-processing of the EEG epochs 60 3.8 Wavelet analyzer unit 61 3.9 Comparator unit 63 3.10 Awake and anesthetized templates. '. 64 3.11 WAV Index derived based on the templates in Fig. 3.5.a: (a) unfiltered index; and (b) filtered (30 s) and scaled index (WAV - black; BIS - gray) 66 3.12 WAV Index: comparison with the bispectral index (WAV - black; BIS - gray) 67 3.13 Zoomed detail: WAV Index (WAV - black; BIS - gray) 68 3.14 Zoomed detail: WAV Index (WAV - black; BIS - gray) 68 4.1 WAV Index: the effect of ocular artifacts 70 4.2 EEG corrupted with ocular artifacts 70 4.3 Influence of various artifacts on the WAV and BIS indexes (WAV-black; BIS-gray) 71 4.4 Uncontaminated baseline EEG and various artifacts: (a) uncontaminated baseline EEG; (b) EEG contaminated with slow blink artifact; (c) EEG contaminated with fast blink artifact; (d) EEG con-taminated with vertical eye movement artifact; (e) EEG contaminated with horizontal eye movement artifact; and (f) EEG contaminated with round eye movement artifact 74 4.5 Power spectra of the uncontaminated baseline EEG and various artifacts: (a) averaged uncontami-nated baseline EEG and individual resting periods; (b) EEG contaminated with slow blink artifact; (c) EEG contaminated with feist blink artifact; (d) EEG contaminated with vertical eye movement arti-fact; (e) EEG contaminated with horizontal eye movement artifact; and (f) EEG contaminated with round eye movement artifact 75 4.6 Fast eye blink artifact: (a) contaminated EEG and corrected EEG; (b) stationary wavelet decompo-sition of contaminated EEG; and (c) power spectra of contaminated EEG and corrected EEG. . . . 79 4.7 Horizontal (left column) and round (right column) eye movement artifacts: (a) contaminated EEG and corrected EEG; (b) stationary wavelet decomposition of contaminated EEG; and (c) power spectra of contaminated EEG and corrected EEG 80 4.8 Slow eye blink (left column) and vertical eye movement (right column) artifacts: (a) contaminated EEG and corrected EEG; (b) power spectra of contaminated EEG and corrected EEG 81 LIST OF FIGURES ix 4.9 Head movement artifacts: contaminated EEG and corrected EEG; and (insert) power spectra of contaminated EEG and corrected EEG 81 4.10 WAV Index: the effect of ocular artifacts and de-noising (WAV prior to de-noising - thin black; WAV after de-noising - thick black; BIS - gray) 83 4.11 WAV Index: the effect of ocular artifacts and de-noising (WAV prior to de-noising - thin black; WAV after de-noising - thick black; BIS - gray) 83 5.1 Time course of both the BIS (gray) and WAV (black) indexes - cases #11 and #9 87 5.2 Time course of both the BIS (gray) and WAV (black) indexes - cases #8 and #22 88 5.3 Induction: 16 cases synchronized at the loss of count 89 5.4 Induction: 4 cases with strong LMA reaction (WAV - black; BIS - gray) 90 5.5 Emergence: 4 cases (WAV - black; BIS - gray) 91 5.6 Emergence: 2 cases (WAV - black; BIS - gray) 91 5.7 Emergence: 2 cases (WAV - black; BIS - gray) 91 5.8 ECT-BIS study worksheet (adapted from [8]) 95 5.9 Intra-individual variability at induction with thiopental - 3 ECT study cases. Patient #3: female, 150 mg thiopental. Patient #4: female, 200 mg thiopental. Patient #5: female, 250 mg thiopental. . 96 B. l Four-channel EEG recording corresponding to a anesthetized patient 118 B.2 WAV Index calculated on a 4-channel recording: comparison between normalized and non-normalized epochs 119 List of Acronyms A C R O N Y M DESCRIPTION A / M Arthroscopy/Menisectomy A C L Anterior Cruciate Ligament A D C Analog Digital Converter A l Artificial Intelligence A R Auto-Regressive A S A American Society of Anesthesiologists A S S R Auditory Steady State Response BIS Bispectral Index Scale B S R Burst Suppression Ratio CNS Central Nervous System C R T Cathode Ray Tube C S A Compressed Spectral Array C T Computed Tomography C W T Continuous Wavelet Transform D O A Depth of Anesthesia D S A Density Spectral Array DSP Digital Signal Processing D W T Discrete Wavelet Transform E C G Electrocardiogram E C o G Electrocorticogram E C T Electro-Convulsive Therapy E E G Electroencephalogram E M G Electromyography E O G Electrooculogram X LIST OF ACRONYMS xi ACRONYM DESCRIPTION ERP Evoked Response Potential FDA Food and Drug Administration FIR Finite Impulse Response fMRI Functional MRI F T Fourier Transform F W T Fast Wavelet Transform HRV Heart Rate Variability ICA Independent Component Analysis ICU Intensive Care Unit IIR Infinite Impulse Response IMP Intermodulation Product LCD Liquid Crystal Display LMA Laryngeal Mask Airway LOC Loss Of Count MAC Minimum Alveolar Concentration M E F Median Frequency M L A E R Mid-Latency Auditory Evoked Response MRA Multiresolution Analysis MRI Magnetic Resonance Imaging NMB Neuromuscular Blocking Agent OA Ocular Artifact OAs Ocular Artifacts OCWT Over-Complete Wavelet Transform ONBs Orthonormal Bases OR Operating Room PCA Principal Component Analysis PDA Personal Digital Assistant PDF Probability Density Function P E T Positron Emission Tomography PSP Postsynaptic Potential LIST OF ACRONYMS xii A C R O N Y M DESCRIPTION QEEG Quantitative E E G QMF Quadrature Mirror Filter R E M Rapid Eye Movement RMS Root Mean Square RSA Respiratory Sinus Arrhythmia RWT Redundant Wavelet Transform SCC Semilinear Canonical Correlation SEF Spectral Edge Frequency SEP Somatosensory Evoked Potential STFT Short Time Fourier Transform SWT Stationary Wavelet Transform UBC The University of British Columbia UBCH The University of British Columbia Hospital VHHSC Vancouver Hospital & Health Sciences Centre WPs Wavelet Packets WT Wavelet Transform ZXF Zero Crossing Frequency Acknowledgements I wish to express my deepest gratitude to Prof. Guy A. Dumont and Prof. Mihai Huzmezan, my supervisors from the Department of Electrical and Computer Engineering, University of British Columbia, for their valuable support, advice and guidance throughout the course of my studies here at the University of British Columbia. I am also indebted to my supervisor Prof. Craig R. Ries from the Department of Anesthesia, University of British Columbia, for his strong support and involvement in the WAV Index project. Thanks are extended to my colleagues in the Division of Control Systems in Departments of Pharma-cology & Therapeutics and Electrical & Computer Engineering, who made the period of my research more pleasant and fruitful. Special thanks are to Stephane Bibian, with whom part of the contributions in the development and validation of the WAV Index are shared, and Lou Ann Mendoza for her efforts invested in interfacing the BIS monitor and help with the real-time implementation of the WAV Index. I also cherish many helpful discussions shared with them. Financial support from the Advanced Systems Institute of British Columbia and Universal Dynamics Ltd. is also gratefully acknowledged. Finally, thanks are extended to Aspect Medical Systems, Inc. for their financial support and equipment loan in the framework of the E C T study. xiii / dedicate this work to the memory of my father Aleksandar, a true scholar. Dragom Taji, s ljubavlju ... xiv C h a p t e r 1 Introduction The goal of anesthesia during surgery is to minimize patient discomfort and risks, while providing surgeons with convenient conditions for performing the procedure. The state of anesthesia is nowadays achieved by administering a combination of drugs to render a patient insensitive to the trauma of surgery, and thereby prevent hazardous hemodynamic and movement responses. Failure to deliver adequate anesthesia may lead to traumatizing experiences of being aware and feeling pain during surgery, whereas excess anesthesia may impair critical organ blood flow. In fact, anesthetizing drugs are among the most dangerous drugs presently used. Therefore, it is crucial to monitor the state of anesthesia. Yet, the search for a reliable method to measure the "depth of anesthesia" is difficult, and compared by many with the search for the Philosopher's Stone. In this introductory chapter, we will first give an insight into the concepts of modern anesthesia and its clinical practice. We will then establish grounds for understanding how to monitor the state of anesthe-sia. Various methods for monitoring the anesthetic state will be introduced, including traditional clinical measures and also those based on the analysis of the Electroencephalogram (EEG). Finally, we will present the first clinically accepted commercial monitor, the Bispectral Index Scale® (BIS®) by Aspect Medical Systems, Inc. This monitor utilizes a sophisticated bispectral algorithm for analysis of the E E G . In spite of its limitations, it is currently the most widely used monitor for measuring the effects of anesthetics on the brain. 1.1 Anesthesia versus Depth of Anesthesia 1.1.1 Early Concepts The science of anesthesiology and its every day practice strikingly evolved since the introduction of anesthesia by inhalation of gases and vapors (nitrous oxide in 1844, ether in 1846 and chloroform in 1847) [9]. Yet, despite this remarkable progress, revealing the essence of general anesthesia and understanding how to measure the depth of anesthesia have moved ahead slowly. While the phrase "depth of anesthesia" is widely 1 CHAPTER 1. INTRODUCTION 2 used, it lacks a strong supporting scientific background. Hence, even nowadays this issue is still hotly debated among academics. John Snow's publication in 1847 was the first systematic observation of the physiological effects of diethyl ether administration [9]. He described progressive levels of ether anesthesia (i.e. first occurs confusion, then excitement and unconsciousness, followed by stillness and absence of reaction to surgical stimuli). He noted that the suppression of more intense stimuli required more ether. In the early days, too "deep" a level of anesthesia, which was needed for performing the surgical procedure, had high risks of morbidity and mortality. (Too "light" a level of anesthesia is also undesirable; the patient being aware of surgical procedure bears tremendous risks of long term psychological damage [10].) In order to improve safety, Snow was the first to attempt delivering a known percentage concentration of anesthetic vapor to the patient based on his or her clinical response. It was obvious that individual anesthetic requirements precluded quantitative assessment of depth of anesthesia when only based on the measurement of the drug concentration. In an attempt to determine a more appropriate way to estimate anesthetic depth, Guedel (1937) [11] proposed his scheme for ether anesthesia. He described four different stages, the third being surgical, which was further subdivided into four planes. He associated different combinations of clinical signs (respiratory patterns, hemodynamic and ocular signs, somatic muscle tone, swallowing and vomiting reflexes) with dif-ferent levels of anesthesia (i.e. stages and planes) in a wide variety of patients. In those days, this scheme had great practical utility for ether, chloroform and cyclopropane anesthesia. As the pace of progress quickened, new anesthetic agents and methods were introduced into clinical practice such as the intravenous induction of anesthesia in adults - a faster technique than inhalational induction. Before the 1930s, clinical anesthesia was employed almost exclusively by the administration of one or two inhaled agents, which yielded unconsciousness, muscle relaxation and analgesia. However, the risks of exposing the patient to a rapid intravenous induction was sadly demonstrated at Pearl Harbor on December 9, 1941 [9]. Many patients with blood loss died due to circulatory collapse after the intravenous administration of thiopentone by inexperienced medical staff under hostile conditions. The concept of "balanced anesthesia", proposed in 1926 by J. S. Lundy of the Mayo Clinic, was then revived. The idea was to wisely combine various agents and techniques so that the anesthetic goals were achieved while reducing the amount of each drug and therefore minimizing the cardiovascular risks. However, Guedel's system started to show serious inconsistencies in estimating the depth of anesthesia administered by different combinations of new methods and agents. In 1942, an intravenous muscle paralyzing drug, curare, was introduced to clinical anesthesia [9]. Al-though Neuromuscular Blocking Agents (NMBs) required artificial rather than spontaneous ventilation, they abolished the need for hazardous deep levels induced by inhaled anesthetics, which were necessary in order to achieve sufficient muscle relaxation (e.g. to facilitate intra-abdominal surgery or tracheal intubation). Soon, the use of NMBs in conjunction with controlled ventilation and lower concentrations of inhaled anesthetics became common practice. However, the probability of both under- and overdosing the patient increased since the significant clinical signs of depth of anesthesia (e.g. respiratory patterns and patient movement in response to surgical stimuli) were no longer available. The anesthesiologist was left only with autonomic CHAPTER 1. INTRODUCTION 3 rather than somatic nervous system responses such as increases in blood pressure and heart rate, sweat-ing and lacrimation. Once again, the previous methods established by Guedel and others were powerless regarding the task of maintaining a proper level of general anesthesia. Early attempts to classify different stages of anesthesia have the inherent misconception that general anesthesia is a multi-featured phenomenon caused by a unitary, non-specific underlying mechanism of anes-thetic action [12]. This implied the following: i) one anesthetic could be freely interchanged by another; ii) the anesthetic effect resulting from the combination of different agents was additive; and iii) in order to hold for all general anesthetics with different chemical structures, the mechanism of anesthesia was non-specific regarding the specific anatomical sites within the Central Nervous System (CNS). Their failure motivated many to question and then revise the basic concepts of anesthesia and the means to measure the anesthetic depth. 1.1.2 Restatement of Concepts In 1950, Gray and Rees [13] divided anesthesia into narcosis, analgesia and relaxation (later renamed "triad of anesthesia"). Woodbridge, in 1957 [14], proposed four components of general anesthesia: i) sensory block of afferent nerve impulses; ii) reflex block of the respiratory, cardiovascular or gastrointestinal tract; iii) motor block of efferent nerve impulses; and iv) mental block, sleep or unconsciousness . He recognized that different drugs could be used to achieve each effect. However, he did not tackle the problem of assessing each of these components. Pinsker [15] confirmed these definitions in 1986; he proposed three components of anesthesia: paralysis, unconsciousness and attenuation of stress response . Although these efforts made the concept of anesthesia more lucid, the understanding of depth of anesthesia did not advance significantly. Prys-Roberts [1] greatly enhanced the understanding of depth of anesthesia by giving a meaningful revision of the anesthesia constituents. Based on the premise that "pain is the conscious perception of noxious stimuli", he defined anesthesia as a state in which "as a result of drug-induced unconsciousness, the patient neither perceives nor recalls noxious stimulation". Thus, he considered unconsciousness as the most important component of anesthesia. He further argued that loss of consciousness is a threshold event, and therefore, anesthesia is an all-or-none phenomenon. Hence, there are no degrees of anesthesia, nor variable depths of anesthesia. In addition, analgesia, muscle relaxation and suppression of autonomic activity are not necessarily components of anesthesia, but important supplements that make the act of surgery possible. They are all discrete desirable pharmacological effects and can be achieved independently by specific drugs. Indeed, many general anesthetic agents1 can exert some or all of these effects up to a certain degree. However, the suppression of sensory perception and the production of unconsciousness are the only properties common to all general anesthetics. In addition, powerful analgesics (opioids) modify the perception of pain without producing unconsciousness. They can achieve a loss of consciousness by doses that completely remove the desire to breathe. Thus, analgesia where consciousness is present cannot be a component of anesthesia. Moreover, in the majority of surgical interventions, which are regularly performed without the use of muscle 1 t h e d r u g s t h a t p r o d u c e t h e s t a t e o f g e n e r a l a n e s t h e s i a - m a i n l y i n t r a v e n o u s a n d i n h a l e d a n e s t h e t i c s CHAPTER 1. INTRODUCTION 4 Noxious stimulus Somatic Sensory Motor I Breathing Anesthesia I (Opioids) I Pain Movement . Muscle . relaxants Autonomic Hemodynamic Sudomotor Hormona l Extradural or Spinal (Opioids) Breathing Arterial pressure/ Sweating Stress Response Heart rate Figure 1.1: Suppression of reflex responses to noxious stimuli by general anesthetics (adapted f r o m [1]). relaxants, sufficient muscle relaxation is achieved by increased concentration of general anesthetic agents. Thus, muscle relaxation is neither a component of anesthesia, nor an alternative to inadequate anesthesia. Finally, Prys-Roberts defined anesthesia in terms of drugs producing unconsciousness (i.e. general anes-thetics) and their modification of noxious stimuli. The act of surgery produces noxious stimuli, which further provoke various reflex responses from the patient. These responses need to be obtunded for the benefit of the patient. The state of anesthesia ensures the suppression of the sensory components, and thus the per-ception of pain. Figure 1.1 shows Prys-Roberts' scheme of reflex responses to noxious stimulation, ranked from left to right in the order in which they are suppressed by general anesthetics. This rank can be used as a surrogate for "depth of anesthesia" [11]. Prys-Roberts' conceptualization and definitions are nowadays widely accepted by academics, although some of the important goals of anesthesia are acquiescently taken for the components of anesthesia, in particularly analgesia2. Stanski [11] agrees that this approach is most appropriate for current clinical anesthesia practice that uses potent inhaled and intravenously administered anesthetics, opioids and muscle relaxants, and their combinations. General anesthesia might be viewed as a "pharmacological intervention" [12] that protects the patient against the trauma of surgery and provides a quiet surgical field. This intervention is a sum of separate pharmacological actions, even if produced by a single drug. Depending on the goal of anesthesia, it can include unconsciousness, anxiolysis, amnesia, analgesia, suppression of somatic motor, cardiovascular, breathing and hormonal responses to the surgery, and muscle relaxation. The achievement of one end-point of anesthesia does not give any assurance regarding the status of another end-point. The concept of anesthesia has obviously overgrown its simple unifying definitions. In modern practice, the term "depth of anesthesia" is not relevant for such a heterogeneous phenomenon. There cannot be a unifying clinical measure of depth of anesthesia: what appears relevant for a single drug may not be relevant 2 Since analgesia by Prys-Rober ts implies conscious perception of pain, to avoid confusion, some authors use term antinoci-ception instead of analgesia to denote the obtunding of the reflex responses to noxious s t imul i . CHAPTER 1. INTRODUCTION 5 for another or a combination of drugs. It is impossible to determine the potency of different actions with a single measure: what is relevant for one anesthesia end-point may not be relevant for another. Thus, anesthesia can only be adequate or inadequate [16]. However, the term "depth of anesthesia" can still be relevant for each of the anesthesia components measured separately [17]. Therefore, the search for a reliable index of depth of anesthesia should be transformed into a search for separate indexes of different components of anesthesia. 1.1.3 A Hypothesis of Anesthesia Ultimately, while there are many desirable goals of anesthesia, modern anesthesia doctrine recognizes at least two as principle components of anesthetic state [12, 2, 18, 19]: 1. unconsciousness (including amnesia); and 2. the obtunding of reflex responses3 to a noxious stimulus (analgesia, antinociception). Indeed, in modern anesthesia practice, the brain is first rendered unconscious, by administering seda-tive/hypnotic drugs (inhaled and intravenous anesthetics). These drugs inhibit the consciousness center, which is within the higher centers of the brain. And then, the transmission of noxious stimuli produced during surgery is blocked from reaching the brain and therefore activating the consciousness center. This can be done at the level of the spinal cord (e.g. by opioids at analgesic levels or by inhaled anesthetics at higher concentrations than those required for unconsciousness) or by the action of local anesthetics at the peripheral nerve, nerve root or spinal cord [2]. Figure 1.2 describes the relationship between unconsciousness and analgesia to produce the state of anesthesia. Conclusively, a monitor for each of these components is necessary in order to assess if the state of anesthesia is adequate. Emotional Stimulus Volatile Anesthetics C O N S C I O U S N E S S C E N T E R Local Anesthetics Peripheral Nerve Noxious 1 Stimulation Figure 1.2: The state of anesthesia: relationship between unconsciousness and analgesia (adapted from [2]). 3somatic motor and autonomic (breathing, cardiovascular, sudomotor and hormonal) responses CHAPTER 1. INTRODUCTION 1.2 Monitoring the State of Anesthesia 6 In order to assess if the state of anesthesia is adequate, we need to monitor various components of anesthesia depending on its goals, including at least unconsciousness and analgesia. However, inadequate analgesia will ultimately trigger consciousness. In the early days of anesthesia practice, the major concern for monitoring the state of anesthesia was to prevent overdosing of the patient. Although overdosing does not imply potential fatalities anymore due to the safer anesthetic drugs and techniques and the existence of antagonists, many anesthesiologists agree that patients usually receive more drugs than needed. Today, cerebral oxygenation and hence anesthetic reversibility is of primary concern. A monitor would therefore provide anesthesiologists with a tool for assessing the CNS well-being as well as a more precise titration of each patient, according to his or her specific needs. This would increase the safety of anesthesia, reduce postoperative side effects and optimize drug consumption in terms of cost. Nowadays, an important concern is also the prevention of under-dosing, which increases the probability of patient awareness. 1.2.1 Prevention of Awareness Awareness during surgery may cause severe physiological and psychological damage to the patient. Vic-tims who suffered this unfortunate experience describe it as the worst of their lives. A group of expert anesthesiologists ranked it just below death as a consequence of inadequate anesthesia [20]. There has been much confusion because of the lack of precision in using the term "awareness". Con-sciousness has been defined in terms of being able to respond to verbal commands [19, 21], but many do not distinguish this ability from being able to form consolidated memory with subsequent recall [22]. Explicit memory is the ability to recall specific events or stimuli. Implicit memory is subconscious, but can have great influence on behavior and feelings. Awareness is the state of being aware (i.e. conscious) [21]. In the context of general anesthesia, this refers to either explicit or implicit memory of intraoperative events, while recall is usually used for explicit memory [22]. A patient might be responding to a verbal command during surgery (so called "wakefulness" or "responsiveness") but have no recall of the event [11]. Furthermore, implicit memory may form even without the patient being responsive during the surgery. It may cause acute or chronic psychosis or nightmares even many years later. However, this is still a controversial topic [21]. The incidence of intraoperative recall varies according to the type of surgery performed and the anesthetic technique used. It has been well described for cesarean delivery procedures under general anesthesia. Recall or intraoperative dreaming occurred in 0.9% and 6.1% of 3076 cases, respectively [20]. Estimations for the general surgical population are from 0.1% to 2% [11, 21, 23, 24, 25]. However, for some major trauma surgeries it was present in as much as 43% of patients or 23% during certain type of cardiac surgeries [2, 20]. Recall with pain is the most worrying complication. Jones [23] has estimated that the incidence of recall with pain during routine general anesthesia is 0.01%, based on a review of the literature. However, the incidence of awareness during surgery without recall is probably much higher. Some estimates are that as much as 80% of general surgical patients (if not paralyzed) can follow commands at some point during apparently CHAPTER 1. INTRODUCTION 7 adequate general anesthesia [19]. This, of course, does not imply that they are aware (i.e. have some form of memory). Therefore, prevention of awareness during general anesthesia is one of the important goals of monitoring the adequacy of the state of anesthesia. So far, none of the clinical measures, methods and strategies has proved to be reliable enough in pre-venting awareness during general anesthesia. Although conscious perception does not imply the ability to remember, nor does the lack of it exclude implicit memory, the probability of recall increases with sustained responsiveness [22]. Also the absence of recall is the only objective criterion for unconsciousness [11]. Thus, if there were means to detect awareness by monitoring the state of anesthesia, we foresee that it would be, most probably, by monitoring the unconsciousness. 1.2.2 Can We Monitor Unconsciousness? The current scientific understanding of the mechanisms of consciousness is still quite limited, which makes the understanding of unconsciousness and anesthesia even more limited. Loss of consciousness is a threshold event. It is an all-or-none response to an anesthetic agent. However, it seems that anesthetics cause a graded progression of the underlying mechanism. Increasing anesthetic effects produce the following progression: disinhibition, sedation, both explicit and implicit amnesia, then loss of conscious awareness and surgical immobility [19]. Stanski defined the depth of anesthesia, in context of its components, as a measurement of the effect of an anesthetic drug [11]. The relationship between the drug dose, plasma concentration and response is determined by pharmacokinetic ("what the body does to the drug") and pharmacodynamic ("what the drug does to the body") concepts. Usually, to measure the drug-effect, one applies stimulus and then observes the response. Frequently, this response is quantal (e.g. unconsciousness, movement). However, it represents only a single point in a continuum of the relationship between the dose of a drug and the effect [17] and does not reflect the depth of the inhibitory effect. Unconsciousness, in the context of general anesthesia, consists of amnesia (absence of recall) and hypnosis (drug-induced sleep) [11]. Therefore, in order to quantify unconsciousness, we can try to measure the effect of hypnotic drugs (i.e. general anesthetics), and thus the depth of hypnosis. Hopefully, we should be able to find a parameter, or a set of parameters, that shows graded changes with increasing hypnotic depth. Unfortunately, at the current state of knowledge, it is not likely that we can quantify amnesia. However, if the depth of hypnosis increases, the probability of recall decreases. Thus, by monitoring hypnosis during surgery, we should be able to assess the likelihood of a patient being aware. In the context of the previously stated anesthesia concepts, by anesthetic depth or Depth Of Anesthesia (DOA) we consider the depth of a particular anesthesia component or end-point, in our case hypnosis (if not stated otherwise). By monitoring hypnosis, we indirectly monitor unconsciousness. By measuring the depth of hypnosis (i.e. the depth of drug-induced sleep), we are measuring the strength of the underlying inhibitory mechanism of consciousness. With deeper levels of hypnosis, the probability of a patient being conscious decreases. CHAPTER 1. INTRODUCTION 8 1.2.3 A n Ideal Monitor It is important to define what criteria should an ideal monitor of DOA meet. We have concluded based on [2], [22], [26] and [27] that an ideal monitor (index) should: 1. show graded changes in order to detect the approach of awareness before its onset; 2. be able to closely reflect changing concentrations of anesthetic agents, as well as the level of hypnosis; 3. have 100% sensitivity and specificity (i.e. the average values yielded by the monitor in two distinct states should be uncorrelated, without overlap of corresponding range of values); 4. be applicable to all anesthetics and their combinations; 5. not be affected by non-anesthetic drugs; 6. be applicable to a wide variety of patients; 7. have a stable baseline value, which should be recovered upon the discontinuation of drugs; 8. have a therapeutic window which coincides with the maximal changes induced by drugs; 9. have automatic rejection of artifacts; 10. have a high temporal resolution with real-time response; 11. be easy to use, with no need for calibration; 12. display results as a scale easy to understand and interpret. As far as we are concerned, these are the guidelines for designing a reliable monitor of anesthesia adequacy. However, to meet all the preconditions is a great challenge, currently still unresolved, which leaves room to an original approach. The task of preventing patient's awareness is the most complicated. In order to detect such an event, which occurs with a very low incidence, high sensitivity (low rate of false-negatives) and specificity (low rate of false-positives) of the device are crucial. 1.3 Intraoperative Assessment of Consciousness In a search for a reliable indicator that would ensure an adequate anesthesia and lack of awareness, various clinical anesthesia methods have been utilized, and different technologies have been introduced, with varying degrees of success. In everyday practice, anesthesiologists administer anesthetic drugs to a perceived level of anesthetic depth that is generally defined by responsiveness to noxious stimulation. Patient responsiveness is usually assessed based on various clinical signs. CHAPTER 1. INTRODUCTION 9 1.3.1 Clinical Signs Two physiological signs have been extensively used in clinical practice for assessment of anesthesia adequacy and as possible predictors of awareness: movement and autonomic/hemodynamic response. Movement For a long time, reflex withdrawal in response to painful stimuli has been one of the most useful signs in clinical practice to assess the adequacy of anesthesia. For years, the concept of Minimum Alveolar Concentration (MAC) of inhaled anesthetic that eliminates a purposeful movement response in 50% of subjects to a skin incision has been used as means of titrating adequate levels of anesthesia induced by inhaled agents. Namely, the concentration needed to prevent a motor response is much higher than that required for loss of consciousness (MAC-awake), but lower than that necessary to eliminate an adrenergic response to stimulus (MAC-BAR) [11]. Movement as a clinical sign is not available for procedures employing neuromuscular agents. In this case, the isolated forearm technique, with a tourniquet preventing paralysis of one limb so that the patient can signal awareness, can be employed [28]. However, this method does not readily identify cognate responses from non-cognate ones. It is used mainly as a research tool, but is not a standard practice in clinical anesthesia. Recent findings have shown that movement responses are mediated by spinal cord, and not by the CNS structures subserving consciousness. In 1993, Rampil showed in rats that the ability of general anesthetic to block a somatic response does not change even after total decerebration. In the same year, Antognini and Schwartz demonstrated in goats that the administration of inhaled anesthetics to the spinal cord alone was necessary and sufficient to eliminate a somatic response [19]. Furthermore, opioids, which are extremely effective in intraoperative movement suppression, do not generally induce unconsciousness. Conclusively, movement does not necessarily imply awareness; it merely suggests a need for analgesia. Hemodynamic and Autonomic responses In today's clinical practice, unexplained changes in blood pressure, heart rate, respiratory parameters, muscle tonus, ocular signs, tears and sweating are considered as a potential signs of awareness. The hemodynamic response to specific intraoperative events, such as laryngoscopy, endotracheal intubation and/or skin incision, has been used to assess the anesthetic adequacy. However, these signs are all indirect, non-specific and unreliable. They heavily depend on the drugs used, the surgical procedure and disease involved. Also, there is a large inter-patient variability. It has been shown that the patients who experienced intraoperative awareness cannot be distinguished from normal patients based on these signs alone [29]. Furthermore, powerful hypnotics produce little suppression of autonomic response, while opioids suppress it without producing sleep (unconsciousness). Like movement response, acute hemodynamic changes are mediated by the spinal cord [19]. CHAPTER 1. INTRODUCTION 10 No set of clinical signs has been proved to be totally specific and sensitive. These signs neither neces-sarily correlate with the concentration of anesthetic agent, nor with the adequacy of anesthesia. Although anesthesiologists extensively use them to monitor the adequacy of anesthesia, they do so qualitatively rather than quantitatively. • 1.3.2 Other Approaches A number of other methods for monitoring anesthetic adequacy have been introduced with different levels of success. The Electromyography (EMG) of facial muscles (FACE®, by Patient Comfort Inc.) analyzes electrical signals originating from muscular tension and movements of facial muscles, especially the frontalis muscle (the one raising the eyebrow), which is particularly resistant to neuromuscular blockade. However, this method seems to be useful mainly for the levels of anesthesia corresponding to the induction and the moment of emergence, but it has no correlation to the drug effect, and has great variability [18, 27, 30]. The lower esophageal contractility method is based on the spontaneous rhythmic activity of the esophagus, which decreases with increasing dose of anesthetics, but is not affected by muscle relaxants. This method was abandoned because patient movement easily obstructed it, it was hard to establish and interpret, and was drug dependent [18, 27]. Recently, Heart Rate Variability (HRV) has been used for monitoring the adequacy of anesthesia. Anemon-I® (by Medical Control SA) is promoted as an analgesia monitor. Fathom® (by Amtec Medi-cal Systems) is another HRV monitor based on a correlation between Respiratory Sinus Arrhythmia (RSA) and awareness [31]. RSA utilizes Electrocardiogram (ECG) and respiratory data that are routinely collected during surgery. It is claimed that this monitor enables the anesthesiologist to evaluate and predict a patient's level of response to anesthesia and painful stimulus. However, it is currently being evaluated and there is not much data available upon its clinical utility regarding the assessment of consciousness. Finally, the brain being the end organ of consciousness, to answer the question of appropriate anesthesia monitoring, it is natural to explore the analysis of the electrical activity of the brain (i.e. EEG). In the following sections, we will introduce various methods of monitoring the anesthetic state that utilize the analysis of the E E G , with special emphasis on the currently most widely used monitor, the BIS by Aspect Medical Systems, Inc. 1.4 The E E G as a Measure of the Anesthesia Adequacy The effect site of anesthetic agents is the central nervous system. Furthermore, the brain is the end organ of consciousness and the target organ for all hypnotic drugs. It is thus reasonable to assume that the cerebral activity contains information concerning the adequacy of anesthesia. Hence, we can try to monitor the depth of the hypnotic component of anesthesia by measuring the effect of anesthetic agents on the brain electrical activity. The electroencephalogram is a non-invasive continuous measure of the electrical activity of the brain. It CHAPTER 1. INTRODUCTION 11 directly correlates to cerebral physiology such as cerebral blood flow and metabolism. Over more than 50 years, the E E G has been thoroughly investigated for monitoring the CNS well-being, and in particular the state of anesthesia, which affects cerebral physiology and alters E E G patterns. In this section, the different EEG-based techniques used to measure the anesthetic depth are briefly reviewed. 1.4.1 Q u a n t i t a t i v e A n a l y s i s o f the Spontaneous E E G The E E G is a result of summated excitatory and inhibitory Postsynaptic Potentials (PSPs), which originate within dendrites of thousands of pyramidal neurons in the cerebral cortex. While awake, PSPs are desyn-chronized since the neurons are firing more independently in the creation of conscious human behavior [32]. Thus, the resulting E E G is a random-appearing signal. The degree of PSP synchronicity determines the magnitude and frequency of the resulting signal. Anesthesia and generally the depression of consciousness increase the cortical synchronicity when highly complex activities present in the awake state are replaced by simpler oscillations. The peak-to-peak amplitude of the E E G typically ranges from 2 to 200 ^/V. The effects of anesthetic drugs on the E E G have been known since the early 1940s [33] when neurophysi-ologists observed different patterns within the E E G corresponding to various levels of anesthesia. The E E G generally changes from a low-amplitude, high-frequency signal in an awake patient, to a high-amplitude, low-frequency signal in deeper hypnotic states. With very high doses of anesthetic agents a burst suppres-sion pattern usually occurs (periods of low amplitude activity or even isoelectricity interrupted by short outbreaks of high amplitude activity), and finally all brain activity may disappear, leading to an isoelectric E E G . However, this is a very simplified classification. Faulconer and Bickford [34] published a review cata-loguing different patterns of the E E G associated with the depth of anesthesia for different anesthetic agents. Bickford has also shown that the E E G measured from the scalp is similar and synchronous to the E E G that can be measured in the depth of the human brain. An extensive review published in 1973 by Clark and Rosner [35] concluded that all anesthetics have an effect on the E E G . Furthermore, they suggested that different anesthetic agents have somewhat different actions on the CNS. Figure 1.3 shows some of the typical E E G patterns occurring during anesthesia that are associated with increasing doses of anesthetics. The effect of opioids on the E E G has also been thoroughly investigated. In 1984, Smith et al. [36] concluded that the E E G probably reflects the depth of anesthesia with high dose of narcotics. Note that at these doses opioids usually introduce the loss of consciousness. Later studies by Scott et al. [37, 38] and Billard and Shafer [39] resulted in similar conclusions. The effect of NMBs has also been investigated. Pichlmayr et al. [40] found no visible change of the E E G . Moreover, the E E G is not influenced by the administration of non-anesthetic agents such as anticholinergic or other cardiovascular drugs. All the above results are strong advocates for the use of E E G as a monitor of DOA. The routine monitoring of the "raw" or unprocessed E E G during general anesthesia is very impractical due to the need for continuous observation and interpretation by an expert, the complexity of instrumentation and the interference sensitivity. All studies have shown that there is a clear correlation between the E E G CHAPTER 1. INTRODUCTION 12 Time (s) Time (s) Figure 1.3: Various EEG patterns (from clinical data) corresponding to increasing doses of anesthetics: (a) awake; (b) general anesthesia; (c) slow delta waves; and (d) burst suppression. and the anesthetic depth. They have also shown that patterns generated through anesthesia were different according to the drugs and their combinations. Intermittent surgical stimuli will further change the E E G patterns. Researchers concluded that neuroelectric recordings couldn't provide a simple and uniform measure of the anesthetic depth. Stanski [11] concluded that in order to use parameters derived from the E E G as a measure of the anesthetic depth, one must first study, for each drug, the correlation between the plasma concentration, the E E G effect and the resulting clinical anesthetic state of the patient. He also stresses that the pharmacokinetic and pharmacodynamic model of the drug must be associated with the parameter to accurately represent and predict the anesthetic depth. Also, when using the E E G to determine the anesthetic depth, great care has to be taken into defining a proper therapeutic window. For instance, in their study of opioids Billard and Shafer [39] observed that the E E G could represent adequately the drug effect, but only for clinically deep anesthesia. Therefore, the use of E E G associated with opioids seems to be irrelevant in preoperative and postoperative pain relief in the Intensive Care Unit (ICU), where only a very shallow anesthetic state is provided. The electroencephalogram monitoring mainly correlates with functions that originate within the cortex, such as consciousness or memory, rather than reflex activity mediated by spinal cord, such as movement or acute hemodynamic response [19]. Only if a drug exerts its effect both on the spinal cord and on the CHAPTER 1. INTRODUCTION 13 brain (e.g. inhaled anesthetics), it is likely that there will be a correlation between the E E G effect and the suppression of movement in response to noxious stimuli. Finally, the E E G only provides information about the current level of suppression of spontaneous activity of cortex, but it will not predict the response of the cortex to noxious stimuli. However, inadequate anesthesia generally causes E E G activation, since during unconsciousness it appears that only noxious stimuli can cause appreciable CNS arousal [41]. Over the last few decades, substantial research has been dedicated to the development of appropriate techniques and methods that can analyze and display E E G information in a quantitative format. Quan-titative E E G (QEEG) techniques and parameters significantly improved since their beginnings due to the progress in understanding the physiology of the brain and technological advances. Nowadays, QEEG meth-ods utilize the computerized techniques for signal analysis of the E E G . The motivation for quantitation of the E E G is to reduce the complexity of analyzing the intraoperative E E G and to provide a parameter for the potential automated closed-loop administration of anesthetic drugs in the future. Time Domain Methods Historically, the first clinical applications utilizing E E G analysis employed time domain methods. In fact Bickford, Faulconer and Soltero [42, 43, 44] developed the first automatic control of anesthesia in the early 1950s, based on the observation that the electrical power of the E E G signal correlates with changes in the rate of drug administration (thiopental or ethyl ether). They quantified the energy of the E E G by rectifying and integrating the signal using an analog technique. When the energy reached some threshold value, a pulse was generated and supplied to an infusion mechanism to deliver a single anesthetic bolus. The system was reported to be successfully employed to control the anesthetic depth in 50 patients undergoing laparatomy [42]. This technique is essentially the same as that of digital total power spectrum analysis. Later on, Bellville [45] developed a method using cyclopropane anesthesia. He observed that increasing concentration of cyclopropane was also increasing the E E G amplitude. He thus used the raw amplified signal directly as a feedback quantity. The above mentioned methods were drug dependent. Also, they were not employed in a balanced anesthesia framework, since the state of anesthesia was achieved by the administration of a single drug, and therefore associated with the risks of too deep levels. Note that they were also insensitive to important frequency shifts in the E E G . Another time domain technique, which tried to estimate the "average" frequency of the E E G signal is the so-called Zero Crossing Frequency (ZXF). This technique measures the intervals between subsequent crossings of the zero voltage level by the E E G signal [46]. Although simple to calculate, this method gives erroneous results for the E E G signal that contains waves dissimilar in frequency and amplitude (e.g. small fast waves superimposed to a slow wave with larger amplitude). Aperiodic analysis [47] is based on a revised ZXF concept, which uses two distinct wave detection schemes: the fast and slow wave algorithms. This technique generates a list of all determined waves with their amplitude (based on a maximum value between two sequential minima) and duration (based on the time interval between two sequential minima), which can be displayed in a spectrum-like fashion. Aperiodic analysis is a mathematically simple method; however, it CHAPTER 1. INTRODUCTION 14 is inferior to frequency domain analysis, except for quantifying the isoelectric E E G . In this case, frequency domain derived parameters become unstable [11]). A very important and widely used time domain E E G parameter is the Burst Suppression Ratio (BSR). It quantifies the burst suppression phenomenon, which carries significant information in the case of head trauma or brain ischemia; however, it is also associated with large doses of general anesthetics and cardiovascular overdose. BSR is calculated as the fraction of the epoch length where the E E G amplitude is suppressed below 5 uV for more than 0.5 s. Power Spectrum Analysis With the development of microprocessors and signal processing tools, researchers have focused their attention on Fourier analysis of the E E G . Power spectrum analysis is used to obtain a frequency distribution of the E E G where the phase spectrum is discarded as uninteresting. Hence, any change in the frequency content of the signal can be visualized. Spectral array data from sequential epochs can be displayed together in stack to show the variations in frequency distribution over time. Typical examples are the Compressed Spectral Array (CSA) and Density Spectral Array (DSA). Pichlmayr et al. [48] have published a thorough review of the effect of the different anesthetic agents on the E E G spectral distribution. The E E G spectral distribution is characterized by different modes. It is common practice to distinguish 4 frequency bands: delta band (0.25 Hz - 3.5 Hz), theta band (3.5 Hz -7.5 Hz), alpha band (7.5 Hz - 12.5 Hz) and beta band (12.5 Hz - 32 Hz). It has been considered for a long time that the E E G signals have a maximal frequency of 30 Hz - 40 Hz. However, very recent findings have highlighted the importance of the gamma band (40 Hz - 60 Hz) in assessing the conscious state [49]. For a normal awake patient, the E E G activity is principally concentrated in the delta and alpha bands. With increasing level of anesthetics, the activity of the alpha band tends to reduce, while the low frequency content of the delta band is increased. However, this is a general observation. For some drugs, other phenomena such as high frequency activity and/or burst suppression appear. Hence, it is reasonable to assume that feature extraction techniques are agent-specific. This is obviously not very practical, when considering that modern anesthesia practice involve many different drugs, each of them interacting with each other. To quantify the effect of anesthetics on the E E G , researchers have tried to derive univariate indexes based on Fourier analysis to characterize its spectral distribution. Among the parameters that have been thoroughly investigated, we can enumerate the following five: - Median Frequency (MEF): the M E F is the frequency that splits the power spectrum distribution into two parts of equal power, see Fig. 1.4. This descriptor is advocated by Schwilden and co-workers [50, 51, 52, 53] who used it for closed-loop anesthesia. However, it has proved to be unreliable. - Spectral Edge Frequency (SEF): the SEF is the frequency at which 95% of the E E G power is present. This descriptor has been proposed in 1980 by Rampil et al. [54]. According to them, the SEF is highly repeatable but presents large inter-subject variations. CHAPTER 1. INTRODUCTION 15 5- power = 6.5 (iV 2.5|-2 1.5 1 0.5 Total power=16.6fiV 2 Theta ratio=3.5 SEF=35 H z Trrrrnin M E F = 4 . 5 H z 8-power = 62.8 uV 2 40 r ' 1 30 25 00 10 Total power=90uV 2 Theta ratio=5.2 o | S E F = 9 H z 5 10 M E F = 1 . 5 H z (b) Figure 1.4: Power spectrum distribution for awake (a) and anesthetized (b) patients. - Total Power: almost all anesthetics tend to increase the overall power of the E E G signal, mostly for the hypnotic levels of general anesthesia. In deep anesthesia, some drugs will generate burst suppression patterns where the signal becomes isoelectric. Therefore, the total power technique can lead to an erroneous conclusion in cases of deep anesthesia. - Delta Power: the E E G of anesthetized patients being characterized by slow waves, one of the easiest methods used to quantify the anesthetic state is to measure the power contained in the low frequency band. However, this method strongly depends on the type of drug used. - Theta Power Ratio: ratio of the power in the band from 2 to 7 Hz to the power in the band from 5 to 7 Hz. This parameter significantly increases in the case of a shift towards lower frequencies. A legitimate question is to determine which of these indexes give the best correlation to the level of hypnosis. Gregg et al. [55] developed a statistical method Semilinear Canonical Correlation (SCC) to determine which parameter offers the best correlation to plasma concentration. They introduced two new parameters: a weighted sum of the logarithms of the powers and a weighted sum of the logarithms of the powers expressed as percentage of the total power. The weights have been optimized in order to reduce a cost function and thus maximize the correlation between the observed and the predicted response. As compared to the ad hoc indexes (e.g. SEF, M E F etc.), the proposed SCC offers less variability in the parameters. Therefore, they better characterize each state. However, it is interesting to note that the SEF also gives fair results. This CHAPTER 1. INTRODUCTION 16 study has only been carried out for the opioid alfentanil, and no further research has been done in order to validate the proposed parameters for other drugs. None of the above mentioned parameters clearly demonstrated its relationship to the DOA. A comparative study by Drummond found none of them to be sufficiently reliable as a sole indicator of DOA [56]. A major limitation of these indexes is that their therapeutic window seems to be limited to general anesthesia. In case of light sedation, they do not drift significantly from their baseline value to allow the anesthesiologist to draw accurate conclusion (see [57] and [58]). Also, there is the ambiguity created by burst suppression and isoelectric E E G . Levy [59] argued that, since E E G activity exhibits multimodal patterns, no univariate descriptor based on processed E E G can serve as a consistent and adequate representation of the power spectrum. Levy thus doubted the reliability of univariate descriptors of DOA based on the E E G power spectrum. Some authors propose the use of these parameters in combination with other parameters such as some of the clinical signs (e.g. SEF and mean blood pressure [60]). Recent findings by Sleigh et al. suggested a new index, the median frequency of the first time derivative of the E E G signal - SE50d, as a fairly accurate parameter in separating the awake from the surgically anesthetized states [49]. E E G modelling The idea of modelling the E E G using an auto-regressive technique dates back to the early 1970s. But its application to anesthesia is more recent [61]. Auto-Regressive (AR) models are based on the assumption that the E E G amplitude is the sum of two components: one is correlated to past values while the other is a random perturbation (white noise). The E E G signal can thus be written as: The parameters a,, i = l,2...n, are obtained by minimizing the sum of the squared error between the real signal and its modelled counterpart. As compared to univariate descriptors, auto-regressive modelling generates a set of parameters that can further be correlated to the anesthetic depth. In order to derive a single index from the AR parameters, a neural network can be trained. Sharma et al. [61] have shown that this technique can lead to accurate results for measuring analgesia (i.e. somatic response in dogs). However, the size of the proposed network is rather large (43 inputs, 43 nodes in the primary hidden layer, 6 nodes in the secondary hidden layer and one output). Also, Sharma uses 4 channel E E G modelled by a 10 t h order AR model and associated with the hemodynamic parameters (e.g. heart beat, blood pressure, etc.) and the MAC of the inhaled anesthetic. Other Methods Recently, various sophisticated mathematical techniques have been investigated in order to extract the relevant information from the E E G , which discriminates between different hypnotic levels. Many of them are based on theory of nonlinear dynamics applied to the E E G signal. They are meant to characterize its chaotic n (1.1) i = l CHAPTER 1. INTRODUCTION 17 behavior and to describe the amount of predictability within. "Chaos" is one of the fundamental properties of the E E G signal. It has been assumed that the E E G signal shows more "order" and less "randomness", and therefore more "regularity" and "predictability", with increasing anesthetic concentrations. Characteristics of disorder, such as fractal spectrum, complexity, spectral entropy and particularly approximate entropy, showed correlation to the level of consciousness, see [62, 63, 64, 65]. However, this is a fairly new and insufficiently exploited field. 1.4.2 Evoked Potentials Evoked potentials are the non-random patterns within the E E G that follow a brief acoustic, visual or somatosensory stimulus. The early cortical response to auditory stimuli, so-called Mid-Latency Auditory Evoked Response (MLAER), has been particularly interesting for the assessment of DOA. M L A E R is seen between 40 and 60 ms after stimulation and consists of positive and negative waves. Several studies have recently demonstrated that increasing anesthetic dose increases the latency of these waves and decreases their amplitude. They correlated the latency of specific waves with increasing levels of hypnosis and the loss of consciousness [21, 22, 66, 67, 68, 69]. However, the utility of M L A E R in providing more appropriate dosing and faster recovery has not yet been evaluated. So-called AEP index has been further developed to incorporate both latencies and amplitudes of M L A E R wave components [70]. Apparently, it is an excellent discriminator between awake and unconscious state. However, it does not seem to be a graded index. Another successful approach has been suggested by Plourde et al., see [71, 72, 73, 74, 75]. They studied the 40 Hz Auditory Steady State Response (ASSR). ASSR is a sinusoidal response of the brain following auditory stimuli, which are repeatedly delivered at sufficient speed to produce overlapping of the responses to individual stimuli. They have found that the unconsciousness caused by general anesthetics is correlated with the profound attenuation of ASSR. It seems that indexes derived from evoked potentials are rather independent on the anesthetic drug, but they do not show graded changes. Their major disadvantage is the need of extensive averaging in order to extract the signal of interest, usually very small, from the background E E G . This leads to a high computational complexity and prolonged measurement time. Nevertheless, a monitor (A-line®) has been developed by Alaris Medical that predicts the shape of the MLAER by an autoregressive adaptive filter, and thus reduces the measurement time to less than 15 s [76]. However, this monitor has not been subjected to extensive clinical validation. For further understanding of monitoring the state of anesthesia and utilizing the E E G analysis for this purpose, the interested reader is referred to excellent reviews by Stanski [11, 77] and Rampil [32]. CHAPTER 1. INTRODUCTION 18 1.5 Bispectral Analysis of the E E G for Monitoring the State of Anesthesia 1.5.1 Bispectral Analysis Historically, most processed E E G parameters have utilized the power spectrum analysis. Their limited success can be explained by the fact that i) they do not quantify all of the information available in the E E G ; and it) they make some unsound assumptions about the nature of the E E G signal. Traditional Fourier analysis assumes that the E E G signal is a stationary signal with a normal (gaussian) amplitude distribution. Its frequency constituents are assumed to be uncorrelated, which implies that the E E G arises from a linear process. Although the statistical properties of the E E G signal may be stable over tens of minutes, they can also significantly change within seconds. The E E G signal rarely has a true gaussian probability distribution [32]. Furthermore, the neuronal activity has a strong nonlinear nature, which results in very complex dynamics. Indeed, the power spectrum does not contain all the information carried by the signal. For instance, phase information is discarded. Thus, phase relationships between signal frequency components cannot be tracked by standard spectral analysis. The physiological background of these phase relationships is not clear. However, it is believed that anesthetic agents tend to synchronize the generation of postsynaptic potentials [32]. Therefore, strong phase relationships may imply a small number of independent E E G generators. This may explain why standard spectral parameters are not able to capture small variations of the hypnotic state such as sedation. For light anesthesia, we can expect small phase shifts of some of the frequency components. These changes are not observable by spectral analysis. Hence, standard parameters fail to characterize sedative states. This reduces their therapeutic window. The most accurate and general technique derived for measuring the level of hypnosis is based on bispectral analysis, a signal processing technique that measures the correlation of phase between different frequency components. In addition, bispectral analysis is capable of tracking changes in signals originating from nonlinear as well as linear underlying generating process [3]. Let us consider a simple nonlinear system, whose output y(k) is equal to the square of its input x(k). Let the input signal be the sum of two sinusoids of frequencies / i and /2, and random, independent phase angles Q\ and #2, respectively; that is, x(k) = cos (f\k + 9\) +cos (/2A; + 82). Then, the output of the system is as following: y(fc) = l+cos[(/i+/2)fc+(0i+02)] + c ^ (1-2) The resulting frequency components 2/i, 2/2, /1+/2 and / 1 — / 2 , so-called Intermodulation Products (IMPs), are obviously dependent on the fundamental frequencies f\ and / 2 (i.e. they are "phase-coupled"). Let us now observe the following signal: z(k) = 1 + cos (fak + 6a) + cos (fbk + 0b) + | cos (/cfc + 6C) + | cos [fdk + 6d), (1.3) CHAPTER 1. INTRODUCTION 19 where fa = /1+/2, /b = / 1 - / 2 , / c = 2/i and / d = 2/ 2, and 0a, t?6, t?c and 0d are random and independent. Although signals y(k) and z(k) have a complete different phase nature, their power spectra are identical. Thus, using this simple example it is obvious that phase coupling, which is typical for nonlinear systems, cannot be tracked by the standard power spectral analysis. The bispectrum of the signal x(k) is computed for all possible frequency triplets fa, fa and fi+fa in the following manner: B(fa,fa) N i = l (1.4) where Xi(f) is the Fourier transform of ith signal epoch of a total of N epochs, and X*(f) is its conjugate. This method is generally used to measure the deviation of a given signal from gaussianity and to detect the presence of quadratic nonlinear properties and phase coupling. Unfortunately, Eq. 1.4 has high complexity. Also, in order to have a reliable measure of phase coupling, it is necessary to examine a large number of epochs. This makes the calculation of the bispectrum quite impractical. However, this problem can be partially alleviated using specific, modified techniques. The amplitude of the bispectrum is influenced by the amplitude of the signal as well as the degree of phase coupling. To quantify the degree of quadratic phase coupling, it is common practice to normalize the bispectrum. The normalized index is called the bicoherence index and is expressed as: where RTP(f\, fa) is the average real triple product, which is an equivalent of the bispectrum of a signal with all the phase angles of the components being equal to zero (maximum phase coupling). The RTP(f\, fa) is calculated as follows: N RTP(fa,fa) = £ |Xi( / i ) | 2 - | X i ( / 2 ) | 2 - \X*\2. (1.6) The bicoherence index is a pure measure of the degree of phase coupling, since both the B(fi,fa) and the square root of the RTP(fa, fa) are equally affected by the amplitude of signal. Figure 1.5 illustrates the bicoherence values for signals y(k) and z(k) from Eqs. 1.2 and 1.3, where fa = 7 Hz and fa — 3 Hz. The bispectrum is non-zero at frequency pairs (6 Hz, 4 Hz) and (10 Hz, 4 Hz), indicating that the components at 4 Hz, 6 Hz, 10 Hz, and 14 Hz are IMPs generated by phase coupling. Generally, a strong phase coupling within the triplet of frequencies fa, fa and fa+fa (i.e. a high bicoherence value at (fa, fa)) means that the components at fa and fa may have a common generator, or that there is some nonlinear interaction between them producing a new dependent component at modulation frequency /1 + / 2 . 1.5.2 Bispectral Index Scale (BIS) Monitor Ning and Bronzino [78] were the first to apply bispectral analysis to E E G in order to characterize sleep patterns in rats. They realized that there was a strong coupling between the frequencies of 6 and 8 Hz during the Rapid Eye Movement (REM) phase of sleep. Sleep patterns being close to patterns obtained CHAPTER 1. INTRODUCTION 20 Figure 1.5: Bicoherence values of signals y(k) and z(k) from Eqs. 1.2 and 1.3 (adapted from [3]). during anesthesia, Ning et al. assumed that this technique might lead to interesting results in monitoring the depth of anesthesia. Their assumption was validated in 1991, when Kearse et al. [79] reported that the bispectral index was more accurate than the spectral edge frequency for anesthesia induced by opioids (alfentanil and sufentanil). These findings were confirmed by Sebel et al. [80] with isoflurane and later on with propofol [81]. A study by Muthuswamy and Roy in 1993 on anesthetized dogs (halothane) also led to similar conclusions [82]. However, they also noticed that the accuracy of the analysis could probably be improved by using bicoherence indexes in association with other parameters such as the MEF, SEF and hemodynamic parameters. Probably the most interesting result was obtained by a research team from Aspect Medical Systems, Inc. Initially, they derived no less than 33 bispectral variables plus a burst suppression ratio that were combined by stepwise discriminant analysis to derive a single index to predict somatic response to surgical incision [83]. The database included data collected during standard surgical procedures from 160 patients. In comparison with other methods (SEF, MEF, etc.), the accuracy of the index was significantly higher. A major advantage of bispectral analysis over more conventional techniques is that it has a wider therapeutic window. Indeed, bispectral analysis can differentiate between different levels of sedation [58]. Based on these findings, Aspect Medical Systems, Inc. further developed a monitor of hypnosis/sedation: the BIS monitor. During the development, this monitor underwent several revisions. A clinical endpoint of a somatic response to surgical incision (which is basically an analgesia endpoint) was substituted with hypnosis/awareness. The BIS monitor was the first monitor of anesthetic effect on the brain approved by U.S. Food and Drug Administration (FDA) in 1996, and for a long time it was the only one. Recently the PSA 4000® monitor by Physiometrix has also been approved. Like the BIS, it is based on a multivariate statistical analysis of various spectral and bispectral discriminants derived from retrospective data. However, there are still not many clinical validation studies published verifying its utility. CHAPTER 1. INTRODUCTION 21 The BIS integrates various EEG-derived bispectral, power spectra and time domain (such as the BSR) parameters into a single value. A proprietary weighted sum of these subparameters was derived empiri-cally from a large database containing clinical data. Approximately 1,500 anesthetic administrations were collected, consisting of 5,000 hours of E E G recordings with corresponding clinical endpoints, and under a va-riety of anesthetic regiments. After artifact detection/removal and subparameter calculations, multivariate statistical data analysis techniques were employed to synthesize a final index that best correlates to clinical endpoints of sedation and hypnosis, and yet is independent of the specific agents used. The BIS monitor is a small compact medical device, easy to use and interpret. A single channel E E G signal is recorded through a simple self-prepping sensor placed on the patient's forehead and connected to the monitor by a single cable. The BIS-sensor uses three or four electrodes, which comprise plastic sponges embedded in conductive gel and placed on an adhesive strip incorporating a printed circuit. Thus, the monitor can automatically check if the electrodes have low-impedance contacts with the scalp (less than 5 kf2). This is readily achieved without previous skin preparation. The raw E E G signal is processed by sophisticated algorithms utilizing bispectral analysis performed in real time. An index representing the level of hypnosis is calculated and displayed on the monitor screen. A value of 100 represents the awake state. With increasing concentration of anesthetics, the index decreases. General anesthesia is obtained for an index between 60 and 40 (see Fig. 1.6). Lower values represent a deep hypnotic state, while values between 90 and 60 generally represent the state of sedation (i.e. light anesthesia). The guidelines for the interpretation of BIS values are illustrated in Fig. 1.6. BIS 100 Light Hypnotic <j State Moderate Hypnotic <J State Deep Hypnotic ^ State Awake Light/Moderate Sedation Deep Sedation General Anesthesia Deep Hypnosis Increasing Burst Suppression Memory Intact Low Probability of Explicit Recall how Probability of Consciousness • Memory Function Lost K | — Q — • Isoelectric EEG Figure 1.6: The Bispectral Index Scale (adapted from [ 4 ] and the web site of Aspect Medical Systems, Inc.). CHAPTER 1. INTRODUCTION 22 Numerous studies have demonstrated the higher accuracy of the BIS in assessing the anesthetic effect on individual patients. According to Aspect Medical Systems, Inc., the use of the BIS monitor, besides helping to prevent awareness during surgery, can also reduce the amount of drugs provided to the patient, which in turn allows a faster wake-up, better recovery and earlier discharge from the ICU. Over the past several years, more than 10,000 cases were included in independent clinical studies that investigated the usefulness of the BIS monitor [84]. The bispectral index has proven to be an acceptably reli-able measure of hypnosis during sedation and general anesthesia with most anesthetic agents including their combinations [16, 58, 85, 86, 87, 88, 89, 90, 91]. Hypnotic agents such as propofol, midazolam or thiopental have a strong depressant effect on the bispectral index; inhaled anesthetics (e.g. isoflurane, sevoflurane, desflurane) depress the BIS to a somewhat lesser degree, whereas opioids have almost no influence on the BIS at clinically relevant concentrations. Finally, nitrous oxide and ketamine appear to have paradoxical effects on the BIS. However, nitrous oxide does not reduce the accuracy of the BIS, but its use lowers the probability of awareness at any given BIS value [92]. Thus, the BIS monitor is most useful during high hypnotic-low opioid anesthesia. BIS guided anesthesia in over 2000 patients has been associated with reduced anesthetic consumption, earlier awakening and faster recovery [92, 93, 94, 95, 96, 97, 98, 99]. BIS is still better than other E E G measures in predicting a somatic response to noxious stimuli; however, its correlation to hypnosis performs significantly better [16, 100] and [101]). Some studies showed that the bispectral index not only predicts a loss of consciousness, but further separates between loss of consciousness and recall [89, 91, 102]. The BIS value of 58 was associated with unconsciousness for propofol and thiopental [103]. For a sufficiently low BIS value (e.g. 40), one can be almost certain that the patient is neither responsive to command nor can have a recall [91, 103, 104, 105]). As the bispectral index rises above 70, the probability of both responsiveness and awareness rises. However, the BIS monitor provides useful trend information in individual patients. However, we have to be careful when adopting specific BIS values as thresholds for various clinically relevant events. It is still uncertain if those thresholds are robust enough for different patients as well as various anesthetic regimens. It has also been shown that the BIS is able to indicate the awakening from anesthesia [94, 106] although some argue that it is not totally reliable when the patient is waking [62]. Also, few cases have been reported where the BIS attributed very low values to patients that could still respond to verbal command (see a review of these cases in [107]). On the other hand, BIS monitoring has been used on more than 2.75 million patients worldwide [84] and so far a total of 63 cases of possible awareness have been reported to Aspect. However, 44% of these were associated with a BIS value greater than 65; in another 15% BIS was not used at the time when awareness occurred, and the remaining cases were also rather inconclusive [108]. A large multi-center clinical study, which is investigating if the BIS monitoring may prevent intraoperative recall, is in progress. It is expected to be completed within this year [19]. Many other applications have been explored. Recent study in 38 psychotically depressed patients revealed a positive correlation between the BIS value immediately before electrical shock and the duration of seizure during electroconvulsive therapy [109]. Positron Emission Tomography (PET) confirmed that the BIS can CHAPTER 1. INTRODUCTION 23 measure the decrease of glucose metabolic rate in the brain associated with the increasing depth of anesthesia [110] 4 . It was also demonstrated that the BIS is consistent across a wide range of age groups [111] and that it reflects the increasing sensitivity of patients to general anesthetics with advancing age [112]. In children [113, 114], it was shown that the bispectral index correlates with sevoflurane concentration, although it was not possible to test its correlation to clinical endpoints. Evaluating the BIS in the ICU, where oversedation is common [115], is an important area of research. Finally, it was reported that the bispectral index can track changes in physiological sleep [116]. Currently, the BIS is the most optimized monitor for assessing the anesthetic drug effect. However, this very useful clinical device has a number of limitations. The bispectral index has been empirically derived from a database that included most of anesthetic agents and their combinations as well as a wide variety of patients. However, for a new anesthetic agent or new patient population, it is necessary to readjust the BIS monitor according to new gathered data. Also, it cannot be used for some anesthetic agents like ketamine, or it has to be used with discretion for nitrous oxide. Consequently, the BIS is somewhat drug-dependent. Also, the BIS monitor can produce erroneous results for patients with neurological disease or under certain physiological conditions (e.g. subdural hematoma, cerebral ischemia, hypothermia), or for those taking psychoactive medications. BIS values were extremely low (BIS = 40 for awake baseline) in one subject who had a genetically determined low-voltage E E G [117]. Although the BIS monitor implements sophisticated algorithms for artifact rejection, false index values may still be produced. Significant E M G activity, particularly excessive facial muscle tone, will increase the value of the bispectral index [118]. Also, from our experience, the BIS can be fooled by intensive ocular activity. Conclusively, reliance on the BIS alone for intraoperative anesthetic management is not recommended. The BIS values should always be interpreted taking into account the accompanying clinical circumstances and anesthetic regimen. Finally, Hagihira et al. thoroughly investigated the practical issues in bispectral analysis of electroen-cephalographic signals [119]. They concluded that at least 3 min are needed for reliable calculations of bicoherence values. There has been a lot of unknown regarding how the BIS is calculated. Also, various versions of the BIS monitor and their specific algorithms contribute to the confusion. Rampil [32] stated that the BIS A-1050 utilizes 2-second epochs and a running average filter over the latest 60 s in order to compute the final index. The monitors A-1050 and A-2000, which we had a chance to use, have the option of 15 s or 30 s smoothing. Thus, according to Hagihira et al., the length of the E E G signal used to produce a single BIS value might not be long enough to give a reliable measurement. However, the BIS value incorporates many other variables beside bispectral ones, which increases accuracy, but adds to a very high computational complexity. A recent study on 19 patients has shown that the actual time delay of the BIS is about 60 s behind the clinical situation; but again, they have used the old version of monitor, BIS A-1000 [120]. Our observations are (using the BIS A-1050 monitor, but also the A-2000 and XP monitors) that this delay is closer to 15 - 30 s. 4 Glucose is a sole sugar used by neurons. CHAPTER 1. INTRODUCTION 1.6 Summary and Project Statement 24 Despite the quantity of work that has been carried out in the search for a monitor of anesthesia adequacy, there does not seem to be a reliable technique that meets the expectations of anesthesiologists. Until recently, the BIS monitor was the only device that had been commercialized. Despite many studies that have confirmed its clinical utility, many practitioners are still reticent of using it since experience has shown that the index gives lots of room for interpretation and, in some cases, erroneous measurements. Furthermore, we can consider the automated control of anesthetic administration as the ultimate goal, which will hopefully improve quality in patient care, increase safety and reduce associated costs. Recently, some studies have demonstrated the potential of closed-loop anesthesia in improving patient outcomes [98, 121]. The pioneering work in this field has been carried out by Bickford and Faulconer back in 1950 [42]. Ever since, closed-loop control of anesthesia has been a challenge waiting for the development of a reliable and robust monitor. New designs of automated anesthesia have been recently proposed using the BIS [121, 122] or M L A E R [123, 124] as a feedback quantity. However, a clinically acceptable solution necessitates significant improvements in sensor technology and artifact removal techniques [125]. The E E G is typically a non-stationary signal. One can thus wonder if standard spectral analysis can capture all the information contained in the signal. Bispectral analysis gives satisfactory results; however, in our opinion, it produces a significant time delay. Wavelets analysis, a powerful emerging signal processing technique, has generated a lot of enthusiasm in the biomedical field. The great majority of signals generated through biological processes tend to be similar to the shapes of wavelet basis functions [126]. Thus, researchers came up with the idea of using wavelets to detect and measure specific patterns, which can be associated to illness and otherwise hidden in a non-stationary signal. The objective of the present work is to analyze the E E G using wavelet transform and feature extraction techniques, and then derive an index correlated to the hypnotic state of the patient. To do so, we use clinical data (raw E E G and the bispectral index) collected during surgical procedures. The main goal of this project is to process the raw E E G using wavelet transform in order to derive an index correlated to the BIS, and hopefully improve upon it. For control, we also want an index that we understand completely and that is transparent to us. In addition, the delay of the bispectral index is a disadvantage for closed-loop control. Wavelet theory will be briefly disclosed in the following chapter. The experimental data sets and data acquisition are fully commented and illustrated in the third chapter. That chapter also includes a math-ematical approach to the feature extraction methodology that we used. Preliminary results obtained by the analysis based on the stationary (redundant) wavelet transform are also included. The fourth chapter is dedicated to an artifact removal technique also based on the stationary wavelet transform. Finally, for validation purposes, the fifth chapter presents more clinical results of the derived index in comparison to the BIS. We conclude this thesis by proposing the work to be carried out in the future in order to move towards a complete solution. Chapter 2 Wavelets Wavelet analysis is a new and still evolving discipline that offers powerful tools for signal analysis. The theory of wavelets was developed by Y. Meyer, I. Daubechies, S. Mallat and other researchers in applied mathematics in the late 1980s. Over the last decade, wavelets overgrown their theoretical framework and became a practical method for solving problems in many fields such as signal processing, numerical analysis and mathematical modeling. Wavelets are used as a basis for a signal or function expansion in the similar way Fourier analysis uses sinusoids. They differ from the traditional Fourier analysis as they offer both time and frequency localization of the information. Many signals of interest comprise non-stationary or transitory features for which Fourier analysis is not always suited. The reason for which is that it focuses on the frequency content of the signal and discards all the time information. The Short Time Fourier Transform (STFT) provides both frequency and time information by using a fixed size of the time window to perform the Fourier Transform (FT). Therefore, for a given STFT analysis we have to make a compromise between the time and frequency resolution (i.e the shorter the time window, the worse the frequency resolution, and vice versa). Wavelet analysis alleviates this problem by introducing a variable window size of its analysis functions. Thus, one of the major advantages of the wavelets is that they can perform local analysis. In particular, they can tradeoff the time resolution for the frequency resolution, and vice versa. This property makes them especially well suited for the analysis of non-stationary signals. Most signals with physiological origin (e.g. E E G , E C G and EMG) have parameters that significantly vary over time and across recording sites, hence exhibiting non-stationary behavior. They are combinations of transient events that are well localized temporally or spatially, and more diffuse rhythms such as a background texture. Both types of signals carry important clinical information and can be handled by wavelets since they capture both large- and small-scale structures in a single efficient analysis protocol. This may explain the success of wavelets in the area of biomedical signal processing. This chapter first presents theoretical basics of wavelets, followed by a brief overview of wavelets in biomedical applications. Special emphasis will be given on E E G analysis. Finally, we will present the rationale behind the analysis of the E E G using wavelets for estimating the anesthetic depth. Before reading 25 CHAPTER 2. WAVELETS 26 h h h t hk h t (c) (d) Figure 2.1: The short-time Fourier and wavelet transforms (adapted from [5]): (a) different modulates of an elemen-tary STFT window; (b) tiling of the time-frequency plane by the STFT; (c) different scales of the analyzing wavelet; and (d) tiling of the time-frequency plane by the WT. the following sections, we refer the reader to the Appendix A, which gives basic overview as well as some definitions and mathematical preliminaries on signal expansion, bases and frames. 2.1 Wavelet Transform: Theory and Implementation A wavelet means a "small wave". A time-domain wavelet is an oscillating amplitude function of time with a small concentrated burst of finite energy [6]. Wavelets are relatively localized in both time and frequency. They have large fluctuating amplitudes during a restricted period of time, and very low or zero amplitude elsewhere. Furthermore, they are band-limited to a relatively restricted range of frequencies. A single wavelet function ^(i) (so-called analyzing or mother wavelet) generates a wavelet family of analysis functions {ijja,b(t)} by scaling (stretching and shrinking) and translating (moving along the time axis) itself over a continuum of scaling values a and translation values b. Figure 2.1.c illustrates several members of one wavelet family at different scales and for the same translation in time. CHAPTER 2. WAVELETS 27 2.1.1 Continuous Wavelet Transform The Continuous Wavelet Transform (CWT) decomposes any finite energy signal f(t) (i.e. f(t) G L2(Tl)) into an infinite set of wavelet coefficients {ca)b} by correlating a given signal with wavelets ipa,b{t) at infinite ranges of scales a and time translations b by computing the following integral: /oo ra,b(t)-f(t)dt = ^atb(t),f(t)), (2.i) -oo where il>Zib(t) is the complex conjugate of V^bM = .7^ (^ r)> a € ft+ and b e H. Note that c a h is an inner product over C2(TZ)The normalization factor ensures equal energies of each wavelet ipa,b(t) with the mother wavelet tpit) (i.e. HVvMH = IIV'WID-Thus, a wavelet coefficient ca,b is a projection of the signal f(t) onto the analyzing function tp(t) at a given scale a and translation in time b. The CWT actually performs correlation analysis, therefore, its output tends to be maximum at those scales and time locations where the signal most resembles the analysis template ipa^b(t). (This is the basic principle for a matched filter, which is the optimum detector of a deterministic signal in the presence of additive noise.) Further, each template ipa,b(t) 1S predominantly located in a certain region of the time-frequency plane with a central frequency inversely proportional to scale a, and a time location centered around b. The scale a is proportional to the duration of the wavelet. For small a (a < 1), ipa,b(t) will be short and of high frequency. Conversely, for a large a (a > 1), ipa,b(t) be long and of low frequency (see Fig. 2.1.c and d). Thus, the CWT acts as a mathematical microscope, zooming in on the small-scale (high-frequency) structures to detect closely spaced events in time, or zooming out on the large-scale (low-frequency) structures to observe global patterns [6]. Unlike the wavelet analysis, for example, the S T F T 1 uses a fixed window size for all scales and therefore can not capture both the small- and large-scale signal structures. Figure 2.1 illustrates the tiling of the time-frequency plane for the STFT and the Wavelet Transform (WT) 2 . The invertibility of the CWT places some restrictions on functions available for use as an analyzing wavelet. In order for such an inverse transform to exist, the so-called admissibility conditions must be satisfied. However, it turns out that these conditions are not too restrictive and that any zero-mean and absolutely integrable function may be an analyzing wavelet (i.e. ip(t) dt = \I>(0) = 0, where ^f(to) is the F T of ip(t)). Thus, we can achieve perfect reconstruction of the original signal / (£) by the inverse CWT: f(t) = Cp rr C a ' b ' t a ' b ( t ) dadb, (2.2) JO J-oa i where C^ is a constant that depends on properties of the mother wavelet, for more details see [5]. 1 T o compute the S T F T , a signal is mult ipl ied by a window function w(t — T) and then the usual F T is performed, that is /oo w*(t-r)f(t) e~iwtdt= (tu(t-T)e-*wt, /(«)). -oo Thus , the S T F T measures the similari ty between f(t) and the analysis functions w(t — r ) e - 1 " t , which are the translates and modulates of an elementary window w(t). Each analysis function has the same time and frequency resolution and a different location in the time-frequency plane. 2 T h e figure shows the Discrete Wavelet Transform ( D W T ) defined in Section 2.1.2, which is the discretized C W T for the dyadic sampling grid (i.e. a = 2J and b = k 2^ , (j, k) e Z2). CHAPTER 2. WAVELETS 28 1 0.5 0 -0 .5 Analyzed signal (length 1000) Al l i i i l l l l l i lM H I 100 200 300 400 500 600 Coefficients Ca,b 700 800 900 1000 Time (samples) " S c a l e o f c o l o r s f r o m M I N t 0 M A X Figure 2.2: The CWT magnitude of a chirp signal. Figure 2.2 shows the C W T 3 of the chirp signal in the form of a time-scale plot, which illustrates graphically the energy distribution in the analyzed signal as a function of both time and scale and is known as scalogram. Since scale is close related to frequency, the scalogram is similar to the time-frequency plot (i.e. spectogram). Note how the change in scale (i.e. frequency) over time is well depicted. The C W T measures analyzed signals with infinite resolution. However, it is highly and unnecessarily redundant and inefficient, and thus is associated with practical limitations. For a signal of the length N characterized by its sample values {/(fc)}fc=o,...,/v-i, the computation of a single coefficient ca^ from Eq. 2.1 (i.e. for a fixed value of a and b) typically requires O(N) operations. However, the number of coefficients that need to be computed is extremely large since the CWT is a two-dimensional function defined continuously over the time-scale plane. Therefore, it cannot be computed using finite precision discrete machines and a more efficient and computationally simpler wavelet analysis is desirable. However, approximations to the CWT can be made to almost arbitrary precision through its samplings. In general, discrete WTs are obtained by sampling the CWT on a two-dimensional grid (a,j,bj^)- In order for the resulting transform to be invertible, it is clear that the choice of time-scale sampling cannot be made arbitrarily. In fact it is required that the underlying family of wavelets form a frame for the space of interest. Discretizations of this sort generally lead to over-complete (redundant) WTs, but a non-redundant transform exists as well. Therefore, the distinction between different types of discrete WTs depends on the way in which the scale and translation parameters from Eq. 2.1 are discretized. 2 . 1 . 2 Discrete Wavelet Transform There is an infinity of possible discretizations of the CWT. It is often convenient to use dyadic (power of two) scales aj = 2j and translations bj>k - k2j, for j, k G Z. This choice of dyadic sampling grid generates 3 A c t u a l l y , a dense sampling of the C W T is shown. CHAPTER 2. WAVELETS 29 the spanning set {ipj,k}- Thus, under appropriate conditions on the function ip(t) for the inverse transform to exists, it is possible to represent a signal /(£) G L2(TZ) through its wavelet expansion /(*)= E dj(fc)-^ -,fe(t), (2-3) j,kez where V'j.fcW = 2~iijj(2~H — k), and factor 2~2 preserves the norm of i/>jtk(t). The wavelet coefficients dj(k) are obtained by the following inner product: dj(k) = (f(t),$Jtk(t)), (2.4) where function tp(t) is the dual wavelet of il>(t)4, and ifij^it) = 2~iip(2~H — k). If the sets{-0j1fc(i)} and {V>j,fe(t)}, j, k & Z, are unconditional bases of L2(TZ) (not necessarily orthogonal), then we say that Eq. 2.4 is a Dyadic Wavelet Transform. This transform is non-redundant, and, therefore, the Parseval's equality holds (i.e. the norm or energy of a given signal can be exactly partitioned in the transformed domain). The tiling of the time-frequency plane for the discrete dyadic WTs is depicted in Fig. 2.1.d. For a signal represented by N samples, the number of coefficients that need to be computed is N. The term Discrete Wavelet Transform (DWT) is commonly reserved for the dyadic W T with orthonormal or biorthogonal wavelet bases. In the case of an orthonormal wavelet basis, ip(t) and -*/>(£) are identical, which significantly simplifies the calculations. Although Orthonormal Bases (ONBs) are convenient, biorthogonal bases offer greater flexibility for some problems, at the expense of increased complexity. The focus of interest in signal processing is the idea of generating wavelet families that form ONBs, and, especially, compactly supported5 ONBs. Wavelet ONBs of compact support allow for DWT implementa-tion with Finite Impulse Response (FIR) filters utilized by an efficient algorithm called the Fast Wavelet Transform (FWT). The FWT has an effective 0(1) complexity per coefficient (i.e. 0(N) for a signal of N samples)6, and, therefore, can be easily implemented on a digital platform. Multiresolution Analysis A major breakthrough of wavelet-based signal processing was due to the discovery of orthonormal wavelet bases [127], and, especially, compactly supported ONBs [128, 129]. Multiresolution Analysis (MRA), pio-neered by Mallatt [127, 130] and Meyer, plays a crucial role in both the theory and the practical imple-mentations of wavelets since it provides the mathematical framework and mechanism for the construction of wavelet ONBs. Discrete WTs associated with MRA-type wavelet bases can be implemented using extremely efficient algorithms. The main idea behind the MRA concept is the specification of the basic resolution structure, which is a sequence of nested subspaces Voo C . . . V 2 C Vi C Vo C V_i C V _ 2 • • • C V_oo, (2.5) 4 N o t e that ip(t) and i/i(t) in Eqs. 2.3 and 2.4 are interchangeable since they are dual to each other. 5 T h i s means that the analyzing wavelet (and hence, the entire family) is compact ly supported (i.e. has a finite durat ion). 6 N o t e that, for example, the Fast Fourier Transform ( F F T ) has the complexity of 0(N logN). CHAPTER 2. WAVELETS 30 where V Q Q = {0} and V _ o o = L^iJVj. Thus, subspaces Vj are uniformly better and better approximations of L2CR) as the index j decreases. We say that the subspace Vj is a jth level approximation of the space L2(K). In addition, the increase in resolution that is gained by stepping up to the higher resolution space Vj_i from Vj should be the same for all j. In order to obtain an uniform increase in resolution across the subspaces, we require each space Vj to be a scaled version of the central space Vo- Precisely, we require for each subspace Vj to have an ONB that can be obtained by scaling the ONB of Vo by a power of 2. Thus, we obtain subspaces {V,-} with the uniform increase in resolution (by a factor of 2) across them. This so-called scale invariance requirement can be formally denoted as: fj(t) G V,- /,_!(*) = fj(2t) G (2.6) or Sj(t) G V, 4 = ^ f0(t) = fj(2H) G Vo, (2.7) where fj(t) is the orthogonal projection7 of the given function /(£) onto Vj and it approximates /(£) to the jth level resolution. The approximation fj(t) approaches the projected function as j decreases. Thus the space that contains high resolution functions will contain those of lower resolution as well. In order to complete the requirements of the MRA, there must exist a scaling Sanction <p(t) G L2(TZ), whose integer translations in time form an ONB for the subspace VQ C L2(1Z), and, therefore Vo = Span{<pk{t)} (2-8) fc where <pk(t) = </>(£ — k), k G Z. In general, we can increase the size of the subspace spanned by changing the time scale of the scaling function. Thus, a two-dimensional family of functions {<Pj,k{t)} is generated by scaling and translating the basic scaling function <p(t), that is <Piik{t) = 2-t<p{2-H-k). (2.9) The span of {<Pj,k{t)~\ o v e r k is Vi = Span{(fJ!k(t)} = Span{(p(2-H - k)} (2.10) fc ' fc for all k G Z. This means that every /(£) G Vj can be written as the following weighted sum: S{t) = YJCj(k)-^k{t). (2.11) fc Thus, for smaller scales when j < 0, we obtain higher resolution subspaces since <fjtk(t) becomes narrower and is translated in smaller steps k2^ in time. For j > 0, <Pj,k(t) is wider and translations are larger time steps, so these scaling functions can only represent coarse (i.e. low resolution) information. 7 If the set {vi,V2, •••} is an orthonormal basis of Vj C L2(Tl) then the orthogonal projection of /(t) £ L2(U) onto Vj is given by CHAPTER 2. WAVELETS 31 The nesting of the spans of {ip(2~H — k)}, as defined by Eq. 2.10 and for all j G Z, is achieved by fulfilling the requirement ip(t) G V_ i . Therefore, since V_i is spanned by {<p(2t — k)}, we conclude that ip(t) can be represented by the following weighted sum: where coefficients h(k) are the so-called scaling filter coefficients and the factor \/2 maintains the norm of the scaling function tp(2t). It is obvious then that subspaces Vj in Eq. 2.10 fulfill the requirements given by Eqs. 2.5 and 2.6 (or 2.7). Actually, the requirements associated with Eqs. 2.5, 2.6 and 2.8 formally define the MRA. Equation 2.12 is called the refinement equation or the MRA equation. It is fundamental for the theory of generating the scaling functions. Figure 2.3.a shows the Haar scaling function </?(£), which is the simple unit-width, unit-height pulse function, while Fig. 2.3.C illustrates two different resolutions in the Haar MRA. For this scaling function, the MRA equation is satisfied for h(0) = h(l) = -A= since ip(t) = (p(2t) + ip(2t — 1). However, representing the function /(£) G L2(TZ) as a weighted sum of tpjtk{t), such as in Eq. 2.11, and with decreasing j for higher resolutions, is not convenient since the subspaces {Vj} are embedded in each other and, thus, are clearly not orthogonal to each other. This means that the direct union of the ONBs that span each Vj is not an overall ONB for the whole space L2(1Z). Therefore, we define a slightly different set of functions i/>j,k(t), so-called wavelet functions, that span the differences between the spaces spanned by the various scales of the scaling function <p(t). We further require that the scaling functions and wavelets are orthogonal to each other. The orthogonal basis are computationally less demanding and they allow the exact partitioning of the signal energy by the expansion coefficients using Parseval's equality. Thus, the orthogonal complement of Vj to V , _ i is defined as Wj by requiring (y>j,fc(£), ipj,i(t)) = 0 for j, k,l G Z. This is denoted by: where Vo is the initial space spanned by the translations of the scaling function ip(t). However, the scale jo of the initial scaling space can be chosen arbitrarily and it usually corresponds to the coarsest signal resolution of interest (e.g. if j0 = 3, we have L2(TZ) = V 3 8 W2 © Wi 0 • • • © W_oo)- Figure 2.4 depicts the nesting of the scaling function spaces Vj, and the wavelet spaces Wj as their orthogonal complements. Since Wo G V _ i , these wavelets can be expressed as a weighted sum of translations of the scaling function f(2t) as defined in Eq. 2.12, that is, (2.12) Vj_i = Vj © Wj with Vj 1 Wj. (2.13) In general, this leads to the following span of L2(TZ) L2(K) = V 0 © Wo © W-i © ... (2.14) TP(t) = ^2hw(k)V2ip(2t-k), keZ (2.15) CHAPTER 2. WAVELETS 32 <p(0 = (p(2r)+ tp(2/ — 1) V(0 = cp(2/)-(p(2?-l) (a) (b) Projection onto V 0 _t I i 1 I I J i J i L -S: -4 -3; -2 -1 0 : : 1' '2 . 3 * S Time (s) Projection onto V.! W_ ± ... 1 W , X W, ! W o i V 0 Figure 2.4: Scaling function and wavelet vector spaces. CHAPTER'2. WAVELETS 33 for some set of so-called wavelet filter coefficients hw(k). It can be shown that the orthogonality requirement poses the following relationship between the wavelet filter coefficients and the scaling filter coefficients: hw(k) = (-l)kh(l-k), (2.16) which, in the case, of a finite even length-TV h(n), is hw(k) = (-l)kh(N-l-k). (2.17) Equation 2.15 generates the analyzing or mother wavelet function ip{t) for the wavelet family {^ •,fc(*) = 2-*^(2-J't-fc)}. (2.18) Figure 2.3.b depicts the simple Haar wavelet for which ^(0) = -~= and hw(l) = Thus, we have constructed a set of functions <Pj,k{t) and ipj,k{t) that can expand any function f(t) G L2(1Z) as following: oo j—jo oo /(*)= E CA>W ¥>*, ,* (*) + E E d i ( * ) ^ , k (2-19) -OO k~— OO or oo j=0 oo /(*)= E c(*0 + E E dj(k)^J,k (for j 0 = 0) (2.20) fc= — o o j = —oofc=—oo , where jo sets the initial space spanned by {fj0,k(t) = 2~i% ip(2~i°t — k)}, and the rest of L^itZ) is spanned by wavelets {ipj,k(t)}. The coefficients of the expansion in Eqs. 2.19 and 2.20 are actually the discrete wavelet transform of the signal f(t). The first summation of the expansion gives a function that is a low resolution or coarse approximation of f(t). Thus, the coefficients c3(k) are called the approximation coefficients or scaling coefficients. The second summation provides additional details of higher or finer resolutions. The coefficients dj(k) are called the detail coefficients or wavelet coefficients. As the detail resolution increases (i.e. j —> —oo), the approximation error goes to zero. For an ONB (or a tight frame), we further have c(k) = (f(t), <pk(t)) and dj{k) = (f(t),i>jtk{t)). (2.21) Then, the Parseval's equality holds and relates the energy of the signal f(t) to the energy in each of the expansion coefficients: oo j=0 oo •^!i/wu2= E iic(fc)ii2+ E E itowii2. (2-22) k= — oo j—~ oo oo -where K is constant due to the redundancy and K = 1 in the case of an ONB (non-redundant expansion). Finally, Daubechies showed [128, 129] that there exist ONBs with the compact support of the scaling function and wavelets. As previously mentioned, this is of a great significance since it allows for DWT computation by an efficient algorithm using FIR filters. CHAPTER 2. WAVELETS 34 Fast Wavelet Transform and Filter Banks Mallatt's fast algorithm [127, 130] for the computation of the DWT coefficients associated with compactly supported ONBs uses discrete convolutions with FIR filters and downsampling, by a factor of two, opera-tions8. Therefore, it is especially well suited for implementation on a digital platform. We then consider only the coefficients h(k) and hw(k) (Eqs. 2.12 and 2.15), and Cj(k) and dj(k) (Eqs. 2.19 and 2.20), that can be viewed as digital filters and digital signals respectively. It can be shown that, in the case of orthonormal wavelet bases9, Eqs. 2.12, 2.15 and 2.19 lead to the Cj{k)= ^ h{2k -m)cj-i(m) (2.23) and dj(k)= hw(2k-m)cj-i(rn) (2.24) where Cj(k) and dj(k), for j, k G Z, are the scaling and wavelet coefficients at the level (i.e. scale) j expressed as a function of the scaling coefficients Cj_i(fc) at the level j — 1. This analysis actually shows that the scaling and wavelet coefficients at the coarser scale j can be obtained by filtering10 the scaling coefficients at the finer scale j - 1 by two digital filters with coefficients ho(n) = h(—n) and h±(n) — hw(—n) followed by down-sampling11. The low-pass filter ho(n) and the high-pass filter hi(n), each followed by a down-sampler, form the so-called two-band analysis bank. The analysis bank splits the spectrum of its input signal (i.e. the vector of coefficients C j _ i ) into two equal parts (i.e. a low and high band). This process is applied iteratively, the so-called iterating the filter bank, starting at the finest resolution level J to obtain the desired low resolution level jo. It can be further shown by combining Eqs. 2.11, 2.9 and 2.15 that a reconstruction of the scaling coeffi-cients at the finer resolution j — 1 can be obtained from the scaling function and wavelet coefficients at the 8 T h e r e are other nonorthonormal (e.g. biorthogonal) M R A - t y p e wavelet families that also have fast algorithms associated wi th them [131]. 9 T h i s is also valid in the case of a tight frame. 1 0 I n the theory of digital signal processing, for an input sequence x(n) (i.e. evenly spaced samples of x(t) in time) and filter coefficients h(n), the output sequence y(n), which is the sequence x(n) "filtered" by a digi ta l filter h(ri), is given by a following convolution sum: J V - l J V - l y(n) = h(k) x(n — fe) = VJ h(n — fc) x(k), k=o k=o where N is the length of h(n). If the number N is finite, the filter is called Finite Impulse Response Filter (FIR) filter. If the number JV is infinite, it is called Infinite Impulse Response Filter (HR) filter. T h e frequency response of a digi ta l filter is the discrete-time F T of its impulse response h(n), that is oo H(J) = VJ h(n)eiun. n= —o 1 1 For an input sequence x(n), we get the down-sampled sequence y(n) by taking every other sample of x(n) (i.e. y(n) = x(1n)) and we denote this operation by y(n) =[2 x(n). In some cases, the down-sampling can be done by keeping the odd samples of x(n) (i.e. y(n) = x(2n+ 1)). Also , the down-sampling can be done by a factor other than two. However, if not stated otherwise, we assume that the down-sampling is done by keeping every other and even sample of the input sequence. CHAPTER 2. WAVELETS 35 Synthesize Figure 2.5: Perfect reconstruction filter bank for the fast wavelet algorithm (adapted from [7]). coarser resolution j, that is, for f(t) £ Vj- i we have cj_i (fc) = c i ( m ) h(k - 2 m ) + Yl di(m) M f c - 2 m ) > ( 2 - 2 5 ) where j, k £ Z. The synthesis operation in Eq. 2.25 actually shows that the scaling coefficients at level j — I can be realized by first up-sampling12 both the scaling and wavelet coefficients at the coarser level j, then by filtering the obtained sequences using two digital filters go(n) = ho(—n) and gi(n) = h\(—n), and finally adding them up. This structure forms the so called two-band synthesis bank and it is applied iteratively starting at the coarsest resolution level jo to obtain the desired high resolution level J. The set of analysis and synthesis operations is known as Mallat !s algorithm or the Fast Wavelet Transform (FWT). We have previously shown that the scaling function ip(t) cannot be chosen arbitrarily in order to generate the MRA. Thus, the appropriate sequences h(n) for the construction of scaling functions according to Eq. 2.12 have to fulfill a set of conditions. This is crucial in the implementation of the DWT associated with multiresolutions and it further poses certain constrains on the filters ho(n), hi(n), go(n) and gi(n) [6]. The up-sampling operation clearly does not lose any information. On the other hand, when down-sampling, there is a possibility of losing information since half of the data is discarded. However, if the filters ho(n), h\(n), go{n) and g\(n) define a perfect reconstruction filter bank13, the aliasing is cancelled out. This is schematically illustrated in Fig. 2.5. In general, the perfect reconstruction filter banks lead to biorthogonal solutions [132]. In the case of orthonormal bases, filters ho{n), hi(n), go{n) and g\(n) satisfy perfect reconstruction requirements, and further, filters ho(n) and h\(n) also satisfy the so called Quadrature Mirror Filter (QMF) 1 2 F o r an input sequence z(n), we get the up-sampled sequence y(n) by inserting a zero between each of the original samples of x(n) (i.e. y(2n) — x(n) and y(2n + 1) = 0) and we denote this operation by y(n) = | 2 x(n). In some cases, the up-sampling can be done by assignments y{2n) = 0 and y(2n + 1) = x(n), or by a factor greater than 2 (i.e. by inserting more than one zero in between the samples of x(n)). However, if not stated otherwise, we assume the first of the above definitions. 1 3 I n general, filters ho{n), h\(n), go(n) and g\{n) define a perfect reconstruction filter bank in z-domain if and only if Ho(z)G0{z) + Hx{z)Gi{z) = 2, H0(z)Go(-z) + H1(z)G1(-z) = 0. CHAPTER 2. WAVELETS 36 condition14. However, for the QMF solution, only IIR filters are allowed (except for the Haar). The alternative solution is to chose the alternating flip filter bank given by Eqs. 2.16 and 2.17. This solution allows the use of FIR filters and the Daubechies filters fit this pattern15 [6, 132]. Figure 2.6 illustrates Daubechies 4 scaling and wavelet function in both time and frequency domains as well as their corresponding discrete filters. A comprehensive elaboration on the necessary and sufficient conditions, both in the time and frequency domain, for the scaling function and the associated filter banks and for various discrete WTs can be found in [5, 7, 132, 6]. In practice, we are dealing with band-limited signals, and, thus, the wavelet coefficients dj(k) (repre-senting signal details) are negligible at scales that are smaller than some scale J. Actually, for some small enough scale J, the scaling coefficients cj(k) are samples of f(t), since the scaling functions ipjtk(t) act as "delta functions" in inner products Cj(fc) = (/(£), <Pj,k(£))16- Thus, we start with the initial sequence cj (or input signal), and we calculate the DWT to as low a resolution j = j 0 as desired through J — jo analysis stages. So, for f(t) G Vj, we have 3=30 /(«) = E cJO v*.,fc(*) + E Edi(fc)^ .*(')' (2-26) kez j=j-ikez which is a finite scale version of Eq. 2.19. Thus the signal S = /(£) is represented by the sum of its coarse approximation Aj0 at the scale jo (i.e. Aj0 = ^ f c Cj0 fj0,k(t)) plus the details Dj at higher resolutions (i.e. Dj = Z)k dj'CO V'j.fcC*) for J — 1 < j < jo). Also, note that for a compactly supported basis, only a finite number of additions over k is necessary. Further, we commonly assume that a signal specified by its sample values corresponds to the resolution 0, that is, the most frequently used initialization rule for the F W T algorithm is c0(k) = f{k), k = 0,...,N-l. (2.27) In order to perform a wavelet decomposition with a certain depth L, the signal length N has to be a multiple of 2L. Using the analysis filter bank (i.e. the decomposition rule specified by Eqs. 2.23 and 2.24), the sequence Co is then split into its half-length subparts {ci(A;)}fc=o,...,JVi-i and {di(fc)}fc=o,...,Ar,-i, (2.28) where A i^ = f • The process is iterated onto the approximation sequence Cj until ones reaches the final level L. The wavelet representation of the signal is provided by the detail (or wavelet) sequences dj, j = 1, ...,L (a total of Y + ••• + = N — ^ coefficients) plus the approximation sequence CL ( ^ T coefficients) at coarse resolution L. This F W T algorithm requires 0(N) operations for L = logN levels of decomposition. The 1 4 T h e high-pass response |i/i(a;)| is a.mirror image of the low-pass response |ifo(^)| with respect to the middle frequency ^ (the quadrature frequency). This can be obtained if the filter hi(n) is of a form ...ho(2), — ho(l), /in(0), — /io(l), ho(2)... (the alternating signs pattern). 1 5 T h e flattest \HO(CJ)\2 leads to the Daubechies wavelets. 1 6 T h e samples of f(t) above the Nyquist rate are good approximation to the scaling coefficients at that scale. CHAPTER 2. WAVELETS 37 Figure 2.6: Daubechies 4 filters: (a) low and high pass filter coefficients and their corresponding scaling and wavelet functions; (b) power spectra of scaling (solid) and wavelet (dashed) functions (adapted from [6]). CHAPTER 2. WAVELETS 38 128 samples Hi(u>)U0-Ho(u)b.(J>-> H,(UJ) w 2 c, (64) H o ( w ) ^ 0 - > > d,(64) •> d 2 (32) c 2(32) H ^ X j h ^ d s (16) Ho(w ) k l > ^ c 3 (16) stage 1 staae 2 stage 3 (a) |H(o.) | BAi BDi BDi V 3 I W 3 I w . 0 j t / 8 7t/4 W , J U / 2 • m rt = 64 H z S1 cz V —j i-» ^ v Q j — A, cz v 2 A,c:v2 D , cz W 3 (b) (c ) Figure 2.7: Three-level DWT (L = 3): (a) two-band filter bank structure; (b) analysis tree; and (c) filter bank frequency responses. input signal can be exactly recovered from these transformed coefficients by applying sequentially the reverse synthesis process described by Eq. 2.25 and for j = L — 1 down to j = 0. The 3-level DWT filter bank structure is illustrated in Fig. 2.7.a for a 1-second signal 5 sampled at 128 Hz. The approximation coefficients are the down-sampled outputs of the low-pass filters, while the detail coefficients are the down-sampled outputs of the high-pass filters. In a dyadic MRA, the decomposition process is iterated so that the approximation coefficients are further decomposed through successive stages. The original signal can be reconstructed from its approximation and details at each stage, for example, for a 3-level signal decomposition, the signal S can be written as S = A3 + D3 + D2 + D\ as illustrated in Fig. 2.7.b (so-called 3-stage analysis tree). Frequency bands corresponding to this analysis tree are depicted in Fig. 2.7.c. The decomposition may proceed until the approximation coefficients consist of a single sample. The nature of the process generates a set of vectors C 3 , 0I3, d 2 , and di , containing the approximation and detail coefficients. These vectors are of different lengths (as indicated in figure), based on powers of 2. They contain information that defines approximation and detail signals (i.e. A3, D3, D2 and D \ , respectively) at different frequency bands (BA3, BD3, BD2, and -Boi) determined by the filter bank frequency response. As expected, these bands are of unequal widths based on powers of 2 (so-called octave bands). CHAPTER 2. WAVELETS 39 2.1.3 Redundant (Over-complete) and Stationary Wavelet Transform There are many possible discretizations of the CWT and few constrains on the spacing between sample points in the time-scale plane17. However, we will discuss a subclass of sampling sets that allows for a fast computation of the transform and its inverse. One important choice of a sampling grid is when the scales are still given by dj — 2J, but the translation steps are chosen to be bk = k, j,k G Z. In this case, the space L^H) is the span of wavelet functions { 2 - £ if}(2~3(t — k))} and they constitute a frame of £ 2 ( ^ ) -This means that any signal /(£) G L2CR-) can be recovered from its transform by a reconstruction algorithm. In general, the resulting transform is redundant since the underlying basis functions are nonorthogonal and over-complete. The terms Redundant Wavelet Transform (RWT) or Over-Complete Wavelet Transform (OCWT) are often used to describe this particular case of a wider class of redundant discrete WTs. One potential limitation of the non-redundant WTs such as the DWT is that they are not translation-invariant. This means that the DWT of a simple translation of the input signal generally results in a nontrivial modification of the wavelet coefficients [133]. Also, at higher levels of decomposition (i.e. lower frequency bands), the DWT has a time sampling rate that might be essentially too low. Although the DWT offers low complexity, it is not necessarily the best solution in applications such as certain pattern recognition tasks, detection of singularities and de-noising (i.e. noise removal or reduction), or, for example, in statistical explorations and local spectra estimations of non-stationary time series. A simple way to alleviate these problems is to compute the DWT for all translated versions of the signal and then average the results [133], or equivalently, to modify the basic DWT algorithm in Fig. 2.5 by removing the down-sampling modules in the analysis and the up-sampling modules in the synthesis. Instead, we modify the low- and high-pass filters at each stage j by padding them out with zeros as follows: (i) the low-pass analysis filter is defined by h$\ri) =T2J'-1 ho(n), (ii) the high-pass analysis filter is defined by h^\n) = t 2 3 - x h\(n), and (iii) the low- and high-pass synthesis filters g^\n) and gf\n) are obtained in a similar manner from filters g0(n) and 51 (n), where ho(n), /ii(n), go(n) and g\(n) are the analysis and synthesis filter banks previously defined for the DWT. The operator "fm denotes the up-sampling by a factor of m (i.e. inserting m — 1 zeros between every adjacent pair of elements of the original sequence). In this way, one obtains a redundant decomposition that is perfectly reversible and closely related to the above defined DWT. Since this transform is translation-invariant, it is commonly denoted as Stationary Wavelet Transform (SWT). Furthermore, if we start with the DWT associated with ONBs, the underlying basis functions form a tight frame with the redundancy factor being equal to 2. Hence, the Parseval equality given by Eq. 2.22 holds, with K — 2. Indeed, we have twice as many basis functions and corresponding expansion coefficients than in the non-redundant case. This redundancy, besides allowing for the translation 1 7 T h e underlying basis functions must constitute a frame of a space of interest. CHAPTER 2. WAVELETS 40 invariance, gives more freedom in design of the filter banks, carries great potential for statistical applications and provides more robustness in the presence of noise. The potential performance benefits from using the SWT come at the cost of increased computation and storage requirements. The approximation and detail sequences at each level of decomposition are of the same length as the original sequence, rather than becoming shorter by a factor of 2 as the level of decomposition increases. The frequency responses of the filter bank are the same as those depicted in Fig. 2.7.c. However, in each of the frequency bands we have as many coefficients as in the original signal, with their time resolution being equal to the time interval between the signal samples. The complexity at level L is increased from 0(2L) to 0(L 2L) for a signal of length N = 2L. 2.1.4 Wavelet Packets The classical wavelet systems realized with 2-band filter banks, such as the standard DWT or SWT, result in a logarithmic frequency resolution. Thus, the low frequencies have narrow bandwidths and the high frequencies have wide bandwidths (see Figs. 2.1.d and 2.7.c). However, there is a variety of generalizations and extensions to the basic wavelet systems, that allow for better control over the tiling of the time-frequency plane. For example, instead of 2-band filter banks, one can use M-band filter banks. Of particular interest are Wavelet Packets (WPs), a simple extension of standard DWT, which allow a finer and adjustable resolution at high frequencies. The WP analysis is obtained by iterating the high-pass branch of the Mallatt's algorithm tree (see Fig. 2.7.a) as well as the low-pass branch. At each stage, details (outputs of high-pass filter \HI(UJ)\) as well as approximations (outputs of low-pass filter \H0(u))\) are further decomposed into low and high frequency signal components. Figure 2.8.a shows the resulting full binary WP tree. Any pruning of this full tree generates a valid WP basis and allows a flexible tiling of the time-scale plane [6]. Accordingly, a given signal can be written in a more flexible way than provided by the standard dyadic decomposition. For example, at level 3 we can write S = Ai+ AD2 + ADD% + DDD3, where DDD3 is the signal component of the narrow high frequency band BpDD3 • WP analysis results in signal decomposition with equal frequency bandwidths at each level of decom-position. This also leads to an equal number of approximation and detail coefficients in each band. This is a desirable feature for data analysis and information extraction. Figure 2.8.b illustrates frequency bands for the 3-level WP decomposition. The cost of this more flexible structure is a computational complexity of O(NlogN). 2.2 Wavelets in Biomedical Appl icat ions Biomedical signals give a fundamental insight into the functioning of a human body, and, therefore, are essential for diagnosing a wide spectrum of diseases. In general, information carried by bioelectric signals are non-stationary, consisting of impulse-like events as well as more diffuse oscillations such as background textures. Also, these signals are usually rather noisy. CHAPTER 2. WAVELETS 41 I 5cz v 0 A AA2 V2 AAA, , V 3 DA4 3 c W 3 r-l AD2a W2 ADA 3 w 2 0 £>ZX43 c W 2 1 DA2 c W 1 0 /MD3 c W100 3 w 1 0 1 r-i DD2cz ADD, W l 1 ( (a) H(u) V 3 A W 3 A ' W 2 o A W 2 1 A W 1 o o A W 1 o i A W 1 1 o A W 1 1 1 - • c o 0 TU/8 T t / 4 3n/8 i t / 2 5ir/8 3TC/4 7n/i W 2 0 u W 2 , = W 2 V0 = V 3 u W 3 u W 2 u ' V V V 3 c V 2 cr V, cr V 0 (b) Figure 2.8: Three-level wavelet packet decomposition: (a) analysis tree (approximations and details); (b) filter bank frequency responses. Wavelets are becoming increasingly popular in the biomedical field since they have a number of fa-vorable properties. They are well suited for the analysis of both small- and large- scale signal structures due to their excellent joint time-frequency localization. Also, wavelets act as a multi-scale matched filter. Thus, since the great majority of signals generated through biological processes comprise of wavelet-like elements [126], wavelets matching certain shapes can be used to detect specific patterns associated with diseases and, otherwise, hidden in a non-stationary signal [134]. The WT can be further viewed as a special spectral analyzer that provides energy estimates in different decomposition bands or other related param-eters. Finally, wavelet-based de-noising techniques, pioneered by Donoho [135], have proven to be simple, and yet very efficient. The range of applications of the wavelet analysis in the field of medicine and biology is extremely large. Two main areas of applicability are medical imaging and biomedical signal processing. Also, the impressive similarity between the WT and the biological information processing, such as the auditory or visual system, has led to a variety of wavelet-based biological models [126]. Medical images, such as those obtained by X-ray Computed Tomography (CT), Magnetic Resonance Imaging (MRI) and mammography, are an invaluable aid for diagnosis. Functional neuroimaging, such as CHAPTER 2. WAVELETS 42 Positron Emission Tomography (PET) and functional MRI (/MRI), is a fascinating and fast developing area that provides a window into the neuronal activity of the brain in vivo. However, medical images consist of data that are indirect, incomplete and generally very noisy projections or observations of the object of interest. They often require statistical analysis for their interpretation. Wavelets have been used in various applications, from image de-noising and enhancement, to statistical analysis of image differences and patterns for diagnosing (e.g. in /MRI). Further, many types of medical images are characterized primarily by their edges, which leads to the extreme sparseness of their wavelet representations. Hence, wavelets have been successfully used for their compression. They have also been applied for feature extraction, for example, in digital mammography. Wavelet based encoding schemes in MRI significantly reduce imaging times thus avoiding hardware costs associated with fast imaging methods. For an excellent overview of wavelets in medical imaging, refer to [7, 126]. The W T has found great applicability in processing of 1-D physiological signals, particularly the E C G and E E G , including the Evoked Response Potentials (ERPs). The W T has been used for the time-frequency analysis and characterization of the heart beat sounds and murmurs [126]. It has been shown that the W T could pick up certain sound components untraceable by other methods. Wavelets have also been utilized for the automatic recognition and timing of the E C G constituent waves (so-called QRS complex and P, T and U waves) [126]. In the U.S. patent No. WO 00/69517, the W T was applied to the E C G signal for discriminating between the normal and abnormal heart rhythms. The algorithm identified the higher amplitude wavelet coefficients, and then compared them with sets of templates derived from the rhythms of known type, in order to reach the decision. Neural networks were also investigated for diagnosing the coronary artery disease. A method, which utilizes the network inputs comprising of parameters obtained by wavelet-based preprocessing of diastolic heart sounds in combination with physical exam variables, has obtained favorable results [136]. More applications of wavelets in biomedical signal processing can be found in [7]. In the case of E E G (or electrocorticogram18 (ECoG)), substantial efforts have been made for the de-tection of epileptic waveforms (spikes, sharp waves and spike-and-wave complexes) and seizures. A reliable short-time prediction of the clinical onset of seizure and localization of the epileptogenic focus are of tremen-dous therapeutic significance for providing a time window during which appropriate intervention can be undertaken to minimize the risks or abate the seizure. Early methods were limited to off-line studies due to excessive computational time needed for performing wavelet analysis. Thus, the development of filter banks for wavelet analysis has been extremely important to make these techniques real-time [134], and hence ex-pand their practical utility. Majority of methods for detecting epileptiform activity are essentially based on the time-frequency localization of signal energy and its relevant changes. Wavelets have proven to be well suited for this task. Moreover promising results of the wavelet-based methods for the real-time detection of epileptic waveforms [137, 138,139, 140] and short-term prediction of seizures [141] has already been achieved. Further, the E E G signal in epileptic patients is commonly monitored over long periods of time in'order, for 1 8 T h e electrical act ivi ty recorded from the cerebral cortex by "in-depth" subdural electrodes, which provides recordings wi th min ima l artifact contamination. CHAPTER 2. WAVELETS 43 example, to localize the seizure focus for its surgical removal. Thus, detecting and keeping only data relevant for characterizing the seizures leads to significant data reduction of long-term E E G recordings [142]. Also, the W T has been applied to the problem of localization of seizures [143, 144]. Particularly scalograms have proven to be a promising tool for both the localization and prediction of seizures [145]. An interesting application of wavelet analysis of the E E G signal by Sartene et.al. [146] offered perspec-tives for automated analysis of polysomnography19 signals. The method utilizes the quantification of the instantaneous energies in conventional E E G frequency bands (S, 9, a and /3) in order to define the sleep stages. It has been demonstrated that the ratio -f^y can be chosen as an index of depth of sleep. Also, the W T has proven to be particularly suitable for studying the interactions between ventilatory and neurological activity. But, as wavelets start to be increasingly used to detect patterns in the E E G , very little literature exists concerning wavelet transform applied for quantifying the anesthetic drug effect. A study of interest has been carried out by R. Roy et al. from the department of biomedical engineering of Rensselaer Polytechnic Institute. They applied wavelet transform to characterize the MLAER of anesthetized dogs [147]. Wavelet coefficients were input to a neural network to determine the anesthetic depth. Fuzzy logic was then used to predict and regulate the drug concentration. Another quite extensive paper was written by Angel et al. [148] on the use of Somatosensory Evoked Potentials (SEPs) to measure the appropriateness of sedation. In this work, wavelets are used in conjunction with Artificial Intelligence (Al) to measure latency changes in the evoked potentials. For more on wavelet analysis of neuroelectic waveforms, see a conceptual tutorial by Samar et al. [149]. The increasing amount of literature on the subject seems to indicate that we are now at the onset of a small technical revolution in signal processing. 2.3 Motivation for Wavelet Analysis of the E E G As presented in Chapter 1, the BIS monitor by Aspect Medical Systems, Inc. is currently the most optimized commercial monitor for measuring the effects of anesthetics on the brain. It is a composite index based both on power spectrum and bispectral analysis of the E E G signal. Thus, unlike parameters derived by power spectrum analysis alone, the BIS contains both power and phase information. Indeed, we have seen that anesthesia increases the cortical synchronicity and, therefore, modifies the latency of some of the frequency components of the E E G . Therefore, for example, for light anesthesia, we can expect small phase shifts of some of the frequency components. For deeper anesthesia, we can expect stronger phase relationships between E E G generators. These changes are not observable by spectral analysis alone, but can be captured by bispectral analysis. Thus, since we want to investigate the use of wavelets for quantifying the anesthetic-induced unconscious-1 9 Po lysomnography is a rather extensive method of investigation of the main physiological functions dur ing sleep, especially neuro-cardio-respiratory ones, and requires a simultaneous recording of a great number of signals ( B E G , Electrooculogram ( E O G ) , E M G , E C G , blood pressure, respiratory variables etc.). CHAPTER 2. WAVELETS 44 Repeating Sequence Sine Wave 8 Hz RI SineWave12 Hz (a) To Workspace Scope 75 150 225 300 375 Time (s) (b) Fi gure 2.9: The three-harmonic signal sout with a drifting phase: (a) Simulink synthesis; (b) phase change over time. ness, we first have to assess if the wavelet decomposition carries similar information as bispectral analysis. To do so, we use a synthetic signal created with Simulink and try to track the changes from the baseline value using wavelets. Figure 2.9.a illustrates the synthesis of the 3-harmonic signal sout. Two of the signal components are fixed both in amplitude and phase. The phase of the third frequency component changes over time as depicted in Fig. 2.9.b. Figure 2.10 shows the power spectrum of the signal sout over time, calculated for 4-second epochs and sampling frequency of 128 Hz. The change of phase over time is not observable. Figure 2.11 shows the change of the absolute value of wavelet coefficients over time, obtained by the wavelet decomposition of the signal sout for 1-second epochs and the same sampling frequency of 128 Hz. The phase drift over time is clearly observable. Thus, the wavelets are capable of capturing the changes in the phase relationships between underlying E E G generators. Hence, wavelets have a potentialfor quantifying the effect of anesthetics on the brain. 2.4 Summary The W T is a natural tool for time-frequency signal analysis since each analysis template ipa^ is well localized in a certain region of the time-frequency plane, with a central frequency inversely proportional to a and a time location centered around b. The multiresolution nature of the analysis enables the W T to zoom in or zoom out on the small scale or large scale structures. This makes this transform particularly suitable for the analysis of non-stationary signals such as the E E G . Indeed, wavelet analysis has achieved great successes in biomedical signal processing over the last decade. The distinction between the various types of WTs depends on the way the scale a and translation b parameters are discretized. The CWT is the most redundant transform, with a and b being continuous variables. The DWT is a non-redundant transform that uses the dyadic scales a = 2J and critically sampled translations b = fc2J. In between these two cases, there are many over-complete (i.e. redundant) transforms CHAPTER 2. WAVELETS Frequency (Hz) Figure 2.10: The power spectrum of the signal sout over time calculated for 4-second epochs. Figure 2.11: The absolute value of wavelet coefficients of the signal sout over time. CHAPTER 2. WAVELETS 46 that use a finer grid for a and b, the SWT being of particular interest, where a — 2J and b — k, (j, k) € Z2. The availability of a fast wavelet algorithm, which uses FIR filters and has an 0(N) complexity, has made the non-redundant WTs with compactly supported ONBs extremely popular. They have found great applicability in data compression and processing (generalized filtering), and when the orthogonality of the transform is of importance. Nevertheless, due to their translation-invariance and robustness, the redundant WTs are generally more suitable for signal analysis, feature detection and extraction tasks, statistical ap-plications and de-noising. However, they are more computationally demanding, and therefore, the choice between these two options is not always obvious. Also, obtaining more general representations such as with WPs, which allow for a more flexible tiling of the time-frequency plane, is straightforward by the simple extension of the standard pyramidal DWT decomposition. Finally, we have successfully tested the ability of the WT to capture phase shifts of underlying signal components. Since the state of anesthesia is associated with the increased synchronicity of the E E G com-ponents, we can conclude that wavelet analysis has potentials for the quantification of the anesthetic effect on the human brain. Such features of W T are to be investigated in the next chapter. Chapter 3 WAV Index: Methodology and Implementation This chapter presents the methodology for the derivation of an index of anesthetic-induced unconsciousness (i.e. hypnosis). The index, further referred to as the WAV Index, is based on the wavelet analysis of the E E G , and the differentiation between two known hypnotic states (awake and anesthetized). In Section 3.1, the experimental data sets that were used to develop the index are introduced. Information extracted from these data relative to the patient's hypnotic state constitutes the foundation on which the feature extraction algorithm of the WAV Index is based. The feature extraction technique and the method to determine the best wavelet suited for the derivation of the index are the subject of Section 3.2. The real-time implementation is further discussed in Section 3.3, where the detailed algorithms and flow charts are provided. Finally, some preliminary results are presented in Section 3.4. Note that the basic ideas behind the proposed technique have already been the subject of a prior publication [150]. 3.1 Experimental Data Sets All data used in this chapter have been acquired using the commercial BIS® 1050A monitor (Aspect Medical Systems, Inc.) updated with the version v.3.4 of the analysis software1. For this purpose, we have used an interface2 which allows the transfer of data from the BIS monitor to the hard drive of a laptop computer through serial communication. Both the processed data (bispectral index, MEF, SEF, etc.) and the raw E E G can thus be accessed and stored. The interface provides new data every second. This setup (i.e laptop and BIS monitor) was approved for use in the Operating Room (OR) by the Department of Biomedical Engineering of the University of British Columbia Hospital. 1 This is actually the analysis software utilized by BIS r 2000 monitor. 2 This interface was developed by Ms. Lou Ann Mendoza, an undergraduate student from the Department of Electrical and Computer Engineering, as part of the requirements for a course project. 47 CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 48 For the purpose of developing the hypnosis index, we first acquired a record of both raw and processed data from a single surgical procedure. 3.1.1 Raw E E G The raw E E G data was obtained from a surgical procedure (bilateral breast reconstruction) that lasted for about one hour and fifty minutes. A single frontal E E G channel was acquired at the fixed sampling rate of 128 Hz using the setup described above. The raw E E G data was transferred via the interface to the hard drive of the laptop computer, and is displayed in Fig. 3.1. a. A number of large spikes within the E E G can CD 100 80 60 40 20 0 0 0.2 0.4 0.6 0.8 1 Time (hours) (a) Event 1 1.2 1.4 1.6 1.8 Event 4 Electrocautery artifacts N Event 2 0.2 0.4 0.6 0.8 1 Time (hours) (b) 1.2 1.4 1.6 1.8 Figure 3.1: Surgical procedure data: (a) raw electroencephalogram; (b) bispectral index. be observed. They are mainly due to surgical manipulations using an electrocautery device which provokes strong interferences that saturates the amplifying unit of the monitor. Note that the E E G signal is indeed quite sensitive to environmental noise. Such heavily corrupted epochs will first need to be detected and excluded from the analysis when deriving the index. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 49 3.1.2 Processed Data The BIS presented in Fig. 3.1.b corresponds to the raw E E G signal in Fig. 3.1.a. This index was computed and displayed in real time by the BIS monitor. Observing the bispectral index, we can distinguish 4 major events that occurred during the procedure: - event 1: the patient arrives in the OR in a mildly sedated state (BIS = 80); - event 2: after injection of different anesthetics, the patient is brought into a deep hypnotic state to facilitate the intubation; - event 3: the patient is moved from the supine position to a sitting position to allow a surgeon to assess the surgical result, which provokes a strong stimulation due to the inherent movement of the tube inserted in the trachea. The arousal in the patient's hypnotic level can be observed; - event 4: end of surgery, the drugs are discontinued, and the patient is reacting. A good index of the level of hypnosis should be able to capture these events as they are representative of the surgical procedure. Note that the large electrocautery interferences are rejected from the analysis. In this case, due to the lack of uncorrupted data, the BIS monitor is not able to provide an index of hypnosis relevant to the patient's state. In exchange, it notifies the user by displaying a negative index. Besides processing the raw E E G in order to derive the bispectral index, the BIS monitor also calculates standard indexes based on spectral analysis, such as the M E F and the SEF indexes. As can be seen in Fig. 3.2.a, the M E F fails to give an accurate representation of the hypnotic state. However, the SEF is more precise (see Fig. 3.2.b). Indeed, the SEF is able to detect the deep hypnotic state (event 2) and the return to the awake baseline value (event 4). However, it fails to detect the third event and does not give a good representation of the sedated states (event 1). It is worth noticing that "raw" indexes derived on short time windows (i.e. few seconds) are usually rather noisy. It is then common practice to filter these indexes to retrieve the main trend. Very few information concerning the BIS monitor filtering technique has been published. However, the user manual refers to a 30-second averaging filter. We will take this value as a basis for filtering our own index for the purpose of a meaningful comparison. 3.2 Feature Extraction Methodology This section investigates the use of wavelet transform to distinguish between two known hypnotic states, that is, the awake and the anesthetized state. Thus, the first task is to determine a feature that characterizes and discriminates between these states. Then, such a feature can be used to determine whether an unknown E E G epoch is more likely to belong to the awake or anesthetized state. Both a general mathematical approach and one illustrative example are presented in this section. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 10 r 50 X LU 2 N X co 0.8 1 Time (hours) Figure 3.2: Median Frequency (MEF) and Spectral Edge Frequency (SEF) indexes based on the EEG signal from Fig. 3.1.a.' 3.2.1 Training Data Sets The E E G data representing the awake state was acquired from a healthy and fully awake subject. A 1-minute recording during which the subject performed a mental task was obtained. To limit artifacts due to eye and muscle movements, the subject was asked to maintain a still position with eyes closed. In addition, Fig. 3.1.b shows that the hypnotic level of the patient corresponding to the BIS values between 30 and 35 was maintained during most of the maintenance phase of the anesthetic procedure. Thus, based upon the values of the bispectral index, and further confirmed by the anesthesiologist's assessment, we extracted from the raw E E G signal depicted in Fig. 3.1.a a set of data which is the representative of the anesthetized state. Care has been taken in selecting only the data that were not corrupted by the electrocautery interferences or the most apparent artifacts. This was achieved by visual inspection of the raw E E G , based on its amplitude and waveform. We have seen in Chapter 1 that the spontaneous E E G signal is a random-appearing signal with the peak-to-peak amplitude typically ranging from 2 to 200 /xV. Therefore, the epochs of too large or too small amplitude as well as those having unusual patterns (e.g. square pulses or spikes) were excluded from the analysis carried on further. Two 1-minute E E G signals corresponding to two known hypnotic states were obtained, which are depicted in Fig. 3.3. These two signals were further divided into epochs of fixed duration. Otherwise indicated, we assume that each epoch is 1 second long. Clearly, the amplitude of each epoch can differentiate between the two states. However, it is a well known fact that the E E G amplitude depends on many different factors such as the electrodes and their location on the scalp. Also, E E G amplitude may vary between patients. Hence, in order to derive a consistent index, we decided to normalize the amplitude of each epoch. By doing so, the most apparent information is lost. On the other hand, any index derived on normalized data will be more robust to electrode placement and CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 51 100 60 epochs (b) Figure 3.3: Training data sets: (a) awake; and (b) anesthetized. impedance (Appendix B illustrates this point). To normalize the amplitude, we calculate the Root Mean Square (RMS) value of each epoch. Let x be an epoch vector that contains the E E G samples. The RMS value of x is defined as: i?M5(x) i N v-£*( f c > 2 N fc=i where N is the number of samples. The normalized epoch x is then defined as: RMS(x)' 3.2.2 M e t h o d o l o g y : A M a t h e m a t i c a l A p p r o a c h Let us now formally define the two previously mentioned training data sets that represent the awake and anesthetized states: ( Tw = {xffl|fc, fc = 1,2, . . .M} (awake), Ta = {xa)fc, fc = 1,2,... M} (anesthetized), where the vectors x#]fc contain N = 128 samples, representing the fcth epoch of either the awake or anes-thetized data set, each containing M — 60 epochs. To characterize the data sets, a particular feature is extracted from each epoch. We define the feature extraction function / as: / : x. ) f e —> /(x. ) f c) = ' f . , f e . (3.1) Each normalized epoch x.^ is thus associated with a feature f . ^ . This feature can be either a scalar or a vector. A particular state can then be characterized by averaging the set {f,,k} over the corresponding training data set. This results in two averaged features f u , and f a , which are representatives of the awake and the anesthetized state, defined as: w M 1_ v ^ M f M ' Z^ fe= l A i w.ki and T — X. f a M ' Z^ fc= l 1a,fc-(3.2) CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 52 In order to assess the hypnotic state of a patient, we claim that it is sufficient to record the patient's E E G and calculate the feature f for each epoch x. Comparing this value to fw and fa, it is possible to calculate the likelihood for the patient to be either awake or anesthetized. Hence, two indexes iw (awake) and ia (anesthetized) are defined such that: { iw = \\f ~ fw\\i, and (3.3) *a = | | f - f a | | l , where the norm || • | | i for a vector v is defined as: N l|v||i = £ K i ) | . j = l The norm 11 • 111 accurately quantifies the difference between f and f, by integrating the distance between these two vectors. Higher degree norms can be used for this analysis as well as the correlation function between two vectors. However, such norms would emphasize larger differences which might be due to outliers, hence leading to noisier indexes. The final index i , so-called WAV Index, is obtained by combining the two indexes as following: i = i a - i w - (3.4) Note that the analysis based on a per epoch basis yields the real-time implementation of the index. The main obstacle is obviously the selection of an appropriate function / . In Chapter 2, we have given the rationale behind the analysis of the E E G using wavelets for the assessment of DOA. Thus, each normalized E E G epoch will be decomposed using the W T into a set of coefficients CL and dj: x—^{{cL;d,}, j = l,2,...L}, (3.5) where L is the level of decomposition. Each vector dj comprises the detail coefficients of the signal in a specific frequency band B^, while the vector ci comprises the approximation coefficients in the lowest frequency band BCL (Figure 2.7 depicts the case where L = 3.). As for the feature used to characterize each E E G epoch x, the Probability Density Function (PDF) of the wavelet coefficients in a chosen detail band Ba. is selected, that is, /:x—>/(*)= f = PDF(B d . ) . (3.6) Such a choice is motivated by the fact that the probability density function emphasizes neither large nor small coefficients. In fact, it tends to focus more on the general content of each wavelet decomposition band. Indeed, this property is interesting when dealing with noise-like signals such as the E E G . Other statistical functions, such as the variance or standard deviation of the wavelet coefficients, can also be considered. As we are interested in the feature extraction for statistical characterization of two distinct hypnotic states, we will use the SWT for the wavelet decomposition of the normalized E E G epochs. Indeed, as mentioned in Chapter 2, the SWT is better suited for this task than, for example, the DWT. This results in vectors cr, and dj, j = 1,2,... L, being of a length equal to the number of samples AT in each epoch. Thus, CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 53 we have more coefficients in the decomposition bands than for the DWT, resulting in more information for statistical feature extraction. Also, having the same number of coefficients in each vector facilitates various statistical comparisons between different bands of decomposition. Now, we are confronted with the task of finding the optimal wavelet and best frequency band for the analysis. 3.2.3 Best Wavelet Selection A difficulty rises when selecting an appropriate wavelet for decomposing the epochs, and choosing the best detail coefficient set dj (i.e. the frequency band Bdd) for carrying out the analysis. To compare the effectiveness of different wavelets, we introduce the discrimination parameter D as follows: D = \%-Ha\\i. (3.7) The discrimination parameter D quantifies the difference between the averaged features fa and fw, which are the PDFs of the wavelet coefficients in a yet to be chosen frequency band Bdr Obviously, to better distinguish between the awake and the anesthetized state, we need to calculate the parameter D through all the bands of decomposition, and then choose the one that maximizes D. Another important requirement is that the averaged features iw and fa must be reliable representatives of the individual epochs from the training data sets Tw and Ta, respectively. Hence, we define two variables, the awake variability vw and the anesthetized variability va as follows: * vw{k) = \\f(kw,k) ~ M i l , and (3.8) Va(k) = ||/(xa,fc) - fa 111, V where x.^ is fcth normalized epoch of the corresponding training data set T., for k = 1,2,... M. Thus, the kth element of the vector v. stores the difference between the feature derived on the individual fcth epoch of the data set T. and the average feature f. that represents all the epochs of that set. For example, in the ideal case, the averaged feature f. should be equal to the features of the individual epochs of the corresponding data set, and thus v. = 0. However, in reality, this is not the case, and therefore, we need to minimize both the mean and variance of the elements stored in the vectors v„, and v a . The best wavelet selection method has been applied to the training data sets presented in Fig. 3.3. An algorithm, which calculates the feature f.^ for a given epoch x.^ through all frequency bands up to a 4-level SWT decomposition, has been coded using standard Matlab wavelet functions. The sets have been processed to derive the averaged features and fa as well as the discrimination parameter D, together with the mean and variance of vectors vw and v a . The program has been iterated for different compactly supported wavelets. Orthogonal families, such as Daubechies, Coiflets and Symmlets, have been investigated, as well as biorthogonal spline and reverse biorthogonal spline wavelets. The 4-level SWT decomposition using Daubechies wavelet filters Dbl - Db303 resulted in the highest values of the discrimination D in comparison 3Daubechies wavelets of order JV, denoted by DbN, have the length of the scaling and the wavelet filter equal to 2 N. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 54 with other wavelet families. All results are displayed in Fig. 3.4. From this analysis, we conclude that the PDF of the band gives the most discriminating power for filters lengths larger than 12 (filters Db6 or higher), see Fig. 3.4.a. The means of v a and vw are minimal (or close to minimal) for the band B^ (see Fig. 3.4.b and c), confirming that this band is the best choice . Finally, Fig. 3.4.d shows the variance of v a and vw for the band B^ • The wavelet Db8 appears to be the best wavelet for the analysis: it has a good discriminating power D and close to the minimal mean and variance of v a and v m , while the corresponding filter is short allowing for faster calculations. Note that the Dbl4 is also suitable for this analysis, having a large discriminating power (close to the maximum value of D) and a smaller mean of v m . However, the filter length and the variance of v a and are also larger. In the following, we will give the advantage to the parameter D over the filter length and chose the wavelet Dbl4 as "optimal" for the analysis. To make an optimal choice a cost-based decision can be employed. In this work, we have limited ourselves to a choice obtained by inspection. The averaged PDFs awake and anesthetized are displayed in Fig. 3.5 for the first four bands of the SWT decomposition using the wavelet filter Dbl4. The awake and anesthetized templates corresponding to the band B^ (Fig. 3.5.a) are thus chosen as representatives of the awake and anesthetized states. Previous analysis singled out the probability density function of the band B^ as the most discriminating. This result is interesting since the B^ decomposition band corresponds to the signal detail in the 32-64 Hz frequency range. In neurophysiology, this particular frequency band, referred to as the 7-band, is often discarded in classical power spectral analysis since it carries a very small amount of the E E G energy. However, as mentioned in Chapter 1, this frequency range carries the information associated with cognitive processing. We appreciate this as a clear indication that the traditional methods can be improved by the wavelet approach. In a similar fashion, wavelet packets can also be used to discriminate between the awake and the anes-thetized state. It has been found that the discriminating power can be increased by performing the WPs analysis using Dbl4 and the bands Baaci3 or Bdad3, see Fig. 2.8. However, the improvement is small while the computational cost associated is very high. Hence, the optimization of the technique using WPs does not seem necessarily useful for this particular application. However, in the cases where a narrower frequency focus is needed, WPs are offering a better solution. Thus, we conclude that the probability density function of the first detail band obtained using the SWT Daubechies 14 filter is the best suited for differentiating the awake state from the anesthetized state. This result has the additional advantage of keeping the necessary decomposition level to 1, which makes the extraction of the feature rather fast. 3.3 Imp lementa t ion This section discusses the implementation of the technique presented in Section 3.2. Figure 3.4: Best wavelet selection for Daubechies filters: (a) discrimination parameter D; (b) mean of the variabil-ity v a; (c) mean of the variability vm; and (d) variance of the variabilities v a and vw for the band Bdx • CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 56 Figure 3.5: Daubechies 14 decomposition - PDFs of the detail coefficients averaged over the awake and anesthetized training data sets: (a) Bd1; (b) Bd2(c) Bd3; and (d) Bd4 • CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 57 3.3.1 Functional Block Diagram Figure 3.6 shows the functional block diagram of the complete system for monitoring the patient's hypnotic state in real time. As previously mentioned, raw EEG data are acquired using the capabilities of the BIS monitor and transferred via the interface to the laptop computer for real-time analysis. The E E G signal is first sensed using the standard BIS sensor (1) placed on the patient's forehead4, and comprising of 3 electrodes forming a single differential channel5. This signal is the input of an amplifier and, further, the Analog/Digital Converter (ADC) unit (2) within the BIS monitor6. This unit amplifies, pre-filters and digitizes the signal. Hence, the very low frequency components (e.g. < 0.5 Hz) as well as higher frequencies (e.g. > 100 Hz) are filtered out. Also, the eventual electromagnetic interferences due to the mains (typically 50 Hz or 60 Hz) are rejected using a notch filter. Epochs of fixed 1-second duration are digitized by the ADC, and acquired at the fixed sampling rate fs of 128 Hz. Digitized epochs containing N = fs samples are then transferred on a per second basis to the Digital Signal Processing (DSP) unit (3) implemented on the laptop computer, where the WAV Index is calculated in real time by means of the afore described wavelet analysis based method. The resulting index is further displayed by the display unit (4) using any standard display device (Cathode Ray Tube (CRT), Liquid Crystal Display (LCD), printer, etc...). In our realization, the display unit is the LCD screen of the laptop computer. The index is displayed both as a trend graph and as a number, and is updated every second. Note that this prototype can conveniently be implemented on a Personal Digital Assistant (PDA) or a handheld PC. The core of our system (i.e. the wavelet analysis based algorithm for calculating the index of hypnosis -the WAV Index) is represented by the DSP unit (3). Its description is the subject of the following subsection. 3.3.2 Analysis Algorithm All parts of the DSP unit, that is the functional blocks of the analysis algorithm, were first developed using standard Matlab functions. The Matlab implementation was used for the off-line calculations of the WAV Index. For the on-line analysis (i.e. real time), the Matlab realization was translated into C functions, which were further integrated into the interface program. Functional blocks of the analysis algorithm are described in detail in the following. 4 T h e forehead is the most convenient region for sensor placement. Also , the frontal and prefrontal lobes (which are at the origin of higher cognitive functions) are located directly behind the forehead. 5 T w o electrodes correspond to the F p z (International 10-20 System) and left mastoid placement, w i t h the common ground electrode placed on F p l . 6 T h e da ta acquisit ion can be performed using other available E E G amplifiers and DSpace platform as an A D C unit. However, the B I S monitor is preferred since it is approved for the use in the O R . CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 58 Figure 3.6: Functional block diagram of the system for estimating the hypnotic state based on wavelet analysis. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 59 Pre-processing Unit Once an epoch is acquired, it is sent to the pre-processing unit (5), see Figure 3.7. It is first stored as a vector x (14) of length N = 128. The mean value x = Yt,k=ix(k) 1S removed (15) to detrend the signal from the offset due to the acquisition process (this offset does not correspond to any physiological activity since the E E G signal is of zero mean). The RMS amplitude (16) of the epoch is then calculated and epochs with amplitudes greater than some maximum value (e.g. 200 uV) and less than some minimum value (e.g. 2 fiV) are rejected7. It is indeed assumed that they contain either artifacts such as ocular and electrocautery artifacts, or isoelectric E E G (i.e. the E E G signal of near to zero amplitude). Chapter 4 discloses and investigates a pre-conditioning technique that removes such artifacts from the E E G (15'). The isoelectric E E G corresponds either to periods of burst suppression patterns or cortical silence. Since these epochs contain noise that might lead to erroneous values of the index, we are excluding them from the analysis, and assign the value 0 to the corresponding index. If the RMS amplitude is within the two bounds (17), a flag (20) takes the value 1, indicating that the current epoch is not corrupted. The algorithm then proceeds to the next stage, that is, the wavelet analyzer unit (6). If an artifact is present (18), it is indicated by the flag value of 0, and the algorithm proceeds to the filtering unit (10). If an isoelectric E E G is detected (19), it is indicative that the patient is in the deepest level of hypnosis. Hence the flag takes the value 1 and the variable W AVun/utered (21), containing the unfiltered WAV Index, takes the value of 0. The algorithm then proceeds to the filtering unit (10). Wavelet Analyzer Unit After the pre-processing stage, the current epoch x (which we assume does not contain any large artifacts), represents the input of the wavelet analyzer unit (6). The wavelet analyzer function calculates the wavelet coefficients in the decomposition band Bd,, which was selected in the previous section as the most discriminating band. This is achieved by applying the SWT and the wavelet filter Dbl4 to the E E G epoch normalized by its RMS value. Further, it is important to note that for obtaining the coefficients in , we actually need to perform only 1-level SWT decomposition. Therefore, the wavelet coefficients are obtained by convolution of the E E G epoch with the wavelet filter Dbl4. Figure 3.8 details the wavelet analyzer unit. In order to compensate for the boundary distortion when calculating the convolution, the current epoch is first extended for I — 1 samples, where I is the length of i i D b i 4 filter (i.e. I = 28). If the previous epoch8 has not been corrupted (22), the current epoch x is extended 7 N o t e that, because of its small amplitude, the E E G signal is very susceptible to various artifacts and noise. T h e detec-t ion/ removal of specific artifacts or interferences is a very complex task. For the t ime being, we w i l l only consider the rejection of the most apparent artifacts characterized by inadmissible amplitudes. Thus , the R M S value of the epoch is ut i l ized as a simple tool for artifact detection. 8 T h e wavelet analyzer unit has access to the stored previous epoch and the corresponding value of flag. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION EEG SIGNAL A M P L I F I E R & ADC UNIT PRE-PROCESSING UNIT /x = [x(\)...x(k)...x(N)],/-~ 14 N = 128 Detrending: X = x - X 15 - 1 5 ' P r e - c o n d i t i o n i n g [..-> , ' 1 16 19 <2jUV 18 >200/iV <200uV and >2pV • • fl^S current 1 r 20 c , - 2 1 WAV„nflllered=0 A flag current 20 flagcurren, = 0 FILTERING UNIT WAVELET ANALYZER UNIT FILTERING UNIT 10 10 Figure 3.7: Pre-processing of the EEG epochs. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION PRE-PROCESSING UNIT 22 WAVELET ANALYZER UNIT 0 , - - 24 ext = [xpmkHB(N-l + l) ... xprevmus(N)] / = length(hM 1 4) = 28 ^t = [-xcurre„,(N) ... -xcurre„,(N-l + \)] / = length(h o m ) = 28 23 Epoch extension: x -> x « , = [ext x] 26 27 -28 Normalization: x„, = -RMS Convolution: ; Cal («) = X hDMA (" - *) „,(*) , « = l,...,A^ + 2-(/-l) 29 C(n) = C e B (« + / - l ) J n = l,...,W 31 pdf = — histogram(C, bins) N bins = 100 25 COMPARATOR UNIT Figure 3.8: Wavelet analyzer unit. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 62 (26) by samples from the previous epoch X p r e v i o u s (23). However, note that if the previous epoch has been corrupted (24), then the current epoch is anti-symmetrically extended (25, 26). The RMS value of the extended epoch x e x t is then calculated (27) and the epoch is normalized (28). By performing the convolution (29) of the extended normalized epoch 5text with the filter hobiA, we obtain the vector Cext- Further, only the undistorted portion of the vector Cext is kept (30) (i.e. N = 128 central samples) and stored in the vector C. The vector pdf storing the probability density function of the wavelet coefficients in B^ band is then obtained. This is achieved by calculating the histogram (31) of the elements of the vector C, and then dividing the resulting vector by N. The histogram is calculated for 100 bins. The number of bins has been empirically chosen. A small number of bins (e.g. 10) does not give a sufficient resolution around zero to differentiate deeper hypnotic levels, thus yielding an "on-off" index. However, a large number of bins (e.g. 1000) increases the computational stress, and results in a noisier index. Comparator The resulting pdf vector is input into comparator unit (7), see Fig. 3.9. This unit compares the pdf vector of a current epoch (11) with two template vectors pdfw (12) and pdfa (13) representing the two known hypnotic states - awake and anesthetized - and stored as look-up tables in the computer memory. The awake template pdf^ , is derived from the awake training data set, which comprises a number of E E G epochs obtained from a group of healthy awake subjects (i.e. population norming). Another possibility for the awake data set is to record the patient's E E G while the patient is still awake (i.e. self-norming). The awake template is the averaged feature f^ extracted from the awake training data (see Section 3.2.2), that is, the PDF of the SWT wavelet coefficients in B^ band averaged over the awake training data set. Each epoch of the awake data set is extended and normalized as in steps 23-28 in Fig. 3.8. Then, the wavelet coefficients of B^ band are obtained by the convolution (29) of the normalized epoch and hrjbi4 filter. The PDF of the epoch is then obtained by calculating the histogram of uncorrupted convolution coefficients for 100 bins (30, 31). Finally, the PDFs of all epochs in the awake training data set are averaged yielding the awake template pdf^. This analysis is performed off-line. Once extracted from the awake training data set, the template pdf^ can be then stored on a mass storage device for future real time comparisons. The anesthetized template pdfa is the PDF of the wavelet coefficients of an isoelectric signal, which corresponds to the deepest level of hypnosis. All coefficients are equal to 0. Hence, this particular PDF is a Dirac function centered at the origin. Therefore, the vector pdfa comprise 100 bins, all containing the zero value, except the central one that contains the value 1. Figure 3.10 shows the awake and anesthetized templates corresponding to the awake baseline and iso-electric E E G , respectively, as well as the PDF corresponding to an intermediate anesthetized state. For example, in Section 3.2.3, we have used such an intermediate PDF (see Fig. 3.5.a) as a representative of the anesthetized state for the best wavelet selection analysis. However, as it will be illustrated later in Sec-tion 3.4, the use of this anesthetized template results in an index that is unable to asses very deep hypnotic CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION u> ro in ro D CC O H < CM S o o 3 2 z 05 U N >H < z < EH H Cd > s Figure 3.9: Comparator unit. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 64 1.0 0.6 \-0- 0.4 Deeply Anesthetized (isoelectric HfiG) 0.2 h - 2 - 1 0 1 2 Amplitude of Wavelet Coefficients Figure 3.10: Awake and anesthetized templates. levels. This is understandable, since this template was derived from the anesthetized training data set Ta in Section 3.2.2 corresponding to the hypnotic levels of 30 < BIS < 35. Deeper levels are unreachable for such an index. The comparison (32) between the pdf (11), which is calculated for the current epoch x in the wavelet analyzer unit (6), and the two templates pdf w (12) and pdf a (13) is achieved using the L l distance metric. This comparison yields two values iw (32) and ia (33) calculated as: Tt • Efc=i \Pdfk ~ pdfw,k|, and N Efe=i \pdfk - pdfa,k\, (3.9) where pdfk, pdfWtk, and pdfa^ denote the k t h elements of the vectors pdf, pdf^, and pdf a respectively. Note again that, instead of the L l norm, any norm of any order can also be used to compare the feature function against both templates. However, since the higher order norms emphasize large differences among the elements of vectors compared, they would lead to a noisier index. An index i (36) is then generated by calculating (35) the difference between iw (33) and ia (34) as: % = (3.10) The output of the comparator unit is then input to the scaling unit (9). Scaling Unit The index i is scaled in order to take values between 0% (corresponding to an isoelectric signal) and approxi-mately 100% (corresponding to the awake baseline) with higher values indicating higher level of consciousness or hypnosis: i = i • scale + of f set, (3-11) where scale and offset are two constants calculated in the off-line analysis, based on the awake and the anesthetized training data sets. Thus, the scale value is computed so that the mean value of index calculated CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 65 for the awake training data set corresponds to « 98%, while the isoelectric E E G gives 0%. The result of the scaling process is further stored into the variable WAVunfiltered (21). Filtering Unit The variable WAVunfiltered (21) contains the unfiltered version of the final WAV Index calculated for the current epoch x. The random character of the E E G dictates that in order to extract a meaningful trend of the patient's hypnotic state it is necessary to smooth this variable using a low-pass filter. A new value WAVunfutered is delivered by the scaling unit (9) for every epoch (i.e. every second). The filtering unit has access to 30 previous values of WAVunfutered and flag stored in the computer memory, including the current ones. However, note that if the current epoch is corrupted with an artifact (i.e. flagCUrrent = 0 (20)), the variable WAVunfutered c a n take any arbitrary value as it will not be used to derive the final value of the index. In this case, the variable WAVunfutered is averaged over the past 30 seconds of data. The result of the averaging filter is stored in the variable. WAV. However, when calculating the average, only uncorrupted epochs are taken into account (by investigating the corresponding flag variable). Furthermore, in order to account for poor signal quality, if more than a certain number of epochs during the last 30 seconds are corrupted (e.g. 20), the monitor is unable to give an accurate estimate of the patient's hypnotic state. In that case, the variable WAV takes the negative value -100%. The output variable WAV of the averaging filter is then sent to the display unit (6 ) . In case of poor signal quality, a message indicating the presence of numerous artifacts is sent to be also displayed by the display unit. 3.4 Preliminary Results In this Section, preliminary results corresponding to the case study presented in Section 3.1 are provided. Note that since the best wavelet analysis showed that the Daubechies 14 wavelet family gives the optimal discriminating power, the analysis and results reported are obtained using this wavelet filter. Figure 3.11.a shows the unfiltered WAV Index derived based on the awake and anesthetized templates in Fig. 3.5.a. It can be observed that this index is rather noisy and thus non-usable for monitoring purposes. Figure 3.11.b shows the filtered WAV Index based on 30-second averaging using the same templates9. This index was also scaled so that it can be compared with the BIS. This procedure can rise questions regarding the benefits of such a scaling. In short, we prefer to use the scale similar to the BIS, since this choice is already familiar to practitioners. The WAV Index captures most of the significant events during the surgical/anesthesia procedure. However, there is a significant difference between indexes in representation of the event 2. Figure 3.11.b reveals the inability of the index to take values below 32 (detail within the circle). This WAV Index fails to quantify deep anesthetic state, which is a serious drawback. However, as 9 N o t e that this first version of the W A V Index has already been published [150]. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 66 Time (hours) Figure 3.11: WAV Index derived based on the templates in Fig. 3.5.a: (a) unfiltered index; and (b) filtered (30 s) and scaled index (WAV - black; BIS - gray). argued previously, this is expected since we derived the anesthetized template based on the hypnotic levels of 30 < BIS < 35, thus limiting the deepest observable hypnotic levels to this range of values. Also, there is some disagreement between indexes in representing the event 1. Figure 3.12 shows the filtered WAV Index based on 30-second averaging using the isoelectric anesthetized template (see Fig. 3.10). This index captures all the significant events during the procedure. Also, the use of the isoelectric anesthetized template enables to capture the deep hypnotic levels significantly better. Furthermore, a strong agreement between the WAV and the bispectral index is obvious. The representations of events 4 and 5 by these two indexes are very close to each other. Some deviation between indexes is evident as well, especially in the intermediate states. The electrocautery artifacts are captured by the two indexes with slight differences. During this period of time, both indexes do not accurately represent the hypnotic state of patient (notified by negative values). Of course, due to different processing techniques, these negative values do not occur at the same time for both indexes. The WAV Index takes values below zero when disturbances and/or artifacts are presented in more than 20 out of 30 consecutive seconds. Also, there is a significant difference between indexes in the representation of the first event (light sedation). This is mainly the result of the presence of ocular and movement artifacts since the patient is still awake at this point. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 67 1 0 0 00 Time (hours) Figure 3.12: WAV Index: comparison with the bispectral index (WAV - black; BIS - gray). Finally, by examining the details in Figs. 3.13 and 3.14, we can observe a significant lead (in the range of 20-25 s) of the WAV Index in comparison to the BIS. This lead indicates that the response of the WAV Index to fast and large changes in the patient's state is faster. This property might be the major advantage of this technique over the BIS analysis. 3.5 Summary In this chapter, we clearly established the potential of wavelet analysis for the extraction of information relative to the hypnotic state of a patient undergoing anesthesia. It has been shown that the PDF of the wavelet coefficients corresponding to the 7-band of normalized E E G epochs exhibit significant difference between different hypnotic states. Based on a remarkably limited amount of data, a methodology that uses the large discrepancy which exists between the PDF of a normal fully awake subject and the PDF of a deeply sedated patient has been presented. These two PDFs were further used as templates against which the PDF derived from an unknown E E G epoch could be properly compared to. The resulting index, so called WAV Index, mirrors very closely the bispectral index. Also, the simplicity of the proposed wavelet-based technique as compared to the bispectral algorithm, resulted in a significant time lead in favor of the WAV Index. This work also carries significant medical interest since the two analyses of the E E G signal based on totally different digital signal processing techniques have produced similar results. CHAPTER 3. WAV INDEX: METHODOLOGY AND IMPLEMENTATION 1.4 1.5 1.6 1.7 Time (hours) Figure 3.13: Zoomed detail: WAV Index (WAV - black; BIS - gray). 100r 9 0 -| J L l I ! I ! I I I I I 0 60 120 180 240 300 360 420 480 540 600 Time (seconds) Figure 3.14: Zoomed detail: WAV Index (WAV - black; BIS - gray). C h a p t e r 4 Pre-Conditioning of the EEG Signal This chapter investigates a wavelet based de-noising technique of the electroencephalogram to correct for the presence of artifacts, and in particular, Ocular Artifacts (OAs). Its application for the enhancement of the WAV Index while patient being awake is also presented. The proposed technique is based on an over-complete wavelet expansion of the E E G using the SWT followed by the thresholding of the coefficients in the lower frequency bands. The reader will be first introduced to the problems caused by OAs for the analysis and interpretation of the E E G signal as well as the prior attempts for their removal. The influence of various artifacts on the WAV Index is illustrated. Because they are the most severe type of artifacts, OAs are thoroughly investigated in Section 4.2. The wavelet-based de-noising method for a single frontal channel is then discussed in Section 4.3. The results are further presented in Section 4.4, demonstrating the potential of the technique for a successful artifact correction. No reference EOG channel is needed. The same approach works both for the eye blink and eye movement artifacts, and, further, for head movements. Finally, the application of the wavelet-based de-noising to the WAV Index is investigated in Section 4.5, along with ways to integrate the method into the real time calculation of the index. The results obtained are also illustrated. Note that the technique presented here can be applied as a separate entity, independent of the WAV Index application. This work has already been the subject of a publication [151]. 4.1 Introduction The electroencephalogram gives researchers a non-invasive insight into the intricacy of the human brain. It is a valuable tool for clinicians in numerous applications, from the diagnosis of neurological disorders, to the clinical monitoring of depth of anesthesia. For an awake healthy subject, the normal E E G amplitude is in the order of 20-50 /xV. This makes the E E G signal very susceptible to various artifacts, causing problems for its analysis and interpretation. In current data acquisition, eye movement and blink related artifacts are often dominant over other electrophysiological contaminating signals (e.g. heart and muscle activity, head and body movements), as well as external interferences due to power sources. Eye movements and blinks 69 CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 70 produce a large electrical signal around the eyes (in the order of mV), known as Electrooculogram (EOG), which spreads across the scalp and contaminates the E E G . These contaminating potentials are commonly referred to as ocular artifact (OA). Figure 4.1 shows the WAV Index calculated based on the E E G acquired at the beginning of a bilateral breast reduction surgery without applied artifact removal. As can be observed, prior to induction, this unmodified WAV Index goes as low as around 55, which corresponds to the hypnotic levels of general anesthesia. This is obviously an erroneous result since, at this point, the patient was awake and free of drugs. The visual inspection of the corresponding E E G signal reveals that it is heavily corrupted by eye blinks and eye movements. However, as can be observed, the bispectral index is not influenced by these artifacts, since the BIS device utilizes an artifact recognition/rejection scheme1. Figure 4.2 depicts 15 s of the E E G signal when the unmodified WAV Index takes one of its minimum values prior to induction. Ann I I Induction 12h57 12h58 12h59 13h00 13h01 13h02 13h03 13h04 13h05 13h06 13h07 13h08 13h09 13h10 13h11 13h12 Time Figure 4.1: WAV Index: the effect of ocular artifacts. 250 r -1501 ' ' ' ' ' ' ' ' ' ' ' ' ' 1 ' 13:02:00 13:02:05 13:02:10 13:02:15 Time Figure 4.2: EEG corrupted with ocular artifacts. The BIS monitor discards corrupted epochs from the derivation of the index [152, 153]. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 71 Therefore, in order for the WAV Index to have any practical value, it is essential to address the issue of artifacts and their removal. In order to assess the influence of artifacts in a more controlled way, we have performed an experiment to record various artifacts and their influence on the WAV and BIS indexes. Thus, the E E G data was acquired from a healthy awake subject. A single frontal channel was recorded. The baseline E E G and eight different artifacts were recorded in the following fashion. For each artifact, the subject was first instructed to keep his eyes shut and still. Sixty seconds of presumably uncontaminated baseline E E G was thus recorded. Then, for the next 60 seconds, the subject was instructed to produce the artifacts, for example, to blink or move his eyes in a predetermined fashion. Finally, another resting period of 60 seconds with no artifacts was recorded. Five ocular artifacts were recorded in this fashion: slow blinks (1 blink per 2 seconds), fast blinks (2 blinks per second), vertical, horizontal and round eye movements, as well as head movements, teeth clenching and chewing. Figure 4.3 shows the influence of various artifacts on the WAV and BIS indexes. As can be observed, the indexes react differently in the presence of artifacts. The WAV Index is extremely sensitive to the OAs and head movements. However, since the head movements are easy to avoid, we are going to concentrate our efforts to the OAs removal. Note that the BIS is also sensitive to the presence of OAs, however, only to a lesser extent. Furthermore, the BIS is able to recognize the artifacts and display the message about their presence thus indicating unreliable results. As mentioned previously, the corrupted epochs are not taken into account when calculating the BIS. Hence, in the case of the extreme presence of artifacts, the bispectral index is not able to assess the patient state and therefore takes a negative value. However, we have noticed that, in the case of fast blinks, the BIS monitor does not recognize the presence of artifacts and does not issue a warning message. This means that, in this case, the bispectral index takes a value of around 40 and displays it as a valid one while the patient is still awake. GO CD I > 1 0 0 Time (min) Slow Blinks Fast Blinks Horizontal Vertical Eye Round Eye H e a d Teeth Chewing and Eye Movement Movement Movement Clenching Swallowing Movement Figure 4.3: Influence of various artifacts on the WAV and BIS indexes (WAV-black; BIS-gray). CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 72 The rejection of epochs contaminated with OAs usually leads to a substantial loss of data. Asking subjects not to blink or move their eyes, or to keep their eyes shut and still, is often unrealistic or inadequate (in contrast to head movements). The fact that the subject is concentrating on fulfilling these requirements might itself influence the E E G . Hence, devising a method for successful removal of OAs from E E G recordings has been (and still is) a major challenge. In the past, significant efforts have been made among researchers to find a way to successfully remove OAs. Widely used time-domain regression methods involve the subtraction of some portion of the recorded E O G from the E E G ([154], [155]). They assume that the propagation of ocular potentials is volume conducted, frequency independent and without any time delay. However, Gasser et al. in [156] argued that the scalp is not a perfect volume conductor, and thus, attenuates some frequencies more than others. Consequently, frequency-domain regression was proposed. In addition, no significant time delay was found, which was consistent with the E O G being volume conducted. In [157] it was reported that, in reality, the frequency dependence does not seem to be very pronounced, while the assumption of no measurable delay was confirmed. Thus, while some researchers support the frequency domain approach for EOG correction ([156], [158]), others dispute its advantages ([157], [159], [160]). However, neither time nor frequency regression techniques take into account the propagation of the brain signals into the recorded EOG. Thus a portion of relevant E E G signal is always cancelled out along with the artifact. Further, these techniques mainly use different correction coefficients for eye blinks versus eye movements. They also heavily depend on the regressing EOG channel. In addition, Croft and Barry [160] demonstrated that the propagation of the E O G across the scalp is constant with respect to OA types and frequencies. They proposed a more sophisticated regression method (the aligned-artifact average solution) that corrects blinks and eye movement artifacts together, and made possible the adequate correction for posterior sites [159]. They claim that the influence of the EEG-to-EOG propagation has been minimized in their method. In an attempt to overcome the problem of the EEG-to-EOG propagation, a multiple source eye correction method has been proposed by Berg and Scherg [161]. In this method, the OA was estimated based on the source eye activity rather than the E O G signal. The method involves obtaining an accurate estimate of the spatial distribution of the eye activity from calibration data, which is a rather difficult task. Due to its de-correlation efficiency, the principal component analysis (PCA) has been applied for OA removal from the multi-channel E E G and it outperformed the previously mentioned methods. However, it has been shown that PCA cannot fully separate OAs from the E E G when comparable amplitudes are encountered [162]. Recently, independent component analysis (ICA) has demonstrated a superior potential for the removal of a wide variety of artifacts from the E E G ([163], [164]), even in a case of comparable amplitudes. ICA simultaneously and linearly unmixes multi-channel scalp recordings into independent components, that are often physiologically plausible. Also, there is no need for a reference channel corresponding to each artifact source. However, ICA artifact removal is not yet fully automated and requires visual inspection of the independent components in order to decide their removal. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 73 Other attempts have been based on different adaptive signal processing techniques ([165], [166], [167], and [168]). The performance of these methods also relies on a cerebral activity to minimally contaminate the E O G reference. Finally, the E E G may contain pathological signals, which resemble OAs. Thus, these signals are most likely to be removed from the E E G recordings, leading to erroneous diagnosis. Therefore, it is often important to distinguish between artifacts and pathological E E G signals prior to artifact removal. In this case, artificial intelligence techniques prove to be somehow helpful [169]. OAs being of a main concern for artifact removal in the application of the WAV Index, we will first investigate their time and frequency characteristics in the next section. 4.2 Ocular Artifacts There are two different originating phenomena for ocular potentials ([154], [170], and [171]). There is a potential difference of about 100 mV between a positively charged cornea and negative retina of the human eye, thus forming an electrical dipole (i.e. the corneo-retinal dipole). Firstly, the rotation of the eyeball results in changes of the electrical field across the skull. Secondly, eye blinks are usually not associated with ocular rotation; however, the eyelids pick up the positive potential as they slide over the cornea. This creates an electrical field that is also propagated through the skull. Hence, ocular potentials spread across the scalp and superimpose on the E E G . The mechanism of origin (eye movements versus eye blinks) and the direction of eye movements determine the resulting E O G waveshape. Vertical, horizontal and round eye movements usually result in square-shaped E O G waveforms, while blinks are spike-like waves. Ocular artifacts decrease rapidly with the distance from the eyes [170]. Therefore, the most severe interference occurs in the E E G recorded by the electrodes placed on the patient's forehead. Yet, this is the most convenient region for their placement. The frontal and prefrontal lobes, which are at the origin of higher cognitive functions, are located directly behind the forehead. Therefore, the task of E O G correction for frontal channels is challenging. As mentioned in the previous section, we have acquired E E G data from a healthy awake subject that comprise five different ocular artifacts (i.e. slow blinks, fast blinks, and vertical, horizontal and round eye movements) as well as the resting, presumably artifact-free, periods. The waveforms in Fig. 4.4 present 10-second signals corresponding to the uncontaminated baseline E E G and the E E G contaminated with different artifacts. The corresponding power spectra are presented in Fig. 4.5. Besides the average power spectrum of the baseline E E G signals, see Fig. 4.5.a, individual power spectra of the E E G during resting (i.e. the baseline) periods are presented as well. Figures 4.4 and 4.5 clearly shows that OAs are large, transient slow waves. They occupy the lower frequency range; from 0 Hz up to 6-7 Hz for the eye movement artifacts, and typically up to the alpha band (8-13 Hz), excluding very low frequencies, for the eye blinks. This is a well-known and documented result [156], which our experiments proved consistent with. Clearly, OA amplitudes are of a much higher order CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 74 150 > 100 (a) 150 > 100 E < (b) 150 100 H 50 E < 0 128 256 384 512 640 768 896 102411521280 Time (samples) (c) E < E < -300 300 150 -150 -300 (e) 0 128 256 384 512 640 768 896 102411521280 Time (samples) (0 Figure 4.4: Uncontaminated baseline EEG and various artifacts: (a) uncontaminated baseline EEG; (b) EEG con-taminated with slow blink artifact; (c) EEG contaminated with fast blink artifact; (d) EEG contaminated with vertical eye movement artifact; (e) EEG contaminated with horizontal eye movement artifact; and (f) EEG contaminated with round eye movement artifact. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 75 Figure 4.5: Power spectra of the uncontaminated baseline EEG and various artifacts: (a) averaged uncontaminated baseline EEG and individual resting periods; (b) EEG contaminated with slow blink artifact; (c) EEG contaminated with fast blink artifact; (d) EEG contaminated with vertical eye movement artifact; (e) EEG contaminated with horizontal eye movement artifact; and (f) EEG contaminated with round eye movement artifact. than those of the uncontaminated E E G and have a characteristic pattern of changes. Vertical eye movement artifacts (Fig. 4.5.d) also seem to produce a rise in the higher frequencies. However, this is most likely due to the increased muscle activity, and it is also present to a lesser degree for the other two eye movements (horizontal - Fig. 4.5.e, and round - Fig. 4.5.f). 4.3 Wavelet-Based De-Noising As previously mentioned, the EOG is the non-cortical activity that contaminates the E E G recordings. Thus, since the brain and eye activities have physiologically separate sources, we will assume that the recorded CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL E E G is a superposition of the true E E G and some portion of the EOG signal. Thus, we have: 76 EEGredt) = EEGtrue(t) + k • EOGit) + DCoffset (4.1) where EEGrec is the recorded contaminated EEG, EEGtrUe is due to the cortical activity, and k • EOG(t) is the propagated ocular artifact at the recording site. The DC0ffset takes into account the zero mean value of the cortical E E G , since this might not be true for the recorded E E G due to the process of data acquisition. We are interested in estimating the ocular artifact based on the analysis of the recorded E E G . By subtracting it from the contaminated EEG, we will then obtain a corrected E E G , which minimizes the effect of the ocular artifact and gives an appropriate representation of the true cortical E E G . The true E E G is a noise-like signal. We cannot observe any clear patterns within it, nor can we simply correlate the particular underlying events with its waveshape [32]. Furthermore, in the awake, conscious state, neurons are firing in a more independent fashion. As a result of this desynchronization, the awake E E G signal is even more random-appearing. The E O G removal can be approached by recovering the regression function k • EOG(t) from the recorded (i.e. noisy) data. For this purpose, in the last decade, wavelet thresholding has emerged as a simple, yet effective technique for de-noising [172]. 4.3.1 Wavelet Thresholding The main statistical application of wavelet thresholding is a nonparametric estimation of the regression function / , based on observations at time points ij. The Si observations are modelled as: Si = f(U) + eu i = l,2,...,N (AT = 2n), (4.2) where £j are independent and identically distributed N(0, a2) random variables at equally spaced time points U. Due to the orthogonality of the wavelet transform, we are allowed to perform filtering in the space of wavelet coefficients. The procedure for suppressing the noise involves: i) finding the coefficients of the wavelet transform of {SJ}; ii) comparing each wavelet coefficient against an appropriate threshold; iii) keeping only those coefficients larger than a threshold; and iv) applying an inverse wavelet transform to obtain an estimate of / . The assumption is that large coefficients kept after thresholding belong to the function to be estimated, and those discarded belong to the noise. This is a fair assumption due to the good energy compaction of the wavelet transform. Some of the function coefficients might eventually be discarded if they are of the same level as the noise coefficients. Thus, this de-noising technique works well for functions whose wavelet transform results in only a few nonzero wavelet coefficients, such as, for example, functions that are smooth almost everywhere, except for only a few abrupt changes [173]. Special care has to be taken when choosing an appropriate threshold, which always involves the estimation of the noise variance a2 based on the data. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 77 4.3.2 Wavelet Estimators and the Choice of Wavelet Transform As discussed in Chapter 2, the discrete wavelet transform is not translation-invariant (i.e. in general, if we apply a DWT to a shifted version of a signal x, we do not get the shifted version of the DWT of a;). Due to this drawback, de-noising with standard DWT often suffers from artificial additional artifacts such as the ringing effects in the vicinity of discontinuity, depending on its actual location [174]. This produces wavelet estimators of irregular visual appearance. In practice, we obtain the DWT of a sequence {x^} by applying a pair of low-pass and high-pass filters, which must comply with certain conditions (see Section 2.1.2). Then, the resulting sequences are decimated (i.e. only every even member of a sequence is kept). By feeding down the low-pass sequence to the next level and repeating the filtering and decimation, we obtain the approximation and the detail of the sequence {x^ at that level. This procedure is repeated until we reach the desired level of decomposition. The decimation steps can be equally carried out by selecting every odd member of each sequence instead of even ones. Furthermore, we can choose to select even members at some levels and the odd members at others. Obviously, the final result will be different, but the orthogonality of the transformation is kept, hence the process can be easily inverted to obtain the perfect reconstruction. Different wavelet transforms corresponding to various selections of decimation steps are referred to as the e-decimated discrete wavelet transforms, and they are all shifted versions of the ordinary DWT applied to the shifted sequence {x^ [133]. Any particular e-decimated DWT defines a particular set of wavelet bases and their time positions, and consequently, the grid of integers for each level at which the wavelet coefficients are localized. Various misalignments between the features in the signal and those of wavelet basis are leading to more or less pronounced artifacts when de-noising. To minimize these artifacts, we apply e-decimated DWTs with different shifts, followed by the averaging over the obtained results ([174], [150]). The stationary wavelet transform of sequence {xi} is, actually, equivalent to applying each possible e-decimated DWT, and then averaging over the results [133]. In practice, the SWT is easily obtained in a similar fashion as the DWT, except that the decimation step is not performed (refer to Section 2.1.3). This leads to the redundant representation of the original signal, with great potential for statistical applications. The approximation and detail sequences at each level L are of the same length as the original sequence, with the complexity increased from 0(2L) to 0(L2L). 4.3.3 Methodology Since OAs occupy lower frequency bands and have significantly larger amplitudes than the noise-like awake E E G , a multi-level wavelet decomposition of EEGrec allows the detection of the presence of artifacts, as they generate much larger coefficients. Thresholding these coefficients (i.e. setting them to zero), and then recomposing the signal will thus correct the E E G . We will use the stationary wavelet transform and its inverse, since it has better sampling rates in the lower frequency bands, hence leading to smoother results. The analysis is based on a 1-second epoch of E E G signal (128 samples). To overcome boundary effects, epochs are extended on both ends, with the samples from the previous epoch at the beginning and flipped CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 78 samples of the current epoch at the end. The length of the epoch extensions has to be greater or equal to the wavelet filter length. For acceptable computational complexity, the analysis has been carried out by performing a 5 level decomposition (frequency bands 0-2 Hz, 2-4 Hz, 4-8 Hz, 8-16 Hz, 16-32 Hz and 32-64 Hz). A Coiflet 3 wavelet filter has been chosen, since the shape of its mother wavelet resembles the shape of the eye blink artifact. This maximizes the amplitude of coefficients corresponding to the eye blink artifacts in the lowest band of decomposition. Furthermore, this choice minimizes the spread of artifacts to higher frequency bands. The filter length keeps the computational complexity low. It has turned out that it works properly for the eye movement artifacts as well. Wavelet coefficients have been thresholded only in the lower frequency bands (i.e. up to 16 Hz). The thresholding procedure sets all coefficients larger than a threshold to zero. This one step procedure is equivalent to estimating the OA with standard wavelet de-noising technique, and then subtracting it from the corrupted E E G . Although in overall the E E G signal is non-stationary, we can assume its stationarity during the period of one epoch. Furthermore, it is true that the variance of the awake E E G signal can significantly change from one second to another, but this variation is small in comparison with the variance of the epochs corrupted by artifacts. Hence, the threshold could be estimated by simple statistical analysis of the baseline E E G , which is presumably artifact-free. Thus, the acquired 60 seconds of the baseline E E G was analyzed and each second was decomposed by 5 level SWT with Coiflet 3 wavelet filter. For each second, the maximum absolute value Mk of wavelet coefficients has been calculated for each band k of decomposition below 16 Hz. The threshold Tk for the band k has been calculated as: Tk=M + 2-std{Mk), (4.3) where Mk and std(Mk) are, respectively, the mean and standard deviation of {Mk}k-i,2,...,60- The proposed methodology enables real time E E G correction in the presence of OAs. 4.4 Results Figure 4.6 shows the results of the de-noising for a single epoch containing the fast blink artifact, while Fig. 4.7 depicts the de-noising of a single epoch containing the horizontal and round eye movement artifacts, respectively. The contaminated and corrected E E G epochs and their power spectra are shown, together with the removed artifacts. In addition, the stationary wavelet coefficients of the contaminated E E G epochs are presented along with the applied thresholds for details at levels 3, 4 and 5 (frequency bands 18-16 Hz, 4-8 Hz and 2-4 Hz), and the approximation at level 5 (frequency band 0-2 Hz). Note the exact preservation of the high frequency content of the original signal, both in amplitude and phase. This is very important in the applications that are exploring various aspects of consciousness, such as the WAV Index, since the 7-band E E G oscillations (40-60 Hz) are preserved. Note also how the DC bias is automatically removed, which is the result of the detrending carried out prior to the analysis. Figure 4.8 shows 10 seconds of contaminated and corrected E E G for the slow blink artifacts as well as the CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL EEG with artrfacl Corrected EEG Removed artifact > 13 Detail coefficients at level 1 at level 2 Detail coefficients I Detail coefficients at level 5 Approximaton coefficients at level 5 16 32 48 64 80 Time (samples) (b) 112 121 EEG with artifact Corrected EEG Frequency (Hz) (c) Figure 4.6: Fast eye blink artifact: (a) contaminated EEG and corrected EEG; (b) stationary wavelet decompositi of contaminated EEG; and (c) power spectra of contaminated EEG and corrected EEG. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 80 Figure 4.7: Horizontal (left column) and round (right column) eye movement artifacts: (a) contaminated EEG and corrected EEG; (b) stationary wavelet decomposition of contaminated EEG; and (c) power spectra of contaminated EEG and corrected EEG. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 81 Slow Eye Blinks -10 s Vertical Eye Movements - 5 s 4 * v1 i t \t > 50 s ° i 256 384 512 640 768 Time (samples) 1024 1152 1280 V (a) Time (samples) Figure 4.8: Slow eye blink (left column) and vertical eye movement (right column) artifacts: (a) contaminated EEG and corrected EEG; (b) power spectra of contaminated EEG and corrected EEG. de-noising result obtained from a 5-second E E G contaminated with the vertical eye movements. Since the per second analysis is based on extended epochs, we are able to cancel out the boundary distortion originating from the convolution of the wavelet filter with the sampled data, thus preserving smooth transitions between the epochs. Figure 4.9: Head movement artifacts: contaminated EEG and corrected EEG; and (insert) power spectra of con-taminated EEG and corrected EEG. CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 82 In addition, since the head movement artifacts resemble somewhat the vertical eye movements (see Figs. 4.8.a and 4.9.a), we have applied the same de-noising technique to these artifacts. Results are shown in Fig. 4.9.b. Finally, note that the signals presented in this section correspond to the eye blinks and movements of extremely high rate (except for the slow blinks), which are rarely encountered in reality. However, these signals were used in order to emphasize the power of the proposed method even in the case of extreme artifacts. 4.5 Application to the WAV Index We have seen that OAs and head movements can make the WAV Index practically unusable in the case of lighter hypnotic levels (awake and light sedation) when the patient is still conscious and can exhibit ocular activity. This activity mainly causes problems at the beginning phase of procedures when patients are typically anxious2. Hence, the E E G signal acquired from the patient prior to induction is subject to the proposed wavelet de-noising technique. Indeed, since the technique is based on the thresholding of wavelet coefficients and thresholds derived from a fully awake subject, we can not apply it for deeper hypnotic states. In this case, the amplitudes of coefficients corresponding to artifact-free E E G are comparable to those due to artifacts in an awake state3. Thus, de-noising would result in thresholding of the uncorrupted coefficients and thus in an incorrect index. Therefore, we have to make sure that the de-noising algorithm is applied only for lighter hypnotic levels. The moment of induction of anesthesia, corresponding to the administration of a large induction dose of an anesthetic, appears suitable for switching off the de-noising technique as the patient rapidly gains unconsciousness. From the perspective of automation in clinical anesthesia (which is the long-term project of our research group), this solution is viable since the computer will be commanding directly the infusion device. Figure 4.10 shows the effect of wavelet-based artifact removal method applied to the E E G corresponding to the WAV Index presented in Figure 4.1. The de-noising technique was applied off-line to the recorded E E G and prior to the induction, and the WAV Index was then calculated. For the derivation of thresholds, the 60-second uncontaminated baseline E E G acquired in Section 4.2 was used. Figure 4.10 shows the WAV Index prior to and after the de-noising as well as the BIS. The de-noising is obviously successful, although the thresholds were derived from a different patient. However, in the ideal case, we would acquire the artifact-free baseline from a respective patient prior to the surgical procedure (this needs patient's cooperation, which might not always be possible). 2 N o t e that dur ing emergence when the patient regains consciousness, these artifacts are not extensively present. In our opinion, it is due to patients being commonly rather t ranqui l at this phase. 3 N o t e that even though the de-noising does not incorporate the artifact recognition scheme, the uncorrupted epochs in the awake state should not be affected (i.e. thresholded). CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 83 12h57 12h58 12h59 13h00 13h01 13h02 13h03 13h04 13h05 13h06 13h07 13h08 13h09 13h10 13h11 13h12 Time Figure 4.10: WAV Index: the effect of ocular artifacts and de-noising (WAV prior to de-noising - thin black; WAV after de-noising - thick black; BIS - gray). Figure 4.11 shows the WAV Index from Fig. 3.12 after de-noising. Note that the WAV Index after de-noising is able to capture the sedated state (i.e. BIS and WAV values of « 80). i Time (min) Fi gure 4.11: WAV Index: the effect of ocular artifacts and de-noising (WAV prior to de-noising - thin black; WAV after de-noising - thick black; BIS - gray). 4.6 Summary Ocular artifact correction can be a challenging task. A multitude of techniques have been proposed in the literature. However, assumptions concerning the propagation of these artifacts across the scalp are still actively discussed. There is no general consensus as of which of these techniques offer the best choice. In this chapter, we have investigated a simple wavelet-based de-noising technique to remove these artifacts from the E E G of awake and lightly sedated patients. Compared with previous methods, this technique neither relies on the availability of reference EOG channel, nor makes the distinction between eyeball movement and eye blink artifact correction. Since in real-life applications both artifacts are present simultaneously, a method based on such discrimination is unreliable. The same method has successfully been tested for head movement artifacts. The nature of E C G artifacts makes them suitable candidates for being removed from CHAPTER 4. PRE-CONDITIONING OF THE EEG SIGNAL 84 the E E G by this method as well. Thus, separate analyses are not required to remove different classes of artifacts. Using the time-frequency localization property of the wavelet transform, and the high sampling feature of SWT in all frequency bands, the proposed technique has clearly shown its potential in correcting for low frequency artifacts, such as OAs, while preserving exactly the high frequency content of the original signal, both in amplitude and phase. This is very important in the applications that are exploring various aspects of consciousness. The selection of thresholds can be further investigated. More conservative thresholds lead to a stronger filtering in lower frequency bands. Thresholds derived based on the E E G of a respective patient (self-norming) are preferred, but the method is robust enough for population-derived thresholds. In the former case, the calibration period is needed to acquire 30-60 seconds of, presumably, artifact-free E E G . If lower accuracy is allowed, the calibration procedure may be discarded by allowing the algorithm to use the parameters derived from the statistical analysis of data from a number of patients. Also, the thresholds can be periodically updated to account for the non-stationarity of the E E G . Adaptive filtering techniques might be investigated in the future to finely tune threshold values in real-time. Finally, the method is simple and computationally efficient. It can also be easily automated and im-plemented in real time. Further, when computational complexity is not of a concern, a higher level of decomposition (e.g. 7 levels) can produce even smoother results. However, there are also some potential limitations of the proposed method to be resolved in the future. The method in its present state is non-specific, that is, some pathological signals of potential interest might be filtered out along with artifacts. The excellent time-frequency localization properties of the SWT are a promising tool that should be investigated in order to separate different artifacts in the wavelet domain. In addition, the method is not well suited for the removal of the E M G artifacts. However, the E M G artifacts are usually filtered out to an acceptable degree by conventional filtering techniques, since their spectra (30-300Hz) does not severely overlap with the E E G spectrum of interest. Furthermore, it has been shown that statistically in about 90% of patients the EMG signal is a minor component of the observed scalp E E G signal [49]. However, we have observed that the teeth-clenching artifact (which can be considered as an E M G artifact) can elevate the WAV Index in deeper hypnotic states. However, since teeth clenching represents a patient reaction, it remains to be answered if the index should reflect it. Thus, the acceptable strategy for dealing with the E M G might be to observe the level of this artifact by calculating the energy conserved in wavelet coefficients in higher frequency bands (e.g. > 100 Hz) and to issue a warning if the presence of this artifact is too strong. C h a p t e r 5 Clinical Validation The promising results based on the first data recorded from a patient undergoing surgery, which were presented in Fig. 3.11 and published in [150], motivated us to obtain more clinical data in order to further improve the index. Therefore, a 16-case audit comprising of various surgical procedures was acquired. Based on these data, the current version of the WAV Index was produced, based on the SWT and isoelectric PDF template, and then implemented in real time. Further, the wavelet-based artifact removal technique was developed. Then, two studies were proposed and carried out with the aim of validating the WAV Index in the clinical environment. These two trials and the results obtained are detailed in the following. 5.1 Arthroscopy/Menisectomy ( A / M ) Observational Study The audit study generated hypotheses for further clinical research testing. Accordingly, an observational study was proposed. The next section presents its clinical research protocol. 5.1.1 Protocol A request for ethical review entitled: "Non-invasive Elecroencephalographic Monitoring During General Anesthesia for Arthroscopic Menisectomy" was submitted by our research group in the Division of Control Systems. After being reviewed and approved by the Clinical Research Ethics Board, the observational study was carried out. The E E G data was acquired from 25 patients undergoing various knee surgeries (arthroscopy, menisectomy or ligament repair) in the Operating Rooms at the University of British Columbia Hospital (UBCH) site of Vancouver Hospital & Health Sciences Centre (VHHSC). Patients unable to communicate orally in English or give a written informed consent were excluded. The study population further included adult patients (i.e. 18 years of age and older) of both male and female gender, with neither allergic history nor contraindication to anesthetics, and of ASA 1 physical status 1, 2 or 3. The BIS monitor (BIS A1050 1 T h e physical status classification system by the American Society of Anesthesiologists (ASA): • A S A 1: a normal healthy patient; • A S A 2: a patient with mild systemic disease; 85 CHAPTER 5. CLINICAL VALIDATION 86 with the analysis software v.3.4) was used as a data acquisition device. The raw E E G data and the BIS response to anesthesia were downloaded on a per second basis from the BIS monitor to the laptop computer via an interface, and the WAV Index was calculated in real time. Both indexes were displayed on the computer's LCD screen. However, the anesthesiologists were blinded to this indication in their assessment of the patient's anesthetic state, as well as the patients. Thus, other than the application of the BIS electrodes to the patient's forehead, patient care during general anesthesia was not changed by this observational study. In addition, various significant events during each procedure have been recorded and time stamped. Therefore, each patient was asked to count down during the administration of the induction bolus, starting from 100. The time corresponding to the Loss Of Count (LOC) was recorded as indicative of the loss of consciousness. The cardiorespiratory (heart rate, blood pressure, etc.) and motor responses during the surgery were also recorded, as well as the drug titration. Finally, the patient reaction to oral command (" squeeze your hand") was recorded every so often during recovery (i.e. emergence) until the patient started to be responsive. 5.1.2 Results Examples of Time Course of the WAV and BIS Indexes In this section we will discuss the results of the A / M observational study. Figure 5.1 shows two cases of this study that illustrate the remarkable correlation between the WAV Index and the bispectral index. At the recovery period (i.e. emergence) of the case #9, the BIS was not able to display a reliable value of index due to the presence of artifacts. However, after a couple of minutes, it reached the same value as the WAV Index. In Fig. 5.2, the good correlation between indexes is still present, however, some differences can be observed, especially during the emergence when the BIS is not able to recover to its baseline value. In addition, the case #22 is a good example that justifies the use of a monitor of hypnosis in clinical anesthesia, and further reveals the benefits which automatic drug delivery might bring to its practice. Namely, the level of patient's hypnotic state during the maintenance phase of a procedure was rather unstable, and, at some points in time, undesirably light. In the next sections, we are going to examine the induction and emergence phases of various cases from the study. • A S A 3: a patient with severe systemic disease; • A S A 4: a patient with severe systemic disease that is a constant threat to life; • A S A 5: a moribund patient who is not expected to survive without the operation; and • A S A 6: a declared brain-dead patient whose organs are being removed for donor purposes. CHAPTER 5. CLINICAL VALIDATION Fi gure 5-1.1 Time course of both the BIS (gray) and WAV (black) indexes - cases ^ 11 and ^ 9. CHAPTER 5. CLINICAL VALIDATION 88 SI8 - A V M SI8 - A V M Figure 5.2: Time course of both the BIS (gray) and WAV (black) indexes - cases #8 and #22. CHAPTER 5. CLINICAL VALIDATION 89 Loss of count 220 Time (s) Figure 5.3: Induction: 16 cases synchronized at the loss of count. Induction The term induction refers to the short period of time that follows the initial administration of anesthetic drugs. Induction is a very important phase in the course of a surgical procedure since during this phase the patient rapidly gains unconsciousness. The period of time closely following induction is usually a very deep hypnotic level which is needed to secure the patient's airway by inserting the endotracheal tube or Laryngeal Mask Airway (LMA). This procedure is very stimulating as it provokes strong gagging reflexes. Figure 5.3 summarizes graphically 16 of 25 cases from the A / M study. All presented cases have been synchronized to the LOC event. It can be observed that the WAV Index consistently starts to drop somewhat before the LOC, which is not the case with the bispectral index that is scattered and delayed. Further, the WAV Index has a clear lead of about 18 s in comparison with the BIS. The remaining 8 cases where excluded from Fig. 5.3, since the insertion of the L M A provoked strong reactions from patients as illustrated in Fig. 5.4. Emergence The term emergence refers to a period of time at the end of a surgical procedure when the patient regains his/her consciousness as a result of drug discontinuation. Figure 5.5 illustrates the course of recovery for 4 cases. The patient reactions to oral command and touch are marked, as well as the moment when the patient starts to be responsive. It is only with this last event that the patient is considered to be conscious. It can be observed that, besides having a faster response to the emergence of the patient, the WAV Index is more consistent in representing this event than the bispectral index (the BIS does not recover its baseline value in cases #8 and #21). CHAPTER 5. CLINICAL VALIDATION 90 Figure 5.4: Induction: 4 cases with strong LMA reaction (WAV - black; BIS - gray). Figure 5.6 illustrates two cases of emergence when the patient is left to regain consciousness without external stimulation by either voice or touch. Finally, Fig. 5.7 illustrates two cases where the emergence is not clearly captured by either the WAV Index or the BIS, since neither index discriminates between the states of patient reacting and patient responding. 5.1.3 Survey A complete survey of the A / M study is presented in Table 5.1 and Table 5.2. The first table summarizes the patient population (gender, age, weight, ASA) and the titration regimen administered during the surgery. In a majority of cases, the anesthesia maintenance was achieved by a combination of propofol (infusion) and inhalational anesthetic. An average MAC value is also given for information. The second table summarized the time shift at induction and emergence of the WAV Index as compared to the BIS. The value of the index at LOC, patient reacting and patient responding is also presented. Both the LOC and the patient reacting events take place during the transition from the awake to unconscious state, or vise versa. Note how the WAV Index captures these transients well as its value rapidly changes at these particular endpoints. This change is indicated by its slope. 5.2 Electro-Convulsive Therapy (ECT) Study Electro-Convulsive Therapy (ECT) can be used for the treatment of acute depression. It involves the delivery of electroshocks that induce seizures. It is advocated by current medical practice that such seizures have the potential for correcting the chemical imbalance within the brain responsible for depression. The therapy CHAPTER 5. CLINICAL VALIDATION 91 Figure 5.5: Emergence: 4 cases (WAV - black; BIS - gray). Figure 5.6: Emergence: 2 cases (WAV - black; BIS - gray). Figure 5.7: Emergence: 2 cases (WAV - black; BIS - gray). CHAPTER 5. CLINICAL VALIDATION 92 z o CD :« .S Q S S o _H pq t i o Z OH + PH < O l — . 0 0 II" CO ^ CO o .3 o - m • l O o Z I 1 Q <D Q •- c .-5 -3 c -.2 o. o '•s S >2 • .5 c < 3 > , CC C SR ce < 2 o 00 o 3 as CO CD a o 3 CC o <i 2 t -O o OQ a g a 3 — .3 pq 8 S OH , PH ft 3 o o Z u O in , 03 OH + J2 o CP O _. •- '3 a 5 .2 3 ' cn ,CD 5 3 .3 w — "o H2 CO O 3 S OQ PH , cti O 3 •a . 0 o cu a a o ca x cc * g Z CD a cc) u 3 CC 3 O 3 2 2 S 1 2 < S < PH 2 PH 2 2 3 to . £ £ 3 a co --H (H — o £ CO O 3 a ° 2 CQ OH , 3 3 3 « ^ o . 3 .O — O -H c2 « , <2 o 3 o tt "3 tt 2 pq 2 PH , PH O O ^ CQ "3 ' C O H< .2 2 -C i—1 +J ... O HJ •9 2 J S .5 _ >s >>-2 CD O 3 HO cti ' Is J o a > 11 "i ce O pq O < s o 3 3 g 2 1 9 CD CD 3 cc , o > 00 0 3 CC o > 3 CC o ce CD .c b .3 ^ _ i 3 ^ cc >8 f 2 < O OH 2 £ l cc . 2 CJ cci o> Z C O . - o • , -° -° a 2 S o 2 2 co Q . tt 3 0 0 — lH kH O OH OH CQ 0 0 0 0 0 _ 1 0 CN 3 3 ce c8 o tt o o tt o OH OH OH OH LO • O o q CN CN 1 0 CN T—i - H T—i —. —. —. —. '3 'a 'a 'a cS cS c3 C8 4J HJ 3 B 3 3 CD <2 <52 CD <HH 3 3 3 3 CO CO CO CO o tt o in OH CN 3 _cS " o o tt o o tt o o tt o o tt o cS "d o o cS HJ a O O o tt o a o a CD HH CD CHH Q HH c ° o Q O -2 o tt o <2 o tt o >S «2 ' o o o p tt tt 0 . 0 o o o <£ H 14 C J OH PH PH 2 + O O S 33 t - t - 3 C N C N S is S § 0 0 0 s 0 f\ ttO 2 < PH PH Pi OH & PH PH OH PH lH &H CD OH OH OH Pi cS T3 cci T3 I O ca P ca *J 3^ *^ a is a RCD o RCD fc N PH + O IO H u 2 ^ 3 -B £ CH ca o ca - - a - < 1 0 N <^ 3 ca „ 3 a ca " o a fa « ca a ca effi T3 T3 3 a co 2; co 3 JS IS) ca -a CN s 2 O ca a 4J rl 3 PH O T4 ca a H J 3 ca CD o OH tsl 1 0 E ca HJ a PH IO ">> a HJ a fa 0 0 1 0 r-H to CN — 1 _ • — 1 3 _ '3 'a 'a j a 'a ca C3 ca 0 ce 4H> HJ a a a a CD CD 3 3 3 3 CO CO CO 2 CO < CO < z j , , w t o ho; T3 a CD O 00 CO ^cf CO O 00 CN O 00 ^ ca 2 ce 2 <D <D Q) 0) <V cd cd P3 cd 2 2 2 22 a s CO < xl o Z 3 o o 3 o o O < CJ co 3 a 3 2 "5 2 PH 3 o 3 3 o o HJ HJ CD CD CD CD CQ CO 2 2 2 22 2 3 o CHAPTER 5. CLINICAL VALIDATION 5 « .2 x i CD tD CD Cfl iw 00 0) 3_= 23 CB 8 3 3 f .§ g H "S pa = CD CD 8 S , > 5 L « 00 i « 6 S 2 ffl <S C C f i i §•2 3 •s a s s CD j?. <D 22 - CD CD u .c ) o £ a " ! 3 s o : co . „ J CD JD _ < I « 5 J cj cfl I > bo CD > C •"•S ° 2 ! t Cfl i o d OT r"S I oi Cfl cc) •g "5 O ra ^ 3 O ho CD e <~ cj CD O CJ . CD ! 3 •§ ra .3 3 3 -O CD -O £ cfl c ^ CD ^ O OT ra CD 3.3 a. ca B CD T 3 > o ra CD ana s o 22 f" 2 pa z 3 1 S 2 f- 9 a5 O G) I> OI ra 3 I \ \ \ \ CD 00 d o o i> t - Oi ' \ \ \ (N H t> q CN o d o o \ \ \ \ CO co d d d d '5 \ \ \ \ \ \ \ 00 CN CS o o CO CN CO , | Oi i o CO Oi to O p^ d_ d^ o_ o d,cT CN l> CO l> <o to o CO CN CN \ \ \ \ \ \ \ \ \ C N Oi rA d CO CN , , , O d , CO O CO Oi CN d d CO 00 O I I > CO eo Oi t -/ \ \ \ / / \ in o i-l co d d \ CN o \ \ \ CN Tf d d 00 CD O 6 O CO ra 3 ^ q q o 6 6 CO / / \ / -cf / / CO o / / / / / O CO IN ^ O 00 CN I O O / / / / / / CO CM / / 00 i> 01 CO H O O rH O -I rH i—I CN d d --5 d_ 0 o iH °? 01 CO to d ^ ra CHAPTER 5. CLINICAL VALIDATION 94 usually consists of a series of repeats. For instance, acute patients can receive a treatment as often as every two days for 2 to 8 consecutive weeks. In some cases, patients must receive one treatment every few weeks to avoid relapse. Prior to the application of the electroshock, patients are induced with an induction dose of thiopental. While the dose might be adjusted according to the patient's response during his/her previous repeat, it is often maintained constant throughout their therapy. Following the anesthetic, a bolus of succinylcholine (muscle relaxant) is administered to avoid muscle movements during the procedure. The seizure can indeed provoke strong muscle contractions that might endanger the patient. No other drug is administered. The isolated forearm technique is used so that the patient can signal his/her discomfort. This technique consists of an inflated cuff that blocks blood circulation. As a result, no muscle relaxant is administered to the isolated limb. This particular procedure is thus a unique opportunity to observe intra-patient variability. By taking recordings of the same patient throughout his/her therapy, and assuming that the drug titration remains fixed, we can observe how consistent the WAV and bispectral indexes are. In an ideal case, each index should respond in a similar fashion under similar circumstances. 5.2.1 Protocol After obtaining patient's consent, a series of 6 repeats is recorded. For this particular clinical trial, we are using the BIS XP Monitor, which is the latest version available. This device is equipped with a serial port for data transfer and storage. We are also using an interface program developed by Aspect Medical Systems, Inc. This interface allows us to download the raw E E G signal as well as the processed data (bispectral index, MEF, SEF, etc.). A number of events are manually recorded during the treatment using the worksheet presented in Fig. 5.8. The Event 1 corresponds to the start of the anesthetic procedure which coincides with the administration of the induction bolus. The patient rapidly becomes unconscious. Any purposeful movement is then recorded. Just prior to the application of the electroshock, a nurse asks the patient to squeeze her. hand in order to check the adequacy of the anesthesia. This stage is marked as the Event 2. The Event 3 corresponds to the moment where the electroshock is applied. Muscle contractions due to the seizure can be best observed on the isolated arm. The end of the motor seizure is marked as the Event 4. Once the treatment is finished, patients are usually transported to a recovery room. The recovery period is also recorded. Every so often, a nurse asks the patient to squeeze her hand. This allows them to assess whether the patient has regained consciousness (marked as the Event 5). Once the patient is awake, nurses ask his/her name (orientation return - considered as the Event 6). Finally, the Event 7 corresponds to the moment when patients are discharged from the ICU. Both pre- and post-ECT hemodynamic parameters (heart rate, blood pressure) are recorded, as well as the drug titration and the placement of electroshock electrodes. The duration of the E E G seizure is also recorded (it is provided by a separate E E G monitor that is used by a psychiatrist who delivers the electroshock). CHAPTER 5. CLINICAL VALIDATION o CO 3 c - a> g ' • ra • ' 8 ro CO o z "D CD *o O CD "D C "D CD CO ro CO S. O tit CO b --. o Sys I — (bpm) IHR 3 E |Dose "3 E |Dose t9> E 1 Dose E I * - CD x £ \ ro co c 0 o CO E lo CO CO E 3 E E oa E -C ' c if . c ra CD c a E o a E V- CO h-g a> z CO L _ 3 N T3 sei ion O O CD UJ UJ Q 25 c p S "c 13 o z CO CD >-I I1 CD 8 o CO ion ( irat Q (Hz) 1 >. c £ LL (SLU) -g a> I Puis P I 1° f o UJ so. o sei Q — o Sys (bpm) IHR I Score | r*. H-O 2 a i-o LU CD 51 >= 1551 M i ICD! 0 -o Figure 5.8: ECT-BIS study worksheet (adapted from [8]). CHAPTER 5. CLINICAL VALIDATION 96 0 20 40 60 80 100 120 140 160 0 20 40 60 80 100 120 140 160 0 20 40 60 80 100 120 140 160 Time (s) Time (s) Time (s) Figure 5.9: Intra-individual variability at induction with thiopental - 3 ECT study cases. Patient #3: female, 150 mg thiopental. Patient #4: female, 200 mg thiopental. Patient #5: female, 250 mg thiopental. 5.2.2 R e s u l t s While this clinical trial started in February 2002, the study is still not complete at this date. Only a few complete cases (6 repeats for the same patient) have been recorded so far due to the limited availability of patients. In Fig. 5.9, the time course of the WAV Index and the bispectral index are plotted for three patients. Time t — 0 corresponds to the Event 1 (i.e. the end of the intravenous injection of thiopental). The time course is plotted up to the application of the electroshock. The comparison with the BIS clearly shows that the WAV Index presents a more consistent behavior and a less erratic time course. Note that the administration of succinylcholine produces fasciculation (muscle twitches) which results in some artifacts in the E E G . The bispectral index is particularly sensitive to these artifacts, while the WAV Index is less so. Note that only 4 time courses are plotted for the patient 4, as the titration was not maintained constant during the trial. This study being still underway, no further analysis of the data recorded so far has been performed. The preliminary results in terms of patient intra-variability still need to be confirmed. However, we feel confident that the cases presented in Fig. 5.9 are representative of the indexes behavior. 5.3 Summary Based on the results of the A / M and E C T studies, we derived a number of conclusions. The WAV Index and the BIS are well-correlated in the majority of cases. In average, a difference of 5.5 % can be observed. However, fast changes in the patient's hypnotic state, such as induction, emergence and intra-operative movements, are usually detected much sooner by the WAV index (i.e. 18.6 s ±9.1 in average during induction and 16.7 s ±16.6 during emergence). Further, the WAV index predicts changes in the patient's hypnotic state better than the BIS. This is particularly visible during induction, see Fig. 5.3, where the WAV index CHAPTER 5. CLINICAL VALIDATION 97 consistently starts to drop somewhat earlier then the patient becomes unconscious. This is understandable when knowing that the WAV Index is based on a per second analysis, and thus reacts instantaneously to the changes in the patient's state. Further, the WAV Index exhibits a smoother behavior than the BIS which tends to react abruptly and with a significant delay. Also, during the 5 cases, the BIS was unable to recover to its baseline value (i.e. 80-100) after emergence. Finally, the preliminary results of the E C T study, which is still in progress, indicate that the performance of the WAV Index exhibits less intra-patient variability. The excellent results obtained so far and the simplicity of the analysis algorithm, have motivated our group to patent this technology. A provisional patent has been submitted in July 2002 [175]. The non-provisional application will be filed within a year to fully protect the invention. C h a p t e r 6 Conclusion This concluding chapter summarizes in Section 6.1 the significance of the research work that has been carried out during the MSc. program here at the University of British Columbia. The main contributions are briefly presented in Section 6.2. Finally, this thesis would not be complete without giving some insights into the research directions that need to be investigated in order to improve our results. This is done in the last section of this chapter. 6.1 Significance of the Work In the modern practice of clinical anesthesia, practitioners rely on a number of physiological signs to assess the state of a patient. Since the first administration of ether as an anesthetic in the mid-nineteenth century, monitoring technology has considerably evolved. Pulse oximeters, capnographs and anesthetic mass spec-trometers are among the most recent advances in that field. However, all these devices quantify a particular parameter indicative of the patient's well-being. The measure of the anesthetic depth by itself has puzzled researchers for many years. New advances were made only with the understanding that anesthesia is a combination of different components, principal ones being hypnosis, analgesia and areflexia. Concerning hypnosis (i.e. anesthetic-induced unconsciousness), it is clear that the cortical activity reflects the patient's state of consciousness. Cortical activity can be assessed by observing the Electroencephalogram (EEG). While fast oscillations usually characterize an awake state, slow and large amplitude oscillations are typical for deeper hypnotic states. Different techniques have been used to derive a univariate index quantifying the changes in E E G activity. Time domain techniques and, later on, spectral and higher order spectral analysis were used. It is only in recent years that this type of technology has found its way into the Operating Room. The BIS® monitor (Aspect Medical Systems, Inc., MA) is probably the most representative. The interest of awareness/consciousness monitors remains to be proven. Anesthesia is a very safe proce-dure, with almost no incidence of intra- and post-operative mortality due to overdosing (in elective surgeries). Also, although bearing risks of serious psychological damage to the patient, the incidence of intra-operative 98 CHAPTER 6. CONCLUSION 99 awareness is extremely low. From the economic point of view, in the cases where law suits were engaged following intra-operative awareness with presence of pain, financial compensations were rather low (an av-erage of 18,000 USD). Thus the question whether such a device would be useful in the practice remains unanswered. With the numerous independent clinical trials assessing the performance of the BIS monitor, new insight has been gained in the advantages of measuring cortical activity during anesthesia. The major advantage of these monitors is that they actually provide a measure of the drug effect on the E E G . There is a clear dependence between the amount of drug administered to a patient and the index value. Hence, using such a device as a titration guide has proven beneficial in improving the anesthesia outcome (recovery time, hemodynamic stability, patient's comfort). In a majority of cases, less drugs than would be otherwise used were administered while maintaining an adequate anesthetic state. Further, such a monitor indirectly provides an insight into a cerebral oxygenation and hence anesthetic reversibility. Conclusively, such a device would give anesthesiologists an insight into the Central Nervous System (CNS) well-being as well as provide them with a tool for a more precise titration of each patient, according to his or her specific needs. The dose-response relationship provided by monitors of consciousness has motivated our group and other researchers to investigate automatic control of clinical anesthesia. While the BIS monitor has been accepted by the medical community, its inherent technology makes it a poor candidate for being used as a feedback quantity for closed-loop control. Indeed, fast transients in the patient's state are captured by the device with a significant delay. This delay is often mentioned by practitioners as one of the main reasons for not using the BIS monitor routinely in their practice. Also, for control purposes, such a delay would seriously undermine the controller performance, in particular its ability to deal with fast unmeasured disturbances caused by external stimuli. The aim of this thesis was to provide an index of hypnosis that reacts instantaneously to changes in the patient's state. This goal was met using wavelet analysis. We believe that the index we developed is superior to the bispectral index. Reactions from the anesthesiologists from the University of British Columbia Hospital were very favorable, and many who once used the BIS monitor would be willing to investigate further the use of our index in their practice. 6.2 Contributions This thesis has clearly established the potential of wavelets for the analysis of the E E G in the context of clinical anesthesia and for the assessment of anesthetic-induced unconsciousness. Wavelet analysis is particularly well-suited for the analysis of non-stationary signals. Its time-frequency localization property makes it particularly useful in a wide variety of applications, such as signal process-ing, de-noising, data compression, etc. While Fourier and short-time Fourier transforms have been used extensively in the past for feature detection and extraction tasks, wavelets provide a more dedicated tool. As far as Depth Of Anesthesia (DOA) is concerned, only very few attempts at using wavelets and wavelet transform to measure hypnosis have been reported in the literature. None of these attempts involved the analysis of passive E E G recordings. They all typically look at the analysis of evoked potentials, which are CHAPTER 6. CONCLUSION 100 particular patterns embedded in the E E G occurring as a result of sensory stimuli. In this work, we have tried to characterize the changes in patients' hypnotic state by quantifying the changes in the Probability Density Function (PDF) distribution of the wavelet coefficients in a particular frequency band. To the best of our knowledge, no prior work has been done in this direction in the past. Our proposed technique yields a highly repeatable index of hypnosis, referred to as the WAV Index. It is based on the comparison of the PDF of the wavelet coefficients extracted from the 7-band (32-64 Hz) with two templates corresponding to the awake state and the isoelectric E E G which occurs at the deepest anesthetized level. We decided to use the Stationary Wavelet Transform (SWT), as it yields more coefficients for statistical analysis. Furthermore, we have developed an original solution for the de-noising of the E E G in order to remove Ocular Artifacts (OAs) that strongly perturb the E E G waveform. Also based on SWT and thresholding of wavelet coefficients, the de-noising technique removes the large spikes and square-like waveshapes that originate from ocular activity, while successfully keeping high frequency information without distortion, neither in amplitude nor in phase. The technique is particularly well-suited for the analysis of the awake and light sedation states, which are typically corrupted by such artifacts. Further, the head movement artifacts are successfully removed as well. Once again, we haven't found any indication in the literature that a similar type of work has been done in the past. While the WAV Index is by itself a novel index of hypnosis, its characterization of the awake and lightly sedated states suffers tremendously from OAs. The association of the de-noising technique with the index thus yields a more robust index of hypnosis. However, note that this system is not yet fully automated, refer to Section 6.3. We envisage that, using the additional information such as E M G , a completely automated solution can be provided. The extensive comparison with the bispectral index carried out in the context of clinical trials have shown that the WAV Index correlates strongly with the BIS. While these two techniques use the same raw signal, the underlying analysis is carried out by two totally different algorithms. The observed similarity between the BIS and the WAV Index tends to prove that both indexes are successful in extracting the relevant information from the E E G . In addition, one remarkable feature of the WAV Index is that its algorithm has been developed based on a very limited amount of data. Two E E G signals acquired from an awake subject and an anesthetized patient were sufficient to provide enough insight to derive the index. Furthermore, the simplicity of the WAV Index yields a computationally efficient algorithm. Since the wavelet analysis is able to capture precisely all of the signal information contained in a very short time window, the index does not depend on an extensive amount of past data. As a result, it exhibits a clear lead time of about 18 s as compared to the BIS. Data acquired during electro-convulsive shock therapy have also indicated that the WAV Index yields more consistent results than the BIS. CHAPTER 6. CONCLUSION 101 In order to summarize the achievements of the WAV Index, we will now argue to which extent we have fulfilled the criteria defined in Section 1.2.3 that an ideal monitor of DOA should meet: 1. The WAV Index shows graded changes of the patient's hypnotic state. Since intra-operative awareness is an event of a very low incidence, extensive clinical testings are needed to assess whether the index detects awareness. 2. Preliminary results of the validation studies suggest that the index is able to reflect changing concen-trations of anesthetic agents (however, further validation is needed), as well as the level of hypnosis. 3. The WAV Index demonstrates a good sensitivity and specificity. The exact values of these parameters are yet to be calculated. 4. The index has been applied for a wide variety of anesthetics used nowadays, and their combinations; it seems to be drug independent. 5. Two validation studies included cases where different non-anesthetic drugs were used; the index has not been influenced by them. 6. The index has been applied to a variety of adult patients, including ASA physical status 1, 2 and 3. 7. The baseline value (the awake state corresponds to index values of 90%-100%) is stable and recovers upon the discontinuation of drugs. 8. The index has a large therapeutic window that coincides well with the changes induced by drugs that were encountered within the cases studied. 9. The WAV index is associated with the method that successfully removes artifacts. This method is close to being fully automated. 10. The WAV index has an excellent temporal resolution (1 s) with instantaneous response to the change in the patient's hypnotic level, since it is based on a per second analysis. 11. There is no need for the index calibration; the index is robust since it is based on the population derived parameters. However, the self-normed parameters can lead to more precise results. In that case, an additional acquisition of 30 to 60-second artifact-free E E G signal from a patient is needed prior to the start of monitoring. 12. The index uses a scale between 0% (isoelectic EEG) and 100% (the awake state) that is easy to interpret and understand. The excellent results obtained so far have motivated our group to patent this technology. More clinical validation is planned in the future. CHAPTER 6. CONCLUSION 6.3 Future Work 102 There are some improvements that would greatly contribute to the proposed WAV Index so that it evolves to a solution comparable to commercially available monitors of hypnosis. In the following, we will briefly discuss the issues that, in our opinion, need to be addressed as in the future. We believe that the core algorithm for assessing in real time the patient's anesthetic state based on the wavelet analysis of his/her electroencephalogram has, in many ways, an advantage over the existing solutions. Also, the proposed de-noising technique offers a lot of potential in performing the complex task of removing typical artifacts from the E E G signal. The solution successfully suppresses OAs as well as head movement artifacts from the E E G in the awake state and light sedation. The association of the de-noising technique with the index yields excellent results in the off-line analysis. Unfortunately, at this moment, we haven't been able to offer an integral solution that is completely automated. Epochs corrupted with ocular and head movement artifacts often greatly resemble, both in the time and frequency domains, to the uncorrupted epochs of deeper hypnotic states. Hence, de-noising of the E E G in deeper hypnotic states removes the relevant information thus leading to incorrect results. Therefore, it is essential to apply the artifact removal method only in lighter hypnotic states, when these artifacts are present and relevant. We initially used the WAV Index values to switch on and off the de-noising. However, results were unsatisfactory. The current solution, which allows for the index to be automated, uses the moment of induction (i.e. administration of the induction anesthetic bolus) to switch off the de-noising. In the context of automated anesthesia this solution is justifiable, since it is a known time stamp. Note that, in this case, the recovery period is not de-noised. Although, during this period, the E E G does not significantly suffer from artifacts, this might be a drawback. Thus, in order for the system to be truly automated, the means of making a more reliable decision on when to apply the de-noising need to be investigated. One potential solution might be to investigate the Electromyography (EMG) signal of facial muscles, especially the frontalis muscle, since the strong presence of this signal is usually indicative of the awake state. However, we see the need for a data acquisition device that is totally transparent1. Another foreseen possibility is the development of an artifact recognition algorithm that would be further associated with the proposed de-noising technique. We feel that wavelets can probably offer the necessary separation of artifacts and uncorrupted E E G epochs in the case of any hypnotic state. Another issue that needs to be resolved is the calculation of the Burst Suppression Ratio (BSR) parameter, which is associated with large doses of general anesthetics, and further carries important information in the case of head trauma or brain ischemia. However, its calculation is straightforward since this time-domain parameter is calculated as the fraction of the defined length of time where the E E G amplitude is suppressed below 5 uV for more than 0.5 s. We did not calculate the BSR parameter since we did not have access to information concerning the electrode impedance values and the electrode ground check periods (during this period the acquired E E G has near to zero amplitude). Also, if we have access to a transparent data acquisition device, the E M G artifacts and their influence x N o t e that the BIS monitor gives the indicat ion of the E M G activity, however, we have no knowledge about its true meaning. CHAPTER 6. CONCLUSION 103 can be thoroughly investigated2. We think that the acceptable strategy for dealing with the E M G might be to observe the level of this artifact by calculating the energy conserved in the wavelet coefficients in the higher frequency bands (> 100 Hz) and to issue a warning if the presence of this artifact is too strong. Investigation of the linearity of the WAV Index in terms of the relationship between the drug dose and the index response would be extremely useful. Defining and testing an intermediate endpoint between the awake and isoelectric E E G , and corresponding values of the WAV Index would further increase the reliability of the index. However, these tasks would require extensive clinical testings. Finally, although in our work we found that Wavelet Packets (WPs) did not improve the results while being more computationally demanding, they should be further investigated (e.g. the 7-band can be split into finer subbands). In particular, if the computational complexity is of no concern, the redundant WPs might be worth further explorations. 2 N o t e that we have not noticed that the influence of these artifacts is significant. Bibliography [1] C. Prys-Roberts, "Anaesthesia: a practical or impractical construct?," British Journal of Anaesthesia, vol. 59, p. 1341, Nov. 1987. [2] P. Glass, "Why and how we will monitor the state of anesthesia in 2010?," Acta Anaesthesiologica Belgica, vol. 50, pp. 35-44, Jan. 1999. [3] J. Sigl and N . Chamoun, "An introduction to bispectral analysis for the EEG," Journal of Clinical Monitoring, vol. 10, pp. 392-404, 1994. [4] J. Johansen and P. Sebel, "Development and clinical application of electroencephalographic bispectrum moni-toring," Anesthesiology, vol. 93, pp. 1336-1344, Nov. 2000. [5] M . Vetterli and J. Kovacevic, Wavelets and Subband Coding. Prentice Hall, 1995. [6] A . Teolis, Computational Signal Processing with Wavelets. Birkhauser, 1998. [7] A . Aldroubi and M . Unser, eds., Wavelets in Medicine and Biology. CRC Press, 1996. [8] C. Ries, UBC CREB #C01-0108. 2002. [9] G. Rushman, N . Davies, and R. Atkinson, A short history of anaesthesia: The first 150 years. Reed Educational and Professional Publishing Ltd, 1 ed., Jan. 1996. [10] A. Stacy, Awareness in the operating room: a patient's view, pp. 10-30. Englewood Cliffs, NJ: Prentice Hall, 1993. [11] D. Stanski, Monitoring Depth of Anesthesia, chapter 30, pp. 1001-1029. Churchill Livingston, third ed., 1990. [12] I. Kissin, "General anesthetic action: an obsolete notion?," Anesthesia Analgesia, vol. 76, no. 2, pp. 215-218, 1993. [13] T. Gray and G. Rees, "The role of apnoea in anaesthesia for major surgery," British Medical Journal, vol. 2, pp. 891-892, 1952. in Johansen: 2000. [14] P. Woodbridge, "Changing concepts concerning depth of anesthesia," Anesthesiology, vol. 18, p. 536, 1957. [15] M . Pinsker, "Anaesthesia: a prapragmatic construct (letter)," Anesthesia Analgesia, vol. 65, p. 819, 1986. [16] P. Sebel, E . Lang, I. Rampil, P. White, R. Cork, M . Jopling, N . Smith, P. Glass, and P. Manberg, "A multicen-ter study of bispectral electroencephalogram analysis for monitoring anesthetic effect," Anesthesia Analgesia, vol. 84, pp. 891-899, 1997. [17] I. Kissin, "Depth of anesthesia and bispectral index monitoring," Anesthesia Analgesia, vol. 90, pp. 1114-1117, May 2000. 104 BIBLIOGRAPHY 105 [18] C. Rosow, "Can we measure depth of anesthesia?," 51st Annual Refresher Course Lectures, vol. 226, pp. 1-7, 2000. Volume 226 corresponds to lecture section. [19 [20; [21 [22 [23 [24] [25 [26 [27 [28 [29 [30 [31 [32 [33; [34; [35; I. Rampil, "Monitoring depth of anesthesia," Current Opinion in Anesthesiology, vol. 14, pp. 649-653, 2001. R. Epstein, "Pro: monitoring the amnestic state during general anesthesia should be a standard of care," Anesthesiology, pp. 506-508, Oct. 2000. G. Schneider and P. Sebel, "Monitoring depth of anesthesia," European Journal of Anesthesiology, vol. 14 (Suppl. 15), pp. 21-28, 1997. J. Drummond, "Monitoring depth of anesthesia," Anesthesiology, vol. 93, pp. 876-882, Sept. 2000. With emphasis on the application of the bispectral index and the Middle Latency Auditory Evoked Response to the prevention of recall. J . Jones, "Perception and memory during general anaesthesia," British Journal of Anaesthesia, vol. 73, pp. 31-37, July 1994. R. Sundlin, G. Enlund, P. Samuelson, and C. Lennmarken, "Awareness during anesthesia: a prospective case study," Lancet, vol. 355, pp. 707-711, 2000. P. Myles, D. Williams, and M . Hendrata, "Patient satisfaction after anaesthesia and surgery: results of a prospective survey of 10,811 patients," British Journal of Anaesthesia, vol. 84, pp. 6-10, 2000. D. Newton, C. Thornton, and C. Jordan, The auditory evoked response as a monitor of anesthetic depth, pp. 274-280. Englewood Cliffs, NJ : Prentice Hall, 1993. C. Pomfrett, "Monitoring depth of anesthesia," The Royal College of Anaesthetists. Bulletin 4, PP- 155-157, Nov. 2000. M . Tunstall, "Detecting wakefulness during general anaesthesia for caesarean section," British Journal of Medecine, vol. 1, p. 1321, 1977. N . Moerman, B. Bonke, and J. Oosting, "Awareness and recall during general anesthesia," Anesthesiology, vol. 79, pp. 454-464, Sept. 1993. M . Struys, L. Versishelen, E. Mortier, D. Ryckaert, J. C. D. May, C. D. Deyne, and G. Roily, "Comparison of spontaneous frontal E M G , E E G power spectrum and bispectral index to monitor propofol drug effect and emergance," Acta Anaesthesiologica Scandinavica, vol. 42, pp. 628-636, July 1998. C. Pomfrett, "Heart rate variability, BIS and depth of anaesthesia," British Journal of Anaesthesia, vol. 82, pp. 659-662, May 1999. I. Rampil, "A primer for E E G signal processing in anesthesia," Anesthesiology, vol. 89, pp. 980-1002, Oct. 1998. M . Rubin and H. Freeman, "Brain potential changes in man during cyclopropane anesthesia," Journal of Neurophysiology, vol. 3, pp. 33-42, Jan. 1940. A. Faulconer and R. Bickford, Electroencephalography in Anesthesiology. Springfield, Illinois: Charles C Thomas, 1960. D. Clark and B. Rosner, "Neurophysiology effects of general anesthetics. I. The electroencephalogram and sensory evoked responses in man," Anesthesiology, vol. 38, pp. 564-582, June 1973. BIBLIOGRAPHY 106 [36] N . Smith, H. Dec-Silver, T. Sanford, C. Westover, M . Quinn, F. Klein, and D. Davis, "EEGs during high-dose fentanyl-, sufentanil-, or morphine-oxygen anesthesia," Anesthesia Analgesia, vol. 63, pp. 386-393, Apr. 1984. [37] J. Scott, K . Ponganis, and D. Stanski, "EEG quantitation of narcotic effect: the comparative pharmacodynam-. ics of fentanyl and alfentanyl," Anesthesiology, vol. 62, pp. 234-241, Mar. 1985. [38] J. Scott, J. Cooke, and D. Stanski, "Electroencephalographic quantitation of opiod effect: comparative phar-macodynamics of fentanyl and sufentanil," Anesthesiology, vol. 74, pp. 34-42, Jan. 1991. [39] V . Billard and S. Shafer, Does the EEG Measure Therapeutic Opioid Drug Effect?, pp. 79-95. Berlin: Springer-Verlag, 1995. [40] I. Pichlmayr, U . Lips, and H. Kiinkel, Fundamentals of Electroencephalographic Analysis, pp. 12-21. Springer-Verlag, 1984. [41] O. Wilder-Smith, "CNS monitoring during general anaesthesia: making pain visible." The 12th World Congress of Anaesthesiologists, Montreal, Canada, June 2000. [42] R. Bickford, "Automatic electroencephalographic control of general anesthesia," Electroencephalography and Clinical Neurophysiology, vol. 2, pp. 93-96, 1950. [43] R. Bickford, "Use of frequency discrimination in the automatic electroencephalographic control of anesthesia (servo-anesthesia)," Electroencephalography and Clinical Neurophysiology, vol. 3, pp. 83-86, 1951. [44] D. Soltero, A. Faulconer, and R. Bickford, "The clinical application of automatic anesthesia," Anesthesiology, vol. 12, pp. 574-582, Sept. 1950. [45] J. Bellville and G. Attura, "Servo control of general anesthesia," Science, vol. 126, pp. 827-830, Oct. 1957. [46] F. Klein, "A waveform analyzer applied to the human EEG," IEEE Transactions on Biomedical Engineering, vol. 23, pp. 246-252, 1976. [47] T. Gregory and D. Pettus, "An electroencephalografic processing algorithm specifically intended for analysis of cerebral electrical activity," Journal Clinical Monitoring, vol. 2, pp. 190-197, 1986. [48] I. Pichlmayr, U . Lips, and H. Kiinkel, The Electroencephalogram in Anesthesia. Berlin: Springer-Verlag, 1984. [49] J. Sleigh, D. A . Steyn-Ross, M . L. Steyn-Ross, M . L. Williams, and P. Smith, "Comparison of changes in electroencephalographic measures during induction of general anaesthesia: influence of the gamma frequency band and electromyogram signal," British Journal of Anaesthesia, vol. 86, pp. 40-58, 2001. [50] H. Schwilden, J. Schiittler, and H. Stoeckel, "Closed-loop feedback control of methohexital anesthesia by quantitative E E G analysis in humans," Anesthesiology, vol. 67, pp. 341-347, Sept. 1987. [51] H. Schwilden and H. Stoeckel, "Quantitative E E G analysis during anaesthesia with isoflurane in nitrous oxide at 1.3 and 1.5 M A C , " British Journal of Anaesthesia, vol. 59, pp. 738-745, June 1987. [52] H. Schwilden, H. Stoeckel, and J. Schiittler, "Closed-loop feedback control of propofol anesthesia by quantitative E E G analysis in humans," British Journal of Anaesthesia, vol. 62, pp. 290-296, Mar. 1989. [53] J. Schiittler and H. Schwilden, Feedback Control of Intravenous Anaesthetics by Quantitative EEG, pp. 194-207. Berlin: Springer-Verlag, 1995. BIBLIOGRAPHY 107 [54] I. Rampil, F. Sasse, N . Smith, B. Hoff, and D. Flemming, "Spectral edge frequency - a new correlate of anesthetic depth," Anesthesiology, vol. 53, p. S12, Sept. 1980. [55] K . Gregg, J. Varvel, and S. Shafer, "Application of semilinear canonical correlation to the measurement of opioid drug effect," Journal of Pharmacokinetics and Biopharmaceutics, vol. 20, pp. 611-635, June 1992. [56] J. Drummond, C. Brann, D. Perkins, and D. Wolfe, "A comparison of median frequency, spectral edge frequency, a frequency band power ratio, total power, and dominance shift in the determination of depth of anesthesia," Acta Anaesthesiologica Scandinavica, vol. 35, no. 8, pp. 693-699, 1991. [57] R. Dwyer, I. Rampil, E. Eger, and H. Bennett, "The E E G does not predict movement in response to surgical incision at 1.0 M A C isoflurane," Anesthesiology, vol. 75, no. 3A, p. A1025, 1991. [58] J. Liu, H . Singh, and P. White, "Electroencephalographic bispectral index correlates with intraoperative recall and depth of propofol-induced sedation," Anesthesia Analgesia, vol. 84, pp. 185-189, 1997. [59] W. J. Levy, "Intraoperative E E G patterns: implications for E E G monitoring," Anesthesiology, vol. 60, pp. 430-434, May 1984. [60] G. Gurman, "Assessment of depth of general anesthesia. Observations on processed E E G and spectral edge frequency," International Journal of Clinical Monitoring and Computing, vol. 11, pp. 185-189, Aug. 1994. [61] A. Sharma and R. Roy, "Design of a recognition system to predict movement during anesthesia," IEEE Trans-actions on Biomedical Engineering, vol. 44, pp. 505-511, June 1997. [62] J. Sleigh and J. Donovan, "Comparison of the bispectral index, 95% spectral edge frequency and approxi-mate entropy of the E E G , with changes in heart rate variability during induction and recovery from general anaesthesia," British Journal of Anaesthesia, vol. 82, no. 5, pp. 666-671, 1999. [63] F.-Z. Shaw, R.-F. Chen, H.-W. Tsao, and C.-T. Yen, "Algorithmic complexity as an index of cortical function in awake and pentobarbital-anesthetized rats," Journal of Neuroscience Methods, vol. 93, pp. 101-110, 1999. [64] H. Viertio-Oja, V . Jantti, P. Talja, H. Tolvanen-Laakso, and A. Yli-Hankala, "Fractal spectrum, bispectrum, complexity, and entropy of E E G during anaesthesia," European Journal of Anaesthesiology, vol. 17 (Suppl.19), p. A271, 2000. [65] J. Bruhn, H. Ropcke, and A. Hoeft, "Approximate entropy as an electroencephalographic measure of anesthetic drug effect during desflurane anesthesia," Anesthesiology, vol. 92, pp. 715-726, 2000. [66] C. Thornton, M . Barrowcliffe, K . Konieczko, P. Ventham, C. Dor, D. N . DE, and J. Jones, "The auditory evoked response as an indicator of awareness," British Journal of Anaesthesia, vol. 63, pp. 113-115, 1989. [67] F. Davies, H. Mantzaridis, and G. Kenny, "Middle latency auditory evoked potentials during repeated transi-tions from consciousness to unconsciousness," Anesthesia, vol. 51, pp. 107-113, 1996. [68] M . Doi, R. Gajraj, and H. Mantzaridis, "Relationship between calculated blood concentration of propofol and electrophysiological variables during emergence from anaesthesia: comparison of bispectral index, spectral edge frequency, median frequency and auditory evoked potential index," British Journal of Anaesthesia, vol. 78, pp. 180-184, 1997. [69] R. Sharpe, D. Nathwani, and S. Pal, "Auditory evoked response, median frequency and 95% spectral edge during anesthesia with desfurane and nitrous oxide," British Journal of Anaesthesia, vol. 78, pp. 282-285, 1997. BIBLIOGRAPHY 108 [70] H. Mantzaridis and G. Kenny, "Auditory evoked potential index: a quantitative measure of changes in auditory evoked potentials during general anaesthesia," Anaesthesia, vol. 52, pp. 1030-1036, 1997. [71] G. Plourde and T. Picton, "Human auditory steady-state response during general anesthesia," Anesthesia Analgesia, vol. 71, pp. 460-468, Nov. 1990. [72] G. Plourde and J. Boylan, "The auditory steady state response during sufentanil anaesthesia," British Journal of Anaesthesia, vol. 66, pp. 683-691, 1991. [73] G. Plourde, "The effects of propofol on the 40-Hz auditory steady-state response and on the electroencephalo-gram in humans," Anesthesia Analgesia, vol. 82, pp. 1015-1022, May 1996. [74] G. Plourde and C. Villemure, "Comparison of the effects of enflurane/N20 on the 40-Hz auditory steady-state response versus the auditory middle-latency response," Anesthesia Analgesia, vol. 82, pp. 75-83, Jan. 1996. [75] G. Plourde, "Auditory evoked potentials and 40-Hz oscillations: an opportunity to study mechanism of action of general anesthetics?," Anesthesiology, vol. 91, pp. 1187-1189, Nov. 1999. [76] W. Jensen, M . Nygaard, and S. Henneberg, "On-line analysis of middle latency auditory evoked potentials (MLAEP) for monitoring depth of anaesthesia in laboratory rats," Med. Eng. Phys., vol. 20, no. 10, pp. 722-728, 1998. [77] D. Stanski, Monitoring Depth of Anesthesia, chapter 29, pp. 1087-1116. Churchill Livingston, fifth ed., 2000. [78] T. Ning and J. Bronzino, "Bispectral analysis of the rat E E G during various vigilance states," IEEE Transac-tions on Biomedical Engineering, vol. 36, pp. 497-499, Apr. 1989. [79] L. Kearse, V . Saini, F. deBros, and N . Chamoun, "Bispectral analysis of E E G may predict anesthetic depth using narcotic induction," Anesthesiology, vol. 3A, p. A175, Sept. 1991. [80] P. Sebel, S. Bowles, V. Saini, and N . Chamoun, "Accuracy of E E G in predicting movement at incision during isoflurane anesthesia," Anesthesiology, vol. 3A, p. A446, Sept. 1991. [81] J. Vernon, S. Bowles, P. Sebel, and N . Chamoun, "EEG bispectrum predicts movement at incision during isoflurane or propofol anesthesia," Anesthesiology, vol. 77, no. 3A, p. A502, 1992. [82] J. Muthuswamy and R. Roy, "Bispectrum analysis of E E G of a dog to determine the depth under halothane anesthesia," in Proceedings of IEEE Bioengineering Conference, pp. 5-6, 1993. [83] S. Bowles, P. Sebel, V. Saini, and N . Chamoun, Effects of anesthesia on the EEG-bispectral analysis correlates with movement, pp. 249-254. Englewood Cliffs, NJ: Prentice Hall, 1993. [84] C. Rosow.and P. Manberg, "Bispectral Index monitoring," Anesthesiology Clinics of North America: Annual of Anesthetic Pharmacology, vol. 19, no. 4, pp. 947-966, 2001. [85] L . Kearse, C. Rosow, and P. Sebel, "The bispectral index correlates with sedation/hypnosis and recall: com-parison using multiple agents," Anesthesiology, vol. 83, p. A507, 1995. [86] R. Flaishon, P. Sebel, and J. Sigl, "Detection of consciousness following thiopental: isolated forearm and bispectral E E G (BIS)," Anesthesiology, vol. 83, no. 3A, p. A515, 1995. [87] K . Leslie, D. Sessler, M . Schroeder, and K . Walters, "Propofol blood concentrations and the bispectral index predict suppression of learning during propofol/epidural anesthesia in volunteers," Anesthesia Analgesia, vol. 81, pp. 1269-1274, 1995. BIBLIOGRAPHY 109 [88] L. Kearse, C. Rosow, P. Glass, and J. Sigl, "Monotonic changes in EEG bispectral index correlate with targeted plasma concentrations of propofol and midazolam," Anesthesia Analgesia, vol. 82, p. S220, 1996. [89] J. Liu, H. Singh, and P. White, "Electroencephalogram bispectral analysis predicts the depth of midazolam-induced sedation," Anesthesiology, vol. 84, pp. 64-69, 1997. [90] V. Billiard, B. Plaud, G. Boulay, R. Bocquet, and B. Debaene, "Monitoring induction and maintenance of sevofiurane anesthesia by bispectral analysis of EEG," Anesthesiology, vol. 85, no. 3A, p. A352, 1996. [91] P. Glass, M. Bloom, L. Kearse, C. Rosow, P. Sebel, and P. Manberg, "Bispectral analysis measures sedation and memory effects of propofol, midazolam, isoflurane, and alfentanil in healthy volunteers," Anesthesiology, vol. 86, pp. 836-847, 1997. [92] L. J. Kearse, C. Rosow, and A. Zaslavsky, "Bispectral analysis of the electroencephalogram predicts conscious processing of information during propofol sedation and hypnosis," Anesthesiology, vol. 88, pp. 25-34, 1998. [93] P. Sebel, F. Payne, T. Gan, C. Rosow, and S. Greenwald, "Bispectral analysis (BIS) monitoring improves PACU recovery from propofol alfentanil/N20 anesthesia," Anesthesiology, vol. 85, no. 3A, p. A468, 1996. [94] F. Payne, P. Sebel, P. Glass, C. Rosow, and J. Sigl, "Bispectral index (BIS) monitoring allows faster emergence from propofol alfentanil/N20 anesthesia," Anesthesiology, vol. 85, no. 3A, p. A1056, 1996. [95] D. Song, G. Joshi, and P. White, "Titration of volatile anesthetics using bispectral index facilitates recovery after ambulatory anesthesia," Anesthesiology, vol. 87, pp. 842-848, 1997. [96] T. Gan, P. Glass, and A. Windsor, "BIS utility study group. Bispectral index monitoring allows faster emergence and improved recovery from propofol, alfentanil, and nitrous oxide anesthesia," Anesthesiology, vol. 87, pp. 808-815, 1997. [97] D. Song, J. van Vlymen, and P. White, "Is the bispectral index useful in predicting fast-track eligibility after ambulatory anesthesia with propofol and desflurane?," Anesthesia Analgesia, vol. 87, pp. 1245-1247, 1998. [98] M. Struys, L. Versichelen, G. Byttebier, E. Mortier, A. Moerman, and G. Roily, "Clinical usefulness of the bispectral index for titrating propofol target effect-site concentrations," Anesthesiology, vol. 53, pp. 4-12, 1998. [99] A. Yli-Hankala, A. Vakkuri, P. Annila, and K. Korttila, "EEG bispectral index monitoring in sevofiurane or propofol anaesthesia: analysis of direct costs and immediate recovery," Acta Anaesthesiologica Scandinavica, vol. 43, pp. 545-549, 1999. [100] J. Vernon, E. Lang, P. Sebel, and P. Manberg, "Prediction of movement using bispectral electroencephalographic analysis during propofol/alfentanil or isoflurane/alfentanil anesthesia," Anesthesia Analgesia, vol. 80, pp. 780-785, 1995. [101] K. Leslie, D. Sessler, and W. Smith, "Prediction of movement during propofol/nitrous oxide anesthesia, per-formance of concentration, electroencephalographic, pupillary, and hemodynamic indicators," Anesthesiology, vol. 84, pp. 52-63, 1996. [102] I. Iselin-Chaves, R. Flaishon, P. Sebel, S. Howell, T. Gan, and J. Sigl, "The effect of the interaction of propofol and alfentanil on recall, loss of consciousness, and the bispectral index," Anesthesia Analgesia, vol. 87, pp. 949-955, Oct 1998. BIBLIOGRAPHY 110 [103] R. Flaishon, A. Windsor, J. Sigl, and P. Sebel, "Recovery of consciousness after thiopental or propofol. Bispec-tral index and isolated forearm technique," Anesthesiology, vol. 86, no. 3, pp. 613-619, 1997. [104] R. Gajraj, M . Doi, H. Mantzaridis, and G. Kenny, "Analysis of the E E G bispectrum, auditory evoked potentials and the E E G power spectrum during repeated transitions from consciousness to unconsciousness," British Journal of Anaesthesia, vol. 80, no. 1, pp. 46-52, 1998. [105] G. Lubke, C. Kerssens, H. Phaf, and P. Sebel, "Dependence of explicit and implicit memory on hypnotic state in trauma patients," Anesthesiology, vol. 90, no. 3, pp. 670-680, 1999. [106] J. Liu, H . Singh, and P. White, "Use of E E G bispectral index to predict awakening from general anesthesia," Anesthesia Analgesia, vol. 78, p. S254, 1994. [107] T. Dushane, "Cons: monitoring the amnestic state during general anesthesia should not be a standard of care," Anesthesiology, pp. 509-512, Oct. 2000. [108] P. Manberg, D. Zraket, L. Kovitch, and L. Christman, "Awareness during BIS monitoring: 2001 update," Anesthesiology, vol. 95, no. 3A, p. A564, 2001. [109] F. Nishihara and S. Saito, "Pre-ictal bispectral index has a positive correlation with seizure duration during electroconvulsive therapy," Anesthesia Analgesia, vol. 94, pp. 1249-1252, 2002. [110] M . Alkire and C. Pomfrett, "Toward the fundamental unit of anesthetic depth: positron emission tomography suggests that bispectral index BIS monitors an important component of anesthetic depth," Anesthesiology, vol. 85, no. 3A, p. A174, 1996. [Ill] T. Katoh, H. Bito, and S. Sato, "Influence of age on hypnotic requirement, bispectral index, and 95% spectral edge frequency associated with sedation induced by sevoflurane," Anesthesiology, vol. 92, pp. 55-61, 2000. [112] C. Bannister, K . Brosius, J. Sigl, B. Meyer, and P. Sebel, "The effect of bispectral index monitoring on anesthetic use and recovery in children anesthetized with sevoflurane in nitrous oxide," Anesthesia Analgesia, vol. 92, no. 4, pp. 887-881, 2001. [113] W. Denman, E. Swanson, D. Rosow, K . Ezbicki, P. C. PD, and C. Rosow, "Pediatric evaluation of the bispectral index (BIS) monitor and correlation of BIS with end-tidal sevoflurane concentration in infants and children," Anesthesia Analgesia, vol. 90, pp. 872-877, 2000. [114] C. Degoute, C. Macabeo, C. Dubreuil, R. Duclaux, and V. Banssillon, "EEG bispectral index and hypnotic component of anesthesia induced by sevoflurane: comparison between children and adults," British Journal of Anaesthesia, vol. 86, pp. 209-212, 2001. [115] C. D. Deyne, M . S. and. J. Decruyenaere, J. Creupelandt, E. Hoste, and F. Colardyn, "Use of continuous bispectral E E G monitoring to assess depth of sedation in ICU patients," Intensive Care Med, vol. 24, pp. 1294-1298, 1998. [116] J. Sleigh, J. Andrzejowski, D. Steyn-Ross, and M . L. Steyn-Ross, "The bispectral index: a measure of depth of sleep?," Anesthesia Analgesia, vol. 88, pp. 659-661, 1999. [117] T. Schnider, M . Luginbuhl, S. Petersen-Felix, and J. Mathis, "Unreasonably low bispectral index values in a volunteer with genetically determined low-voltage electroencephalographic signal," Anesthesiology, vol. 89, pp. 160-161, 1998. BIBLIOGRAPHY 111 [118] J. Bruhn, T. Bouillon, and S. Shafer, "Electromyographic activity falsely elevates the bispectral index," Anes-thesiology, vol. 92, pp. 1485-1487, 2000. [119] S. Hagihira, M . Takashina, T. Mori, T. Mashimo, and I. Yoshiya, "Practical issues in bispectral analysis of electroencephalographic signals," Anesthesia Analgesia, vol. 93, pp. 966-970, 2001. [120] G. Baker, J. Sleigh, and P. Smith, "Electroencephalographic indices related to hypnosis and amnesia during propofol anaesthesia for cardioversion," Anaesth Intensive Care, vol. 28, pp. 386-391, 2000. [121] M . Struys, T. D. Smet, L. Versichelen, S. V . de Velde, R. V. den Broecke, and E. Mortier, "Comparison of closed-loop controlled administration of propofol using BIS as the controlled variable versus "standard practice" controlled administration," Anesthesiology, vol. 95, pp. 6-17, 2001. [122] T. Sakai, A. Matsuki, P. White, and A. Giesecke, "Use of an EEG-bispectral closed-loop delivery system for administering propofol," Acta Anaesthesiologica Scandinavica, vol. 44, pp. 1007-1010, 2000. [123] A. Nayak and R. Roy, "Anesthesia control using midlatency auditory evoked potentials," IEEE Transactions on Biomedical Engineering, vol. 45, pp. 409-421, Apr. 1998. [124] G. Kenny and H. Mantzaridis, "Closed-loop control of propofol anaesthesia," British Journal of Anaesthesia, vol. 83, pp. 223-228, 1999. [125] J. Schiittler and H. Schwilden, "Present state of closed-loop delivery in anesthesia and intensive care," Acta Anaesthesiologica Belgica, vol. 50, pp. 187-191, 1999. [126] M . Unser and A. Aldroubi, "A review of wavelets in biomedical applications," in Proceedings of the IEEE, vol. 84, pp. 626-638, Apr. 1996. [127] S. Mallat, "Multiresolution approximations and wavelet orthonormal bases of L^iTZ)," Trans. Amer. Math. Soc, vol. 315, pp. 69-87, 1989. [128] I. Daubechies, "Orthonormal bases of compactly supported wavelets," Comm. Pure Appl. Math., vol. 41, pp. 909-996, 1988. [129] I. Daubechies, Ten lectures on wavelets, vol. 61 of CBMS-NSF Regional Conference Series in Applied Mathe-matics. SIAM, Philadelphia, 1992. [130] S. Mallat, "A theory for multiresolution signal decomposition: the wavelet representation," IEEE Trans. Pattern Anal. Machine Intell, vol. 11, pp. 674-693, 1989. [131] P. Abry and A. Aldroubi, "Designing multiresolution analysis-type wavelets and their fast algorithms," J. Fourier Anal. Appl, vol. 2, pp. 135-159, 1995. [132] G. Strang and T. Nguyen, Wavelets and filter banks. Wellesley-Cambridge Press, 1996. [133] G. Nason and B. Silverman, The stationary wavelet transform and some statistical applications, pp. 281-299. [134] S. Schiff, A. Aldroubi, M . Unser, and S. Sato, "Fast wavelet transformation of E E G , " Electroencephalography and Clinical Neurophysiology, vol. 91, pp. 442-455, Dec. 1994. [135] D. Donoho, "De-noising by soft-threshholding," IEEE Trans. Information Theory, vol. 41, no. 3, pp. 613-627, 1995. BIBLIOGRAPHY 112 [136] M . Akay, Diagnosis of coronary artery disease using wavelet-based neural networks, ch. 19, pp. 513-524. CRC Press, 1996. [137] L. Senhadji and F. Wendling, "Epileptic transient detection: wavelets and time-frequency approaches," Journal of Clinical Neurophysiology, vol. 32, pp. 175-192, 2002. [138] F. Sartoretto and M . Ermani, "Automatic detection of epileptiform activity by single-level wavelet analysis," Clinical Neurophysiology, vol. 110, pp. 239-249, Feb. 1999. [139] C. D'Attellis, S. Isaacson, and R. Sirne, "Detection of epileptic events in electroencephalograms using wavelet analysis," Annal of Biomedical Engineering, vol. 25, no. 2, pp. 286-293, 1997. [140] I. Clark, R. Biscay, M . Echeverria, and T. Virues, "Multiresolution decomposition of non-stationary E E G signals: a preliminary study," Computer in Biology and Medicine, vol. 25, no. 4, pp. 373-382, 1995. [141] I. Osorio, M . Frei, and S. Wilkinson, "Real-time automated detection and quantitative analysis of seizures and short-time prediction of clinical onset," Epilepsia, vol. 39, no. 6, pp. 615-627, 1998. [142] R. Sirne, S. Isaacson, and C. D'Attellis, "A data reduction process for long-term EEGs," IEEE Magazine in Engineering in Medicine and Biology, vol. 18, pp. 56-61, Jan-Feb 1999. [143] G. Benke, M . Bozek-Kuzmicki, D. Colella, G. Jacyna, and J. Benedetto, "Wavelet-based analysis of E E G signals for detection and localization of epileptic seizure," Proceedings of SPIE Conference on Wavelet Applications, vol. 2491, pp. 760-769, Apr 1995. [144] S. SchifF, J. Milton, and S. Weinstein, "Wavelet transforms and surogate data for E E G spike and seizure localization," Optical Eng., vol. 33, no. 3, pp. 2162-2169, 1994. [145] J. Benedetto and D. Colella, "Wavelet analysis of spectogram seizure chirps," in Proceedings of SPIE Confer-ence, vol. 2569b, pp. 512-521, 1995. [146] R. Sartene, L. Poupard, J.-L. Bernard, and J.-C. Wallet, Sleep images using wavelet transform to process polysomnographic signals, ch. 13, pp. 355-382. CRC Press, 1996. [147] J. Huang, Y . - Y . Lu, A. Nayak, and R. Roy, "Depth of anesthesia estimation and control," IEEE Transactions on Biomedical Engineering, vol. 46, pp. 71-81, Jan. 1999. [148] A. Angel, D. Linkens, and C. Ting, "Estimation of latency changes and relative amplitudes in somatosensory evoked potentials using wavelets and regression," Computers and Biomedical Research, vol. 32, pp. 209-251, 1999. [149] V . Samar, A . Bopardikar, R. Vao, and K. Swartz, "Wavelet analysis of neuroelectric waveforms: a conceptual tutorial," Brain and Language, vol. 66, pp. 7-60, 1999. [150] S. Bibian, T. Zikov, G. Dumont, C. Ries, E . Puil, H . Ahmadi, M . Huzmezan, and B. MacLeod, "Estimation of the anesthetic depth using wavelet analysis of electroencephalogram," in Proceedings of the 23rd Annual International Conference of the IEEE Engineering in Medicine and Biology Society, vol. 1, pp. 951-955, 2001. [151] T. Zikov, S. Bibian, G. A. Dumont, M . Huzmezan, and C. Ries, "A wavelet based de-noising technique for ocular artifact correction of the electroencephalogram," in Proceedings of the 24th Annual International Conference in IEEE Engineering in Medicine and Biology Society, 2002. BIBLIOGRAPHY 113 [152] Aspect Medical Systems, A-2000 BIS Monitor: Serial Port Technical Specification, 1998. [153] N . Chamoun, J. Sigl, and C. Smith, "US Patent #5,458,117," Oct 17 1995. [154] G. Gratton, M . Coles, and E. Donchin, "A new method for off-line removal of ocular artifact," Electroen-cephalography and Clinical Neurophysiology, vol. 55, no. 4, pp. 468-484, 1983. [155] R. Vergler, T. Gasser, and J. Mocks, "Correction of EOG artifacts in event-related potentials of the E E G : aspects of reliability and validity," Psychophysiology, vol. 19, pp. 472-480, 1982. [156] T. Gasser, L. Sroka, and J. Mocks, "The transfer of EOG activity into the E E G for eyes open and closed," Electroencephalography and Clinical Neurophysiology, vol. 61, pp. 181-193, 1985. [157] J. Kenemans, P. Molenaar, M . Verbaten, and J. Slangen, "Removal of the ocular artifact from the E E G : a comparison of time and frequency domain methods with simulated and real data," Psychophysiology, vol. 28, no. 1, pp. 114-121, 1991. [158] J. W. M . Verbaten and J. Slangen, "The removal of the eye-movement artifact from the E E G by regression analysis in the frequency domain," Biology and Psychology, vol. 16, no. 1-2, pp. 127-147, 1983. [159] R. Croft and R. Barry, "Removal of ocular artifact from the E E G : a review," Neurophysiology Clinics, vol. 30, no. 1, pp. 5-19, 2000. [160] R. Croft and R. Barry, "EOG correction: a new aligned-artifact average solution," Electroencephalography and Clinical Neurophysiology, vol. 107, no. 6, pp. 395-401, 1998. [161] P. Berg and M . Scherg, "A multiple source approach to the correction of eye artifacts," Electroencephalography and Clinical Neurophysiology, vol. 90, no. 229-241, 1994. [162] T. Lagerlund, F. Sharbrough, and N . Busacker, "Spatial filtering of multichannel electroencephalographic recordings through principal component analysis by singular value decomposition," Journal of Clinical Neuro-physiology, vol. 14, no. 1, pp. 73-82, 1997. [163] R. Vigario, "Extraction of ocular artifacts from E E G using independent component analysis," Electroen-cephalography and Clinical Neurophysiology, vol. 103, pp. 395-404, 1997. [164] T.-P. Jung, S. Makeig, C. Humphries, T.-W. Lee, M . McKeown, V . Iragui, and T. Sejnowski, "Removing electroencephalographic artifacts by blind source separation," Psychophysiology, vol. 37, pp. 163-178, 2000. [165] E. Ifeachor, B. Jervis, E. Morris, E. Allen, and N . Hudson, "A new computer-based online ocular artifact removal (OAR) system," IEE Proceedings A, vol. 133, pp. 291-300, 1986. [166] K . Rao and D. Reddy, "On-line method for enhancement of electroencephalogram signals in presence of electro-oculogram artefacts using non-linear recursive least squares technique," Medicine and Biology Engineering Coraput., vol. 33, pp. 488-491, 1995. [167] S. Selvan and R. Srinivasan, "Removal of ocular artifacts from E E G using an efficient neural network based adaptive filtering technique," IEEE Signal Processing Letters, vol. 6, no. 12, pp. 330-332, 1999. [168] R. Ksiezyk, K . Blinowska, P. Durka, W. Szelenberger, and W. Androsiuk, "Neural network with wavelet preprocessing in E E G artifact recognition," in Medicon'98, 1998. [169] E. Ifeachor, M . Hellyar, D. Mapps, and E. Allen, "Knowledge-based enhancement of human E E G signals," in IEE Proceedings in Radar and Signal Processing, 1990. BIBLIOGRAPHY 114 [170] D. Overton and C. Shagass, "Distribution of eye movement and eye blink potentials over the scalp," Electroen-cephalography and Clinical Neurophysiology, vol. 27, p. 546, 1969. [171] T. Eelbert, W. Lutzenberger, B. Rockstroh, and N. Birbaumer, "Removal of ocular artifacts from the EEG - a biophysical approach to the EOG," Electroencephalography and Clinical Neurophysiology, vol. 60, no. 5, pp. 455-463, 1985. [172] D. Donoho and I. Johnstone, "Ideal spatial adaptation by wavelet shrinkage," Biometrika, vol. 81, pp. 425-455, 1994. [173] The Mathworks Inc., MA, Matlab user's guide. Wavelet toolbox, 1997. [174] R. Coifman and D. Donoho, "Translation invariant de-noising," Lecture Notes in Statistics, vol. 103, pp. 125-150, 1995. [175] S. Bibian, T. Zikov, G. Dumont, C. Ries, E. Puil, H. Ahmadi,»M. Huzmezan, and B. MacLeod, "Method and apparatus for the estimation of the anesthetic depth using wavelet analysis of the electroencephalogram," Patent Pending, #60/395,313, July 12th 2002. [176] C. Burrus, A. Gopinath, and H. Guo, Introduction to wavelets and wavelet transforms: a primer. Prentice Hall, 1998. Appendix A Bases and Frames This Appendix has been adapted from [176, 6, 5]. For precise mathematical details or an overview, please refer to these references. A . l Bases A set of functions or vectors vk(t) spans a vector space V (or V = Span{vk}) if any element of that space can be fc expressed as a linear combination of members of that set, that is, all elements in V are of the form where k e Z and ak € C or TZ. An inner product1 {f(t),g(t)) and a norm \\f\\ = ^(f, f) are usually defined for such 1 A n inner product on a vector space V over C or TZ, is a complex valued function {•, • ) , defined on V x V wi th the following properties: (i) (x + y, z) = (x,z) + (y,z), (ii) (x,ay) = a{x,y), (iii) {x,y)* = {y,x), and (iv) {x,x) > 0, and (x,x) = 0 if and only if x = 0; where x, y and z € V , a e C or TZ, and the operator * denotes a complex conjugate. 2 O f part icular interest for signal processing is a space L^CR-) that consist of complex-valued finite energy signals defined on TZ. A function f(t) £ Li{1Z) satisfies f(t) = Y ak' Vk^ (A.l) fc The inner product between 2 functions on L<2(1Z) is given by ten where f*(t) denotes a complex conjugate of f(t). The norm of / is defined as: 115 APPENDIX A. BASES AND FRAMES 116 If the space contains not only all the functions that can be expressed by a linear combination of members of the spanning set, but also the functions which are the limits of these infinite expansions, then we say that the space V is the closure of the spanning set {vk(t)}, and we denote it by V = Span{vk}-k A set {vk(t)} is a basis set or basis for a given space V if the set of {ak} in Eq. A. l is unique for any particular }{t) £ V. This implies that the basis functions are linearly independent (i.e. J2k a k ' vk(t) = 0 only if all ak = 0). A basis set is orthogonal if (vk{t),vi(t)) = 0 for all k ^ I. For example, in three dimensional Euclidean space, orthogonal basis vectors are coordinate vectors that are at right (90°) angles to each other. The basis set is orthonor-mal if (vk(t),vi(t)) = 5(k — I). This means that in addition to being orthogonal, the basis vectors are normalized to unity norm, that is, ||^ fc(t)|| = 1 for all k. If we have an orthonormal basis, we can express any f(t) g V based on an inner product as: /(*) = £</(*).«*(*)> (A.2) k since by taking the inner product of Vk{t) with both sides of Eq. A. l , we get ofc = (/(t),ufc(t)>. (A.3) The calculation of the expansion coefficients such as in Eq. A.3 is called the analysis part, while the computation of the signal from the coefficients and expansion functions such as in Eq. A. l is called the synthesis part. Further, for an orthonormal basis, the Parseval's equality holds. This means that the norm or energy of a signal can be partitioned in terms of the expansion coefficients ak, that is, for every /(£) e V ll/WII2 = £K/ (0 ,M t )> | 2 = £|a*|2- (A.4) 12 \ak k k Although orthonormal bases are very convenient, the more general case of nonorthogonal bases is more flexible for some problems. In the case of a so-called biorthogonal basis, the Eq. A. l is still valid. However, a dual basis set {vk(t)} exists, whose elements are not orthogonal to each other, but to the corresponding element of the basis set {vk(t)}, that is, <«j(t). «*(*)> =S(k-l). (A.5) Thus, /(*) = £</(*).«*(*)> ' Mt) = J2&k' Vk^< (A-6) where ak = (f(t),vk(t)). The Parseval's equality becomes \\f(t)\\2 = Y^(Mt),f(t)r-(Mt)j{t))- (A-7) Although a biorthogonal system requires an additional or dual basis, it is more general and allows for a larger class of expansions. A .2 Frames The necessary condition for a set of functions to be a basis is the linear independence of that set, that is, the uniqueness of its expansion. This means that no element of the set can be written as a linear combination of the APPENDIX A. BASES AND FRAMES 117 others. If the set of functions or vectors is dependent, and yet allows the expansion given by Eq. A.6, which is not unique anymore, then the set is called a frame. Thus, a frame is, a spanning set. To be a frame, an expansion set {ipk(t)} must satisfy k for some 0 < A and B < oo and for all functions f(t) in the space. Dividing Eq. A.8 by !|/(i)||2 shows that A and B are bounds that "frame" the normalized energy of the expansion coefficients. If A = B then the expansion set is called a tight frame. In this case we have A||/W||2 = £ | ( ^ ( t ) , / ( t ) > | 2 , (A.9) fc which is the generalized Parseval's equality for tight frames. If A = B = 1, the tight frame becomes an orthogonal basis. For frames, there is no strict Parseval's theorem and Eq. A.8 shows that the energy in the transform domain cannot be exactly partitioned. In the case of a tight frame, the partitioning can be done exactly with Eq. A.9 and the signal can be expanded as follows: f(t) = A-1J2(Mt)J(t))-Mt), (A.10) k which is the same as the expansion using an orthonormal basis except for the A"1 factor. This factor is a measure of the redundancy in the expansion set (e.g. A = 2 means there are twice as many vectors as needed to cover the space). Thus, frames are an over-complete (redundant) version of a basis set, and tight frames are an over-complete version of an orthogonal basis set. In the case of a frame that is neither a basis nor a tight frame, a dual frame set exists so that the analysis and synthesis can be done in a similar way as for a nonorthogonal (biorthogonal) basis. In the case of a tight frame, calculations are very similar to those for an orthogonal basis. The use of frames and tight frames implies a certain amount of redundancy, which gives robustness to the expansion so that it is more error insensitive. However, redundancy also means computational inefficiency. A.3 Unconditional Bases A basis {vk(t)} is an unconditional basis for a space V if every convergent series of the form ^2kak ' vk(t) is uncon-ditionally convergent; that is, if every arrangement of its terms converges to the same element. An unconditional basis is a very robust basis where the coefficients drop off rapidly for all members of the space. This basis is further a bounded unconditional basis if there are constants 0 < A < B < oo such that A < \\vk\\ < B. The fundamental idea of bases or frames is to represent a continuous function by a sequence of expansion coeffi-cients. In the case of orthogonal bases and tight frames, the Parseval's equality relates the norm of the function to the norm of the coefficients (see Eqs. A.4 and A.9). For an unconditional basis, not only the norm of the function in the space can be related to some norm of the coefficients in the expansion, but the absolute values of the coefficients carry the sufficient information to establish the relation. Therefore, if only the norm of the function is of interest, there is no condition on the sign or phase information of the expansion coefficients. Appendix B Effect of Amplitude Normalization In this very short appendix, we wish to motivate our choice for normalizing the amplitude of EEG epochs prior to analysis. To do so, we are using a raw 4-channel EEG recording provided by Aspect Medical Inc. The recording is 30 minutes long. No information concerning the patient is available. The raw EEGs are plotted in Fig. B.l. 0 5 10 15 20 25 30 Time (min) Figure B. l : Four-channel EEG recording corresponding to a anesthetized patient. These data have been obtained using electrodes located at different position on the patient's scalp. It appears clearly that the channels 1 and 2 have larger amplitudes than the channels 3 and 4. Using the extracting technique described in Chapter 3, the WAV Index on normalized and non-normalized epochs has been obtained and is displayed in Fig. B.2. Clearly, since the four EEG signals have been recorded from the same patient under the same conditions, we should expect similar indexes. The normalization of the amplitude of each epoch gives a more consistent index. Another positive aspect is that the discriminating power of the PDF at level d\ is increased. As a result, the index is less noisy. 118 APPENDIX B. EFFECT OF AMPLITUDE NORMALIZATION 119 100 80 60 40 20 0 Normalized epochs 10 15 20 Time (min) 25 30 100 80 60 40 20 0 Non-normalized epochs 10 15 20 Time (min) 25 30 (a) (b) Figure B.2: WAV Index calculated on a 4-channel recording: comparison between normalized and non-normalized epochs.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Monitoring the anesthetic-induced unconsciousness (hypnosis)...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Monitoring the anesthetic-induced unconsciousness (hypnosis) using wavelet analysis of the electroencephalogram Zikov, Tatjana 2002
pdf
Page Metadata
Item Metadata
Title | Monitoring the anesthetic-induced unconsciousness (hypnosis) using wavelet analysis of the electroencephalogram |
Creator |
Zikov, Tatjana |
Date Issued | 2002 |
Description | This thesis presents a novel method for the discrimination between various levels of anesthetic-induced unconsciousness (i.e. hypnosis) in the context of intra-operative anesthesia. This is achieved by means of the wavelet analysis of the patient's Electroencephalogram (EEG). While monitors of hypnosis are already available, they rely on complex signal processing algorithms, and in some instances, on multivariate statistical analysis of a large quantity of data. To date, the BIS® monitor (Aspect Medical Systems, Inc., MA) provides the most accepted solution. It displays an index calculated based on bispectral analysis, and has been the subject of many clinical trials over the past 5 years. As such, it is widely known and accepted by the medical community, particularly in the United States. Wavelet analysis can be viewed as a generalization of Fourier analysis that introduces time localization in addition to frequency decomposition of a signal. This major advantage of wavelets makes them particularly suitable for the analysis of non-stationary signals such as the EEG. In this thesis, we present a simple wavelet-based technique that calculates, in real time, an index of hypnosis (so-called WAV Index) based on a patient's electroencephalogram. Wavelet analysis significantly reduces the computational complexity to perform this task. Further, neither a large number of patients nor an extensive amount of clinical EEG data is needed for the index derivation. In addition, a technique for the de-noising of the electroencephalogram to remove ocular artifacts is proposed. These artifacts, which result from ocular activity (eye blinks and movements), strongly perturb the EEG signal by superimposing large spikes and step-like patterns. Hypnotic states such as the awake state and light sedation are typically corrupted by these waveforms. The solution proposed successfully removes these artifacts from the EEG, thus providing a clean signal for the wavelet analysis. This technique further works in the case of head movement artifacts. Finally, clinical validation has been carried out to assess the performance of the proposed index. While being particularly well-correlated to the bispectral index, the proposed index has a faster response to fast transients and it is more consistent in terms of patient intra-variability. |
Extent | 12147208 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-09-22 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
IsShownAt | 10.14288/1.0065055 |
URI | http://hdl.handle.net/2429/13058 |
Degree |
Master of Applied Science - MASc |
Program |
Electrical and Computer Engineering |
Affiliation |
Applied Science, Faculty of Electrical and Computer Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2002-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_2002-0619.pdf [ 11.58MB ]
- Metadata
- JSON: 831-1.0065055.json
- JSON-LD: 831-1.0065055-ld.json
- RDF/XML (Pretty): 831-1.0065055-rdf.xml
- RDF/JSON: 831-1.0065055-rdf.json
- Turtle: 831-1.0065055-turtle.txt
- N-Triples: 831-1.0065055-rdf-ntriples.txt
- Original Record: 831-1.0065055-source.json
- Full Text
- 831-1.0065055-fulltext.txt
- Citation
- 831-1.0065055.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0065055/manifest