METHOD AND DEVICE FOR CALCULATING AT LEAST ONE HEART RATE VARIABILITY INDEX FROM AN RR SERIES

Abstract

The calculation process is implemented by a processing unit starting from an initial RR series or an RR series derived from the initial RR series. The processing unit is configured to detect local minima (Pji) and / or local maxima (Mji) in a sliding main window (F), to divide the main window (F) into n (n≥3) sub-windows (Fj), to calculate a lower curve (Cmin) and / or an upper curve (Cmax), and to calculate at least one heart rate variability index (IV) from the portion of the lower curve (Cmin) in the central window(s) (F2,…,Fn-1), excluding the first sub-window (F1) and the last sub-window (Fn), and / or from the portion of the upper curve (Cmax) in the central window(s), excluding the first sub-window (F1) and the last sub-window (Fn). Figure 5

Description

Title of the invention: METHOD AND DEVICE FOR CALCULATING AT LEAST ONE HEART RATE VARIABILITY INDEX FROM AN RR SERIES Technical field

[0001] The present invention relates to the calculation of at least one characteristic index of heart rate variability from samples of an RR series constructed from a cardiac signal. This variability index can in particular be used to quantify the activity of the Autonomic Nervous System (ANS) Prior art

[0002] From a physiological point of view, the heart of a living being, isolated from any external influence, contracts automatically in a very regular way like a metronome, under the action of the sinus node which generates an independent nerve impulse, and causes a spontaneous contraction of the cardiac muscle. The heart is not, however, isolated, but is connected to the Autonomic Nervous System (ANS).

[0003] The Autonomic Nervous System (ANS) is a part of the nervous system that regulates certain physiological processes such as, for example, blood pressure, heart rate, body temperature or respiratory rate. It is composed of two branches: a "sympathetic" branch, accelerating, and a "para-sympathetic" (or vagal) branch, inhibiting. Its functioning is "autonomous", that is to say, without conscious effort from the person. The ANS is in connection with the organs and sensory receptors of the organism from which it receives information on the internal and external environment of the subject. In response to this information, it stimulates physiological processes, through its sympathetic branch, or inhibits them, through its para-sympathetic branch. This regulation is necessary in order to maintain the internal homeostasis of the organism.

[0004] One of the important functions of the ANS concerns cardioregulation, that is, the regulation of heart rate and blood pressure to meet the oxygen needs of tissues, for example, during physical exertion, stress, pain, emotion or other physiological events. This is a constant, continuous regulation, ordered by the two branches of the ANS, sympathetic (stimulation) and vagal (inhibition).

[0005] This phenomenon is particularly visible on the heart rate, which undergoes these influences, acceleration and deceleration, causing variability of the instantaneous heart rate. However, outside of these influences, the heart rate is under the control of the "sinus node", a nerve node located in the tissues of the right atrium of the myocardium, acting like a biological clock at a very regular rhythm. This variability in heart rate is therefore a direct reflection of the activity of the ANS.

[0006] Furthermore, heart rate is a physiological parameter that is simple to implement and that can be measured non-invasively, for example by means of electrodes for measuring the electrocardiogram (ECG), a plethismographic pulse sensor (PPG) placed on the finger, ear, forehead or wrist depending on the use of the device, or any other sensor such as ultrasound, microphones, piezoelectrics, smartphone cameras, etc.

[0007] In many situations, it is interesting to know the activity of the ANS. In particular in medicine, during surgical interventions, in resuscitation or intensive care situations, or even in outpatient monitoring, but also in non-medical fields such as sport, physical activities in general or relaxation.

[0008] In particular, knowledge of the activity of the ANS can be of valuable assistance in developing a diagnosis of many clinical situations. On this subject, one can refer for example to the following publication: Lacroix D, Logier R., Kacet S., Hazard JR, Dagano J. (1992): “Effects of consecutive administration of central and peripheral anticholinergic agents on respiratory sinus arrhythmia in normal subjects, J. of the Autonomie Nervous System”, Vol 39, pages 211-218.

[0009] To study these fluctuations in heart rate, various techniques for spectral analysis of a signal representing the evolution over time of the instantaneous heart rate (or frequency) have been developed since 1970, for example, and which is obtained after sampling an analog bioelectric signal, characteristic of the heart rate of a living being, and subsequently called a cardiac signal.

[0010] An example of a method for analyzing heart rate variability known to date consists of: - to acquire a cardiac signal, by any invasive or non-invasive means [for example, and without limitation, acquisition of an electrocardiographic (ECG) signal by means of an electrocardiograph, or use of a blood pressure sensor connected to a catheter introduced into an artery, or use of an infrared pulse sensor], - to construct from this signal a series called RR which is made up of a plurality of samples (RR;) representing the time intervals which separate two successive heartbeats, - and to carry out a spectral analysis of the RR series.

[0011] More particularly, the spectral analysis of an RR series from a cardiac signal is usually carried out in two main stages.

[0012] In a first step, we calculate the spectral density curve of the RR series, by example between 0 and 2Hz, using different known methods. The most commonly used method is to calculate the fast discrete Fourier transform of the RR series, in predefined time windows, weighted by means of a predefined weighting window. Depending on the intended implementation, this can be a rectangular weighting window, or for example a Kaiser, Hamming, Hanning or Bartlett weighting window. The computational time windows can also be predefined and fixed, or it can be a computational time window, of predetermined size, which is slid over time. For example, the Fourier transform is performed in a sliding time window of 256 seconds, applied to the RR series, and subjected to Kaiser weighting to limit the edge effects due to windowing.

[0013] In a second step, from the spectral density curve obtained at the end of the first step, the spectral powers (areas under the spectral density curve) are automatically calculated between predetermined frequency limits, and possibly adjustable by a user. These spectral power calculations make it possible to obtain quantitative information, which is characteristic of the activity of the Autonomic Nervous System (ANS), and thus constitute a means of investigating and analyzing cardiac regulation by the ANS.

[0014] The aforementioned spectral analysis method, however, has several drawbacks.

[0015] An RR series is non-stationary. Consequently, applying a Fourier transform to this type of non-stationary series gives results that are imprecise, or even erroneous, and that cannot be interpreted with certainty. The calculation of the spectral density curve by fast Fourier transform (or equivalent) is relatively heavy in terms of computing power and / or computing time, which currently makes this spectral analysis method unsuitable and difficult to implement in real time, particularly in embedded systems.

[0016] To obtain an acceptable frequency resolution, the fast Fourier transform must also be calculated in relatively large time windows (for example 256s), which corresponds to a large number of samples of the RR series. It follows that this method of spectral analysis is accompanied by a memory effect which delays the taking into account of a change occurring in the RR series.

[0017] Furthermore, a method for calculating a heart rate variability index has already been proposed, in particular in international patent application WO2006 / 032739, which can more particularly be used to quantify the activity of the ANS. This method makes it possible to overcome the aforementioned drawbacks inherent in prior art methods using spectral analysis and is based on calculating a variability index from local minima or local maxima of a resampled, filtered and normalized RR series, which are detected in a sliding main window. More part Specifically, in this method, the main window is divided into several sub-windows, a curve is determined by linear interpolation, connecting all the local minima (when the samples of the RR series are respectively equal to the time intervals ôti which separate two successive heartbeats) or connecting all the detected local maxima (when the samples of the RR series are respectively equal to the inverse 1 / ôti of the time intervals which separate two successive heartbeats), and heart rate variability indices are calculated from this curve.

[0018] In this type of method, the sliding of the main window is typically carried out with an overlap, the sliding step of the main window having a duration less than the duration of the main window. This overlap avoids the loss of information and promotes the real-time calculation of the variability index.

[0019] By implementing this type of method, it has nevertheless been observed at certain times that oscillations or untimely variations over more or less long periods of time, in the value of the variability index. These oscillations or untimely and temporary variations do not correspond to physiological reality and constitute artifacts which can harm the reliability of this index. Presentation of the invention

[0020] The main objective of the present invention is to propose a method for calculating at least one heart rate variability index from local minima and / or local maxima of an RR series which are detected in a sliding main window with overlap, which method has been improved so as to reduce, and preferably to eliminate, the aforementioned untimely and temporary oscillations or variations of the variability index which do not correspond to physiological reality. Summary of the invention

[0021] The invention thus has as its first object a method for calculating at least one heart rate variability index, implemented automatically by a processing unit, from an initial RR series consisting of a plurality of samples (RRi) which are a function respectively of the time intervals (ôti) which separate two successive heartbeats or from an RR series derived from the initial RR series, in particular a resampled and / or normalized and / or filtered RR series.

[0022] This method of the invention comprises a step of calculating (a) at least one heart rate variability index in a sliding main window (F) of duration (T), with a sliding step of at least one sample and of a duration (d) less than the duration (T) of the main window; the processing unit is configured to, during said calculation step (a), detect local minima (Pji) and / or local maxima (Mji) in the main window (F), to divide the main window (F) into n (n>3) sub-windows (Fj) consisting of a first sub-window (FJ), one or more central windows (F2,...,Fn_i) and a last sub-window (Fn), to calculate a lower curve (Cmin) connecting at least the last local minimum of the first sub-window (FJ, at least each local minimum (Pj;) included in the central sub-window(s) (F2,...,Fn_i), and at least the first local minimum of the last sub-window (Fn) and / or an upper curve (Cmax) connecting at least the last local maximum of the first sub-window (FJ, at least each local maximum (Mjj) included in the central sub-window(s) (F2,...,Fn_i), and at least the first local maximum of the last sub-window (Fn), and to calculate said at least one index (IV) of heart rate variability from the portion of the curve lower (Cmin) in the central window(s) (F2,...,Fn_i), excluding the first sub-window (FJ and the last sub-window (Fn) and / or from the portion of the upper curve (Cmax) in the central window(s), excluding the first sub-window (Fi) and the last sub-window (Fn). .

[0023] To carry out the invention, it has been demonstrated that the aforementioned untimely and temporary oscillations or variations of the variability index could be linked to an uncertainty in the calculation of the aforementioned curve at the extreme edges of the main window, this uncertainty in calculation being able to be amplified by the overlap in the sliding of the main window. In the invention, the exclusion of the first and last sub-windows makes it possible to overcome this drawback.

[0024] In the present text, the term "cardiac signal" refers to any physical signal characteristic of the instantaneous heart rate of the living being. For the implementation of the invention, different invasive or non-invasive techniques can be used to acquire this cardiac signal. A known invasive technique consists, for example, of using a blood pressure sensor connected to a catheter introduced into an artery. Among the known non-invasive methods (and which will be preferably chosen), we can cite, for example, the use of an infrared pulse sensor, the use of an ultrasound sensor allowing the detection of cardiac cycles, of the type of sensor implemented in a cardiotocograph, or the acquisition of an electrocardiographic signal (ECG).The acquisition of an electrocardiographic signal (ECG) is in practice the most commonly used method, because in addition to its non-invasive nature, it allows obtaining a more precise cardiac signal than that obtained, for example, by means of an infrared pulse sensor.

[0025] In the present text, the term "initial RR series" generally designates a series of several successive samples (RRi), obtained after sampling an analog cardiac signal characteristic of the heart rate of the living being, each sample (RRi) being generally a function of the time interval (ôti) between two successive heartbeats. More particularly Each sample (RRi) can be equal or proportional to the time interval (ôti) between two successive heartbeats or can be equal or proportional to the inverse (1 / ôti) of this time interval.

[0026] In the preferred embodiment described below with reference to the appended figures, this initial RR series is more particularly constructed from the R waves of an ECG signal. This is not, however, limiting of the invention. In the case of an ECG-type cardiac signal, the so-called “RR” series can be constructed by using the other depolarization waves (P, Q, S or T) of the ECG signal to construct the RR series, the precision being, however, less good than by using the R waves of the ECG signal. Also, when the cardiac signal is not an ECG signal, the samples of the initial RR series are not calculated by determining the time interval separating two successive R waves of the ECG signal, but are more generally determined by detecting in the cardiac signal the time interval between two successive heartbeats.

[0027] The variability indices calculated using the method of the invention can advantageously be used to quantify the activity of the ANS or to characterize or detect any stimulus having an effect on the activity of the ANS, and resulting in a variation in the variability of the heart rate.

[0028] More particularly, the method of the invention may include the following additional and optional characteristics, each additional and optional characteristic being able to be combined with any of the other technical characteristics of the method: - The processing unit is configured to determine, in each central sub-window (F2,...,Fn.i), an area (Aj) which is delimited at least by the lower curve (Cmin) and is above the lower curve (Cmin), and which extends over the entire duration (D) of the sub-window and / or an area (Aj) which is delimited at least by the upper curve (Cmax) and is under the upper curve (Cmax), and which extends over the entire duration (D) of the sub-window and to calculate said at least one index (IV) of heart rate variability from said areas (Aj). - The processing unit is configured to determine, in each central sub-window (F2,...,Fn.i), an area (Aj) which is delimited by the lower curve (Cmin) and by the upper curve (Cmax) and which extends over the entire duration (D) of the sub-window, and to calculate said at least one index (IV) of variability of the heart rate from said areas (Aj) - The main window (F) is preferably divided into at least four sub-windows (n>4) and the processing unit is configured to select from among said areas (Aj), the one with the lowest value, and to calculate said at least one index (IV) of heart rate variability from this minimum area (Smin). - The processing unit is configured to calculate the ratio (Smin / S) between said minimum area (Smin) and the total area (S) of the main window (F) and to calculate said at least one index (IV) of heart rate variability from this ratio (Smin / S). - The duration (T) of the main window (F) is preferably between 32s and 256s, and more preferably is equal to 64s. - The number (n) of sub-windows is preferably greater than or equal to 4. - Sub-windows can have identical durations. - The duration of the first sub-window (FJ) is preferably equal to the duration of the last sub-window (Fn). - the calculation step (a) is carried out on an RR series which is obtained at least by resampling the initial RR series or an RR series derived from the initial RR series, at a sampling frequency (Fe) greater than or equal to 4Hz, and preferably equal to 8Hz. - The duration (d) of the sliding step is greater than or equal to the inverse (1 / Fe) of the sampling frequency (Fe), and preferably is equal to 1s. - the calculation step (a) is carried out on an RR series which is obtained at least by normalizing the samples (RRi) of the initial RR series or the samples (RRi) of an RR series derived from the initial RR series, over the entire duration (T) of the main window (F). - the calculation step (a) is carried out on an RR series which is obtained at least by filtering the initial RR series or an RR series derived from the initial RR series in at least one predefined frequency band [fl - f2]. - Filtering allows filtering in at least one high frequency band [0.15Hz-0.4Hz]. - Filtering allows filtering in at least one low frequency band [0.04Hz; 0.15Hz], - The processing unit is configured to display the heart rate variability index(es) (VI) in real time on a screen.

[0029] Another subject of the invention is a device for calculating at least one index (IV) of heart rate variability, comprising at least one processing unit configured to automatically implement the above-mentioned calculation method.

[0030] Preferably, the device further comprises means for acquiring a cardiac signal and means for sampling this analog cardiac signal, and the processing unit is configured to also construct said initial RR series from the signal sampled by the sampling means.

[0031] Preferably, the processing unit is configured to automatically implement the above-mentioned calculation method in real time, as the cardiac signal is acquired and the initial RR series is constructed.

[0032] Another subject of the invention is a computer program product comprising program code instructions and making it possible, when executed by at least one processing unit, to carry out the above-mentioned calculation method. Description of figures

[0033] Other characteristics and advantages of the invention will appear more clearly on reading the detailed description which is given by way of non-limiting and non-exhaustive example, and with reference to the appended drawings in which: - [Fig.l] [Fig.l] schematically represents the main elements of an example of a device for acquiring and processing a cardiac signal, which can be configured to calculate indices of variability of the heart rate in accordance with the method of the invention, - [Fig.2] [Fig.2] represents the set of waves (PQRST) characteristic of a heartbeat in an ECG signal, - [Fig.3] [Fig.3] represents an example of a digital ECG signal, obtained after sampling an analog ECG signal, - [Fig.4] [Fig.4] represents an example of an initial RR series (also designated initial RR signal) constructed from the signal in [Fig.3], and - [Fig.5] [Fig.5] represents an example of calculation of areas on an RR series after normalization and filtering of the signal and after centering and reduction of amplitude in a range of values ​​without unit between -0.1 and +0.1. - [Fig.6] [Fig.6] represents another example of calculating areas on an RR series after normalization and filtering of the signal and after centering and reduction of amplitude in a range of values ​​without units between -0.1 and +0.1. Detailed description

[0034] Device for acquiring and processing a cardiac signal

[0035] [Fig.l] shows an example of a device which allows the acquisition and processing of the cardiac signal of a living being (hereinafter referred to as “patient”) in accordance with the invention.

[0036] This system includes:

[0037] - usual means for acquiring an ECG signal, comprising several electrodes measuring 1 connected as input to an electrocardiographic monitor (ECG) 2,

[0038] - means 3 for processing the ECG signal delivered as output by the ECG monitor 2.

[0039] The ECG signal processing means 3 comprise a converter analog / digital 4, and a processing unit 5. The input of converter 4 is connected to the output of the ECG monitor 2, and the output of the converter 4 is connected to an input port of the processing unit 5.

[0040] The processing unit 5 can be implemented with various technologies and is preferably a programmable electronic processing unit, for example of the microcomputer type, programmable automaton, electronic card with microprocessor(s) or microcontroller(s) or electronic card with programmable electronic circuit(s) of the FPGA type. The processing unit 5 can also be implemented with specific integrated electronic circuit(s) of the ASIC type.

[0041] In the particular embodiment illustrated in [Fig.l], and not limiting the invention, the processing unit 5 is constituted by a microcomputer, the converter 4 being connected to a communication port of this microcomputer (for example RS232 serial port or USB port).

[0042] In other embodiments, the analog / digital converter 4 and the processing unit 5 can be integrated into an on-board system, for example of the “holster” type and capable of being worn by a patient.

[0043] In operation, the electrodes 1 are applied to the body of a patient, and the ECG monitor outputs in the usual manner an analog electrical signal, called the ECG signal, which for each heartbeat, has for example the form of the signal shown in [Fig.2].

[0044] With reference to [Fig.2], for each heartbeat, this electrocardiographic signal (ECG) consists of a set of electrical waves: - P-shaped, which corresponds to the depolarization of the atria, and which has a low amplitude and a dome shape; - the PQ space which reflects the atrioventricular conduction time; - R fund considered in practice as a marker of ventricular systole, or heartbeat, the QRS complex reflecting ventricular contraction, and - T background which reflects ventricular repolarization.

[0045] This analog ECG signal is digitized by the converter 4, with a predetermined sampling frequency (Fe), for example 256 Hz.

[0046] The sampled signal delivered at the output of the converter 4, (signal shown in [Fig. 3] is processed by the processing unit 5, by means of specific calculation software which is described in detail later. This calculation software is stored in the memory of the processing unit 5, and allows, when executed, to automatically calculate, from the digital signal delivered by the analog / digital converter 4, one or more indices (IV) of variability of the heart rate.

[0047] A preferred variant of this calculation software will now be detailed.

[0048] Example of algorithm of the software for calculating the frequency variability index (fs) cardiac

[0049] In a preferred embodiment of the invention, the main successive steps of the algorithm of the software for calculating heart rate variability index(es) are as follows: 1. Acquisition of RRi samples from the initial RR series and resampling 2. Selection of RRk samples from the resampled RR series included in a main time window F of duration T [T (in seconds) > 1 / f] 3. Signal normalization. 4. Filtering 5. Division of the main time window F into n [n >3] sub-windows Fj 6. Detection of local minima Pj; and local maxima Mj; in the main window F 7. Calculation of Cmin and Cmax curves 8. Calculation of areas Aj 9. Calculation of one or more IV indices of heart rate variability 10. Sliding (with overlap) of the main window F and repeating the calculation from step 2

[0050] In practice, the system can be programmed to be used in real time or in delayed time.

[0051] When the system is used in deferred time, step 1 is first carried out in real time so as to acquire all the RRi samples over the entire desired analysis period; all of these successive RRi samples are stored in an acquisition file in the memory of the processing unit. In a second step, steps 2 to 10 are carried out in a loop, in deferred time, on the RRi interval values ​​stored in the acquisition file.

[0052] When the system operates in real time, the acquisition step 1 on the one hand, and the other processing steps 2 to 10 on the other hand are executed by two separate software modules operating in parallel, the first acquisition module (step 1) feeding the second processing and calculation module (steps 2 to 10) by means of a buffer file or equivalent.

[0053] Steps 1 to 10 will now be detailed.

[0054] Step 1: Acquisition of RRi samples from the initial RR series and resampling

[0055] The acquisition of the RRi samples of the initial RR series is carried out by a first software sub-module which is supplied as input with the successive digital data constituting the digitized ECG signal (signal of [Fig.3]) delivered by the analog-to-digital converter 4. Each data (or point) of the ECG signal is defined by the instantaneous amplitude ECG; of the ECG signal, and by the instant h sampling rate (f = n; / Fe, with n; sample number and Fe representing the sampling frequency of converter 4).

[0056] The first RR sample acquisition sub-module is designed to automatically detect each successive peak Ri in the digital signal delivered by the converter 4, and to automatically construct an RR series ([Fig.4]) consisting of a succession of RR samples. Each RRi sample is defined by the pair of coordinates: ti [a sampling instant (or number)]; time interval ôti (expressed for example as a multiple of the sampling period 1 / Fe) separating a peak Ri from the following peak Ri+i (in another variant it could be the preceding peak R, J.

[0057] Usually in itself, Fonde R being most often the finest and widest part of the QRS, it is preferably used to detect the heartbeat with very good precision, the time interval oti corresponding in practice to the time separating two successive heartbeats. Nevertheless, in another variant, one could envisage using other waves (for example Q wave or S wave) of a heartbeat of the ECG signal to detect and construct the RR series.

[0058] The initial RR series ([Fig.4]) delivered by the aforementioned first sub-module is resampled by a second software sub-module at a predefined frequency f. The objective of this resampling is to obtain at output an RR series whose RRk samples are equidistant from a temporal point of view, that is to say in other words an RR series whose sampling instants are equidistant temporally.

[0059] This resampling can be carried out in a manner known per se by interpolation, and for example by linear interpolation, in order to constitute a series of samples which are temporally equidistant.

[0060] The resampling frequency f is preferably greater than twice the patient's physiological maximum heart rate.

[0061] In order to take into account all possible physiological situations, the resampling frequency f is preferably set to a value greater than 4Hz, and is for example 8Hz.

[0062] Step 2: Selection of RRi samples included in a main time window of duration T [T (in seconds) > 1 / f]

[0063] This step amounts to isolating a number N of successive RRi samples (N=Tf).

[0064] The duration T of the main window F is preferably greater

[0065] or equal to between 32s. Below 32s, there is a risk that the variability index of heart rate is not physiologically significant.

[0066] The duration T of the main window F is preferably less than or equal to 256s, in order to reduce calculation times and allow implementation substantially in “real time”.

[0067] As an indication, we choose for example and preferably a main window F having a duration T of 64 seconds, which corresponds to 512 successive RRk samples (N=512), for a resampling frequency f of 8Hz.

[0068] The following steps 3 to 9 are applied to the RRk samples of the resampled RR series included in this main window.

[0069] Step 3: Normalization of the signal

[0070] This step is carried out by means of a software sub-module which firstly calculates the norm S of the signal from step 2 in accordance with the usual formula below: 100711 S =

[0072] Where p is an integer greater than 1, and S; are the discrete values ​​of the signal.

[0073] In another variant, the standard can also be calculated using the following formula:

[0074] Preferably, for the implementation of the invention, p will be chosen equal to 2.

[0075] In a second step, the software sub-module performs a normalization of the signal by dividing each value Si of the signal by the previously calculated norm S.

[0076] This normalization makes it possible to overcome the influence of the average value of the heart rate.

[0077] It is essential that this normalization be carried out over the entire width of the main window, i.e. taking into account in the normalization the (N) RRk samples of the main window. In contrast, a normalization carried out on each sub-window (Fj) would be totally ineffective.

[0078] Preferably, during this step, the software sub-module also performs a centering and a reduction of the amplitude of the RRk samples in a unitless range, for example between -0.1 and +0.1, such as in the example of [Fig.5]. This operation makes it possible to overcome the influence of inter-individual amplitude variations.

[0079] Step 4: Filtering

[0080] This step consists of applying a band-pass filter to the RRk samples of the resampled and normalized RR series included in the main window so as to retain only the frequencies included in a predefined frequency band [fl; f2].

[0081] More particularly, the frequency band [fl; f2] is equal to or is included in the frequency band [0.04Hz; 5Hz].

[0082] To carry out this band-pass filtering step, for example, a high-pass digital filter having a cut-off frequency at the frequency fl will be used, in series with a low-pass digital filter having a cut-off frequency at the frequency f2. It is also possible to use a recursive selective filter with infinite impulse response (IIR filter) centered on the frequency band [fl; f2].

[0083] The high-pass filter (cutoff frequency fl) aims to filter very low frequencies, and incidentally to remove many artifacts in the signal. In practice, the cutoff frequency fl is greater than or equal to 0.04 Hz.

[0084] The low-pass filter (cut-off frequency f2) aims to filter high frequencies, typically higher than 5Hz and preferably higher than 1Hz, because in practice they do not contain any interesting information. The implementation of this low-pass filtering, although preferential, is however optional and not essential to the realization of the invention.

[0085] In a preferred embodiment, the frequency band [fl; f2] is equal to [0.15Hz; 0.4Hz] in order to retain only the high frequency variations of the RR series (respiratory variations) relating to the parasympathetic tone.

[0086] In another preferred embodiment, the frequency band [fl; f2] may be equal to [0.04Hz; 0.15Hz]. This low frequency range is known to be influenced by both the sympathetic system and the parasympathetic system.

[0087] Step 5: Splitting the main time window F into n [n >3] sub-windows F I

[0088] The software sub-module corresponding to step 5 automatically divides the main time window F into n [n >3] sub-windows Fj.

[0089] In an alternative embodiment, the duration of the first sub-window Fi is preferably (but not necessarily) equal to the duration of the last sub-window Fn.

[0090] In a simplified embodiment, all the sub-windows Fj can have identical durations.

[0091] As an indication and not a limitation of the invention, in the particular example of [Fig.5], the main window has been divided into four (n=4) sub-windows Fh F2, F3 and F4 of identical durations D (duration D of each sub-window equal to 16s for a main window F of 64s), namely a first sub-window FH two central sub-windows F2 and F3 and a last sub-window F4.

[0092] In this particular example of [Fig.5], the two windows Fi and F4 are the two extreme windows which are excluded and are not taken into account for the calculation of the heart rate variability indices IV in step 9.

[0093] Step 6: detection of local minima Pji and local maxima Mji in the main window

[0094] Step 6 is carried out by means of a software sub-module which detects local minima Pji and local maxima Mj; of this signal.

[0095] This detection is preferably carried out on the principle of measuring the maximum amplitude between two successive zero crossings of the RR series. Thus, a local minimum is Pji or a local maximum Mj; is detected as being the point having the maximum amplitude relative to zero between two successive zero crossings of the RR series.

[0096] Detection can also be achieved by means of an algorithm detecting a slope inversion (or equivalently a change in the sign of the signal derivative).

[0097] The detection of local minima and local maxima can also be optimized to retain only the minima Pji whose value is lower than a predefined or calculated threshold and to retain only the maxima Mj; whose value is higher than a predefined or calculated threshold, in order not to retain the local minima or local maxima which are not significant and linked to non-physiological rebounds.

[0098] In the example of [Fig.5], the local minima that were detected are as follows: - First sub-window Fi: Pu, P[2, Pu - Central sub-window F2: P2i, P22 - Central sub-window F3: P3[, P32, P33 - Last sub-window F4: P4h P42

[0099] The local maxima that were detected are as follows: - First sub-window Fi: Mu, Mn, Mn - Central sub-window F2: M2b M22, M23 - Central sub-window F3: M3b M32 - Last sub-window F4: M4b M42, M43

[0100] Step 7: Calculation of the Cmin and Cmax curves

[0101] In step 7, the processing unit 5 automatically calculates: - the lower curve Cmin which connects the local minima Pji detected in the previous step and - the upper curve Cmax which connects the local maxima Mj; detected in the previous step

[0102] These curves Cmin and Cmax are determined by interpolation, and preferably by linear interpolation, that is to say by calculating the straight line segments connecting the local minima (Cmin curve) and the straight line segments connecting the local maxima.

[0103] These curves Cmin and Cmax can also be determined by any other type of mathematical interpolation and in particular cosine, polynomial or polynomial by parts interpolations.

[0104] Step 8: Calculation of the areas Aj

[0105] The processing unit 5 automatically calculates for each central sub-window F2 to Fn_i [i.e. windows F2 and F3 in the example of [Fig.5]], excluding the two extreme sub-windows [first sub-window Fi and last sub-window Fn, i.e. excluding windows Fi and F4 in the particular example of [Fig.5] with four sub-windows), Do Aj which extends, over the entire duration D of the central sub-window between the lower curve portion Cmin and the upper curve portion Cmax in the central sub-window

[0106] Step 9: Calculation of one or more IV indices of heart rate variability

[0107] The processing unit 5 automatically calculates at least one variability index of the heart rate IV from the areas Aj calculated in the previous step.

[0108] In a preferred embodiment, the processing unit 5 automatically calculates at least one heart rate variability index IV from the smallest (Smin) of the areas calculated in the previous step, i.e. in the particular example of [Fig.5], from the smallest of the areas A2 and A3 [Smin= min (A2; A3)]. In this particular example of [Fig.5], the minimum area Smin is equal to A2.

[0109] More particularly and in a non-limiting manner of the invention, the index IV is a function of the minimum area Smin, and for example is calculated by means of an affine function with Smin as a variable, and can more particularly be calculated by means of the following formula:

[0110] IV =100.(a.Smin +b) / S [YES] with: - S representing the total surface area of ​​the main window F - a and b are coefficients that can be determined empirically.

[0112] In the particular example of [Fig.5] (main window F of 64 seconds between -0.1 and +0.1), the total surface area S of the main window F is 12.8.

[0113] This IV index advantageously makes it possible to quantify the activity of the parasympathetic tone.

[0114] Preferably, the value of this index IV can be visualized by being displayed on a screen by the processing unit 5

[0115] Step 10: Sliding (with overlap) of the main window F and repeating the calculation from step 2.

[0116] The processing unit 5 shifts the main time window of duration n seconds, by a time step equal to d seconds less than the duration T of the main window F (d < T), and loops back to step 2.

[0117] The step d for sliding the main calculation window influences the sensitivity of the index.

[0118] As an indication, and in practice, very satisfactory results have been obtained in choosing a step d worth 1s. Other variants

[0119] In another embodiment variant, the aforementioned filtering step (Step 4) could be carried out before the selection step 2, for example continuously as the RRi samples of the initial RR series are acquired. Also, in another embodiment variant, this filtering step could be carried out before the normalization step (abovementioned step 3).

[0120] In another variant, an index IV making it possible to quantify the activity of the parasympathetic tone can also be calculated from only the curve Cmin connecting the local minima Pj;, for example by calculating for each central sub-window an area Aj, which is delimited at least by the lower curve Cmin and is above the lower curve Cmin, and which extends over the entire duration D of the central sub-window. This is for example an area Aj which extends between this curve Cmin and a predefined horizontal line, which may more particularly but not necessarily be the abscissa axis, as illustrated for example in [Fig.6].

[0121] In another variant, an index IV can also be calculated from only the curve Cmax connecting the local maxima Mj;, for example by calculating for each central sub-window an area Aj, which is delimited at least by the upper curve Cmax and is under the upper curve Cmax, and which extends over the entire duration D of the central sub-window. This is for example an area Aj which extends between this curve Cmax and a predefined horizontal line, which may more particularly but not necessarily be the abscissa axis.

[0122] In another variant, the Cmin curve can be calculated by taking into account only: - the last local minimum of the first sub-window (i.e. for [Fig.5], the point Pn of the sub-window FJ, the other local minima of the first sub-window not being taken into account, - the first local minimum of the last sub-window (i.e. for [Fig.5], point P4i of sub-window F4), the other local minima of the last sub-window not being taken into account, - the local minima of the other central sub-windows (i.e. sub-windows F2 and F3 for [Fig.5]).

[0123] In another variant, the Cmax curve can be calculated by taking into account only: - the last local maximum of the first sub-window (i.e. for [Fig.5], the point Mn of the sub-window FJ, the other local maxima of the first sub-window not being taken into account, - the first local maximum of the last sub-window (i.e. for the [Fig.5], the point of the sub-window F4), the other local minima of the last sub-window not being taken into account, - the local maxima of the other central sub-windows (i.e. sub-windows F2 and F3 for [Fig.5]).

[0124] In another variant, a heart rate variability index IV can be calculated in the same way as the AUCmoy parameter or the AUCmax parameter described in international patent application WO2006 / 032739, but excluding in the calculation the areas calculated in the first sub-window and in the last sub-window of the main window.

[0125] In the embodiments which have been described with reference to the appended figures, the calculation algorithms are applied to an initial RR series whose samples (RRi) are equal or proportional respectively to the time intervals (ôti) which separate two successive heartbeats. In another variant embodiment, the invention can be implemented from an initial RR series whose samples (RRi) are equal or proportional respectively to the inverse (1 / ôti) of the time intervals between two successive heartbeats. In this case, by applying the aforementioned steps 1 to 7, curves are obtained which are the inverse of those of [Fig.5].

[0126] The above-mentioned IV variability indices can generally be used to measure the variability of the heart rate in a living being (conscious or unconscious), and therefore to study the effects on the ANS (Autonomic Nervous System) of any stimulus likely to modify the activity of the ANS.

[0127] As non-limiting and non-exhaustive examples of the invention, the IV indices can be used in the medical or surgical field for: - the assessment of the pain felt by a conscious living being that is not mechanically ventilated (i.e. a living being whose respiratory rate is arbitrary and variable, and is not imposed by a controlled ventilation device, in contrast in particular with a patient under general anesthesia), for example a conscious living being under local anesthesia of the epidural type, - the assessment of the pain felt by an unconscious living being, and in particular a living being under general anesthesia; in the latter case, the assessment of the level of pain makes it possible to indirectly know the level of analgesia during general anesthesia.

[0128] As non-limiting and non-exhaustive examples of the invention, the IV indices can also be used, outside the medical or surgical field, for example for the evaluation of the stress felt by a living being.

Claims

Claims

1. Method for calculating at least one index (IV) of heart rate variability, implemented automatically by a processing unit (5), from an initial RR series consisting of a plurality of samples (RRi) which are respectively a function of the time intervals (ôti) which separate two successive heartbeats or from an RR series derived from the initial RR series, which method comprises a step of calculating (a) at least one index of heart rate variability in a sliding main window (F) of duration (T), with a sliding step of at least one sample and of a duration (d) less than the duration (T) of the main window, and in which the processing unit (5) is configured to, during said calculation step (a), detect local minima (Pj;) and / or local maxima (Mjj) in the main window (F), to divide the main window (F) into n (n>3) sub-windows (Fj) consisting of a first sub-window (FJ), one or more central windows (F2,...,Fn_i) and a last sub-window (Fn), to calculate a lower curve (Cmin) connecting at least the last local minimum of the first sub-window (Fi), at least each local minimum (Pj;) included in the central sub-window(s) (F2,.. .,Fn_i), and at least the first local minimum of the last sub-window (Fn) and / or an upper curve (Cmax) connecting at least the last local maximum of the first sub-window (FJ, at least each local maximum (Mj;) included in the central sub-window(s) (F2,...,Fn_i), and at least the first local maximum of the last sub-window (Fn), and to calculate said at least one index (IV) of heart rate variability from the portion of the lower curve (Cmin) in the central window(s) (F2,...,Fn i), excluding the first sub-window (Fi) and the last sub-window (F n) and / or from the portion of the upper curve (Cmax) in the central window(s), excluding the first sub-window (FJ and the last sub-window (Fn).;

2. Method according to claim 1, in which the processing unit (5) is configured to determine, in each central sub-window (F2,...,Fn_i), an area (Aj) which is delimited at least by the lower curve (Cmin) and is above the lower curve (Cmin), and which extends over the entire duration (D) of the sub-window and / or an area (Aj) which is delimited at least by the upper curve (Cmax) and is under the curve upper (Cmax), and which extends over the entire duration (D) of the sub-window and to calculate said at least one index (IV) of heart rate variability from said areas (Aj).

3. Method according to claim 1, in which the processing unit (5) is configured to determine, in each central sub-window (F2,...,Fn_i), an area (Aj) which is delimited by the lower curve (Cmin) and by the upper curve (Cmax) and which extends over the entire duration (D) of the sub-window, and to calculate said at least one index (IV) of variability of the heart rate from said areas (Aj).

4. Method according to claim 2 or 3, in which the main window (F) is divided into at least four sub-windows (n>4) and in which the processing unit (5) is configured to select from among said areas (Aj), the one whose value is the lowest, and to calculate said at least one index (IV) of variability of the heart rate from this minimum area (Smin).

5. Method according to claim 4, wherein the processing unit (5) is configured to calculate the ratio (Smin / S) between said minimum area (Smin) and the total surface area (S) of the main window (F) and to calculate said at least one index (IV) of heart rate variability from this ratio (Smin / S).

6. Method according to any one of claims 1 to 5, in which the duration (T) of the main window (F) is between 32s and 256s, and preferably is equal to 64s.

7. Method according to any one of claims 1 to 6, wherein the number (n) of sub-windows is greater than or equal to 4.

8. A method according to any one of claims 1 to 7, wherein the sub-windows have identical durations.

9. Method according to any one of claims 1 to 8, wherein the duration of the first sub-window (FJ is equal to the duration of the last sub-window (Fn).

10. Method according to any one of claims 1 to 9, in which the calculation step (a) is carried out on an RR series which is obtained at least by resampling the initial RR series or an RR series derived from the initial RR series, at a sampling frequency (Fe) greater than or equal to 4Hz, and preferably equal to 8Hz.

11. Method according to claim 10, wherein the duration (d) of the sliding step is greater than or equal to the inverse (1 / Fe) of the sampling frequency (Fe), and preferably is equal to 1s.

12. Method according to any one of claims 1 to 11, in which the calculation step (a) is carried out on an RR series which is obtained at least by normalizing the samples (RRi) of the initial RR series or the samples (RRi) of an RR series derived from the initial RR series, over the entire duration (T) of the main window (F).

13. Method according to any one of claims 1 to 12, wherein the calculation step (a) is carried out on a series (RR) which is obtained at least by filtering the initial RR series or an RR series derived from the initial RR series, in at least one predefined frequency band [fl - f2].

14. The method of claim 13, wherein the filtering provides filtering in at least one high frequency band [0.15Hz-0.4Hz].

15. A method according to claim 13 or 14, wherein the filtering allows filtering in at least one low frequency band [0.04Hz; 0.15Hz],

16. Method according to any one of claims 1 to 15, wherein the processing unit (5) is configured to display the heart rate variability index(es) (IV) in real time on a screen.

17. Device for calculating at least one index (IV) of heart rate variability, comprising at least one processing unit (5) configured to automatically implement the calculation method of any one of claims 1 to 16.

18. Device according to claim 17, further comprising means (1, 2) for acquiring an analog cardiac signal and means (4) for sampling this analog cardiac signal, and in which the processing unit (5) is configured to also construct said initial RR series from the signal sampled by the sampling means (4).

19. Device according to claim 18, wherein the processing unit (5) is configured to automatically implement the calculation method of any one of claims 1 to 16 in real time, as the cardiac signal is acquired and the initial RR series is constructed.

20. Computer program product comprising program code instructions and allowing, when executed by at least one processing unit, to carry out the calculation method of any one of claims 1 to 16.

Images