Myoelectric Pattern Identification of Stroke Survivors Using Multivariate Empirical Mode Decomposition

This study presents a novel feature extraction method for myoelectric pattern recognition using a multivariate extension of empirical mode decomposition (EMD), namely multivariate EMD (MEMD). The method processes multiple surface electromyogram (EMG) channels simultaneously rather than in a channel-by-channel manner. From mode-aligned intrinsic mode functions (IMFs, representing signal components over multiple scales) derived from the MEMD analysis, normalized amplitude distributions of the same-mode/scale IMFs across different channels were calculated as features, which serve to reveal the underlying relationship in the aligned intrinsic scales across multiple muscles. The proposed method was assessed for identification of 18 different functional movement patterns via 27-channel surface EMG signals recorded from the paretic forearm muscles of 12 subjects with hemiparetic stroke. With a linear discriminant classifier, the proposed MEMD based feature set resulted in an average error rate of 4.61 ± 4.70% for classification of all the different movements, significantly lower than that of the conventional time-domain feature set (7.14 ± 6.15%, p < 0.05). The results indicate that the MEMD based feature extraction of multi-channel surface EMG data provides a promising approach to modeling of muscle couplings and identification of different myoelectric patterns.


INTRODUCTION
The electromyogram (EMG) signal is the electrical manifestation of muscle activation.Surface EMG recordings play an important role for various biomedical applications such as identification of muscle activity patterns for diagnosing abnormal movements [1][2][3] and control of assistive devices for restoration of limb function [4][5][6][7][8][9][10].In these applications, multi-channel surface EMG signals recorded via a number of electrodes or high-density electrode arrays have been utilized to provide more detailed information than conventional mono or bipolar recordings [11][12][13][14][15][16][17].The currently available methods for identification of different myoelectric patterns routinely describe surface EMG patterns by concatenating all the features extracted separately from individual channels [5][6][7][11][12][13].Such an approach ignores co-activation or coordination of multiple muscles, which may be useful for identification of different myoelectric patterns [18].In contrast, a simultaneous processing of multiple surface EMG channels to extract cross-channel information may be helpful for effective interpretation of muscle couplings and identification of various myoelectric patterns.
Multivariate empirical mode decomposition (MEMD) is a recently developed approach for simultaneous processing of multiple data channels.It is a multivariate extension of original empirical mode decomposition (EMD) that acts as a fully datadriven algorithm adaptively decomposing a given time series into a number of narrowband oscillatory components, namely intrinsic mode functions (IMFs) [19,20].As the EMD is a greedy and complete transformation, decomposition of time series into IMFs can preserve the relevant information contained in different frequency components, and allows an efficient method for multi-scale analysis to search discriminating features [21].Previous efforts have shown that distinct patterns of electrophysiological data can be successfully characterized by features extracted across different IMFs [21][22][23][24][25][26].However, applying the standard EMD approach on each individual channel of multivariate time series data may suffer from a mode-misalignment problem, which is manifested by the same common oscillation modes across multivariate data appearing in the different-index IMFs.This would lead to the fact that the IMFs of different channels are not matched either in number or in scale [18,19,23].Such misaligned-IMFs derived from standard EMD can impact the distinctiveness of the resultant features, thereby affecting the performance in pattern recognition [23].This problem can be mitigated by employing the MEMD, which is designed to simultaneously process multivariate data to result in scale-aligned IMFs [19].Thus, the aligned IMFs extracted by MEMD facilitate the extraction of cross-information regarding the associations among multiple channels for better characterization of distinct patterns in multivariate data [19,23].
In this study, we investigated the feasibility of extracting cross channel information for myoelectric pattern identification by simultaneously processing a number of EMG channels using MEMD.Taking advantage of the mode alignment property of the MEMD, a feature extraction method was proposed for calculating amplitude distributions of the same-scale IMFs across multiple channels.The performance of the proposed method was evaluated with multi-channel EMG data recorded from hemiparetic stroke patients during performance of different intended wrist and hand movements of the affected limb.The advantage of the MEMD based feature extraction or modeling of muscle couplings was illustrated by enhanced identification of different myoelectric patterns of the stroke subjects.

262
Myoelectric Pattern Identification of Stroke Survivors Using Multivariate Empirical Mode Decomposition

BACKGROUND
The EMD acts as a fully adaptive data-driven method that is able to decompose a timeseries into a finite set of IMFs, denoted as ), by means of the so-called sifting algorithm.The EMD decomposition of a time series can be described as: ( where represents IMFs with order number of i ranging from 1 to N, and is a residual usually considered to be the (N+1)-th IMF.Therefore, can be used for concise description of the EMD method: (2) These resultant IMFs represent the inherent oscillation modes (also considered as frequency components over different scales of that time series) [20].Although the EMD method is an investigative approach to the characterization of univariate time series, the mode-misalignment introduced by standard EMD limits its further application on the multivariate time series data.Recently, MEMD has been proposed to extend the application of EMD to multivariate data [19], thus inspiring the improved analysis of signal structure couplings among multiple data channels based on MEMD.
The key of MEMD algorithm is the calculation of local mean.Unlike the EMD in which the local mean is computed by taking an average of upper and lower envelops obtained by interpolating the local maxima and minima of the signal, MEMD directly works with multivariate data (representing an n-variable time series, here denoted as ) of which the local maxima and minima are not well defined.To deal with this problem, multiple n-dimensional envelopes are generated by taking signal projections along different directions in n-dimensional spaces.These envelops are then averaged to obtain the local mean.Suppose that is a set of vectors (indexed by k) along the directions represented by angles on a (n -1)-sphere.The algorithm is briefly summarized as follows [19]: 1) Generate a point set based on the Hammersley sequence for sampling on an (n -1)-sphere.2) Calculate a projection of the multivariate input data along a direction vector , for all k, thus giving 3) Locate the time points according to maxima of the set of projected for all values of k, to obtain multivariate envelope curves .
5) For a set of K direction vectors, calculate the mean m(t) of the envelop curves: (3) 6) Iterate on the detail d(t) = v(t)m(t) until it becomes an IMF in a multivariate form.The above procedure is further applied to v(t)d(t).The stoppage criterion for multivariate IMF is similar to that for the univariate IMFs except that the equality constraints for the number of extrema and zero crossings are not imposed, as extrema cannot be properly defined for the multivariate signal.By projection, MEMD directly processes multivariate signals to produce the aligned IMFs.

METHODS
This section presents the MEMD-based feature extraction from multi-channel EMG signals for enhanced myoelectric pattern identification, as well as the data collection and signal processing framework for performance evaluation.

Dataset Description
The multi-channel data used in this study were part of a database collected in our previous study [12], including high density surface EMG recordings from stroke patients during their performance of different movements involving the affected arm.This dataset was selected for testing the proposed method.Twelve stroke patients with hemiparesis (Subjects 1-12) participated in this study and their physical characteristics and clinical assessment are presented in Table 1.The study was approved by the Institutional Review Board of Northwestern University (Chicago, USA).All participants gave their informed consent before the experiment.
During the experiment, each subject was asked to perform (or intend to perform) different elbow, hand and wrist functional movements of the affected arm following a video demonstration of a healthy person.In this study, only hand and wrist functional movements were analyzed (Table 2).The experiment protocol comprised 18 trials, with each trial consisting of 5 repetitions of the same movement pattern.For each repetition, the subject completed the task over 3 s and then relaxed for a period of 5-20 seconds.They were also allowed to rest for 3-5 minutes between trials to avoid muscular and/or mental fatigue.
During the performance of each movement, high density surface EMG signals of 89 electrodes were collected above the upper arm, forearm and hand muscles in the affected side of each subject.The Refa128 EMG system (TMS International BV, Netherlands) was used for the recording, and the sampling rate was 2000 Hz per channel.Since the 18 movement patterns investigated in this study involved wrist and hand motions, only the electrodes placed on the forearm and hand muscles were selected for analysis.It has been reported [27,28] that the use of bipolar (i.e., single Myoelectric Pattern Identification of Stroke Survivors Using Multivariate Empirical Mode Decomposition differential) EMG recordings helps to improve EMG classification performance, and it is indeed more clinically relevant.Hence, 27-channel bipolar surface EMG data were extracted from the high density surface EMG recordings to assess the performance of the proposed feature extraction method.Figure 1 illustrates the electrode positions of the 27 bipolar channels, including 24 channels arranged in a 3 ¥ 8 grid formation around the forearm and the other 3 channels placed to target the first dorsal interosseous, thenar and hypothenar muscles, respectively.

Signal Segmentation and Decomposition using MEMD
For the multi-channel surface EMG signals recorded from each trial, the onset and offset of each contraction were obtained from the dataset.Such information was determined manually when the dataset was created [12].Based on the determined onset Journal of Healthcare Engineering • Vol.  and offset, the multi-channel signal segment corresponding to each contraction was reshaped into a multivariate form as a 27-dimensional time series for MEMD analysis.
The resultant scale-aligned IMFs for each channel were obtained at the same time, thus facilitating the extraction of features reflecting the couplings among these channels.The MEMD was implemented in MATLAB environment (version 2012a, The Mathworks, Inc., MA, USA), using available source code from the webpage of MEMD algorithm's proposers [29].For each voluntary muscle contraction, the multi-channel surface EMG signal segment, along with their IMFs derived from MEMD analysis, was further segmented into a series of analysis windows with a window length of 256 ms and a window increment of 128 ms.A larger window length would carry more stable neural control information facilitating improvement of classification accuracy.Meanwhile, however, it leads to a longer system delay.The window of 256 ms was selected considering the trade-off between the classification performance and the imposed delay [28].The following EMG feature extraction and pattern classification were then performed on these analysis windows.

Feature Extraction
For each analysis window, features regarding the couplings among multiple channels can be extracted from the corresponding scale-aligned IMFs derived from different channels.Suppose that the IMFs corresponding to an analysis window is denoted as IMF i Locations of electrodes for 27-channel bipolar EMG signal recordings.Each pair of the hollow and solid circles represents a bipolar EMG channel.
absolute value (MAV) and the waveform length (WL), are often used to account for the amplitude or oscillatory intensity of surface EMG time series, both the MAV and WL were therefore calculated for the i-th IMF from the channel m.These MAV or WL parameters, derived from the same-scale IMFs of different channels, were then normalized as a feature vector describing the amplitude distribution of a specific scale across channels.Multiple normalized feature vectors could be obtained by varying the IMF scale order i.
The effect of the scale combination for feature extraction on the myoelectric pattern recognition performance was also investigated.It was found that IMFs with scale order higher than 6 were very weak as compared to those with scale order not higher than 6, indicating that low-frequency signal components located at high-order IMFs did not carry much information (as shown in the following experimental results).In this case, all IMFs with scale order higher than 6 were summed up as a single low-frequency scale, denoted as 7-th scale for convenience.Consequently, four different scale combinations for feature extraction were examined in this study, including scales from 1 to 6 (denoted as {1~6}), scales from 2 to 6 (denoted as {2~6}), scales from 2 to 7 (denoted as {2~7}), and scales from 2 to 4 (denoted as {2~4}).These scale combinations were selected to assess the contribution of the IMFs with the first and highest scale orders to surface EMG classification, with a basic assumption that the majority of signal information was carried within the IMFs with scale order ranging from 2 to 6 (this was confirmed in the following experimental results).With each scale combination, the feature vectors normalized across channels, for both the MAV and WL derived from all possible scales, were concatenated to form a single feature vector for each analysis window.
In this study, the routine TD feature set including four statistics, namely the MAV, the number of zero crossings (ZC), the number of slope sign changes (SSC), and the WL [5,7,[11][12][13], was also used to classify the different movements of each stroke subject for comparison purposes.For each analysis window, all TD features extracted from 27 channels were concatenated to form a 108-dimensional feature vector.

Classification
User-specific classification was conducted to train and evaluate linear discriminant classifiers (LDCs) using the data collected from the same subject.The LDC used in this study worked through two steps.The first step represents linear discriminant analysis (LDA, or more strictly, multivariate extension of LDA), which is a method for finding a linear transformation/projection of features that maximized the ratio of between-class variability to within-class variability for the optimal separation of two or multiple classes (also known as Fisher linear transformation [30]).The second step involves a linear classifier using the maximum a-posteriori probability (MAP) rule and Bayesian principles [30], where the within-class density of each class is usually modeled as a multi-variant Gaussian distribution.With a well-trained classifier, the final decision of an unknown sample can thus be made based on the maximum likelihood (ML) estimation.The use of LDC is due to its simple structure, low computational complexity, and high efficiency [7,11].Note that for a C-class problem, the LDA is able to produce feature transformation to at most C-1 dimensions [30].Thus, with C = 18 in this study,

DISCUSSION
Due to recent advances in amplifier design and electrode manufacture techniques, multi-channel surface EMG consisting of a number of electrodes has been widely applied for myoelectric pattern recognition [11][12][13][14][15]. Effective modeling of muscle couplings by processing a number of EMG signals simultaneously (rather than in a manner of channel-by-channel processing) might facilitate decoding neural control information from multi-channel surface EMG recordings.MEMD is developed to overcome the mode-misalignment problem introduced by the standard EMD [19], and is able to simultaneously decompose data from multiple surface EMG channels into their scale-aligned components (e.g., IMFs), thereby enabling cross-information to be extracted as features describing couplings across channels [18].In this study, an MEMD-based feature extraction method was proposed to calculate amplitude distribution of the same-scale IMFs across multiple channels.These normalized amplitude distribution vectors can reveal the underlying relationship (the quantitative description of the multi-channel signal couplings) in the aligned intrinsic scales across multiple surface EMG channels, thus facilitating identification of the distinct coactivation and coordination patterns of multiple muscles.The proposed MEMD-based feature extraction method was applied to multi-channel EMG pattern recognition in this study.A dataset including 27-channel surface EMG signals acquired from 12 stroke subjects during their performance of 18 different movement patterns was used to assess the performance of the proposed method.The improved classification performance using the MEMD-based feature extraction compared with the conventional TD feature set demonstrates its advantages for multichannel EMG pattern recognition.
Among the examined four MEMD-based feature sets, {2~6} outperformed the others with statistical significance (p < 0.05).Compared with feature set {2~6}, either adding the first-order scale ({1~6}) or adding the last-order scale ({2~7}) for feature extraction leads to reduced classification performance (p = 0.026 and p = 0.049, respectively).This may be due to the fact that the first-order scale primarily carried high-frequency noise, while the last-order scale mainly contained low-frequency motion artifacts of the surface EMG signals (see Figure 2).These findings were consistent with a previous study [18].It follows that both these scales lack discriminative neural control information, and should be discarded from the scale combination during feature extraction to facilitate multi-channel myoelectric pattern recognition in stroke subjects.In this case, considering that MEMD is a recursive approach, it is unnecessary to produce IMFs with scale orders higher than 6, thus greatly saving computational power.In contrast, the average error rate of feature set {2~4} was significantly lower than that of {2~6} (p = 0.019), but was comparable with that of the TD feature set (p value close to 1).This also demonstrated that {2~6} was the most appropriate MEMD-based feature set in this study.
The MEMD algorithm was employed in this study to provide an adaptive framework to decompose multi-channel myoelectric signals into scale-aligned IMFs, thus offering practical tools for modeling muscle couplings to facilitate EMG pattern identification.Stroke patients were studied as a test case of pattern recognition of multi-channel 270 Myoelectric Pattern Identification of Stroke Survivors Using Multivariate Empirical Mode Decomposition surface EMG data.Stroke, however, is not a critical factor in the methodology.In fact, the proposed MEMD-based feature extraction method can be applied to EMG pattern recognition involving other population as well.In addition to the myoelectric control applications, the MEMD analysis of multi-channel surface EMG signals can also be applied to identify abnormal muscle activation patterns in patients with neurological disorders.One limitation of the proposed method is the increased computational burden imposed by the MEMD.Consequently, to process the same amount of data, the MEMD-based feature extraction may take several minutes, as opposed to several seconds required for conventional TD feature extraction.The difficulty in real time implementation of the MEMD limits its utility in the development of a myoelectric pattern recognition control system.In contrast, for neurophysiological investigations which do not require real time implementation (such as identification of abnormal muscle activation patterns in patients with neurological disorders), the computational complexity imposed by the MEMD is less critical for offline analysis.

CONCLUSION
This study presents a novel feature extraction method using multivariate empirical mode decomposition (MEMD) to simultaneously process multiple surface EMG channels.From scale-aligned intrinsic mode functions (IMFs) derived from the MEMD analysis, normalized amplitude distributions of the same-scale IMFs across different channels were calculated as features which serve to reveal the underlying relationship in the aligned intrinsic scales across multiple muscles.The proposed method was assessed for identification of 18 different functional movement patterns via 27-channel surface EMG signals recorded from the paretic forearm muscles of 12 subjects with hemiparetic stroke.The MEMD-based feature set achieved better performance than the conventional time domain (TD) feature set in classifying different movements.The current results suggest that the MEMD-based feature extraction of multi-channel surface EMG data can provide a promising approach to modeling of muscle couplings and identification of different myoelectric patterns.
Figure 1.Locations of electrodes for 27-channel bipolar EMG signal recordings.Each pair of the hollow and solid circles represents a bipolar EMG channel.

Figure 2 .
Figure 2. Illustration of selected two channels of surface EMG signals recorded during wrist extension performed by Subject 9, and two corresponding sets of scale-aligned IMFs after MEMD analysis.Two vertical dashed lines in each channel represent determined onset and offset of the muscle contraction, respectively; between them, the multiple channels of surface EMG signals are decomposed simultaneously.