WO2012076923A1

WO2012076923A1 - Method and system for automatic scoring of sleep stages

Info

Publication number: WO2012076923A1
Application number: PCT/IB2010/003410
Authority: WO
Inventors: Christopher Letellier; Claudia Lainscsek
Original assignee: Universite De Rouen
Priority date: 2010-12-08
Filing date: 2010-12-08
Publication date: 2012-06-14

Abstract

A method to produce a scoring for sleep stages comprises: • acquiring an electroencephalogram (EEG) time series; • defining a sliding time windows having a predetermined time duration of the time series; • applying a Delay Differential Equations model onto the time series, said model having at least one parameter varying according to sleep stages; • scoring each sliding time windows based on the value of said at least one parameter.

Description

METHOD AND SYSTEM FOR AUTOMATIC SCORING OF SLEEP STAGES.

Technical Field

The invention relates to the field of automatic scoring of sleep stages. Background Art

Traditionally, the polysomnographic recordings are categorized into stages according to Rechtschaffen and Kales (R&K) (A. Rechtschaffen and A. Kales, editors. A manual of standardized terminology, techniques and scoring system for sleep stages of human subject. Washington DC: US Government Printing Office, National Institute of Health Publication, 1968.). Their rules were introduced for improving the reliability of sleep scorings from different laboratories. The R&K manual involves parameters, techniques and wave patterns of polysomnographic recordings mainly based on a spectral analysis. The scoring, usually accomplished by well-trained neurologist, consists in classifying all 20-30 s pieces of an approximately 8 h recording into one of the next six stages of vigilance:

But interrater variability of sleep stages scorings is known to have an impact on the evaluation of patients' sleep. The old scoring R&K rules were recently modified to overcome this weakness. For instance, for avoiding large interrater departures, stages 3 and 4 merged into a single stage, named SWS (Slow Wave Sleep). These rules are still mostly based on a frequency analysis (spectral analysis). In spite of that, it is yet concluded that no visual based scoring system will ever be perfect, as all methods are limited by the physiology of the human eye and visual cortex, individual differences in scoring experience, and the ability to detect events viewed using a 30-second epoch. And recent studies only showed slight improvements with the new rules and the interrater agreement in sleep stages scorings is slighly greater than 72%.

It is therefore needed to develop automatic sleep scoring techniques. Most of the computer-assisted identification of sleep stages were based on the R&K manual, that is, mainly based on a spectral analysis. In fact, automatic sleep scoring techniques try to reproduce what is done by neurologists with the help of R&K rules: in other words, these latter rules were translated into algorithmic procedures. They are using a multivariate time series made, for instance, of EEG, EOG and EMG signals. But if some of these techniques can lead to an overall epoch-by-epoch agreement of 80%, they require a quite complex decisional tree.

Since the emergence of "chaos theory", it is known that nonlinear dynamical systems may often produce behaviors characterized by a broadband spectra, leaving useless spectral analysis. Since there is no doubt that the dynamics underlying EEG is non linear, dynamical invariants as correlation dimension, Lyapunov exponents, Recurrence Plots Quantifiers were thus also used but the overall agreement was rarely greater than 60 or 70 %. Clearly none of these dynamical invariants was enough for a reliable identification of all sleep stages and combination of few of them applied to different signals was required to reach an overall agreement around 80 %. The size of the data set to safely compute dynamical invariants is often too large to be of a practical use in sleep scoring.

Some techniques based on neural networks, that is, on modelling, were also developed. But the neural networks were used to distinguished different features exhibited in the spectral domain. One of these techniques were not able to distinguish more than the REM sleep from non-REM sleep. Another technique was correctly scoring sleep stages but required two EEG channels, one horizontal electro-oculogram channel and one chin electromyogram channel. An adaptative database method— still based on spectral analysis— was also developed but needed strong preprocessing and the validation with manual scoring was not provided. Another portable one- channel EEG system — BioSomnia — based on criteria defined by Rechstschaffen and Kales revealed a lack of REM sleep detection. A last automatic sleep classification was recently proposed: based on a linear discriminant analysis, this technique was able to distinguish awake, slow-wave sleep and rapid eye movements sleep stages. But a specific sensor, a head accelerometer, is required and must be added to conventional sensors.

Therefore, there is a need for an automatic method to classify sleep stages with reliability at least equivalent to a practitioner and/or without any particular device other than classical EEG system.

Disclosure of Invention

To better address one or more concerns, in a first aspect of the invention a method to produce a scoring for sleep stages comprises:

· acquiring an electroencephalogram (EEG) time series;

• defining a sliding time windows having a predetermined time duration of the time series;

• applying a Delay Differential Equations model onto the time series, said model having at least one parameter varying according to sleep stages;

• scoring each sliding time windows based on the value of said at least one parameter.

The technology here proposed only uses a single scalar time series— one EEG channel — and is able to return a hypnogram that presents an epoch-by-epoch agreement around 80%. It does not use any prior knowledge on sleep classification like the R&K manual. The technique presented is based on dynamical classification using Delay Differential Equations (DDE). The use of DDE is justified by three arguments: i) a spectral analysis seems not to be appropriate since brain dynamics is very likely non linear by nature, ii) modeling techniques are known to be efficient from short data sets, and iii) DDE techniques are known to be very robust against noise and work for relatively low sampling rates.

In particular embodiments:

• the Delay Differential Equations model is selected by a genetic algorithm;

• parameters of the Delay Differential Equations model are estimated using Singular Value Decomposition; • the scoring parameter is selected from the model parameters as having the largest cross-correlation with a corresponding hypnogram;

• the scoring parameter is selected from the model parameters as having the lowest means of variance computed with a sliding variance technique;

• the time series is normalized in the range [0, 1 ] before the definition of the sliding time windows; and/or

• the same model is applied to different time series coming from a same type of sleep electroencephalogram data.

The method proposed here is an automatic classification scheme. The technique is built only upon the data by minimizing a model error without any prior knowledge on sleep stages as defined in the sometimes mathematically inconsistent R&K manual. It only uses dynamical information on the dynamics underlying the data retrieved from a global model. Global modelling techniques are known to be reliable even with short time series. As a consequence, the method disclosed here may be used with window sizes that can be as short as 2 s instead of 30 s windows as commonly used for visual scoring.

Advantageously, it is a one-step procedure allowing very fast scoring

(the scoring of a 8 hour recording is completed in around 1/2 minute with an average laptop).

It is noise insensitive and works for relatively low sampling rates. No smoothing procedure is required. It may thus be applied to the time series recorded at a sampling rate of 64 Hz as commonly used at hospital. This opens advantageously the possibility for real time scoring, which could be very crucial for some clinical settings.

The disclosed technique only needs a single scalar time series, one EEG channel. EMG and EOG measurements are not required, thus reducing the equipment and the number of electrodes to apply to patients.

In a second aspect of the invention a system to produce a scoring for sleep stages comprises: • an interface to an electroencephalogram system for acquiring an electroencephalogram (EEG) time series;

• a calculator for:

• scoring each sliding time windows based on the value of said at least one parameter;

• a human interface to display a hypnogram based on said scoring.

In a third aspect of the invention a computer program product is able to implement any of the method steps described here above when loaded and run on computer means of an EEG system or connected to an EEG system Depending on the type of signal, a particular embodiment may be preferred as easier to adapt or as giving a better result. Aspects of these particular embodiments may be combined or modified as appropriate or desired, however.

Brief Description of Drawings

These and other aspects of the invention will be apparent from and elucidated with reference to the embodiment described hereafter where:

- Figure 1 is a flowchart of a first embodiment of the invention;

- Figure 2A is an hypnogram of a subject scored visually by a neurologist;

- Figure 2B is an hypnogram of the same subject as Figure 2A scored automatically with a 25s window;

- Figure 2C is an hypnogram of the same subject as Figure 2A scored automatically with a 2s window;

- Figures 3A and 3B are hypnograms for the same subject as Figure 2A with, in black line, the visually scored hypnogram as reference and in grey line the automatic scored hypnograms with or without low-pass filtered parameter estimation; - Figure 4 is a cross-correlation between automatically scored and visually scored hypnograms for 35 different subjects;

Mode(s) for Carrying Out the Invention

In reference to Figure 1 , the parts of the processing scheme are:

· 1 . Data Acquisition, step 1

The data can be recorded with any EEG data acquisition system even with low sampling rates. For instance, the sampling rate for the data tested was 64 Hz.

• 2. Data Pre-Processing, step 3.

In order to optimize the delay differential model, the whole data set is shifted into the positive domain and normalized to be within the range [0; 1]. Normalizing the data allows to avoid singularities during parameter estimation of the model. Normalizing the data is not a necessary step but improves the adequacy between dynamics changes and parameter evolution.

· 3. Derivatives, step 5.

The derivatives are computed numerically with an algorithm that should not increase the noise level at any derivation, thus destroying the dynamics underlying the data. Good candidates that are known to be reliable are:

· centre difference schemes such as

dx _ 1 x(t - n) - x(t + n) . . .

dt ^~ N j-{ n

with N = 4, and/or

where N is the number of forward and backwards points for computing the derivative.

N should not be too large for the usual centered difference scheme (1 ); N = 4 is an example giving good results. In the second case, N should be greater; N = 9 is an example giving good results.

• fitting a polynomial on a sliding window, then analytically deriving the polynomial.

• 4. Structure selection of the model using a genetic algorithm, step 7 Delay Differential Equations (DDE) is used in a polynomial form.

The terms depend on the degree of non-linearity D and the number of delays Ντ. For instance, depending on these modeling parameters, the DDE would have the form

and so on.

It is known that optimal models with as few terms as possible better identify distinguishing features of the underlying dynamics. Therefore, it is better to not use full term models.

The time t and the different delays η are multiples of the sampling time 5t. The selected delays are not necessarily multiples of a single delay η as used in the common reconstruction procedure. Indeed, it is known that nonuniform delay embeddings can provide better results to represent the dynamical evolution of a system with different timescales. Such a non-uniform reconstructed space helps to improve the quality of the obtained model. These delays η are thus selected independently by a genetic algorithm.

The genetic algorithm (GA) automatically selects all modeling parameters from the data. These modeling parameters seem to be more sensitive to the measurements performed— the chosen EEG channel— than the dynamics (the patients). When different patients are evaluated using the same EEG channel, the same model form may advantageously thus be applied, avoiding a new structure selection procedure.

These parameters are:

• (a) the parameter space dimension Dp, that is, the number of terms in the model is selected by the GA for a given number of delays, degree of nonlinearity and maximal number of model terms. Models with a smaller number of terms are preferred by the GA in order to avoid a reduction of the model quality as commonly observed with models made of too many terms.

• (b) the degree of non-linearity D: it should be large enough to capture the relevant couplings required to characterize the underlying dynamics.

• (c) the number of delays Ντ : it should be large enough to capture all the relevant timescales.

• (d) the values of the delays η.

· 5. Parameter identification, step 9.

Parameters of the DDE model selected by the GA in the previous step are estimated using Singular Value Decomposition (SVD) or any parameter estimation technique. This is done with a sliding window technique.

Parameter with the largest cross-correlation between one model parameter evolution and the hypnogram is selected as the scoring parameter. Another way to select the scoring parameter would be to compute the variances of each output model parameter with a sliding window technique (variancogram).

The length of these windows would have to be about half of the length of typically expected length of the sleep stages. The model parameter that has the lowest mean of these variances would then be a good candidate for the scoring parameter.

Note that structure selection as well as the selection of the scoring parameter has advantageously to be done only once as long as the data are typical sleep EEG data.

These steps may probably become necessary again for EEG data in very special settings such as e.g. babies, special diseases, different sampling rates, different electrode placements, different sleep parameter identification among others.

· 6. Filtering Parameter Evolution, step 1 1 .

The next step for getting an hypnogram, is to filter outputs of the selected scoring parameter with a lowpass Butterworth filter.

• 7. Scoring into Different Sleep Stages, step 13. To get integer value as used for plotting a hypnogram, the outputs of the previous step are rounded to integer values within the range 1 -6 as shown in right column of the following table.

Details of a genetic algorithm used at step 7 for selecting a model are now given.

A genetic algorithm is a search algorithm that is based on natural genetics.

A given problem is encoded as an array (population) of artificial strings (chromosomes).

In the cases considered here, where an optimization problem has to be solved, the guesses for possible solutions are encoded. The GA is split in two parts, the first one is devoted to the estimation of the time delays and the second part is used for selecting the structure of the DDE models. In the delay-selection part, possible delay-combinations are encoded into binary strings, while in the model-selection part different guesses for models are encoded.

These chromosomes can be strings of 1 's and 0's. The GA will then manipulate this representation of the solution, but not the solution itself.

A GA also must have a criterion for discriminating good from bad solutions according to the relative fitness of these solutions. This will be used to guide the evolution for future generations. In the case considered here, discrimination criterion uses the minimum of the mean square residual without using the visual scoring. This way a model is obtained that reflects distinguishing properties of the data and verifies the manually obtained scoring. Other cost functions such as best separation of the data classes may also be a good choice.

After having encoded the problem in a chromosomal manner and having found a discrimination strategy for good solutions, an initial population of encoded solutions has to be created. This can be done by using a random number generator without any prior knowledge of possibly good solutions. The initial population is in fact a set of random guesses for the optimal solution. For the model-selection part it is a set of different DDE models and for the delay-selection part it is a set of possible delay-combinations.

The evolution of this initial population towards later generations will be done by applying genetic operators in an iterative process. The most common genetic operators are (1 ) selection, (2) recombination, (3) mutation, and (4) junk.

Selection allocates greater survival to better individuals. Better solutions are preferred to worse ones. Additional new, possibly better individuals not present in the original population have to be created. This is done via recombination and mutation. Recombination combines bits of parental solutions to form a better offspring. It combines parental traits in a novel manner. Mutation on the other hand modifies a single individual. It is a random walk in the neighborhood of a particular solution.

In the GA used here, the starting population size may be 200 individuals for the delay selection part and between 10 and 50 for the model structure selection part depending on the order of nonlinearity and the number of delays wanted. The population size changes during the run of the GA: If the residual is constant for around five generations, the population size is increased in order to escape possible local minima of the residual. To further avoid getting trapped in local minima the best five individuals are taken out for five iterations. During these five iterations the mutation rate is increased by 5% and more junk is added. After those five generations the five best that had been removed are added again and also crossed over with the existing generation. This procedure is done twice for each run.

The GA is stopped when the residual does not change for 7 generations. The GA is started with models up to five terms and two delays. More parsimonious models even with a slightly larger residual are preferred to models with more terms. If this does not give satisfying results, the number of terms in the model and/or the number of delays have to be increased.

As an example, the described procedure is now applied to an EEG time series from the data base of the sleep laboratory of the Bois-Guillaume Hospital (Rouen University Hospital, France).

It was a patient sleeping under non invasive mechanical ventilation. Since these patients are not necessarily sleeping well, this test case cannot be considered as a simplest test case. The hypnogram is shown in Fig. 2a.

The GA returned a model with two terms, two delays and D = 2. The selected model

JC— CI γ JC^ 2 & 2 i ( ^ )

has therefore two parameters, ai and a₂. The delays are τι = 18t = 16 ms and τ₂ = 35t = 47 ms.

In the case of model (3) parameters, ai and a₂, show about the same cross-correlation with the hypnogram for subject 15. For instance, parameter ai can be selected as the scoring parameter. While the visual scoring was performed using 30 s windows, the automatically scored hypnogram shows a good agreement for 2 s as well as for 25 s windows (Fig. 2b and 2c). The mean overall agreement for all of the 35 subjects of our data base between the visually and the automatically scored hypnograms is 84% using 25 s windows and 80% using 2 s windows. This leads to 95% of overlap between the two different scorings obtained with these two window sizes (Fig. 4).

Raw parameter estimation is shown in Fig. 3a. Once the Butterworth filter is applied, the parameter evolution fits rather well the visually scored hypnogram (Fig. 3b). In order to investigate the interplay among sleep stages and other quantifiers, it may be useful to use the filtered signal shown in Fig. 3b. But for a comparison with the hypnogram provided by visual scoring it is better to use a discrete signal, that is, the filtered parameter evolution encoded with integers as shown in Fig. 3b.

The method may be implemented by a computer program product that is able to implement any of the method steps as described above when loaded and run on computer means of an EEG system or connected to an EEG system. The computer program may be stored/distributed on a suitable medium supplied together with or as a part of other hardware, but may also be distributed in other forms, such as via the Internet or other wired or wireless telecommunication systems.

For instance, Fig. 5, a system 50 to produce a scoring for sleep stages comprises:

• an interface 52 to an electroencephalogram system 54 for acquiring an electroencephalogram (EEG) time series;

· a calculator 56 for:

• a human interface 58 to display a hypnogram based on the scoring. The calculator 54 may be a general purpose computer running a computer program to run the disclosed steps or include an integrated circuit that may be arranged to perform any of the method steps in accordance with the disclosed embodiments.

Industrial Applicability

The disclosed method and system are useful to develop automatic supervision of sleep periods.

While the invention has been illustrated and described in details in the drawings and foregoing description, such illustration and description are to be considered illustrative or exemplary and not restrictive; the invention is not limited to the disclosed embodiment.

Other variations to the disclosed embodiments can be understood and effected by those skilled on the art in practicing the claimed invention, from a study of the drawings, the disclosure and the appended claims. In the claims, the word "comprising" does not exclude other elements and the indefinite article "a" or "an" does not exclude a plurality.

Claims

Method to produce a scoring for sleep stages comprising:

• acquiring (1 ) an electroencephalogram (EEG) time series;

• defining (3) a sliding time windows having a predetermined time duration of the time series;

• applying (7) a Delay Differential Equations model onto the time series, said model having at least one parameter varying according to sleep stages;

• scoring (13) each sliding time windows based on the value of said at least one parameter.

Method according to claim 1 , wherein the Delay Differential Equations model is selected by a genetic algorithm.

Method according to claim 2, wherein parameters of the Delay Differential Equations model are estimated using Singular Value Decomposition.

Method according to claim 2 or 3, wherein the scoring parameter is selected from the model parameters as having the largest cross-correlation with a corresponding hypnogram.

Method according to claim 2 or 3, wherein the scoring parameter is selected from the model parameters as having the lowest means of variance computed with a sliding variance technique.

Method according to claim 1 , wherein the time series is normalized in the range [0, 1 ] before the definition of the sliding time windows.

Method according to claim 1 , wherein the same model is applied to different time series coming from a same type of sleep electroencephalogram data.

8. System to produce a scoring for sleep stages comprising: • an interface (52) to an electroencephalogram system (54) for acquiring an electroencephalogram (EEG) time series;

• a calculator (56) for:

• a human interface (58) to display a hypnogram based on said scoring.

A computer program product that is able to implement any of the method steps according to claim 1 to 7 when loaded and run on computer means of an EEG system or connected to an EEG system.