Method and apparatus for determining sleep states

An apparatus is provided for detecting Macro Sleep Architecture states of a subject such as WAKE, NREM and REM sleep from a subject's EEG. The apparatus includes an EEG digital signal assembly of modules arranged to convert analogue EEG signals into digital EEG signals. A bispectrum assembly is responsive to the EEG digital signal assembly and converts the digital EEG signals into signals representing corresponding bispectrum values. A bispectrum time series assembly, in electrical communication with an output side of the bispectrum assembly, generates at least one bispectrum time series for a predetermined frequency. A macro-sleep architecture (MSA) assembly is responsive to the bispectrum time series assembly and is arranged to produce classification signals indicating classification of segments of the EEG signals into macro-sleep states of the subject.

This application is a National Phase Application of International Application No. PCT/AU2009/001554, filed on Nov. 27, 2009, which designated the U.S., and claims priority to Australia Application No. 2008906167, filed on Nov. 28, 2008, the entire contents of each of which are hereby incorporated by reference.

FIELD OF THE INVENTION

The present invention is concerned with a method and apparatus for automatically determining states of sleep of a subject. Embodiments of the invention find particular application in aiding the diagnosis of sleep disorders and assessing daytime sleepiness.

BACKGROUND

The reference to the prior art in the following discussion is not to be taken as any representation or admission that such art forms part of the common general knowledge. The disclosures of each of the publications referred to herein are hereby incorporated by reference in their entireties and for all purposes.

Obstructive sleep apnea hypopnea syndrome (OSAHS) is a serious sleep disorder with high prevalence among the population [1, 2]. In the USA, about 24% of men and 9% of women fall within Medicare guidelines [1] for treatment. In Singapore 15% of the total population is at risk [3]. Over 1.2 million Australians experience sleep disorders costing the country $10.3 billion (in 2004) [4]; OSAHS is the commonest disorder (66% of the total).

OSAHS is characterized by breathing interruption during sleep. Full closure of the airways is known as obstructive Apnea, and a partial closure is defined as obstructive Hypopnea (See Appendix A for technical definitions). The common symptoms of OSAHS are excessive daytime sleepiness and intermittent snoring [5, 6].

OSAHS is a major risk factor for downstream complications such as stroke, diabetes and cardiovascular disease [6, 7]. It is also known to be associated with cognitive deficiencies, low IQ in children, fatigue and accidents. It is responsible [7] for 11,000-43,000 traffic accidents per year in NSW. Untreated patients are known to utilize twice the national health resources prior to diagnosis [8]. OSAHS is treatable. If diagnosed early, its devastating secondary complications can be thwarted. However, over 90% individuals with OSAHS are estimated to be undiagnosed at present [2].

The standard test for OSAHS diagnosis is Polysomnography (PSG) [9]. PSG is a technique to monitor multiple neuro-physiological and cardio respiratory signals, over the course of night. It requires a full-night sleep-laboratory stay in a specifically equipped sleep-suite, connected to over 15 channels of measurements. The 6-8 hours of sleep data is then subjected to a complex and time-consuming manual process (Sleep Scoring) to identify events of Apnea/Hypopnea and a type of sleep a disturbance known as EEG-arousals (EEGA). The outcomes of PSG test are summary measures of OSAHS severity such as the Respiratory Disturbance Index (RDI) and the EEG Arousal index (ArI) etc (please see appendix B for details).

EEG in the Diagnosis of Sleep Disorders

Sleep Scoring and Macro-Sleep-Architecture

Sleep is essentially a neuropsychological phenomenon; EEG still remains the cheapest and the most portable technique for the functional imaging of the brain during sleep. It is also the technique with the highest temporal resolution available. In the current practice of PSG Scoring, EEG is regarded as an indispensable signal when a definitive diagnosis is desired. Thus, in-facility diagnostic PSG tests always include EEG. Electromyography (EMG) and Electroocculography (EOG) signals are also needed for the correct EEG-centred interpretation of sleep states.

In diagnostic PSG tests EEGs are essential for the following tasks:(i) to define EEG-arousals and identify sleep fragmentation. EEG-arousals are also used as one parameter in defining Hypopneas (see Appendix A).(ii) to score the Macro Sleep Architecture (MSA) of sleep. It is a process in which sleep is classified into three macro states: (1) Wake State (SW), (2) Rapid Eye Movement (REM) Sleep State (SR), and (3) Non-REM Sleep State (SN). The MSA is extremely important in the diagnosis of OSAH. In addition, Sleep MSA may be used in the diagnosis/monitoring of a range of sleep disorders including Narcolepsy, Insomnia, Sudden Infant Death Syndrome and Depression etc.Some important uses of MSA in PSG includes:(a) Estimating the Total Sleep Time, TST, which is needed for the computation of the RDI index. The TST is also useful as a summary indicator of the quality of sleep during a PSG test. Note that the TST, which is defined using EEG, can be significantly different from the Total Time in Bed (TTB) measured with a clock.(b) Estimating clinically important descriptors of sleep such as the Sleep Efficiency (SE), Sleep Latency (SL), REM Latency (RSL), the total time spent in REM, and the percentage of time spent in REM.(c) Expressing most of the clinically relevant sleep parameters separately for REM and NREM sleep, before providing an overall number. Some Examples are: the Arousal Index in REM sleep, The Arousal Index in NREM sleep, RDI in REM sleep and RDI index in NREM sleep. The reason behind this is that the REM/NREM classification provides fundamental information about sleep and its diagnostic characteristics.

An EEG may be broadly divided into four major frequency bands [10], Delta (δ, 0.1-4 Hz), Theta (θ, 4.1-8 Hz), Alpha (α, 8.1-12 Hz), and Beta (β, >12.1 Hz).FIG. 1shows the EEG activity at different states, (a) awake drowsy state, (b) light sleep (NREM sleep Stage 1 and Stage 2), (c) deep sleep (NREM sleep Stage 3 and Stage 4) and (d) REM sleep. These frequency bands are heavily used in sleep scoring.

The scoring of MSA is done manually using the rules laid down by Rechtschaffen and Kales (R&K, 1968) (see appendix C for a summary) [11]. Manual scoring relies on visual extraction of specific features in two EEG channels (usually C3-A2 and C4-A1 of the International 10/20 system), two channels each of EMG and EOG. Thus, six channels of electrophysiological data have to be visually interpreted, simultaneously taking care of difficulties such as measurement artifacts.

This process is time consuming (typically 1-2 hours per patient), costly (hundred of dollar per recording) and prone to inter and intra scorer variability [12-15]. The scorers from different laboratories tend to agree less than scorers from the same laboratories, due to differences in interpretation and subjective implementations. For example, the mean epoch by epoch agreement between the scorers from three sleep laboratories in the USA for healthy subjects is 76% (range 65-85%) which decreases to 71% (range 65-78%) in the OSAHS cases [16]. A similar result (76.8%) has been reported by European laboratories based on a large database of 196 recordings from 98 patients [13].

The accurate computation of sleep parameters such as the RDI, ArI and REM Latency is important for the clinical diagnosis of a range of sleep disorders. Thus, the final diagnostic accuracy heavily depends on the precise scoring of MSA. However, due to the subjectivity associated with the scoring process there exists a significant variability in the PSG results between the technicians of the same laboratory and across the different sleep laboratories [14, 15]. For example, [12] reported that a patient can get two different diagnoses in two different laboratories which might range from as low as RDI=4.9 to as high as RDI=79.

In order to overcome the problems associated with manual scoring and cater to the ever increasing demand for PSG testing, several researchers have proposed automatic sleep scoring systems [17-22]. After publication of R&K's rules in 1968, several authors tried to automate them and achieved various degrees of agreement with human scorers. With the advancement in digital signal processing techniques, several other used frequency spectral analysis [22], neural network analysis [18, 20], multidimensional scaling and wavelets techniques or expert system approaches [21] to develop automatic sleep staging systems. However, a reliable and accurate method with sufficient precision suitable for in-facility PSG as well as other take-home OSAHS screening devices does not exist yet. Despite the inherent subjective nature, human scoring is still considered the golden method of MSA scoring. Existing methods for automatic MSA scoring have the following shortcomings:1. The R&K rules depend on visual features in sleep EEG and were originally proposed specifically for manual scoring. Most of the automated techniques try to implement R&K rules and depend on morphological features like k-complexes, vertex-waves and spindles. These characteristics are severely altered in disease states such as OSAHS [6, 10, 23, 24]. Consequently their performances decreases in OSAHS [19, 25]. In addition, the detection of visual features is a highly subjective process.2. The agreement between automatic and human classifications is smaller than the agreement between human scorers [25]. This result is to be expected because the R&K criteria are based on visual features and humans are better than machines in visual pattern recognition.3. The differentiation of REM/NREM/WAKE is critical in OSAHS diagnosis. However, automated methods had difficulties in distinguishing wake state from Stage 1 of NREM sleep and REM sleep.4. Automated scoring techniques currently available in PSG equipment need expert human intervention for manually editing the outcomes, and thus are not truly automated systems [25]. The human intervention makes them subjective and time consuming.5. Existing methods have not been tested under disease conditions such as OSAHS, Periodic Leg Movement Syndrome (PLMS) or upper airway respiratory syndrome (UARS), where sleep is corrupted with frequent EEG arousals, apnea events, and recording artifacts.6. All existing techniques depend on recording multiple physiological signals, making them unsuitable for portable monitors used for OSAHS screening.

The AASM definition [26] of micro-sleep is: “ . . . an episode lasting up to 30 seconds during which external stimuli are not perceived. The PSG suddenly shifts from waking characteristics to sleep”. It is generally believed that micro-sleep is closely associated with excessive diurnal sleepiness. Excessive daytime sleepiness and spontaneous micro-sleep are two major consequences of OSAHS [27], contributing to motor vehicle and work related accidents. It is estimated to affect 12% of the adult population [27].

Clinically, sleepiness is commonly expressed by the measure Sleep Latency (SL), which is the length of time required to fall asleep. The common tests for measuring SL are Multiple Sleep Latency Test (MSLT) and Maintenance of Wakefulness Test (MWT), technical details for which are provided in Appendix E. In these tests SL is computed as the time from the start of recording to the sleep onset. To technically identify sleep onset, sleep technician have to simultaneously look at multiple signals. It is a tedious and a subjective process resulting in high inter-rater, as well as intra-rater, variability [12]. The SL also provides valuable information in the diagnosis of other widespread diseases such as insomnia.

Even though diurnal micro-sleep is an important phenomenon related to sleep disturbances, there is no objective system of measurements to detect micro-sleep or express its severity. In routine PSG tests targeted for OSAHS diagnosis, episodes of micro-sleeps are not scored.

It is an object of the present invention to provide a method and apparatus that addresses one or more of the various problems discussed above in relation to prior art methods for determining sleep-related parameters.

SUMMARY OF THE INVENTION

According to a first aspect of the present invention there is provided a method of determining sleep states from an EEG signal of a subject, the method comprising the steps of:electronically processing the EEG signal to generate a third or higher order spectrum of said signal;electronically processing said spectrum to produce at least one spectrum time series for a predetermined frequency; andelectronically processing said spectrum time series for compliance with predetermined criteria to tangibly classify segments of the EEG signal as corresponding to particular macro-sleep states of the subject.

In one embodiment the segments of the EEG signal are tangibly classified by indicating their classification on an electronic display for viewing by a user.

Preferably the step of electronically processing the spectrum to produce at least one spectrum time series produces a first spectrum time series and a second spectrum time series for corresponding first and second frequencies.

Preferably the spectrum comprises a bispectrum and the spectrum time series comprises a bispectrum time series.

The step of electronically processing said bispectrum time series for compliance with predetermined criteria will preferably include processing the first bispectrum time series to classify a corresponding one of said segments as Wake or Sleep.

In a preferred embodiment the method involves a step of electronically processing the second bispectrum slice time series of said segment classified as Sleep to further classify said segment as NREM or REM sleep.

Preferably the step of electronically processing the bispectrum to produce at least one bispectrum time series produces a first bispectrum time series and a second bispectrum time series for corresponding first and second frequencies.

The step of electronically processing said bispectrum time series against predetermined criteria will preferably include processing the first bispectrum time series to classify a corresponding one of said segments as Wake or Sleep.

In a preferred embodiment the method involves a step of electronically processing the second bispectrum time series of said segment classified as Sleep to further classify said segment as NREM or REM sleep.

According to the preferred embodiment, the first frequency falls within a range of 1-8 Hz and the second frequency is greater than 9 Hz.

The method may include a step of electronically smoothing the second bispectrum time series prior to processing said second bispectrum time series against predetermined criteria.

The step of electronically processing the first bispectrum time series against predetermined criteria will preferably include comparing values of a segment of the first bispectrum time series to a threshold value computed from the first bispectrum time series.

Preferably the step of electronically processing the second bispectrum time series against predetermined criteria includes comparing values of a segment of the second bispectrum time series to a threshold value computed from the second bispectrum time series.

The method may include a step of indicating the particular macro-sleep states of the subject with an electronic display or another sensory signal modality such as a sound or tactile indication.

Preferably the step of generating a bispectrum by electronically processing the EEG signal comprises an indirect estimation method

The indirect method may include applying a Fourier transform to cumulant values corresponding to the EEG signal.

Alternatively, the step of generating a bispectrum by electronically processing the EEG signal may comprise a direct estimation method.

The method may include a step of electronically processing the first bispectrum time series to produce an index indicating sleepiness of the subject.

Preferably the step of electronically processing the first bispectrum time series to produce said index comprises determining a fraction of time over a predetermined period that said series approaches a predetermined threshold value.

According to a further aspect of the present invention there is provided an apparatus for detecting sleep states of a subject including:an EEG digital signal assembly of modules arranged to convert analogue EEG signals into digital EEG signals;a spectrum assembly responsive to the EEG digital signal assembly and arranged to convert the digital EEG signals into signals representing spectrum values;a spectrum time series assembly in electrical communication with an output side of the spectrum assembly and arranged to generate at least one spectrum time series for a predetermined frequency; anda macro-sleep architecture (MSA) assembly responsive to the spectrum time series assembly and arranged to produce classification signals indicating classification of segments of the EEG signals into macro-sleep states of the subject.

Preferably the spectrum assembly is a bispectrum assembly that is arranged to convert the digital EEG signals into signals representing bispectrum values. Alternatively, other higher order spectrum assemblies might also be used, such as a trispectrum.

Preferably the spectrum time series assembly comprises a bispectrum time series assembly arranged to generate at least one bispectrum time series for a predetermined frequency.

In one embodiment, the macro-sleep architecture (MSA) assembly is responsive to the bispectrum time series assembly.

Preferably the EEG digital signal assembly includes:one or more EEG electrode ports;an analogue signal conditioning module in electrical communication with the EEG electrode ports;an analogue to digital converter in electrical communication with the analogue signal conditioning module; anda bandpass digital filter arranged to process output from the analogue to digital converter to produce the digital EEG signals.Preferably the bispectrum assembly comprises:a cumulant calculator;a Fourier transform module responsive to the cumulant calculator to produce a signal representing values of a bispectogram; anda bispectrum time series estimator arranged to produce at least one bispectrum time series signal representing values at a predetermined frequency of a slice of the bispectogram.

The MSA assembly preferably includes a comparator to determine if the value of the at least one bispectrum time series signal exceeds a predetermined value over a segment.

Preferably the apparatus includes a sleepiness index calculator responsive to the bispectrum time series estimator and arranged to generate a signal indicating sleepiness of the subject based on the bispectrum time series signal for the predetermined frequency falling within the EEG Delta band.

In a preferred embodiment the apparatus has a display to display a macro sleep status determined by the MSA.

The user interface may include a display in communication with the sleepiness index calculator to display the sleepiness index.

Preferably the apparatus is arranged as a driver sleep management system, wherein the sleepiness index calculator is incorporated into a sleepiness index alert unit, said alert unit controlling one or more of:a driver intervention module for increasing the wakefulness of the driver;a wireless communications module for communicating with a base station;a vehicle control interface for appropriately and safely immobilising the vehicle; anda data logger for recording sleepiness index values of the driver.

According to a further aspect of the present invention there is provided a media readable by machine, tangibly embodying a program of instructions executable by the machine to cause the machine to perform the previously described method to determine sleep states from an EEG of a subject.

DETAILED DESCRIPTION OF PREFERRED EMBODIMENTS

Method

EEG Data Acquisition

The clinical data acquisition environment for this work is the Sleep Diagnostic Laboratory of The Prince Alexandra Hospital, Australia. The data was recorded using clinical Polysomnography (PSG) equipment (Siesta, Compumedics®, Sydney, Australia). Patient preparation, electrode placement and instrumental set-up were done by an experienced sleep technician according to AASM guidelines [9]. Table 3 describes the demographic details of the subjects studied.

Database A: From each subject in this database, routine PSG data was collected. In a typical PSG test, signals such as ECG, EEG, EMG, EOG, nasal/oral airflow, respiratory effort, body positions, body movements and blood oxygen saturation are monitored. EEG data was recorded from both hemispheres using electrode positions C4, C3, A2, and A1, based on the standard international 10-20 system of electrode placement [10] illustrated inFIG. 2A. The subject population includes individuals with symptoms such as daytime sleepiness, snoring, tiredness, lethargy etc and are suspected of sleep apnoea. They were referred to the hospital for a routine PSG test.

Database B: In this database, EEG data from 19 different positions (seeFIG. 13) on the skull were recorded in addition to the conventional PSG channels described under Database A above. Electrode locations were as specified in the Standard 10-20 system of EEG electrode placement [10]. The data recorded in this database were from healthy volunteers, defined as individuals without any symptoms of the OSAHS syndrome.

Database C: This database contains data from the subjects referred for the Multiple Sleep Latency Test (MSLT) [28] or Maintenance of Wakefulness Test (MWT) [6]. Typical physiological signals which are recorded in MSLT and MWT are (i) 4 channels of EEG (C3, C4, O1, O2 with reference to A1 and A2). (ii) left and right eye Electro-oculograms (EOGs), and (iii) submentalis electromyogram (EMG).

The collected data was segmented into sub-records of length M samples for further analysis. Note: In the field of sleep medicine these segments are conventionally termed as epochs. So from here onward they will be referred to interchangeably as epochs or segments.

HOS-Based Analysis of EEG

Let xi(k), i=1, 2, 3 . . . , N, denote the k-th sample of the i-th segment of digitized EEG data where N is the total number of segments in a recording. The samples xi(k) are modelled as:
xi(k)=hi(k)=ei(k)+wi(k)  (1)
where ei(k) is a white non-Gaussian process, and hi(k) is a stable, possibly non-minimum phase kernel representing the underlying system generating the EEG segment xi(k). The term wi(k) represents measurement noise within the frequency band of interest, which is traditionally modelled as a white Gaussian process. The inventors have relaxed this constraint and allow the measurement noise to be either a white or colored Gaussian noise, or, any noise process with a symmetrical probability density function.

A method according to an embodiment of the present invention will now be described that discounts noise processes ei(k) and wi(k) and capture the features of the underlying EEG generating mechanism via the kernel hi(k) thereby enabling classification of EEG into REM/NREM/WAKE stages, as well as enable us to define a measure of sleepiness. The inventors achieve that in the 3rd-order spectra domain in order to keep the phase information in hi(k) intact. While the third order spectra is preferred by the inventors and an implementation discussed on the third order spectra is described in detail herein other spectra of order greater than two may also be used. For example the fourth order spectrum, i.e. the trispectrum might be employed.

The Bispectrum Estimation

The third order spectrum of a signal is known as the bispectrum. The bispectrum can be estimated as the 2D Fourier transform of the third order cumulant sequence Cxi(τ1,τ2) [29]. Transforming (1) to the Cumulant domain, we obtain:
Cxi(τ1,τ2)=E{xi(k)·xi(k+τ1)·xi(k+τ2)}  (2)
where the sequence xi(k) is assumed to be zero mean. (If this condition is not satisfied, we can easily subtract the mean of the signal to make it zero mean). The bispectrum Bxi(ω1,ω2) of xi(k) is given by,
Bxi(ω1,ω2)=γ·H(ω1)·H(ω2)·H*(ω1+ω2)+Bwi(ω1+ω2)  (3)
where H(⋅) is the one dimensional Fourier transform of xi(k), γ is the skewness of ei(k), and Bwi(ω1,ω2) is the bispectrum of the noise process wi. Since the bispectrum of a Gaussian process is zero, Bwi(ω1,ω2)=0, leaving a high SNR signal:
Bxi(ω1,ω2)=γ·H(ω1)·H(ω2)·H*(ω1+ω2)  (4)
which depends on the EEG system response hi(k).

The bispectrum obtained using (4) will be a complex number. Unlike the power spectrum (2ndorder statistics) based on the autocorrelation, bispectrum preserves Fourier phase information. Thus, the inventors have conceived that the EEG System Response hi(k) can be estimated from the bispectrum keeping the true Fourier phase intact. In contrast, power spectrum (or autocorrelation based) techniques lose phase information, and the EEG System Response estimated from it will be the minimum-phase equivalent of the original response.

Equation (4) points out that Bxi(ω1,ω2) carries information on the EEG system response which, the inventors have realised, provide diagnostic features to score sleep. Note that the response H(ω) appears in (4) as a non-linearly transformed quantity. Bxi(ω1,ω2) is a multi-dimensional signal, and using it to derive features for staging sleep can be difficult.

In [30] the problem of signal reconstruction from the bispectrum was considered, and it was proved that any 1-dimensional slice of the bispectrum carries sufficient information to estimate a system response within a time-shift, as long as the chosen slice is not parallel to any one of the frequency axes or to the diagonal at 135 degrees. The inventors have understood that such a result may be utilized, in a preferred embodiment of the invention, to reduce the computational complexity of the HOS techniques and to identify features easily for sleep staging. The flexibility offered by the choosing of arbitrarily oriented and shifted oblique slices is also advantageous in avoiding unfavourable regions in the Bispectrum.

In the frequency domain, a quantity Pi(ω;ϕ,ρ) can be defined for each data segment xi(k) such that Pi(ω;ϕ,ρ)=Bxi(ω, ϕω+ρ) describes a one-dimensional slice inclined to the ω1-axis at an angle tan−1ϕ and shifted from the origin along the ω2-axis by the amount ρ, (−π<ρ<π) [30].

2.22 The Slice-Bispectrogram and the Bispectrgram Time Series (BTS)

The slice Pi(ω;ϕ,ρ) carries complete information on the EEG system response (i.e. the underlying EEG generating system) according to the model we have adopted. Thus, in order to describe the overnight data in a graphical way, we define a matrix SB(ω;ϕ,ρ) such that:
SB(ω;ϕ,ρ)=[P0(ω;ϕ,ρ)|P1(ω;ϕ,ρ)| . . .Pi(ω;ϕ,ρ), . . .PN−1(ω;ϕ,ρ)]  (5)
where the ithcolumn of SBrepresents a vector [Pi(ω;ϕ,ρ), −π<ω<π]T. We call this matrix the Slice-Bispectrogram. When displayed as an image, it illustrates the time evolution of the slice Pi(ω;ϕ,ρ) during the course of the night.

The inventors have conceived that features to stage sleep into REM/NREM/WAKE stages can be derived from SB(ω;ϕ,ρ). A set of times series, called Bispectrum-Time-Series (BTS, ξf), is formed by considering a set of fixed ω, i.e. ω=ω0, ω=ω1, ω2. . . ωN−1as follows:

Symbol ξfrepresents the Bispectrum-Time-Series at frequency f. The inventors have found that it is possible to choose particular values for ‘f’ such that the Bispectrum-Time-Series, ξfcarries sufficient information to characterize macro-sleep-stages. In Section 3, the implementation details to estimate MSA will be explained with derived results.

Results and Discussion

Implementation details of a method according to a preferred embodiment of the invention will now be described. The results of the method will then be compared to a set of clinical data with a view towards exploring the method's performance.

FIG. 2shows a block diagram according to an embodiment of the present invention. The following discussion will describe the method which entails processing EEG signals, or EEG signal data, to produce indications as to which MSA states segments of the EEG signals correspond to and in one embodiment a sleepiness index. The various steps of the method are programmed as tangible instructions in a computer readable media for execution by one or more processors of the computer. The computer may include a suitable analog to digital converter to receive EEG signals from electrodes that are in contact with a patient. Alternatively, the computer may process previously logged EEG signals stored on a magnetic, electronic memory or optical medium for example.

Implementation Details

EEG is a low-frequency signal and the frequency band of present interest is contained within 1-45 Hz. Thus xi(k) is filtered using a 5thorder, zero-phase digital Butterworth bandpass filter f(k) with lower cut-off frequency fi=1 Hz and higher cut-off frequency fh=45 Hz to remove out-of-band noise, including the ubiquitous power line interference at 50 Hz. All results shown herein have been derived with filtered EEG data. Let the filtered segments xi(k) be denoted by yi(k).

B). Estimation of the HOS:

The bispectrum can be estimated via estimating the 3rdorder cumulant (see (2)) and then taking a 2D-Fourier transform (see (3)), this method, known as the indirect method of estimating the bispectrum, is used in the preferred embodiment of the invention discussed herein. The procedure to form the estimate Cyi(τ1,τ2) of the cumulant Cyi(τ1,τ2) of yi(k) is outlined in steps (S1)-(S4) below.

(S1) Segment yi(k) into J records of length L samples each. Subtract the mean of each record to form the zero-mean sub-segments yij(k). Estimate the 3rdorder cumulantsCyij(τ1,τ2) of each sub-segment yij(k) using (6) as defined in [30].

where Q is the length of third-order correlation lags considered in the computation.

(S3) Apply a bispectrum window function to the overall cumulant estimateCyi(τ1,τ2) to obtain the windowed cumulant function Cwyi(τ1,τ2). We used the minimum bispectrum-bias supremum window described in [29] for the purpose.

(S4) The bispectrum Byi(ω1,ω2) of the segment yi(k) was estimated as the 2-D Fourier transform of the cumulant estimateCwyi(τ1,τ2).
Byi(ω1,ω2)=ΣT1=−cT1=+cΣT2=−cT2=+cCyi(τ1,τ2)e−j(τ1e1+τ2e2)(8)

It should be noted however that the bispectrum slice may also be estimated directly in the frequency domain as
Byi(ω1,ω2)=γ·Yi(ω1)·Yi(ω2)·Y1*(ω1+ω2)  (8A)
where Yi(ω1) is the 1-Dimensional Fourier transform (1D FFT) of the measured, filtered EEG signal and y is the skewness of yi.

FIG. 3shows typical mesh and contour plots of the bispectrum magnitude Abs(Byi(ω1,ω2)) during different macro sleep states, SW, SNand SR. Abs (⋅) stands for the computation of the modulii of each and every entry. The bispectra shown inFIG. 3were computed from short EEG segment of length 10 s, randomly taken from the EEG data of subject ID 5 (Table 3). The salient features seen inFIG. 3were consistently observed across the subjects in Databases A, B and C.

According to mesh plots in theFIG. 3, the bispectrum magnitude varies during different sleep states. It is low during the SWand SRstates and increases considerably during the SNstates. A hypothesis test (two tailed; student-t statistic) showed that the mean bispectrum magnitude during SNsleep state was different from mean bispectrum magnitude during SWor SRstates at the level of significance σ=0.001.

The contour plots of the bispectrum magnitude in theFIG. 3shows the location of the peaks at different sleep states. The comparison of contour plot for SWor SRwith that for SNstate, shows that, bispectrum during SWstate have multiple peaks at different frequencies (ω≈1-20 Hz) whereas peak in bispectrum magnitude during NREM state is just concentrated in narrow frequency band (ω<5 Hz). Bispectrum during SRshows multiple peaks in low frequency region (ω<10 Hz), similar to SWstate however, the bispectrum peak at high frequency (ω>15 Hz) completely diminishes during REM states. These features in the EEG bispectrum were consistently seen over the whole night sleep EEG data.

C). Estimation of the Slice-Bispectogram:

Further implementation details of the method described in Section 2.22 will now be discussed. Without a loss of generality (see [30]), set ϕ=1 and set ρ=0 so that the slice of the bispectrum considered is inclined to the ω1-axis by 45 degrees and passes through the origin (i.e. the line described by ω1=ω2in the (ω1,ω2)-plane) as symbolized by SB(2πf;1,0).

FIG. 4(a)shows the slice bispectogram magnitude (Abs(SB(2πf;1,0)).FIG. 4(b)shows corresponding sleep states scored manually by a sleep technician. Y-axis ofFIG. 4(a)is the frequency axis, f=1-40 Hz. X-axis ofFIG. 4is the epoch number. InFIG. 4we display the first contiguous 300 epochs out of the total recorded sleep epochs (N=810) for that particular subject, for the sake of display clarity.

FIG. 4(a)clearly illustrates the variation in the bispectogram magnitude in the different frequency bands with the change in sleep states, consistently over the night. During the SWstates bispectogram shows comparably high magnitude in the higher frequencies. With the appearance of SNstates magnitude in the higher frequencies (f>10) decreased whereas that of the lower frequencies (f<10 Hz) increased. Again as the sleep state changed to SRstates the magnitude in all the frequencies decreased.

D). Estimation of the Bispectrum-Time-Series:

From the Slice-Bispectogram (SB) the, Bispectrum-Time-Series (BTS, ξf) was estimated at two frequencies SB(2πf0;1,0) and SB(2πf1;1,0). For the results reported herein, the inventors set f0=2 Hz and f1=20 Hz for their ability to discriminate REM/NREM/WAKE states. Note that the entries of BTS are complex valued. The following definitions are used from hereon: ξ2=Abs(SB(2πf0;1,0)), and ξ20=Abs(SB(2πf1;1,0)). Once again Abs(⋅) stands for the computation of the modulii of each and every entry in BTS. The inventors have found that in a preferred embodiment, ξ2and ξ20can be used efficiently and accurately to identify different macro-sleep states (in MSA) and also to define a sleepiness index via the time-density of micro-sleep events. To make the time-series a dimensionless quantity, normalisation is applied using (8B) before using them for classification.

In the case of ξ2, pre-processing (smoothing, de-trending and equalization) is applied before using the time series in classification work in order to improve the performance; ξ20, however did not require such processing. To reveal the slow changes and remove the outliers from ξ2, ‘Loess Smoothening Method [31, 32] is used, which is based on local regression using weighted linear least squares and a 2nddegree polynomial model. De-trending with a least-squares-fit straight line was used to remove the trend in ξ2. Histogram equalization is often used in image processing to increase the contrast of the image. In the presently described embodiment, histogram equalization is applied to ξ2after the steps of smoothening and de-trending to obtain the time series ξ′2.FIG. 5shows the effect the pre-processing on ξ2.

FIG. 6(a)andFIG. 6(b)show the ξf, at frequencies 2 Hz and 20 Hz respectively.FIG. 6(c)shows technician scored macro sleep states. EEG data xi(k) for this figure is taken from subject ID 6 of Database A, Table 3.FIGS. 6(a), 6(b) and 6(c)illustrate explicitly the characteristics of Slice-Bispectogram. The magnitude of 2 Hz component varies in a cyclical fashion, synchronously with the macro sleep states. Magnitude is high during the SNsleep states and it decreased during SRand SWstates. The magnitude of the 20 Hz component consistently remained low during the sleep; however it increased considerably with the episodes of SWstates during the night. Thus, ξ20is a preferred platform for separating wake from other sleep states.

The steps for estimating the MSA from the two bispectrum time series ξ20and ξ′2according to the preferred embodiment, will now be described. The following method is also illustrated in the flowchart ofFIG. 7.

E). MSA Scoring Algorithm Based on Bispectrum Time Series ξ′2and ξ20

Bispectrum time series, ξ20and ξ′2were used to classify each data segment xi(k) (see (1)) into the classes WAKE (SW)/REM (SR)/NREM (SN) sleep. The length of data segments (M) were set to the standard ‘epoch length’ of M=30 s as used in routine sleep scoring.

The classification of the segment xi(k) into the two categories WAKE/SLEEP can be done quite easily based on ξ20(seeFIG. 6(b)). However, it cannot be used to differentiate between SLEEP states NREM and REM. The series ξ′2can be used easily to identify NREM segments, but cannot be used to separate WAKE from REM stages. According to the preferred embodiment of the invention, both ξ20and ξ′2are used together in a sequential decision process to classify xi(k) into the three categories SW/SR/SNsleep.

The identification of sleep states was done based on direct threshold operations on the ξ′2and ξ20as described in steps (T1)-(T3) below:(T1) WAKE/SLEEP Classification: For classifying each i=1, 2, 3 . . . N, segments of EEG data xi(k) into SWand SLEEP (SR/SN) states we used ξ20. We set a threshold ε20computed from ξ20. The condition ξ20≥ε20was tested for all the N segments. If condition ξ20≥ε20is true in a segment then that segment was classified as SWelse as SLEEP.(T2) NREM/REM Classification: After the step (T1), the EEG data segments are separated into two states WAKE and SLEEP. To further classify sleep segments into NREM and REM sleep, used ξ′2is used. The threshold ε2is previously computed from ξ2. The condition ξ2≥ε2is tested in each segment previously classified as Sleep in step (T1). If the condition is satisfied then that segment is classified into NREM sleep else into REM sleep. The method for calculating proper thresholds (ε2and ε20) is described in Appendix F.(T3) Estimation of Sleep Parameters: From the estimated sleep states of each patient, we computed TST, SE, RSL and percentage of SW, SRand SNsleep, using the standard definition given in Appendix B.
3.2 Performance Evaluation of the MSA Scoring Method on Clinical Data

For the purposes of this specification, macro sleep states estimated using the preferred embodiment of the present invention are termed the HOS based Estimated Sleep States (HESS), and the reference states estimated by the expert human scorer is called the Technician Scored Sleep States (TSSS). The TSSS will be considered the ground truth, against which the performance of the HESS will be compared.

The performance of the HESS was evaluated on a clinical database of EEG/PSG data (see Table 3.1, Database-A), and compared with that of the TSSS estimate. In order to test the stability and reliability of the algorithm, subjects with varying degrees of OSAHS severity were included. The test subjects population (23 subjects) had the mean RDI of 27.23±26.24, ranging from as low as 0.6 and as high as 94.4. The mean arousal index (ArI) for the subject population was 29.05±19.23.

To compute the HESS, EEG data from only one channel is required. In Database A, EEG data were recorded with the electrode positions C4, C3, A1 and A2 based on the standard international 10-20 system of electrode placement. The results presented in this section are based on the EEG data from the electrode position ‘C4’ with reference to ‘A2’.

To compute the TSSS, an expert human Sleep Scientist manually scored PSG data based on R&K criteria [7] and clinical practice parameters formulated by the AASM. As known from previous studies, the intra-scorer variability for the Sleep Scientist was 84%. According to the TSSS scoring, of the total of 20,714 epochs from the 23 subjects, 3508 were in SW, 14098 were in SNand 3108 were in SR.

The performance of the HESS in scoring macro sleep states was evaluated in several different ways. Using TSSS as the ground truth, the following performance measures were computed:Sensitivity, Accuracy and the Positive Predictive Value (PPV).The agreement between TSSS and HESS using Cohen's Kappa statistic [31] (please see Appendix D on how to interpret Kappa values).Pearson's correlation coefficient and Altman-Bland plots [34] were used for the descriptive analysis of TST, SE, RSL and percentage of SW, SN, SRsleep computed from HESS.Hypothesis testing to establish the statistical significance of the results at the level of significance σ=0.001.

FIG. 9 (a)-(b)andFIG. 10 (a)-(b)shows the PSG scored sleep states and automated estimated sleep states after Step (T2) of section 3.1 (E), for two subjects (Subject ID=6 and 19). According to these figures there is high resemblance between the technician scored sleep states and automated sleep states. The agreement between TSSS and HESS in classifying N epochs (N=810 for subject 6 and N=1018 for subject 19) into SW, SNand SRwas 78% (kappa 0.61) and 80% (kappa=0.69) respectively forFIG. 9 (b)andFIG. 10 (b).

InFIG. 9(b)andFIG. 10(b)a few brief episodes of REM are seen, which are not REM sleep episodes according to TSSS. From the literature it is known that the duration of a typical REM episode is between 10-20 minutes in first sleep cycle which increases in the latter half of the night. Hence, in order to decrease False positive prediction of the REM episodes and to increase the accuracy of HESS a REM-continuity rule was used. According to this rule those periods of REM sleep were targeted which were for less than D mins. These periods of SRsleep were reclassified into either SNstate or SWstates. For the reclassification a new threshold ε′20was defined from the ξ20time series. Appendix F shows the method to compute D and ε′20. In the preferred embodiment of the present invention D=8 min and ε′20=min(ξ20). If condition ξ20=ε′20was satisfied then that epoch was classified as SWelse as SN.

FIG. 9(c)andFIG. 10(c)shows the HESS results after applying REM combining rule. After the rule agreement between TSSS and HESS increased to 81% (kappa=0.66, σ=0.001) and 86% (kappa=0.71, σ=0.001) respectively forFIG. 9andFIG. 10. REM combining rule was used in all the further classification.

Table 5 shows the contingency table, comparing the scoring performance of the HESS with that of TSSS. Table 6 gives the sensitivity, accuracy and PPV statistics for HESS computed from the Table 5. The overall agreement between TSSS and HESS in separating the 20,714 epochs into the two states Sleep (SNand SRpooled together) and Wake is 83.78% (kappa=0.48, σ=0.001, z=61). Also, there was a high agreement between TSSS and HESS in classifying SW, SNand SRstates, of 75.17% (kappa=0.54, σ=0.001, z=94). Table 7 lists the agreement and kappa values between TSSS and HESS for all the 23 subjects. The mean agreement between TSSS and HESS for 23 subjects was 77.6%±7.2% and mean kappa value was 0.58±0.08.

In the PSG test several sleep parameters such as, Total sleep time (TST), Sleep efficiency (SE), REM sleep latency (RSL) and percentage of SW, SN, SRsleep are directly computed from the scored MSA. Accurate estimation of these parameters substantially depends on the precise classification of macro sleep states. These parameters were computed for each patient from the HESS and compared with those computed using TSSS.

FIG. 11consists of scatter plots, showing the relationship between the parameters TST, SE, RSL and percentage sleep computed using TSSS and HESS. The correlation for TST, SE and RSL was 0.84(σ=0.001), 0.62(σ=0.001) and 0.79(σ=0.001) respectively.

FIG. 12show the Altman-Bland plot for TST, SE, RSL and percentage sleep. These plots show that there is a good agreement between TSSS and HESS in estimating MSA parameters with very small insignificant bias. The bias for TST was −16 minutes (σ=0.32) and for SE was only −3.7% (σ=0.24). The bias for RSL and percentage of SW, SN, SRsleep were also low and insignificant, 14 min (σ=0.26), 1.3% (σ=0.34), −4% (σ=0.13) and 2.8% (σ=0.10) respectively.

3.3 Performance of HESS Algorithm with Other Electrode Position

The performance of the previously described HESS method according to a preferred embodiment of the invention will now be discussed, when EEG from electrode locations other than the traditional (R&K-recommended, 1968) C3, C4, A1 and A2 positions are used. For this data from the subjects in database B, table 3.2 was used. The International 10/20 electrode locations (seeFIG. 13), were investigated with a view towards identifying the best electrode locations suitable for single-channel EEG-based sleep scoring.

Following the customary practice in sleep EEG, site Fpz was used for the ground electrode. The reference electrode for measurement sites located on the left hemisphere of the scalp was A1; similarly, measurements on the right hemisphere used A2. This choice of the ground and reference electrodes comply with the R&K recommended measurements on C3 and C4. This makes a direct comparison between the performances of traditional and the HESS-based techniques possible.

FIG. 14shows the estimated sleep states using HESS algorithm using EEG data from the electrode locations Fp2, F3, Fz and P3over the scalp. Table 8 summarizes the statistical results for all electrodes that were used. A statistically significant (σ=0.0001) agreement between TSSS and HESS for all the electrode positions was found.

The time series ξ20proved highly successful in classifying EEG epochs into the two groups Sleep and Wake. Wake states correspond to high values of ξ20; sleep states are associated with a consistently low magnitude (s0) approaching zero for all practical purposes. It was noted that the magnitude of ξ20gradually moves from a high value towards s0as a person is falling asleep. Before finally settling down to s0corresponding to the state Sleep, the magnitude of ξ20briefly touches s0several times (seeFIG. 15 (a)). Such events are associated with EEGs that are characteristic of Sleep, and the inventors define them as micro-sleep events.

According to a preferred embodiment of the present invention, the stability of the ξ20time series to identify micro-sleep events and then use them to form a measure of sleepiness, termed the Sleepiness Index (SI). We define the SI as the fraction of the time the magnitude of ξ20maintained its value corresponding to sleep, i.e., s0, computed over a time frame of 10 seconds. Thus, SI can vary from 0 (no episodes of micro-sleep during the current 10 s period) to 1 (micro-sleep/sleep events completely covers the current 10 s period).

This definition allows the SI index to be computed real time, making it a useful tool to monitor the sleepiness of an individual. Sleep Onset can be defined as the instant at which SI reaches 1, and the Sleep Latency can be defined as the time duration from the ‘lights out’ to the Sleep Onset. The SI can be applied to automate the scoring of MSLT and the MWT tests, which are two important diagnostic tools used in general sleep medicine. In addition to the use in sleep medicine, the SI-index has the potential to be used in situations such as the sleepiness monitoring of long-range truck drivers; to facilitate this use, SI has been defined in such a way that it can be computed in real-time.

InFIG. 15(a)we show the time series ξ20where a person is falling sleep. The red bars under the series ξ20curve indicate micro-sleep events. The TSSS sleep states are shown inFIG. 15(b). Note that the sleep technologist has not attempted to identify states of micro-sleep.FIG. 15(a)andFIG. 15(b)graphically illustrate the close correspondence of the ξ20and the sleep/wake states of an individual. Moreover, it is seen that the gradual development of sleep over time is captured well by the gradual change in ξ20.

To test the capability of the SI to estimate the sleep onset and the sleep latency (SL), we computed SI for the subjects in Database C. In all the computations, the segment length was M=30 s, and the segment overlap was set to 29 s. Thus, the sleep state was assessed every second, based on the data of duration 30 s. The SI index was computed using the sleep state information covering a period of 10 s.FIGS. 17 and 18shows the SI and TSSS for the subject ID 28. Data inFIG. 17andFIG. 18is from the 2ndand 3rdnap of the MSLT test (see Appendix E for details of the MSLT test). Table 9 gives the SL as computed from the SI and TSSS for all the 4 naps in the MSLT for all the subjects in database C.

According to Table 9, the SL computed via the SI closely matches the Sleep Technologist identified Sleep Latency. The SI based approach, however, has been fully automated and provides consistent results. It does not depend on multiple physiological signals; only one channel of EEG is sufficient for the purpose. On the contrary, Sleep Technologist requires multiple signals and depends on subjective methods to estimate the SL.

Finally, it is of great interest to explore how the micro-sleep events defined via ξ20correspond with micro-sleep defined by the AASM. As an illustration, inFIG. 16we show the EEG (C4-A1, C3-A2), EOG (left and right eye) and EMG data from the epochs35,36,37, and38as marked onFIG. 15(a). Epochs35and36ofFIG. 16clearly follow the micro-sleep characteristics as defined by AASM, whereas epochs37and38do not. The micro-sleep defined by AASM agrees with the micro-sleep events defined via ξ20according to the preferred embodiment of the present invention.

The SI index that has been described can be used to identify the level of sleepiness of an individual, on a real-time basis. This is expected to contribute significantly to the monitoring of sleepiness of truck drivers etc. with a view towards timely intervention to prevent events of micro-sleep/sleep and accidents.

4.1 A Dedicated Sleep Architecture Estimator

Referring now toFIG. 19there is depicted a block diagram of a Sleep Architecture Estimator1according to an embodiment of the present invention. The following embodiment makes use of circuit blocks that are arranged to analyse bispectra, however as previously discussed, other higher order spectra might also be used such as the trispectrum. It will be realised that the blocks that are depicted are provided as separate interconnected circuit modules, as shown. Alternatively in another embodiment they may be implemented as virtual modules by a machine such as a suitably programmed high speed computer, preferably including one or more dedicated digital signal processing integrated circuits. In the latter alternative an aspect of the invention encompasses a program storage device, for example an optical or magnetic disk or memory circuit. The program storage device is readable by the machine and tangibly embodies a program of instructions executable by the machine to cause the machine to perform a method according to the present invention, for example preferably the method described above in sections 2. and 3. The computer preferably includes a suitable analog to digital converter with ports to receive EEG electrodes for monitoring a subject. The computer also includes a graphical display for displaying classifications of sleep that are determined during execution of the program for viewing by a human operator. A keyboard and mouse are also included in order that the operator can vary parameters of the program, such as thresholds for example.

The Sleep Architecture Estimator1includes ports for connection of EEG electrodes2. The ports are connected to an analogue signal conditioning module2, which contains circuitry receiving and conditioning the EEG signals. This circuitry typically includes AC coupling to high input impedance differential amplifiers, noise filtering and limiting. The conditioned analogue signal is passed to an analogue to digital converter module6which includes an antialising filter. The digital signal from ADC6is passed to a bandpass digital filter8to remove out-of-band noise, including mains power distribution frequency hum.

Modules2to8, as described above may collectively be referred to as an EEG digital signal assembly3, for converting analogue EEG signals into digital EEG signals.

A Segmentation Module10is coupled to the output of the bandpass digital filter8in order to receive the filtered digital signal. Segmentation Module10is arranged to segment the filtered digital signal from the bandpass digital filter8into a number of records, each comprising a plurality of samples. The Segmentation Module10is further arranged to subtract the mean of each record to form a digital signal representing zero-mean sub-segments. A 3rdOrder Cumulant Calculator12is coupled to the output of segmentation module10. The Cumulant Calculator12is arranged to produce a signal representing an estimate of the 3rdorder cumulants of each sub-segment. An Averaging Module14receives the 3rdorder cumulant signal and processes it to produce an overall cumulant estimate signal.

A Bispectrum Window Processor16receives the overall cumulant estimate signal and is arranged to produce a windowed cumulant signal in accordance with a minimum bispectrum-bias supremum window.

The output of the Bispectrum Window Processor16is coupled to a 2-D Fourier Transform Module18, which is arranged to process the windowed cumulant signal in order to produce a bispectrum signal that represents an estimate of the bispectrum of the corresponding segment.

Modules10to18, as described above, may collectively be referred to as a bispectrum assembly13to convert the digital EEG signals from bandpass digital filter8into signals representing bispectrum values.

In an alternative embodiment the bispectrum assembly is arranged to directly determine the bispectrum estimate according to equation (8A) as previously described. According to the direct calculation embodiment of the invention the Cumulant Calculator12is not required. Direct method estimation is computationally much less expensive, but has somewhat higher estimation variance (leading to poorer performance).

A Slice Bispectogram Estimator20receives the output from the Fourier Transform Module18and is arranged process the bi-spectrum signal to calculate values for a slice of the bispectrum. As previously described, the particular slice of the bispectrum that is favoured is the line described by ω1=ω2in the (ω1,ω2)-plane. The Slice Bispectogram Estimator calculates the values for the slice and generates a representative digital signal that is output to the Bispectogram Time Series Estimator22.

The Bispectrum Time Series Estimator22is arranged to process the output from the Slice Bispectogram Estimator in order to produce two Bispectrum Time Series data signals at first and second frequencies f0and f1. These two frequencies are preferably set at 2 Hz and 20 Hz respectively although they may be user adjusted within limits.

Modules20and22may be collectively referred to as a bispectrum time series assembly21for generating bispectrum time series at the first and second frequencies.

A pre-processing module24receives the output from the Bispectogram Estimator and is arranged to operate on the 2 Hz time series data signal to apply smoothing, de-trending and histogram equalization.

The pre-processed 2 Hz time series data signal, and the 20 Hz time series data signal, are passed to an MSA estimator module28. The MSA estimator is arranged to processes the time series data signals and includes comparators to compare them to predetermined threshold values as described in steps T1, T2 and T3 of section 3.1 E above. On the basis of those comparisons the MSA Estimator28produces a sleep state signal that indicates any one of Waking, Non-REM and REM sleep states. The sleep state signal is passed to a User Interface32, as will be discussed further below. In another embodiment of the invention further modules may be incorporated to provided a method of diagnosing obstructive sleep based on the determined macro sleep states.

A Sleepiness index Calculator26processes the 20 Hz time series data signal from the Bispectogram Time Series Estimator22. The Sleepiness Index Calculator26is arranged to determine the fraction of the time that the magnitude of the 20 Hz time, in each segment, that the series data signal maintains a value approaching zero (i.e. at or less than (so) as discussed in section 3.4 above).

Modules24to28may be collectively referred to as a macro-sleep architecture assembly25for producing classification signals indicating classification of segments of the EEG signals into macro-sleep states such as Sleep, Wake and NREM sleep and REM sleep.

A User Interface32is provided which communicates with the MSA Estimator module28, the Sleepiness Index Calculator26, the Bandpass Digital Filter8and the Bispectogram Time Series Estimator22. The User Interface32includes a visual display36to display information to a user such as the EEG waveform, sleep state, total sleep time, current f0and f1values and current sleepiness index value. Operator buttons38are disposed about the screen to allow an operator to enter control parameters via an interactive menu system which is displayed upon the screen in use.

FIG. 20shows the Sleep Architecture Estimator in use connected, via electrodes40to a patient42with display36showing the patient's EEG signal and particular information in regard to the patient's sleep state.

The Sleepiness Index may be used for in-facility MSLT and MWT tests, i.e. automated, objective estimation of sleep onset, sleep latencies etc, using just one EEG channel. Furthermore, the described apparatus may be used to automatically score macro sleep states in overnight PSG tests.

4.2 A Driver Sleepiness Monitoring System

Referring now toFIG. 21, there is depicted a block diagram of a Driver Sleepiness Monitoring System44, according to an embodiment of the present invention in use.

The sleepiness monitoring system includes a Sleepiness Index Alert Unit46which incorporates modules2to22and26of the Sleep Architecture Estimator1ofFIG. 19. The Sleepiness Index Alert Unit46receives EEG signals from driver48via a headband mounted set of electrodes50and in responses produces an ongoing sleepiness index. Although a cable connection is shown to the set of electrodes50, a wireless connection may be used instead.

In the event of the sleepiness index rising above a predetermined threshold for longer than a minimum time frame, or more than a minimum number of times in a predetermined time period, then the Sleepiness Index Alert Unit46is arranged to activate alarm50, which requires manual switch-off by the driver within a specified time period, in order to increase the wakefulness of driver48.

In the event of the sleepiness index falling below a predetermined threshold for longer than a minimum time frame, or more than a minimum number of times in a predetermined time period, then the Sleepiness Index Alert Unit46is arranged to activate the driver intervention module (eg. alarm)50in order to increase the wakefulness of driver48.

The Sleepiness Index Alert Unit46also communicates with a Vehicle Control Interface54that in turn is arranged to control the vehicle's engine via the engine management system56.

A data logger52is provided, that is coupled to the Sleepiness Index Alert Unit46and which is arranged to record EEG and sleepiness index data of the driver48. A wireless communications module58, capable of communicating with a home base via a cellular phone network, is coupled to the Sleepiness Index Alert Unit. Consequently the Sleepiness Index Alert Unit46is able to send and receive data relating to the driver's state of wakefulness to the home base. This information may be used to organise a replacement driver for example.

The Driver Sleepiness Monitoring System44may also be equipped with a biofeedback assembly containing an audio stimulus source, e.g. MP3 recordings or radio receiver, to provide stimulus such as music, sports commentaries, etc. to the driver in order to lower the driver's sleepiness index. A continual display of the sleepiness index may also be displayed to the driver including prompts, for example requests for the Driver to stop for a cup of coffee, upon the sleepiness index climbing above a predetermined threshold.

A method and apparatus for estimating Macro Sleep Architecture (MSA) and Sleepiness

Index (SI) as the indicator of sleepiness, using single channel of EEG data, has been described.

The inventors' intensive statistical analysis of the 20,714 epochs from 23 subjects showed a significant agreement between the proposed Higher-Order-Spectra Based technique (HESS) and Technician Scored Sleep States (TSSS). On an average there was an agreement of 77.63% (±9.9%) between TSSS and HESS in classifying sleep into three macro sleep states, State Wake (SW), State REM (SR) and State NREM (SN). This level of agreement is considered an excellent result in clinical sleep scoring, observing that it is on a par with intra-scorer agreement [12] of expert human technologists. Note that the agreement reported on the HESS technique was computed on subjects spanning a large range of RDI, whereas the intra-scorer agreement of humans deteriorates further in the presence of OSAHS cases [13]. In addition, the HESS requires only one channel of EEG, whereas the TSSS uses multiple modes of signals.

The MSA parameters such as NREM/REM stages, Total Sleep Time, Sleep efficiency and REM sleep latency were computed by HESS with significantly high accuracy and correlation, using fully automated algorithms. These are the most important EEG-based parameters resulting from a PSG test. Our ability to compute the parameters reliably, objectively and using a single channel of EEG should make a dramatic impact on the diagnosis and treatment of OSAHS as well a range of other sleep disease such as insomnia.

MSLT and MWT are the main measures of sleepiness used in the clinical environment. They are, however, complicated tests requiring access to sleep laboratories and the services of experienced sleep technologists. Even then, the computation of parameters such as the SL and Sleep Onset are fraught with subjective elements. The Sleepiness Index provides an elegant solution to these problems. It can be automatically computed from a single channel of EEG data at real-time.

SI index is reliable and it is expected that it can be also used as a tool to monitor the sleepiness of humans in hazardous environments such as driving or operating mining equipment.

Appendix A

In Appendix A we provide the standard definitions of Apnea, Hypopnea and Arousals as formulated by professional sleep disorders organizations such as the AASM [9] and ASDA [33]. These definitions are routinely used by sleep physicians around the world to diagnose and treat OSAHS.

Definition of Apnea (AASM):

Cessation of airflow ≥10 s (with oxygen desaturation undefined) ORCessation of airflow <10 s (but at least one respiratory cycle) with oxygen desaturation ≥3%.These events can occur with or without arousal

Definition of Hypopnea (AASM): an event to be defined as hypopnea it should fulfil criterion (a) or (b), plus (c) of the following:(a) A clear decrease (≥50%) in amplitude from base line of a valid measure of breathing during sleep. Baseline is defined as the mean amplitude of stable breathing or the mean of three largest breaths in the two minutes preceding onset of the event.(b) A clear amplitude reduction but does not reach criteria (a) but is associated with oxygen de-saturation >3% or an arousal.(c) Event lasts ≥10 seconds.

Definition of Arousals (ASDA) [33]: American Sleep Disorder Association defines EEG arousals (EEGA) as abrupt shift in EEG frequency, which may include theta, alpha activity and/or frequencies greater than 16 Hz (but not sleep spindles) subjected to following scoring rules:the subject must be asleep for a minimum period of 10 s before declaring an Arousal event,EEG frequency shift must be sustained for a 3 s duration or more, and,EEG arousal from REM sleep requires presence of simultaneous increase in the sub mental EMG amplitude.
Appendix B [10]

Respiratory Disturbance Index (RDI) is defined as the average number of apnea/hypopnea events per hour of sleep.

Arousal Index (ArI) is defined as the average number of arousal events per hour of sleep.

Total time in Bed (TTB) is the time spent on the bed from the ‘lights off’ when the recording starts in the night to the ‘lights on’ in the morning when the recoding ends.

Total Sleep Time (TST) is defined as the actual sleep time (total time spend in REM and NREM sleep).

Sleep Efficiency (SE) is defined as the percentage ratio 100×(TST/TTB) %.

REM sleep Latency (RSL) is defined as the time between the sleep on-set and the first occurrence of REM sleep.

Appendix C [6]

Rectschaffen and Kales (R&K) in 1968 [11] laid the rules for manual scoring of sleep stages. These rules are very broad and complex, in here we are just giving a short summary of these rules. R&K rules require EEG, EMG and EOG channels for implementation.

Macro Sleep State Wake (SW): When a person is awake and active, EEG shows high beta activity (frequency >16 Hz) with low voltage. As the person starts relaxing (with eyes closed) and gets drowsy EEG shows abundance of alpha (8-12 Hz) activity. Voluntary random slow eye movements can be seen in EOG. EMG is tonic with relatively high voluntary activity.

Macro Sleep State NREM (SN): after the sleep onset, EEG shows low voltage and mixed frequency activity. Theta (4-7 Hz) activity increases. The NREM state is further subdivided into four stages; in Stage 1 features such as sharp vertex waves appear; EOG shows slow eye movement; EMG activity decreases. In Stage 2, EEG activity is full of features such as spindles and k-complexes; EOG activity disappears; EMG shows very low tonic activity. In Stage 3 and Stage 4 EEG shows an abundance of delta activity (high voltage (>75 microV), low-frequency (<3 Hz) waves); no EOG activity and very low tonic EMG activity.

Macro Sleep State REM (SR): EEG shows low voltage mixed frequency (range of high theta and slow alpha) sawtooth wave like activity. EOG shows phasic random eye movement activity. In EMG tonic activity gets suppressed and phasic twitches appears.

Kappa statistic is widely used in situations where the agreement between two techniques should be compared. In Appendix D we provide a guideline for interpreting the Kappa values.

Multiple Sleep Latency Test: In MSLT a series of nap opportunities (4-6) are presented to the subject undergoing test at 2 hour intervals beginning approximately 2 hour after morning awakening [6]. Subjects in MSLT test are instructed not to resist themselves from falling sleep. Electrophysiological signals which are recorded to detect sleep onset in MSLT are (i) 4 channel EEGs (C3, C4, O1, O2 with reference to A1 and A2). (ii) left and right eye Electro-oculograms (EOGs), and (iii) submentalis electromyogram (EMG). Each MSLT recording goes for at least 20 minutes. If no sleep onset occurs within 20 minutes than nap opportunity is terminated. If sleep onset occurs before 20 minutes than test session is continued for 15 minutes after sleep onset. Sleep Latency in MSLT is defined as the elapsed time from the start of the test to the first 30-sec epoch scored as sleep. Pathological sleepiness is defined if mean SL is less than or equal to 5 minutes.

Maintenance of Wakefulness Test: The test procedures are similar to MSLT. The only difference is the objective of the test and instructions given to the test subject. The subject is instructed to attempt to remain awake and objective of the test is to assess the function of the underlying wakefulness system.

Appendix F

The method of computing values for the thresholds, ε2, ε20, ε′20and D used in the REM combining rule, used in Section 3.1, (E) and section 3.2 respectively will now be described.

(a) Thresholds (Δ), ε2and ε20: Values for the thresholds (Δ) ε20and ε2were decided after a careful search process, such that the Agreement between HESS classification (WAKE/NREM/REM) and TSSS was optimized. Initially ε2and ε20were set at high values. Both the thresholds were then decreased consecutively in several steps and Δ chosen as the combination of ε2and ε20, please see (9)-(11). For each combination (300) of thresholds set the agreement between TSSS and HESS was computed in estimating, SW, SN, SR. Those values for ε2and ε20for which agreement was optimised were selected.FIGS. 13A and 13Bgraphically illustrate the results of this search. From theFIG. 7it was possible to obtain the absolute values for ε20and ε2at which agreement between HESS and TSSS is optimised. However, it was observed that by using varying threshold (threshold varies for every patient depending on the BTS) such that ε20=median(ξ20)+c and ε2=min(ξ′2), the agreement between TSSS and HESS further increases,FIG. 13A, 13B, where ‘c’ is a constant, whose value was set at 0.03 after a search process. These values of Δ were used in step (T1) and (T2) of section 3.1. At these values the best overall agreement between TSSS and HESS of 75.17% was achieved. The agreement was 74 and 73 when c=0.02 and 0.04 respectively.
ε2ϵ{1,0.67,0.44,0.3,0.20,0.13,0.09,0.058,0.04,0.026,0.017,0.01,0.007,0.005, 0.003,0.0022,0.0015,0.001,0.00067,0.0004}  (9)
ε20ϵ{10,5,2.5,1.25,0.625,0.312,0.156,0.078,0.039,0.019,0.009,0.0048,0.0024, 0.0012,0.0006}  (10)
Δ={ε2,ε20}  (11)

(b) D for REM-combining Rule and Threshold ε′20: To get the best estimation of D, we varied D from 0 to 10 minutes and applied HESS algorithm on the subjects from the database A. Again D was fixed after a careful search process, at which agreement between TSSS and HESS was optimised.FIG. 8(AT the end of this document) shows the variation of agreement with D. From theFIG. 8, D=8 min was set. Epochs excluded from the REM sleep state due to REM rule were re-classified into NREM sleep or WAKE state. We used BTS at f=20 Hz for the re-classification as it has shown the excellent ability to distinguish sleep from wake (seeFIG. 5). We set the threshold constant at ε′20=median(ξ20). We checked for the condition ξ20>ε′20. If condition was satisfied then epoch was reclassified as SWelse as SN.

TABLE 2List of Symbols used hereinSWMacro Sleep State “wake”SNMacro Sleep State “NREM”SRMacro Sleep State “REM”ωDiscrete Frequency (radians per second)rPearson's Correlation CoefficientσSignificance level in statistical testszZ statistics in hypothesis testtt statistic in hypothesis testδDelta frequency range of EEGθTheta frequency range of EEGαAlpha frequency range of EEGβBeta frequency range of EEGNTotal Number of EEG EpochsξfBispectogram Time Series at frequency f

TABLE 3Demographic details of the subjects studied.Sub. IDAgeSexBMIRDI(3.1) Database A: Subjects with PSG recordings145F32.80.6234F23.51.5362F28.683.5463F20.64.1536F34.14.7652M24.54.9744M38.85.6844F49.26.8940F27.68.51044F38.113.81135M28.716.81250F20.118.31350M31.820.81469F23.823.81561M28.133.51671F29.733.61762M30.137.51856F34.440.11963F41.0545.82029M36.848.22147M31.860.12233F39.783.12327M45.594.4(3.2) Database B: Subjects with PSG + Full EEG recording2427M23.77.12532M30.218.5(3.3) Database C: Subject with MSLT recordingSub. IDAgeSexBMISL2634M26.65.52757F19.615.6PSG—Polysomnography, EEG—Electroencephalogram, MSLT—Multiple Sleep Latency Test.

TABLE 4Variation in Accuracy and PPV value of HESS with change in Threshold, ε2andThreshold, ε20used in the section 3.1 and block diagram of FIG. 7.MinMin + (μ − min)/2μμ + (max − μ)/2max(SN+(SN+(SN+(SN+(SN+ε2SWSNSRSR)SWSNSRSR)SWSNSRSR)SWSNSRSR)SWSNSRSR)ε20Accuracy961008961008627474873674888728749287MinPPV1888—881888—87508852925986437158864069Accuracy9780697807646777694067907232659471{Min +PPV1788—881788—88488845745387396851883666(μ − min)/2}Accuracy9770697706656178644361906634589464μPPV1788—881788—88458841704887366546893361Accuracy98403984—3723361385530793948278638{μ +PPV1788—881788—88288431572886274827892645(max − μ)/2}Accuracy100———100———100———100———100———MaxPPV17———17———17———17———17———In our classification algorithm we used two Bispectrum time series, ζ′2and ζ20. ε2was computed from ζ′2and ε20from ζ20. To find out the best value for ε2and ε20we varied both and computed accuracy and PPV of HESS. As we can see from the table at ε20= mean(ζ20) and ε2= minimum(ζ′2), HESS gives the best results. Min is the Threshold set at the minimum value in the time series, Max is the threshold set at maximum value of time series and mean is the mean value of the time series.

TABLE 5Contingency table for TSSS vs. HESS sleep states scoring.TSSSAwakeSleepSWSNSRTotalHESSAwakeSW222317614014385SleepSN8681095946712294SR289134123403970Total338014061320820714

TABLE 7Comparison between percentage of SW, SN, SRsleep computedusing TSSS and HESS, across the 23 subjects.SubSW-SW-SN-SN-SR-SR-Agree-IDTSSSHESSTSSSHESSTSSSHESSmentkappa1171570711324880.64271477761610810.613372559673726740.5948772712022880.655141465602126810.666142167582021830.69730246162914760.56813965532237760.449192165661613640.4910121667602124830.6711101876651417800.5712311958661115820.661371271782220830.651427305963147820.6815725857085720.5516151772651318760.5117191864381744610.4818202266621416650.4719242666611013860.7120161966631818750.652113207571129700.5222192168631316890.552391574581727710.43Last two column shows the % agreement and kappa statistic between TSSS and HESS.

TABLE 8Comparison of estimated sleep states computed usingsingle channel of EEG at different location overthe skull with technician scored sleep states.SensitivityPPVTSTSESleepWakeSleepWakeAgreementKappaTSSS75573.6C3739729386958187.50.78C47507392879585.589.50.81Fp179277.286.649.582.757.375.20.56Fp279277.287.451.383.359.476.40.58F374172.289.77691.472.684.50.73F477875.8937290.378.7840.71Fz7707594.278.292.382.887.20.78Cz78276.286.953.583.959.469.40.44T387585.39541.78274.864.60.41T483881.795.556.88682820.68P375573.692.779.792.779.786.80.77P477175.294.478.692.583.584.60.74T581679.590.9528467.155.10.31T686083.895.24883.678.3570.35O178976.990.761.686.870.567.20.47O285483.294.648.383.676.254.60.27TST—total sleep time, SE—sleep efficiency, PPV—positive predicted value.

REFERENCES

1. Young T, et al.,The Occurrence of Sleep-Disordered Breathing among Middle-Aged Adults. N Engl J Med, 1993. 328(17): p. 1230-1235.2. Young T, et al.,Estimation of clinically diagnosed proportion of sleep apnea syndrome in middle aged men and women. Sleep, 1997. 20(9): p. 705-706.3. Puvanendran, K. and K. L. Goh,From snoring to sleep apnea in a Singapore population. Sleep Res Online, 1999. 2(1): p. 11-14.4. “Wake up Australia”, Access Economics, for the Steering Committee of Sleep Health Australia. Available from the Australian Sleep Association Website, 2005. www.sleepaus.on.net/wakeup.pdf.5. Abeyratne, U. R., A. S. Wakwella, and C. Hukins,Pitch jump probability measures for the analysis of snoring sounds in apnea. Physiological Measurement, 2005. 26(5): p. 779-798.6. Kryger, M. H., T. Roth, and W. C. Dement,Principles and Practice of Sleep Medicine.3 ed. 2000: W. B. Saunders. 1336.7. Australian Sleep Association, Submission to the Parliament of the Commonwealth of Australia. 1999.8. Ronald, J., et al.,Obstructive Sleep Apnea Patients Use More Health Care Resources Ten Years Prior to Diagnosis. Sleep Research Online, 1998. 1(1): p. 71-74.9. Kushida, C., et al.,Practice parameters for the indications for polysomnography and related procedures. Polysomnography Task Force, American Sleep Disorders Association Standards of Practice Committee. SLEEP, 1997. 6(20): p. 406-422.10. Niedermeyer, E. and F. L. d. Silva,Electroencephalography: Basic principles, clinical applications, and related fields.5 ed. 2005: Lippincott Williams & Wilkins. 1309.11. Rechtschaffen, A. and A. Kales,A manual of standardized terminology and scoring system for sleep stages of human subjects. Brain Information Service/Brain Research Institute University of California at Los Angles, 1968 (204).12. Collop N A,Scoring variability between polysomnography technologists in different sleep laboratories. Sleep medicine, 2002. 3(1): p. 43-47.13. Danker H, et al.,Interrater reliability between scorers from eight European sleep laboratories in subjects with different sleep disorders. Journal of Sleep Research, 2004. 13: p. 63-69.14. Suzuki M, et al.,Discrepancy in polysomnography scoring for a patient with obstructive sleep apnea hypopnea syndrome. he Tohoku journal of experimental medicine, 2005. 206(4): p. 353-360.15. Whitney C W, Gottlieb D J, and E. al,Reliability of scoring respiratory disturbance indices and sleep staging. Sleep, 1998. 21(7): p. 749-757.16. Robert G N, et al.,Interobserver Agreement Among Sleep Scorers From Different Centers in a Large Dataset. Sleep, 2000. 23(7).17. B. Kemp, et al.,A model-based monitor of human sleep stages. Biological Cybernetics, 1987. 57(6).18. B. Kim, K. P.,Automatic sleep stage scoring system using genetic algorithms and neural network. Engineering in Medicine and Biology Society, Proceedings of the 22nd Annual International Conference of the IEEE, 2000. 2: p. 849-850.19. Djordje P, et al.,Accuracy of Automated Sleep Staging Using Signals from a Single Forehead Site.20. M. Qianli, et al.,Sleep-stage Characterization by Nonlinear EEG Analysis using Wavelet-based Multifractal Formalism. Engineering in Medicine and Biology Society, Proceedings of the 27th Annual International Conference of the IEEE, 2005: p. 4526-4529.21. O. Pacheco and F. Vaz,Integrated system for analysis and automatic classification of sleep EEG. Engineering in Medicine and Biology Society, Proceedings of the 20th Annual International Conference of the IEEE, 1998. 4: p. 2062-2065.22. Van Hese P, Philips W, and et al,Automatic Detection of Sleep Stages using EEG.23rd Annual Conference—IEEE/EMBS Oct. 25-28, 2001, Istanbul, TURKEY, 2001. 2: p. 1944-1947.23. Chervin R D, Burns J W, and R. D L,Electroencephalographic Changes during Respiratory Cycles Predict Sleepiness in Sleep Apnea. Am. J. Respir. Crit. Care Med., 2005. 171(6): p. 652-658.24. Dingli, K., et al.,Electroencephalographic spectral analysis: detection of cortical activity changes in sleep apnoea patients. European Respiratory Journal, 2002. 20(5): p. 1246-1253.25. V. Svetnik, J. Ma, and e. al.,Evaluation of automated and semi-automated scoring of polysomnographic recordings from a clinical trial using zolpidem in the treatment of insomnia. Sleep, 2007. 30(11).26. American Sleep Disorders Association and D.C.S. Committee,International Classification of Sleep Disorders: Diagnostic and Coding Manual. American Academy of Sleep Medicine, 2005.27. Klink, M. and S. F. Quan,Prevalence of reported sleep disturbances in a general adult population and their relationship to obstructive airways diseases. Chest, 1987. 91(4): p. 540-546.28. CARSKADON M A, et al.,Guidelines for the multiple sleep latency test(MSLT):a standard measure of sleepiness. Sleep, 1986. 9(4): p. 519-524.29. Mendel, J. M.,Tutorial on higher-order statistics(spectra)in signal processing and system theory: theoretical results and some applications. Proceedings of IEEE, 1991. 79(3): p. 278-305.30. Abeyratne, U. R.,Blind reconstruction of non-minimum-phase systems from1-D oblique slices of bispectrum. Vision, Image and Signal Processing, IEE Proceedings, 1999. 146(5): p. 253-264.31. Loader, C.,Smoothing: Local Regression Techniques. Handbook of Computational Statistics, 2004.32. Loader, W. S. C. a. C.,Smoothing by Local Regression: Principle and Methods.33. Michael Bonnet, D. C., Mary Carskadon (Consultant), Paul Easton (Consultant), Christian Guilleminault (Chairman), Ronald Harper, Boyd Hayes, Max Hirshkowitz, Periklis Ktonas, Sharon Keenan, Mark Pressman (Consultant), Timothy Roehrs, Jack Smith, Jim Walsh, Steven Weber, Philip Westbrook,ASDA Report EEG Arousals: Scoring Rules And Examples: A Preliminary Report From The Sleep Disorders Atlas Task Force Of The American Sleep Disorders Association. Sleep, 1992. 2(15): p. 173-184.34. Feinstein A R and Cicchetti D V,High agreement but low kappa: I. The problems of two paradoxes. Journal of Clinical Epidemiology, 1990. 43(6): p. 543-549.35. Feinstein A R and Cicchetti D V,High agreement but low kappa: II. The problems of two paradoxes. Journal of Clinical Epidemiology, 1990. 43(6): p. 551-558.

In compliance with the statute, the invention has been described in language more or less specific to structural or methodical features. The term “comprises” and its variations, such as “comprising” and “comprised of” is used throughout in an inclusive sense and not to the exclusion of any additional features. It is to be understood that the invention is not limited to specific features shown or described since the means herein described comprises preferred forms of putting the invention into effect. The invention is, therefore, claimed in any of its forms or modifications within the proper scope of the appended claims appropriately interpreted by those skilled in the art.