EP1588280A2

EP1588280A2 - System and method of processing electromyographic signals for the diagnosis of parkinson's disease

Info

Publication number: EP1588280A2
Application number: EP03815575A
Authority: EP
Inventors: Maria Chiara Carboncini; Gennaro De Michele; Bruno Rossi; Stefano Sello; Soo-Kyung Strambi
Original assignee: Enel Produzione SpA
Current assignee: Enel Produzione SpA
Priority date: 2003-01-31
Filing date: 2003-07-31
Publication date: 2005-10-26
Also published as: WO2004066832A2; AU2003253298A8; WO2004066832A3; AU2003253298A1; ITFI20030024A1

Abstract

A method of processing electromyographic signals (f(t), g(t)) for the diagnosis of Parkison’s disease in a patient, comprising a processing step based on a Wavelet analysis of said electromyographic signals, said analysis being capable of determining a value (P.D.I) of global cross-correlation energy, said value being indicative of the level of gravity of the disease in said patient.

Description

TITLE

SYSTEM AND METHOD OF PROCESSING ELECTROMYOGRAPHIC SIGNALS

FOR THE DIAGNOSIS OF PARKINSON'S DISEASE

DESCRIPTION Field of the Invention

The present invention relates to a system and a method of processing electromyographic signals particularly useful for the diagnosis of Parkinson's disease and, more particularly, for an objective evaluation of the level of gravity of the illness.

Backgrounds of the invention

The methods currently employed for assessing the level of gravity of Parkinson's disease and monitoring its time evolution are essentially of a clinical and/or an instrumental type.

The semi-quantitative methods, such as the Unified Parkinson's Disease Rating Scale (U.P.D.R. S . ) , form part of the first type.

This method consists of the filling in by the patient of questionnaires in which the symptoms are attributed numerical values that increase with the gravity of the symptom. The specific scale used for this purpose depends on the particular purpose that the examiner sets himself: Unified Parkinson's Disease Rating Scale (U.P.D.R.S. ) , Hoen Yar (H.Y.), Fatigue scale, etc., Since the semi-quantitative evaluation consists of a questionnaire that has to be completed by the patient under examination, it can easily be altered by subjective factors, among them the emotive state Of the patient, which also exerts a fundamental influence on his motorial performance .

The second type comprises C.A.T. and N.M.R. analyses of the brain and the analysis of the kinetic and electromyographic (EMG) parameters recorded during specific movements and on particular muscles .

The instrumental Evaluations (N.M.R., C.A.T.) are reliable examinations that detect in a detailed manner the areas of the brain involved in the worsening of the pathology.

However, these examinations do not provide any information about the clinical characteristics of the subject and some days have to pass before the results become available.

As far as the patient is concerned, moreover, they involve a massive exposure to strong radiations for the entire duration of the examination, which therefore cannot be repeated at brief intervals.

Examinations based on the analysis of kinematics and electromyographic (MG) parameters are generally founded on: a) quantification of the ^' area of the electromyographic trace .

However, the evaluation of the area of the trace does not always represent a significant index, and this not only because equal areas can correspond to different traces, but also because it is possible for significant variations of the patient's clinical state can give rise to area variations of little statistical significance; b) visual description of the trace.

However, the evaluation of the gravity of the illness based on a direct examination of the EMG graphs can hardly be univocal, because there can exist seemingly similar graphs that correspond to clinically different situations: this is a serious drawback when situations of incipient illness are being considered.

Furthermore, like the greater part of biomedical signals, even the EMG signals are of the non-stationary type with highly complex time-frequency characteristics. They generally consist of brief high-frequency components closely , spaced in time, followed by low-frequency components differing only very slightly in frequency.

The mathematical instrument most commonly employed in biomedical applications is the Fourier transform, which decomposes a given signal into its frequency components.

However, this technique requires the signal that is being examined to be stationary, that is to say, its frequency content must not undergo variations in time.

But, as a general rule, biomedical signals are non- stationary and have an intermittent and variable frequency content .

The limits associated with Fourier analysis can be partly overcome by using a mobile-window Fourier transform, also known as Gabor transform, which supposes the signal be almost stationary in a brief interval of time. Nevertheless, even this procedure is associated with limitations and inaccuracies that become particularly evident during the signal partition phases, so that the Gabor transform becomes de facto disadvantageous for the analysis of signals that involve different frequency ranges .

Objects and summary of the Invention

The object of the present invention is to solve the aforesaid problems associated with prior art techniques and to provide a method as described in claim 1.

The present invention also concerns a system as described in the independent claim 12, specifically _.suitable for implementing the aforesaid processing method. Another object of the present invention is to provide a computer product, particularly a software, capable of implementing a method in accordance with the present invention.

The principal aim of the present invention is to provide a methodology suitable for determining efficient qualitative and quantitative parameters useful in the assessment of the level of gravity of Parkinson's disease. This is done by means of a particular analysis of the EMG signals: the Wavelet cross-correlation analysis.

The principal advantage of the method in accordance with the present invention is that when the Wavelet cross- correlation analysis is carried out on a pair of EMG signals, it makes it possible to find characteristic Wavelet Analysis Maps .(W.A.M.) and a numerical index, known as the Parkinson Disease Index (P.D.I.) , correlated in an unequivocal manner with the level of gravity of the disease. A second advantage derives from the fact that the processing of the maps and the calculation of the P.D.I, can be carried out also by means of a dedicated hardware system connected with the data acquisition equipment used for detecting the EMG signals, thus making it possible for data about the patient's state to be obtained in real time and to make a quick assessment of the effect of pharmaceuticals and therapies .

Brief description of the drawings

Further advantages, characteristics and use of the present invention will be brought out by the following detailed description of a preferred embodiment thereof, which serves purely as an example ■ and is not to be considered limitative in any way, the description making reference to the attached drawing in which:

Figure 1 shows a block diagram of a system in accordance wit the present invention; - Figure 2 shows a flow diagram of the processing method in accordance with the present invention; and

Figures 3 to 9 show schemes and graphs of the tests carried out for the validation of the method in accordance with the present invention. Detailed description of the invention

Referring to Figure 1, a system in accordance with the present invention comprises means for the acquisition of data corresponding to electromyographic signals (EMG) .

A central control unit, comprising means for processing the data, receives these data as input, processes them and, as output, provides corresponding information relating to predetermined parameters in a format suitable for being displayed on a monitor, printed or reproduced in some other manner. According to a first embodiment of the system in accordance with the present invention, the analysis and processing of the input data (identifying the EMG signals) is carried out off line.

The acquisition of the EMG signals is carried out with the help of detection equipment of a conventional type, capable of memorizing the acquired signal in digital format in the form of files on any kind of support (floppy disk, CD or other) .

Each file, in ASCII format for example, contains some general information for the processing (type of signal, number of data items, acquisition rate) and a series of number pairs indicating respectively the time variable and the value of the selected EMG channel .

The files are then read by the control unit, which thereupon commences the processing procedure for the purpose of obtaining the previously indicated maps and gravity indexes .

All the processing steps of the method in accordance with the present invention are implemented by means of software and specific procedures designed specifically for that purpose . The remainder of the present description will be dedicated to a compete treatment of the method from the analytical point of view, but will not provide a detailed description of any particular implementation software of the method, since such implementation can depend on the employed hardware platform, the chosen programming language and also on the different approaches in the numerical implementation of the algorithms, all parameter within the knowledge of a person skilled in the art.

In a second preferred embodiment, on the other hand, a system in accordance with the present invention comprises the instrumentation for acquiring the EMG signals, their transfer in real time to the central control unit, which processes them in real time and instantaneously provides the results o the analysis. The subsequent Figure 2 shows a flow diagram of the processing method in accordance with the present invention.

The macro-operations on the input data are represented in sequence in the diagram of Figure 2. At a fist step SI, the cross-correlation Wavelet transform of the EMG signals is calculated.

Preferably, the Wavelet cross-relation analysis is applied in an iterative manner to a sample pair of EMG signals acquired during the signal measurement session. As a general rule, the analysis of 10 pairs of EMG signals is sufficient for an adequate statistical valuation (Step S2) of the results that will bring out the most probable values and the respective indeterminacies of the calculated integral and local Wavelet quantities.

This statistical analysis leads to the determination of the Wavelet parameters and the P.D.I. A subsequent step S3 of the method is provided for the calculation of a Wavelet Average map (W.A.M.) useful for obtaining local graphical representation of all the Wavelet cross-correlation maps for a given patient.

This makes it possible to display (step S4) the value of the P.D.I, and the morphology of the energy content of the most representative cross-correlation for a given patient both in time and in frequency and to carry out qualitative comparisons between different patients or for a given patient in the course of time. In particular, it makes it possible to assess with considerable sensitivity the effect of a given pharmacological therapy on Parkinson patients by making direct use of the EMG signals.

In the following a detailed mathematical treatment is given to describe the calculation algorithm that has been used.

As is known, the continuous wavelet transform represents an optimal localized decomposition of time series, x (t) , as a function of both time t and frequency (scale) a:

where ψ is an analyzing wavelet verifying the following admissibility condition:

where: Cψ and ^(ω) are the admissibility constant and the Fourier transform of the wavelet function, respectively.

The parameters a and r denote the dilation (scale factor) and translation (time shift parameter) , respectively.

The local wavelet spectrum is then defined as :

where k₀ denotes the peak frequency of the analyzing wavelet ψ .

From the local wavelet spectrum we can derive a mean or global wavelet spectrum P (k) :

P(Jfc)= f P(k, t)dt

—00

which is related to the total energy E of the signal x (t) by:

+∞

E = j P(k)dk

0

The relationship between the ordinary Fourier spectrum P_F(ω) and the global wavelet spectrum P (k) is given by:

^' indicating that the global wavelet spectrum is the average, of the Fourier spectrum weighted by the square of the Fourier transform of the analysing wavelet ψ shifted at frequency k. To the end of the present invention the Morlet wavelet function has been advantageously selected which allows to get amplitude and phase information on the signals under examination at a given scale a-.

Let W_f ( -t, T) and W_g [a, τ) be respectively the wavelet transforms of two real signals: f (t) and g (t) .

We define the wavelet cross-scalogram, (also cross- correlation wavelet transform) as :

W_fg(a, r) = W_f*(a, τ)W_g(a, τ)₁

where the symbol *" indicates the "complex conjugate" operator. When the analyzing wavelet is complex, the wavelet cross-scalogram W_fg(fl, r) is also complex and can be written as:

W_fg(a, T) — CoWfg(a, r) — iQuadWf_g(a, r)

If f (t) and g(t) e ^~ (<R) the following important relation holds :

J f(t)g(t)dt = l/c j ] CoW_fβ(a, τ)dτda,

0 -co

which relates the Co-spectrum to the correlation integral of the signals. The local wavelet cross-correlation spectrum is given by:

I W_fg(a, τ) |² = 1 CoW_fg(a, r) + | QuadW_fa(a, τ)

The integration of the local wavelet cross- correlation spectrum over r gives the global wavelet cross-spectrum.

The integration of the global wavelet spectrum over CL is a scalar quantity related to the whole cross correlation energy i.e. the index that has been indicated as P.D.I, heretofore. Furthermore, wavelet analysis allows for a precise evaluation of the intermittency property which is a well-established feature for many complex physical systems.

The energy at a given scale is not evenly distributed in time and/or in space and any variable alternates high activity with quiescence.

Mathematically, we quantify the intermittency level, F, through the deviation of the related distribution function with respect to a Gaussian behavior, by the kurtosis parameter, k:

U; ₌ ,C_3 = -^--3

where μ (for i=2 , 4) denotes the i th central moment. The correction constant 3 is related to the kurtosis of the Gaussian distribution.

If K>0 , i.e. the distribution has sizeable tails which extend much further from the mean than standard deviation, the variable is said to be "intermittent" .

To the ends of the present invention, as local variable we consider the time evolution of the wavelet power spectrum coefficients, at a given scale.

With this quantity we can obtain important information on intermittent behavior at different frequencies, of the cross-correlation energy related to a couple of EMG signals.

With a view to validating the method described hereinabove, specific tests have been carried out both on patients (previously diagnosed in accordance with conventional methods) and on healthy subjects.

The experiments were carried in the manner about to be described, making reference also to Figure 3.

All the subjects were asked to sit on a chair with their right arm unsupported, unrestrained, and free to move . They were then instructed to keep their right arm straight, with the elbow extended and the forearm in the mid-prone position, while grasping the handle of a strap attached to a strain gange.

They were then asked to carry out an isometric flexion at about 70% of the maximum voluntary effort (isometric phase of the task) for about 1 second.

The force thus generated could be monitored on an oscilloscope screen.

Following an acoustic "go" signal, the subjects then had to extend the arm as quickly as possible, reaching the furthest extension they could.

The final position was maintained for 1 second. In this way it was possible to study several aspects of muscle activation: the relaxation of the antagonist muscle before the beginning of the movement, the three-phase pattern associated with the fast limb displacement and, finally, the postural activity necessary to maintain the arm extended in the horizontal plane.

This sequence was repeated ten times at regular intervals of 30 - 120 seconds in order to extract statistically significant results related to the recorded signals.

The recording of the electromyographic signals (EMG) was obtained by means of an integrated system, comprising an 8 -channel electromyograph and two telecameras connected with a main processor (an ELITE system, BTS for example) . The electromyographic activity was measured by means of surface electrodes spaced about 1 centimeter apart and arranged longitudinally with respect to the fibres of the muscle.

The recordings of the EMG signals were differentially amplified, analog-to-digital converted (at 500 Hertz, and anti-aliasing filtered with 12-bit resolution) and stored in a computer. For the analysis the activity of the major pectoral muscle indicated as p in figure 3, and of the posterior detailed muscle indicated at pd in figure 3, were recorded.

Furthermore, for all ill subjects there was made a diagnosis based on a qualitative scale to evaluate their motorial functions. This scale (UPDRS) consists of six sections that define motorial, cognitive, therapeutic and relational functions. Each function is assigned a score ranging from 0 (least ill) to 4 (most ill) and this value provides information regarding the progress of the illness.

Figure 4 show an example of the results obtained when the processing method in accordance with the present invention is applied to a pair of EMG signals recorded in the case of a diseased patient during the previously described ballistic movement.

Figure 5 provides the same information in the case of a healthy subject.

In both these figures the two uppermost graphs represent the respective recordings of the analyzed EMG signals .

The central panel shows the Wavelet map of the power cross-correlation obtained by using a chromatic scale for correlating the energy intensity of the signal . The horizontal axis represents the tine variable expressed in seconds, while the vertical axis shows the frequency values in accordance with a logarithmic scale.

The panel on the right shows the global power spectrum obtained by a time integration of the corresponding local power spectrum. The selected frequency range is comprised between 7.4 and 165 Hz, because this interval contains all the significant components of the cross-correlation energy.

At first sight, using qualitative criteria only, one can already deduce some important characteristics typical of the two different classes of cases: for the Parkinson patient (diseased) the power distribution seems very broad both in time and in frequency, with high values of the Wavelet power integral . In the case of a healthy subject, on the other hand, the power content is well localized and involves only a small range of frequencies and therefore implies a lower value of the Wavelet power integral .

According to the present invention, a more irregular and dispersed movement, well visualized and highlighted by means of Wavelet analysis, is related to the lower muscle synchronization and coordination. The unsteady and intermittent character (more extensive tails than those associated with a Gaussian distribution) quickly appears from a purely qualitative evaluation of the Wavelet maps. More precisely, Figure 6 shows the measure of intermittence as a function of the frequency calculated in both cases .

The average intermittency in the case of a diseased patient (full line) is generally greater than that of a healthy subject (broken line) .

In particular, around 18 and 35 Hz the Gaussian behavior of the cross-correlation energy evidenced by a healthy subject appears clearly intermittent for a diseased patient. Significant differences can also be noted around 60, 80 and 140 Hz.

Expressed as a function of the frequency, the intermittency derived from the Wavelet analysis represents a typical quantity of immediate interest for diagnostic purposes .

With a view to obtaining more reliable results, a statistical analysis can be advantageously carried out on the set of Wavelet maps, mainly by means of the integral quantities obtainable from the cross-correlation results, especially on the P.D.I, parameter.

The statistical procedure used for this purpose is as follows: as set of possible realizations there is selected the set of test sequences repeated for each subj ect ; the sample parameter is assumed to have a normal distribution; the confidence level (95%) is evaluated for the first two distribution moments; an F-test is carried out (using a Fischer-

Snedecor distribution) in order to assess the compatibility of the groups, i.e. to estimate the probability of having the same second order moment for different subjects.

As already mentioned, the test sequence was carried out and recorded for ten consecutive times in the case of each patient.

It is thus possible to represent the results of this statistical analysis in a graph (Figure 7) of the global

Wavelet power (P.D.I.) for each of the subjects examined.

The full squares in Figure 7 represent values relating to healthy subjects, while the empty squares represent values relating to diseased patients.

Figure 8 shows the global Wavelet power (P.D.I.) together with the associated statistical uncertainty.

Using this graph, it therefore becomes possible to correctly discriminate the presence of the disease and to obtain a clear and reliable classification as compared with the previous clinical evaluations.

The experimental data show that there exists a threshold between the two classes of subjects. This threshold corresponds to a P.D.I, value of about 3500 (a.u. ) .

By way of confirmation of the results obtained, the table reproduced below (Table 1) shows the results with the quantitative U.P.D.R.S. classification. Table 1

It can be noted that the condition of Patients 7, 9, 11, 12, 14 and 15 can be better evidenced and distinguished by making use of the P.D.I, parameter rather than the U.P.R.D.S. scale.

Use of the Wavelet technique as described

- ereinabove therefore makes it possible for the time evolution of the level of the illness in a given patient or the efficacy of a therapeutic treatment to be realistically followed on the basis of the variations of the P.D.I, parameter, which - by way of example - are shown for two patients in Figure 9. The present invention has so far been described in the form of its preferred embodiments, which constitute examples and are not to be considered limitative.

It is to be understood that other embodiments are possible, and all of these have to be considered as coming within the scope of the invention as defined in the claims attached hereto.

Claims

CLAIMS 1. A method of processing electromyographic signals (f (t) , g (t) ) for the diagnosis of Parkinson's disease, characterized in that it comprises a processing step based on a Wavelet analysis of said electromyographic signals, said analysis being capable of determining a ^"value (P- D.I.) of global cross-correlation energy, said value being indicative of the level of gravity of the disease in said patient .

2. A method in accordance with claim 1, comprising a step of providing at least one pair of electromyographic signals (f (t) , g (t) ) , each signal being relative to a predetermined muscle ( p, pd) of said patient.

3. A method in accordance with claim 2, wherein said at least one pair of signals (f (t) , g (t) ) is acquired during a predetermined sequence of movements of said muscles (mp, pd) .

4. A method in accordance with any one of the preceding claims, wherein said analysis comprises the calculation of a wavelet transform (W_f(-i_/τ), Wg-( ,τ)) for each of said electromyographic signals fff , g (t) ) .

5. A method in accordance with claim 4 , wherein each of said Wavelet transforms W_f( ,τ), W_g(α,τ) is based on the Morlet function.

6. A method in accordance with claim 4 or claim 5, wherein said analysis comprises the calculation fir each pair of signals of a cross-correlation Wavelet transform in accordance with the relationship:

W_fg(a, T) = W_f*{a, τ)W._g(a, r),

7. A method in accordance with claim 6, wherein said analysis comprises the calculation of a local Wavelet cross-correlation spectrum in accordance with the relationship :

I W , r) I² = I CoW_fg(a, r) \² ÷ \ QuadW_Sg{a,r) |²

8. A method in accordance with claim 7, wherein said analysis comprises the calculation of a global Wavelet cross-spectrum by means of an integration of the local Wavelet cross-correlation spectrum with respect to the variable t .

9. A method in accordance with claim 8 , wherein said analysis comprises the calculation of a value (P.D.I.) of global cross-correlation energy obtained by means of an integration of the local Wavelet cross-correlation spectrum with respect to the variable .

10. A method in accordance with any one of the preceding claims, comprising a step of applying a first statistical model to the set of electromyographic signal pairs (f(t), g(t)) acquired for the same patient.

11. A method in accordance with claim 7, comprising a step of applying a second statistical model to the set of calculated local Wavelet cross-correlation spectrums .

12. A system of processing electromyographic signals (f(t), g(t)) for the diagnosis of Parkinson's disease, in a patient, comprising means for the digital acquisition of said electromyographic signals (f(t) g (t) ) , characterized in that it comprises means of processing said signals capable of carrying out a Wavelet analysis of said electromyographic signals, said analysis being capable of determining a value (P.D.I.) of global cross-correlation energy, said value being indicative of the level of gravity of the disease of said patient.

13, A system in accordance with claim 12, wherein said acquisition means comprise an apparatus for detecting said electromyographic signals.

14. A system in accordance with claim 12 or claim 13, comprising also means for the analog-digital conversion of said acquired signals.

15. A system in accordance with any one of claims 12 to

14, wherein said analog-digital conversion is effected in real time .

16. A system in accordance with any one of claims 12 to

15, wherein said converted signals are acquired in real time by said acquisition means.

17. A system in accordance with any one of claims 12 to 15, wherein said acquisition means acquire at least one pair of electromyographic signals (f (t) g (t) ) , each signal being relative to a predetermined muscle (mp, pd) of said patient .

18. A system in accordance with claim 17, wherein said at least one pair of • signals (f (t) g (t) ) is acquired during a predetermined sequence of movements of said muscles (mp, pd) .

19. A system in accordance with any one of claims 12 to 18, wherein said processing means are capable of carrying out the calculation of a Wavelet transform (W_f(β,τ), Wg(fl,τ)) for each of said electromyographic signals (f (t) g (t) ) .

20. A system in accordance with claim 19, wherein said Wavelet transform (W_f(α,τ), W_g-(Λ_/τ)) is based on a Morlet function.

21. A system in accordance with claim 19 or claim 20, wherein said processing means are capable of carrying out the calculation of a Wavelet cross-correlation transform for each pair of signals in accordance with the relationshi :

W_f3(a,τ) = W (a, r)W_g(a, τ),

22. A system in accordance with claim 21, wherein said processing means are capable of carrying out the calculation of a Wavelet cross-correlation transform for each pair of signals in accordance with the relationship:

! W_fa( , r) |² = j CoW_fg(a, r) f + \ QuadW_fa(a, r) |²

23. A system in accordance with claim 22, wherein said processing means are capable of carrying out the calculation of a global Wavelet cross-spectrum by means of the integration of the local Wavelet spectrum with respect to the variable t .

24. A system in accordance with claim 23, wherein said processing means are capable of carrying out the calculation of a value (P.D.I.) of global cross- correlation energy obtained by means of an integration of the global Wavelet cross-spectrum with respect to the variable a.

25. A system in accordance with any one of claims 12 to 24, wherein said processing means are capable of applying a first statistical model to the set electromyographic signal pairs acquired for one and the same patient .

26. A system in accordance with claim 25, wherein said processing means are capable of applying a second statistical model to the set of the calculated local Wavelet cross-correlation spectrums .

27. A system in accordance with any one of claims 12 to 26, wherein said processing means comprise one or more software moduli.

28. A computer product, in particular a program, stored on a storing support or functioning on a data processing system, capable of implementing a method in accordance with any one of claims 1 to 11.