Microsoft Word - 10ejtmMyoRev1_3_ 131-138 Graupe E-pub Apr08.doc Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 2010; 1 (3): 131-138 - 131 - Recognition and prediction of individual and combined muscular activation modes via surface EMG analysis Daniel Graupe Department of Electrical & Computer Engineering and Dept. of Neurology & Rehabilitation Medicine, University of Illinois, Chicago, IL, USA Abstract The paper discusses how recognition of individual and combined muscular activation modes (functions) and the prediction of intended such modes can be accomplished by identifying parameters of noninvasive surface EMG signals. It outlines the mathematical analysis of surface EMG signal to facilitate such recognition and related prediction, including recognition of intention (in terms of attempts) to activate motor functions from the EMG, without accessing the CNS itself, in cases where a patient, say, a high-level amputee does not have the final-activation muscles and joints. The EMG activity thus allows to interpret and recognize CNS commands from minute variations in the parameters of surface EMG signals that record changes in the firing of motor neurons triggering contractions in related muscle fibers. We note that although in popular media this is sometimes referred to as detection of “thoughts”, no thoughts are detected, but only motor-outcomes of thoughts as found in the EMG signal. Examples of concrete cases where such recognition or prediction were accomplished in the author’s lab and in devices that came out of that lab, are given as are references to these in the literature over the last 35 years. Key Words: prediction, combined muscular activation modes, surface EMG European Journal Translational Myology - Myology Reviews 2010; 1 (3): 131-138 1. The role of the surface EMG signal 2. Analisis 2.1 The Basic Temporal Model of the Surface EMG Signal 2.2 Non-Stationarity: Stochastic Time Series Approach 2.3 Non-Stationarity of Model: The Wavelet Approach 2.4 Neural Network for Clusteing 2.5 Discrimination between Different Motor Functions from EMG Parameters 3. Applications 3.1 Discrimination of Six Contraction Modes in EMG Recorded at Shoulder Electrode: To Control Multi- Degree-of-Freedom Above-Elbow Prosthesis 3.2 Discrimination of Several Muscle Contraction Modes in EMG from Electrodes - To Control Activation of Steps for Walking - via Electrical Stimulation of Paraplegia 3.3 Fixed-Time-Interval Prediction of Onset of Sleep Apnea via Submental EMG 4. Conclusions and perspectives The paper discusses how recognition of individual and combined muscular activation modes (functions) and the prediction of intended such modes can be accomplished by identifying parameters of noninvasive surface EMG signals. The EMG activity thus allows to interpret and recognize CNS commands from minute variations in the parameters of surface EMG signals that record changes in the firing of motor neurons that trigger contractions in related muscle fibers. Recognition and discrimination between several modes of muscular contractions were first reported in 1975 [7], as was the related mathematical derivation. It was subsequently employed by others, such as in [3,18,20]. Such recognition and prediction is important for controlling prosthetic and orthotic devices [7,14] and for controlling electrical stimulation on paraplegic patients [9,10,12]. It is also important for detecting onset of tremors in patients suffering from neurological diseases such as Parkinson Disease (PD) and Essential Tremors [1], in order to adaptively control deep brain electrical stimulation (DBS). Other applications are to predict and detect onset of sleep apnea and hypopnea events [24], to allow early warning of such events, to allow high-level above-elbow [13] and shoulder- disarticulation [14] amputees to activate at will (via a single pair of shoulder-EMG electrodes to activate three below-amputation joints of a prosthesis (six functions), namely, elbow-flexing/extension, wrist rotation (wrist pronation/supination) and hand’s grasp- opening/closing [13,14], or to activate the same six Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 132 - functions on a robotic arm, as was performed and filmed in the author’s lab since the early 1970’s. The paper outlines the mathematical analysis of surface EMG signal to facilitate such recognition and related prediction, including recognition of intention (in terms of attempts) to activate motor functions from the EMG, without accessing the CNS itself, in cases where a patient, say, a high-level amputee does not have the final-activation muscles and joints. We note that although in popular media this is sometimes referred to as detection of “thoughts”, no thoughts are detected, but only motor-outcomes of thoughts as found in the EMG signal. In the coming sections, the paper will discuss the role of the surface EMG signal, analysis methods for retrieving complex information from the surface-EMG signal and specific examples of concrete applications, mainly from the author’s lab. 1. The role of the surface EMG signal We note that parameters of surface EMG signals for a given limb function or combination of such functions (which is executed by contractions of certain related muscle fibers that are activated by their appropriate motor neurons), start to change, at least slightly, even before the actual desired function is executed. Furthermore, in order to execute such functions, the motor cortex must first balance the contractions of other muscle fibers upstream from the joint involved to balance the body for balanced execution of the desired end function. Hence, changes in shoulder and upper arm muscle take place well above an elbow that is to be bent or a wrist that is to be rotated. The above serve to allow a person amputated at the shoulder to be able to control his/her missing elbow (actually, an elbow prosthesis) or even his/her missing hand’s grasp function. It also allows to discriminate between different arm functions from the same surface EMG signal, since the surface EMG is a spatial integration of many hundreds or more motor neurons that may activate different muscle fibers. Furthermore, not all these fibers contract in the same manner for different major functions (elbow or grasp, elbow bending or extension). The model in Fig. 1 is a schematic for the mathematical input-output model. The filter block (denoted: Muscle and skin tissue) is the transfer function between input and output. It encompasses the propagation of the motor unit’s propagation-response through the muscle fibers [21] and the effects of skin tissue all the way to the skin electrode where the surface EMG is measured. Through the recording of surface EMG at the quadriceps muscles of complete thoracic-level paraplegics during electrical stimulation at that same location, we verified this model. The recording of this signal, which we termed as Response- EMG, was possible since, in complete paraplegics the only EMG that can be measured at this stimulated site is due to the stimulation. By separating the resulting EMG from the artifact of the stimulus, the resulting EMG was nothing but a very strong motor-unit action potential, which was in fact the sum of action potentials from many motor-units firing simultaneously due the stimulation pulse [6,8,10]. Fig. 2 shows the response EMG at a paraplegic’s quadriceps during FES stimulation immediately after the patient stood up (via FES, of course) and then after prolonged stading, when muscle fatigue started to take its toll. At extreme left the stimulus’ artifact is visible. We note that similar results are given in [17]. The EMG is a neuromuscular signal that is, of course, controlled by the CNS. It allows indirect reading of certain “thoughts”, namely those related to muscular activity, without having to record it at the CNS itself. All such applications require adequate placements of the EMG electrodes so that they access the muscle fibers whose contractions are of interest and hence, indirectly, the motor neurons of interest. They also require utilization of a broad frequency band in the recording and the subsequent signal processing of the recorded EMG signals. Finally, and not least, they require utilization of rigorous advanced signal processing methods to deal with stochastic time varying signals and to be sensitive to detect changes in the mathematical parameters involved, such as AR or ARMA analysis [4], wavelet transform (WT) analysis [16] and blind signal separation (BSS) algorithms, such as ICA [2] and AMUSE [22], in presence of environmental and recording noise. With the advent of extremely fast (GHz rates) and miniature computer hardware, the practicality of these applications and their cost are no more prohibitive, as processing speed is well above what is needed even for the most complex algorithms that are and will be involved, noting that all biological signals and processes lie in or below the KHz range (a million times below the hardware rate) and that micro- miniaturization is widely available. 2. Analisis 2.1 The Basic Temporal Model of the Surface EMG Signal In Fig. 1 we show surface EMG signals as an output y(t) of a system whose input in an array of (seemingly – for the outside observer) random trains of firings (impulses) of individual motoneurons denoted as xi (t), where t denotes time and where i = 1, 2, . N, N being the number of neurons considered (recorded) at a given pair of surface EMG electrodes. The muscle and skin fibers between the motor neurons above and the EMG skin electrode (that is non-invasive) is considered as a low-pass filter that transforms the train of pulses into the surface EMG signal y(t) as above. The input array {xi} of impulses of random arrival time can be considered to be a Poisson process. Furthermore, in [15] it is shown that if a Poisson input is applied to a low pass filter having a dominant time Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 133 - Fig. 1 The Surface EMG Signal (a) Schematic of the signal (b) The surface EMG signal as a stochastic process (time series)                            Fig. 2 Response EMG during stimulation in paraplegia constant that is considerably longer than the mean of the Poisson process then the filter’s output converges to a Gaussian process. In [7] we have shown that the surface EMG signal satisfies or approximates a Gaussian process, in [19] it is experimentally shown that indeed the surface EMG is Gaussian. Since a Gaussian process is linear, we have shown in [7] that it can be modeled by a linear AR (pure autoregressive), as in eqn (1) or a mixed ARMA (AR and moving average model), as in eqn (2), namely: y(k) = a(1) y(k-1) + a(2) y(k-2) + .... + a(m) y(k-m) + w(k) (1) y(k) = a(1) y(k-1) + a(p) y(k-p) + w(k) + b(1) w(k-1) + .. + b(q) w(k-q) (2) that are equivalent representations [4], k denoting discrete time (k = 0, 1, 2,....), and where w(k) is discrete white noise, satisfying: E[y(k-i) y(k-j] =0 if i differs from j for all i and j (3) and E [y 2 (k-i)]= W < Infnite (4) I and j being integers and E[...] denoting Expectation. We note that w(k) is the innovations of y(t), being the true residual to which y(t) converges, namely, the prediction error when the model orders m, p and q of eqn. (1) or (2) are correct [4,7]. It is shown in [1] that whereas theoretically m should be infinity, in practice, some low value (between 4 and 10) suffices. Orders p, q of eqn (2) can be computed via a(i) of eqn (1), as shown in [4]. 2.2 Non-Stationarity: Stochastic Time Series Approach The problem of distinguishing between modes or combinations of modes of muscle contractions, of predicting the onset of certain oncoming types of contraction or of activating certain limb or other motor functions, as in robotics, prostheses, orthotics, or of control of electrical stimulation, including deep brain stimulation, reduces to identifying the parameters {ai} of eqns, (1) or (2) or their equivalents. In doing so, one must keep in mind that the parameters {ai} of eqns. (1) and (2) are not stationary, but always change over time. Parameter identification is in itself a stochastic process that requires a minimal number of samples of the EMG data to be at all possible, the number being at least as high as the number of parameters that need be identified. We cannot artificially sample at a high rate to get more samples. Sampling must be at least at the Nyquist rate, namely twice the highest frequency at which there is any information in the signal [4]. But over-sampling is wasteful and may yield non- information. Furthermore we must guarantee that we use an algorithm whose identification results are mathematically rigorously proven to converge (convergence “in probability” or “with probability one”) and unbiased (not converging to a constant error). Algorithms not proven to converge “in probability” are useless. If the conditions of convergence in probability and unbiasedness are satisfied, and given enough samples, we shall converge to the true but unknown parameters, if the number of parameters is adequate [4]. Ref. [4] guides us also to finding an adequate number of parameters. The data length must be such that the parameter estimates seems to reach steady state values. A rule of thumb is that we need at least 100-120 samples per parameter if a statistically efficient [4] algorithm is employed. In cases of signal tendency towards instability (it never reaches infinity but may start to grow rapidly for some time), one has to follow the discussion in [4] on unstable processes, but not to avoid identification. Still, even at 120 samples per 4 parameters (480 samples total) at a sampling rate of 350-400/sec (adequate for surface EMG), we need approx 1 second to complete the identification (computer time is no detriment with today’s technology). This may be just about the limit to what one may require for EMG control of fast movements. When movements change, then parameters change too and we must re-identify. Now, there are algorithms that require many more samples than others. One must therefore select an algorithm that can be mathematically proven to be statistically efficient (namely, satisfying the Cramer- Rao bound [4] with an equality). Fortunately there is a rather simple algorithm that is convergent in probability unbiased and that can be proven to be statistically efficient. This is the sequential - (recursive-) least squares algorithm [14]. Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 134 - Fig. 4 Single EMG channel controlling 6 functions (three joints) of above elbow prosthesis. Top: Schematic of Prosthesis. Bottom: System’s block d Fig. 3 EMG signal from single source forms 2 different 2-dimensional vectors to detect two fuctions Doing so will rigorously satisfy that we can cope with non-stationarities in {ai}. To avoid getting stack it is best to restart identification from 0 with each new identification run. Even better is to run N parallel identification algorithms (N being the number of samples per one complete run, say, 480 in our example) and to consider only the identified parameters from the last run that was completed. Thus identification is updated at each sample. 2.3 Non-Stationarity of Model: The Wavelet Approach An alternative to the algorithm of Section 3.2 is to perform a Wavelet Transform (WT) of the data y(t) above [19]. The WT is a time-frequency transformation of data, such that it yields a time- varying frequency representation of that data. It can be viewed as an extension of frequency transform analysis (the Fourier Transform) for processes whose parameters (say, spectrum) varies over time. A parallelism between the WT and the discrete stochastic time-series analysis of Sect. 2.2 is discussed in [23]. The continuous time WT (CWT) is given by: where XWT ( τ, s) is the CWT of a time-domain signal x(t) , τ denoting the dilation, namely, the location of the WT as it is shifted in time (or in space) and s denoting scale and corresponds to the frequency information, and the asterisk denoting complex conjugate. Hence, XWT (τ, s) is clearly a two dimensional function, of time and frequency. The dicrete form of eqn. (5) yields the discrete WT (DWT), which is easily computed. It is conveniently presented in terms of an array of freqency bands (band-pass filters) where the spectrum in each band changes over time. 2.4 Neural Network for Clusteing Since the parameters as obtained via 2.2 or 2.3 are never exactly repeatable since the EMG signal itself is never repeatable [even an opera singer can never duplicate the exact signature of his/her voice, namley the voice signal, which rather resembles Fig. 1 b each time (s)he repeats any given aria], one must cluster parameters that represent any given muscular or joint contraction/movement [5,6]. An afficient way to perform such clustering is doing it with artificial neural newtworks (NN). Differnt NNs can be used, such as ART, Counter-propagation, LAMSTAR and others [5]. In [11] we presented an ART algorithm for clustering EMG signals for controlling electrical stimulation for walking by thorachic-level paraplegic patient. 2.5 Discrimination between Different Motor Functions from EMG Parameters Discrimination between different motor functions or modes of muscle contraction is accomplished by detecting differences in one or more parameters that are consistent for a given mode of contraction relative to others. The parameters considered are either time series parameters (Sect. 2.1), or changes over time in Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 135 - (a) Elbow flexion (b) Elbow extension           (c) Grasp closing (d) Grasp opening Fig. 5 High above-elbow amputee operating 4 limb functions by EMG from shoulder electrode pair Fig. 6 EMG electrode placement to activate above-elbow prosthesis (patient as in Fig. 5) one or more frequency DWT bands (Sect. 2.2) that are different at certain times than those at same band for a different mode of contraction. Neural networks can be employed also to aid in such discrimination, especially if they are somewhat hidden from direct observation. Center of clusters can serve as parameters for discrimination. Fig. 3 illustrates the basic mathematical problem of discrimination in terms of m-dimensional parameter vector-space [in case of Fig. 3, we employ a 2-D space, describing identification via only 2 parameters, namely, a1 and a2 – see eqn. (1) above]. Obviously, discrimination requires that the sub-space of one function (in case of artificial limb control by EMG as in Fig. 5c – Grasp Closing) will not overlap another function to be discriminated. If overlap occurs over a very small percentage of a subspace, training of the user (patient, in this case) must be carried out to reduce the probability of discrimination error. Else, either finer analysis must be done (more sensitive and more Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 136 -               (a)-schematic (b) T-9 complete Fig. 7 EMG electrodes at above spinal cord lesion of thoracic-level complete paraplegic patient, to control activation of Electrical Stimulation for ambulation repeatable identification algorithms to reduce the size of the subspaces, if at all possible), improved NN clustering methods, or electrode placement must be changed. 3. Applications 3.1 Discrimination of Six Contraction Modes in EMG Recorded at Shoulder Electrode: To Control Multi- Degree-of-Freedom Above-Elbow Prosthesis Since our initial study of [7], we have applied the signal discrimination via the mathematical approaches of Sect. 2 to control up to six limb functions (Elbow extension, elbow flexion, wrist pronation, wrist supination, grasp opening, and grasp closing) to several high above–elbow and shoulder disarticulation amputees [13,14]. See Figs. 4, 5 and 6. In the above, we have used a single pairs of EMG electrodes, located on biceps, between biceps and triceps (this gave often a “richer” mixture of muscle fibers involved in the artificial-limb’s functions that we wanted to activate by EMG), or shoulder (if the stump did not allow a lower placement). We requested the patient to pretend and try to use his/her remaining muscles as if he could a ctually perform the intended function. Still, it is possible, if necessary to have patient try to perform other task, as long as those do not activate other limbs or joinys. Response was of the order of a second or more. However, this was all done from 1961 to 1973 (then our lab switched to work on paraplegia - FES for walking - and other topics) and computer processing time was very slow. With today’s (2010), the delay will be below a second, and this mainly since a reasonable number of samples must be collected for parameter identification to converge. 3.2 Discrimination of Several Muscle Contraction Modes in EMG from Electrodes - To Control Activation of Steps for Walking - via Electrical Stimulation of Paraplegia Since early 1982 [10] we fisrt succeded to have thoracic-level paraplegic patients (T-1 to T-12) activate steps via functional electrical simulation (FES), where the command for taking a left step or a right step came from processing EMG signals taking from EMG electrodes placed on the patient’s chest as shown in Fig. 7) as we reported in [10] and in later publications [8,9,11,12]. All patients had a complete SCI lesion, and hence no motor function and no sensation below the lesion (Fig. 7). As in the cases discussed in Sect. 3.1, we requested the patients to contract soulder, chest and arm muscles as if they were healthy and could walk. 3.3 Fixed-Time-Interval Prediction of Onset of Sleep Apnea via Submental EMG In our recent work on predicting the onset of sleep apnea and hypopnea events [24], we found that Submental EMG, we developed an artificial neural- Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 137 - network-based system that employs 5 sources of continuously monitored data for that prediction. The neural network (NN), namely, the LAMSTAR network also ranks the source signals concerning their significance to successful prediction and computed the Submental EMG to be the most significant of all sources. Prediction is of an event to occur within 30 seconds, between 30 and 60 seconds and between 60 and 120 seconds. Success rate in our study was 81% [24]. 4. Conclusions and perspectives We have shown that surface EMG is a powerful and totally noninvasive tool for medical diagnosis and prediction in many fields of medicine, provided that they are processed with adequate mathematical extraction and parameter identification algorithms. Furthermore, the progress in ultra-high speed computation, the most powerful mathematical procedures are already fast enough for real-time control of patient functions, be they prosthetics, orthoses, elecrical stimulation, (includiding deep brain simulation, paraplegia, etc.), or other body functions and for early warning of events (sleep apnea, epilepsy and beyond). Also, technology allows this to be achieved compactly, to fit the needs and comfort of patients. Corresponding Author Prof. Daniel Graupe, Dept. of Electrical & Computer Engineering and Dept. of Neurology & Rehabilitation Medicine, University of Illinois, Chicago, IL 60607, USA. E-mail: dangraupe@gmail.com References [1] Basu I, Graupe D, Tuninetti D, Slavin K. Stochastic modeling of the neural activity in the subthalamic nucleus and model parameter identification from Parkinson Disease patient data. 2010, ECE Report, Univ of Illinois at Chicago, February 2010. [2] Cichocki A, Amari SI. Adaptive Blind Signal and Image Processing, John Wiley, NY. WT 2002. [3] Doerschuk PC, Gustafson DE, Wilsky AS (), Multifunctional upper-extremity prosthesis control-signal generation using EMG signal processing, Proc. JACC (Joint Automatic Control Conf.), Denver, CO., 1979, pp. 414-419. [4] Graupe D, (1989), Time Series analysis, Identification and Adaptive Filtering, 2nd Edition, Krieger Publishing Co., Malabar, FL. 1989. [5] Graupe D. Principles of Artificial Neural Networks, Second Edition, Krieger Publishing Co., Malabar, FL. 2008. [6] Graupe D, Basseas S, Deng J. EMG response to FES at quadriceps of paraplegics and its use to control FES levels to combat muscle fatigue, Symp. International College of Surgeons on Spinal Cord Injuries, Las Vegas, NV, 1986. [7] Graupe D, Cline WK. Functional separation of EMG signals via ARMA identification methods for prosthesis control purposes, IEEE Trans Sys Man and Cyber 1975; 5: 252-256. [8] Graupe D, Kohn KH. A critical review of EMG control of functional electrical stimulation in paraplegics. CRC Critical Rev in Biomedical Eng 1988; 15. [9] Graupe D, Kohn KH. Functional Electrical Stimulation for Ambulation by Paraplegics, Krieger Publishing Co, Malabar, FL. 1994. [10] Graupe D, Kralj A, Kohn KH. Computerized signature discrimination of above-lesion EMG for stimulating peripheral nerves of complete paraplegics, Proc IFAC Symp Prosthetics Control, Columbus, OH, 1982. [11] Graupe D, Kordylewski H. Artificial Neural Network Control of FES in Paraplegics for patient responsive ambulation. Trans IEEE on Biomed Eng 1995; 42: 699-707. [12] Graupe D, Kohn KH. Transcutaneous functional neuromuscular stimulation of certain traumatic complete thoracic paraplegics for independent short-distance ambulation. Nerol Res. 1997; 19: 323-333. [13] Graupe D, Magnusen J, Beex A,. Microprocessor system for multifunctional control of upper-limb prostheses via myoelectric signal identification, IEEE Trans. Automatic Control, 1975, 23(4):538- 544. [14] Graupe D, Salahi J, Kohn KH. Multifunctional prosthesis and orthosis control via microcomputer identification of temporal pattern differences in single-site myoelectric signals, J Biomed Eng 1982, 4(1): 17-22. [15] Papoulis A. Probability, Random Variables and Stochastic Processes, McGraw Hill, NY, 1965. [16] Rao RM, Bopardikar AS. Wavelet Transforms, Addison Wesley, reading, MA. 1998. [17] Sandercock T, Faulkner JA, Alberts TW, Abbrect PH. Single motor unit and fiber action potentials during fatigue. J Appl Physiol 1985; 1073-1079. [18] Saridis GN, Newman MA, Gelfand S, Upper- limb EMG statistical analysis. Proc. JACC (Joint Automatic Control Conf.), Denver, CO., 1979; 399-404. [19] Scott RN. Myoelectric energy spectra. Med Biol Eng 1967, 303-305. [20] Sherif MH, Gregor RJ, Lyman J. The effect of load on myoelectric signals – The ARMA representation. Proc 2nd Ann Conf of IEEE Eng In Med & Biol Soc 1980. Prediction of muscular activation modes via surface EMG European Journal Translational Myology - Myology Reviews 1 (3): 131-138, 2010 - 138 - [21] Stalberg E, Tronely J. Single Fiber Electromyography, Mirvalle Press, Old Woking, Surrey, UK. 1979. [22] Tong L, Liu RW, Soon VC, Huang YF. Indeterminacy and identifiability of blind identification. IEEE Trans Circ & Sys 1991; 38: 499-509. [23] Veselinovic, D, Graupe D. On the parallelism in the enemble sense between time series models and disctrete wavelet transforms of stochastic signals. IEEE Trans Circ & Sys 2000; 47: 485- 489. [24] Waxman JA, Graupe D, Carley DW. Automated prediction of Apnea and Hypopnea using LAMSTAR artificial neural networks. Amer J Respiratory and Crit Care Med, 2010; 181, in press. Also web edition of the journal at: http://intl-ajrccm.atsjournals.org.