New Methods for Analyzing Complex Biomedical Systems and Signals
View this Special IssueResearch Article  Open Access
Vector Autoregressive Hierarchical Hidden Markov Models for Extracting Finger Movements Using Multichannel Surface EMG Signals
Abstract
We present a novel computational technique intended for the robust and adaptable control of a multifunctional prosthetic hand using multichannel surface electromyography. The initial processing of the input data was oriented towards extracting relevant time domain features of the EMG signal. Following the feature calculation, a piecewise modeling of the multidimensional EMG feature dynamics using vector autoregressive models was performed. The next step included the implementation of hierarchical hidden semiMarkov models to capture transitions between piecewise segments of movements and between different movements. Lastly, inversion of the model using an approximate Bayesian inference scheme served as the classifier. The effectiveness of the novel algorithms was assessed versus methods commonly used for realtime classification of EMGs in a prosthesis control application. The obtained results show that using hidden semiMarkov models as the top layer, instead of the hidden Markov models, ranks top in all the relevant metrics among the tested combinations. The choice of the presented methodology for the control of prosthetic hand is also supported by the equal or lower computational complexity required, compared to other algorithms, which enables the implementation on lowpower microcontrollers, and the ability to adapt to user preferences of executing individual movements during activities of daily living.
1. Introduction
The methods for estimating the intention of an amputee to control the movement of a prosthetic hand often relies on the interpretation/decoding of the electrical activity of muscles (electromyograms, EMG) recorded on the skinsurface of the residual limb [1, 2]. Even though information regarding the muscle activity could be obtained in various ways, commercially available prosthetic hands commonly use only a few surface EMG channels. Furthermore, the classification algorithms implemented in such prostheses generally comprise only calculation of an amplitude based EMG feature (e.g., root mean square) that is thresholded to obtain binary control of a single hand function [3]. Even in the case of multifunction prostheses with a plurality of independent joints, the same control paradigm is still employed where switching between different grasps can depend on nonEMG inputs [4], such as smartphone interfaces, specific movements measured using inertial sensors, object recognition algorithms based on camera systems [5], or EMG inputs such as cocontractions (contractions of multiple muscles). Relatively low percentage of users is capable of fully exploiting the capabilities of multifunction prosthetic hands. This suggests that the control strategies need to be improved in order to achieve an intuitive link between intention and the resulting artificial hand function.
With the ongoing development of dexterous prosthetic hands [6–10], the gap between the potential dexterity offered by the new hands and the actual capabilities of humanmachine interfaces is even more evident. In contrast to the commercial devices, a variety of EMG recording techniques [11, 12], signal preprocessing [13], and classifiers have been proposed in the academia [14]. One of the approaches for increasing selectivity and sensitivity of the EMG signals includes surgical intervention in order to place recording electrodes on top of or inside targeted muscles [15]. This technique was shown to provide superior inputs for any type of EMG based controller; however, it is limited by need of surgical intervention and may not be widely accessible. As an alternative, EMG picked up at the surface of the skin provides a safe, noninvasive interface, at the cost of selectivity. Considering that the underlying mechanism of a classification algorithm is to localize the source of the electrical activity (active muscle), the methodology employed for similar problems is often based on multiple recording sites. Therefore, using a multichannel EMG approach provides an interesting perspective for improving the performance of surface EMG controllers. Increasing the number of the input signals also enables advanced data processing and classification techniques. Indeed, computational methods such as support vector regression [16], treestructured neural network [17], Bayesian inference [18], ICA clustering [19], hidden Markov models [20], nonnegative matrix factorization [21], and various pattern recognition approaches [22–25], demonstrated promising classification results while using multiple discrete EMG channels or highdensity surface EMG electrodes.
As one of the approaches for decoding finger movements using only multiple EMG signals, in this study we propose a novel classification algorithms based on the combination of EMG feature extraction and piecewise modeling of the feature temporal dynamics, incorporated in hierarchical hidden semiMarkov models (HHSMM) [26–28] and an online Bayesian classifier implemented through model inversion. We developed and assessed two variations of the aforementioned method and included special cases of these algorithms with reduced computational requirements. In order to assess our algorithm but also to provide an extensive evaluation of the impact of different EMG features on the classification accuracy, we compared its performance with some of the states of the art classifiers. This work is an extension of the study presented in [29].
2. Materials and Methods
2.1. Subjects
Five able bodied subjects participated in the experiment. All the subjects provided informed consent, and the study was approved by the Regional Ethical Review Board in Lund, Sweden.
2.2. Protocol
The methodology presented in this research relies on multielectrode EMG signals for movement classification. During the measurement procedure, a subject was comfortably seated with the right hand resting in a neutral position. The subject was asked to perform a movement to match a hand image shown on the screen in front of him/her. The requested movements were flexion/extension of all five fingers and adduction/abduction of the thumb (12 different movements). Movement and resting periods between movements were of equal length (5 seconds) and were timed by a LabVIEW (National Instruments, Austin, TX) custom application. The participant was visually prompted by the program, which synchronously acquired the EMG signal and the current cue annotation. Two different datasets (one for training and one for testing) each consisting of five repetitions of each movement, totaling 60 movements, and the rest states, were stored with the intended class information, for offline analysis.
2.3. EMG Recording
The dataset used in this paper was previously recorded and used in the publication by Huang et al. [30]. The EMG signals were acquired using a modified version of the amplification and acquisition system [31] and a LabVIEW program. The system acquired 16 channels of EMG, sampled at 1.6 kHz per channel and with a bandpass filter between 0.5 Hz and 800 Hz with 16bit resolution and a gain of 56 dB per channel. The Red Dot Ag/AgCl (3M Health Care, Germany) electrodes that were used in the study were placed on the forearm of the participants as shown in Figure 1. For the 16 channel recordings, twelve electrodes were placed on the superficial flexor muscles on the volar side of the forearm and four electrodes were placed on the superficial extensor muscles on the dorsal side of the forearm. Electrodes were placed in order to cover as many independent muscles as possible. A sample of recorded EMG signals is provided in Figure 2.
(a)
(b)
2.4. Feature Extraction
A set of commonly used EMG features in the time domain was selected as an additional comparison criterion in this study. The choice of time domain features, versus spectral or timefrequency domain features, is in line with the tendency of deriving a control chain that could be implemented in an embedded system with limited processing power/speed. Among the relatively high number of reported features, we chose threeamplituderelated and threefrequencyrelated timebased EMG features. We implemented functions for the following EMG features as reported in [13, 32, 33].
2.4.1. Mean Absolute Value
Mean Absolute Value (MAV) is the favored EMG feature in many myoelectric control applications. It is calculated within a moving average window of the rectified EMG signal:where denotes window length and the th EMG sample within current window.
2.4.2. Root Mean Square
Root Mean Square (RMS) calculation is, similar to MAV, also windowed function:In the case of the signal with mean value close to zero, such as surface EMG (sEMG) recorded with relatively large transcutaneous electrodes, the RMS becomes standard deviation of the input time series.
2.4.3. Variance
Variance (VAR) in the case of sEMG signal could be simplified enabling fast implementation in an embedded system. The resulting formula is the following:
2.4.4. Slope Sign Change
Slope Sign Change (SSC) is a time domain method used to estimate frequency feature if the EMG signal. The calculation of the SSC relies on counting changes between positive and negative slope among three consecutive samples. To limit SSC calculation only to periods with EMG activity, the threshold function is imposed in the feature extraction method:
2.4.5. Zero Crossing
Zero Crossing (ZC) is the function that counts the number of consecutive EMG samples that change sign within the sliding window. Similar to SSC calculation, the threshold function is imposed to remove calculation of the ZC during periods without pronounced EMG activity:
2.4.6. Willison Amplitude
Willison amplitude (WA) is a measure related to superimposed action potentials that make the EMG signal. The WA is the number of consecutive differences between consecutive samples that exceed set threshold:To produce comparable input conditions, the sliding window for all feature calculations was the same (250 ms width and 25 ms overlapping).
In the case of SSC, ZC, and WA calculations, thresholds were set based on the estimated white noise level during recording. As each recording session started with a rest period, it was convenient to use first 0.5 s to calculate thresholds that are used throughout the same recording. A sample of calculated features is shown in Figure 3.
All feature extraction calculations were conducted using MATLAB, with moving window size of 250 ms and a displacement of 25 ms.
2.5. Classification Algorithm
Here, we will extend our previous work [29] in which we derived Bayesian online classifier using vector autoregressive hierarchical hidden Markov models (VARHHMM), with a classifier based on vector autoregressive hierarchical hidden semiMarkov models (VARHHSMM) [26–28]. The flowchart for the classifier’s components is shown in Figure 4. The basic presumption of this approach is that the feature dynamics in time domain during a movement could be estimated by a number of VAR(1) segments defined via a unique set of parameters. A movement denotes a finger movement represented in the EMG feature space.
From the basic flowchart, it can be noted that the main difference between the two algorithms is in the top layer of the hierarchy. In the case of the VARHHMM, this layer only takes signal likelihood as the input to for the layer that identifies the current movement. In the case of the VARHHSMM, the current movement label will also depend on the expected duration of individual movements when identifying the active moment label. In what follows, we will refer to these two classifiers as hierarchical hidden Markov model (HMM) and hierarchical hidden semiMarkov model (HSMM), where we dropped the VAR and hierarchical labels for readability.
We assumed here that the selected time domain EMG features could be represented by a number of segments. To accommodate signal dynamics in realtime and enable creating a classifier suitable for embedded implementation, we defined individual segments as vector autoregressive (VAR) models of the firstorder VAR(1). Using the VAR process for approximating the dynamics of a segment of an EMG time series allows for efficient prediction of feature dynamics and implicitly accounts for the noise presented in the recorded EMG signals. Formally, the VAR(1) process is defined aswhere is the latest acquired signal feature in vector form, is the mean value of the signal, denotes the segment within the th finger movement, is the transition matrix of the VAR(1) model, and is the random normal variable with zero mean. Note that the dynamics within each segment and movement are completely defined by the following set of parameters , where denotes covariance of .
Based on predicted and measured sample vectors we can calculate which of the possible (fitted) VAR(1) models provides the best description of the recorded time series. The transition between individual segments (between different VAR(1) models) is considered as stochastic; thus, in our realization it is modeled as a HMM (as in [26]). The transitions between segments of a single movement are constrained by additional hidden states which define the currently active movement. For this reason, we introduced another HMM chain that acts as the second level of the hierarchy. The transitions between possible VAR(1) models of a single movement and different movements are defined through the transition matrices for each HMM chain. The transition matrices on both levels of the hierarchy are estimated during the training procedure that relies on Viterbi expectation maximization algorithm [34].
By applying an approximate Bayesian inference procedure to the hierarchically arranged HMM and the emission distributions defined with VAR(1) models, we can classify each of the measured data points in a time series. The Bayesian inference for realtime signal classification relies on estimating the posterior probability for each movement. The posterior probability is calculated by combining expectations over movement probabilities with the evidence provided by the emission distributions (observation likelihoods). The expectations (predictions) are estimated using the posterior distribution at previous time sample and the transition matrices of movements and associated segments. Using the values of the posterior probability at the latest acquired signal and the current movement , we can classify the current sample aswhere is the classification label obtained by our algorithm.
The HSMM algorithm is similar to the HMM algorithm with the only difference being that each movement includes a sequence of durations associated with it. This approach increases the number of free parameters with an additional vector comprising prior duration probabilities of individual movements. In other words, the probability of switching between movements is also influenced by the model of movement duration. An illustration of a specific realization of the movement switching sequence of the HSMM is shown in Figure 5.
The implementation of the HSHMM algorithm relies on a timer, or counter (Cnt) in the case of equidistant sampling. Upon entering a state/movement, the counter is set to predefined value related to that particular state duration and with every new sample is decremented. When the counter reaches 0, the state/movement is free to change. Upon switching to a new state, the counter is reset to zero. This additional condition is defined aswhere denotes movement at time sample and denotes transition matrix between movements.
In the case of HSMM, the duration of a movement is provided in form of a prior distribution. As the measurement protocol was executed in automated manner, all performed movements are of the same duration (~5 s) with some variability related to subject’s reaction times to a visual cue. In order to simulate more realistic enduser scenario, we defined prior distribution to have a duration range with 2second variability ( s). As the simplest scenario derived a priori, the imposed distribution was defined as uniform within the whole range. To estimate influence of the range of the imposed distribution, the accuracy score was also calculated for ranges from s to s.
As the number of segments per movement is the free parameter of both algorithms, it was of specific interest to evaluate performance of the classification with respect to a movement division by the automated optimization procedure. Specifically, the case of a single segment per movement was analyzed independently. The rationale for this lies in the extracted features dynamics in the time domain. When visualizing the recorded signals, a movement in a feature space reassembles a step response of a dynamic system, with a relatively short transient compared to the constant part. When ignoring the noise part of the vector autoregressive model (see above) we can rewrite the update equation aswhere the current sample is a function of previous sample and the distance of the previous value from the terminal value in the limit. The constant matrix defines the rate of change along each dimension of the recorded signal. This form of the basic VAR representation is convenient for illustrating the signal trajectory after the onset of a movement (see Figure 6). The advantages of using single segments per movements area significantly shorter optimization time (explained in the following section) and a simplified realtime processing implementation as the lower HMM layer is removed.
2.6. Optimization
The estimation of the free parameters of all models presents the only computationally demanding procedure of the presented methodology. During the optimization, the free parameters of the VAR models (number of segments per movement, mean signal values of individual segments , covariance matrix , transition matrix , and between segment’s transition matrix which defines the HMM structure) are evaluated using a Viterbi algorithm coupled with an expectation maximization algorithm (EM) [35]. The flowchart of the optimization procedure is shown in Figure 7.
The optimization procedure is an iterative process where from the initial values of VAR parameters the Viterbi algorithm determines the most likely sequence of segments and estimates transition matrices. The following step includes reestimation of VAR parameters based on a maximum likelihood calculation and the most likely separation of training data suggested by the Viterbi algorithm. The procedure shown in Figure 7 is iterated with the loglikelihood set as the convergence criterion. When the convergence criterion reaches a value of , the iterations are stopped. The whole optimization method was based on the implementation of the Viterbi EM algorithm provided in the pyhsmmautoregressive Python library [36].
Bearing in mind that the optimization method could converge towards a local optimum, we repeated the whole optimization procedure for 1000 different initial values of free parameters and only the solution with the lowest value of Bayesian Information Criterion (BIC) [37] was kept as final. This iterative extension of the Viterbi EM algorithm helped us to improve the convergence of the parameter estimation to optimal parameters values. The Python 3.6 codes for the parameters optimization and classification simulation are available upon request sent to the corresponding author.
2.7. Classifiers Used for Benchmarking
To evaluate the effectiveness of our algorithm, we used some of the most commonly used classification methods including(i)Linear Discriminant Analysis (LDA),(ii)Quadratic Discriminant Analysis (QDA),(iii)nearest neighbors (knn), with ,(iv)Support Vector Machine with the ErrorCorrecting Output Codes (SVM ECOC),(v)Linear Classifier with the ErrorCorrecting Output Codes (LC ECOC),(vi)Linear Discriminant Analysis with the ErrorCorrecting Output Codes (LDA ECOC),(vii)Naive Bayes (NB),(viii)Random Forest (RF),(ix)Decision Tree (DT).
The implementation of selected classifiers was done using function in the MATLAB 2016b (The MathWorks Inc., Natick, MA). Similar to the VARHHMM algorithms, the classifiers were trained on only the training data set and tested on the test set and no postprocessing of classannotations was performed.
3. Results
The ground truth for classification in this paper is set to be the visual cue presented to the subject. Although there is the noticeable latency between the visual cue and the muscle activity in both movement onset and cessation as shown in Figure 2, we decided to use this as a reference to provide the same basis for all of the classification methods. This was done to force classifiers to adapt to the transient processes between consecutive states. As expected, this approach results in lower overall classification accuracy but provides interesting insights into classifiers performances during dynamic, continuous changes of input parameters.
The main metric for comparing features and performance of the classifiers was accuracy. For the accuracy calculation, each feature sample (one observation every 25 ms) was treated individually. Moreover, as the rest state is directly involved in active motor driving by prosthesis control, it was also considered equal to the other movement classes. With this evaluation paradigm, each classification point was categorized as a true positive (TP) if the predicted class matches the visual cue label, or as a false positive (FP) if there is mismatch between predicted class and the visual cue. This approach generates a large number of data points, which was in our case approximately 20000 per dataset that could be used to investigate the performances of various features and classifiers.
The cumulative results for all subjects are presented in Table 1. As the accuracy results for the subjects do not fit into normal distribution, the accuracy values presented in the table are median values among all subjects. The differences between classifier and features are emphasized by ranking top three features and classifiers overall but also top three results for each EMG feature. From Table 1, it could be noted that HSHMM outperforms other classifiers for the majority of tested feature extraction algorithms. It is also interesting that the HSMM with only one segment per movement closely follows performance of the HSMM with multiple segments per movement. Among feature extraction algorithms, the MAV features resulted in the best performance across different classifiers and subjects; RMS and WA lead to slightly lower results compared to MAV, while using VAR, SSC and ZC failed to produce similar results. Among classifiers, only NB underperformed, while other classifiers produced mutually comparable results.

As the average accuracy across all subjects and all classifiers was the greatest in the case of MAV feature, the results obtained using this feature were analyzed in more detail. Table 2 shows accuracy results for each classifier and each subject when the MAV feature was used. The results show that all VARHHMM variations are topranked, with the HSMM algorithm having the highest median accuracy score. Although at the top ranks, the difference of median values is not prominent.

With the focus on VARHHMM methodology, to illustrate ability to separate individual finger movements and resting periods, the confusion matrix is shown in Figure 8. For this example, we selected the best performing feature method among all classifiers (MAV) and the subject with the median results (Subject 4). This was done to analyze in more detail the case close to the enduser scenario. The common indicator of the classification effectiveness is the confusion matrix which is used to reveal weak points of a classifier. Figure 8 is showing the confusion matrix for the HSMM algorithm. It could be noted that majority of movements (2, 4, 5, 6, 7, 11, 12, and 13) are classified with just sporadic misclassifications. For these movements, most of the errors are the consequence of the latency between the visual cue (used as the golden standard in this paper) and the actual muscle activity of the targeted movement. This effect is responsible for the movementtorest and resttomovement misclassifications. The rest class, except the errors related to delayed transitions, is generally well classified. The movements (3, 8, 9, and 10) have more samples that are labeled incorrectly as some other movement. These errors, if consistent, could significantly impede control of a prosthesis. Thus, to examine the behavior of the VARHHMM implementations in a realtime application, the classification outputs were further analyzed.
Figure 9 contains some of the most severe examples of the misclassification errors produced by the HSMM algorithm. Figure 9(c) is an illustration of the proper movement classification where only errors occur during transition between rest and the movement. Compared to this sample, Figure 9(a) depicts an error of the labeling and keeping the incorrect movement active almost from the beginning of the movement onset, Figure 9(b) depicts misclassification at the middle of the movement period, and Figure 9(d) depicts the indecisiveness of the algorithm to label samples at the transient period. Out of these examples, the most severe in terms of realtime application is the first case that could lead to the execution of the false actuation. The other cases are relatively easy to correct using a state machine and/or majority voting algorithm that prevents situations that are less likely to occur by the user’s intent.
Other metrics that were implemented to calculate classifier behavior in the realtime are motion selection (MS) and motion completion (MC) times. The calculation of MS and MC is derived from the papers of Li et al. and OrtizCatalan et al. [38, 39]. MS was computed as the delay between cue onset and the first correctly classified sample. Due to the difference in features frequency, the condition for the MC was slightly modified compared to other papers. In our algorithm, a new EMG feature was produced every 25 ms, while in [38] every 100 ms, so instead of 10 correct classifications, we waited for 40 correct classifications until the MC condition was fulfilled. This change was done in order to have comparable results with the relevant publications. The same metric for MS and MC was employed for all subjects and all features. Figure 10 is showing MS and MC times mean values across features and classifiers. It could be noted that HSMM has the lowest delay between cue and the first correctly classified sample, indicating that the method has a relatively low detection threshold compared to the other classifiers. It is also clear that the HSMM works the best with amplitude based features (MAV, RMS, and WA). Another interesting outcome is that HSMM with only one segment per movement achieves similar results as HSMM with multiple segments per movement when fed with MAV or RMS features of the EMG signal. On the other hand, the regular HHMM performed slightly worse on average than the best classifiers, such as LDA, SVM ECOC, L ECOC, and LDA ECOC. Similar to HSMM, reducing number of segments per movement for the HHMM does not result in significant drop in performance related to MS and MC times.
We also analyzed another parameter, namely, the optimal number of segments per movement that results from the optimization procedure. For both HHMM and HSMM algorithms, two segments per movement most frequently produced the highest accuracy score (Figure 11). Another interesting parameter that was optimized and analyzed is the order of the vector autoregressive process () in the HHMM and HSMM models with one segment per movement. The highest accuracy was achieved with 3 successive points, but what was not expected was relatively good accuracy with only one autoregressive sample (prediction relies only on the previous observation).
As the performance of HSMM algorithm depends on the predefined distribution of possible movement durations, analyzing this parameter was of great importance. During the measurements, the duration of movements was governed by the visual cues that were presented in automated manner. This way, the resulting movements have roughly the same duration of 5 s. To artificially introduce expected movement duration variability, we expanded the uniform duration distribution around central point of 5 s. With the extension of the expected movement duration, we evaluated classification accuracy. As presented in Figure 12, the dependency between classification accuracy and the range of permitted movement duration is almost linear in the midrange. This finding is valuable for defining expected movement duration, while taking into account the tradeoff between accuracy and flexibility of the classifier to adapt to free movements by a user.
4. Discussion
In this paper, we presented novel algorithms for classifying features from surface EMG signals. Both of the proposed algorithms (HHMM and HSMM) are variations of the VARHHMM algorithm which combines vector autoregression, hidden Markov models, and Bayesian inference. The main focus of the presented results is the comparison between different EMG classifiers, including some of the most commonly used ones and VARHHMM variants. The results presented in Table 1 indicate that the HSMM algorithm with a priori fixed movement duration outperforms other classifiers. It should be also noted that the results of HSMM algorithm presented in this table are with the ad hoc uniform movement durations distribution between 3 s and 5 s. It is envisioned that, in the prolonged prosthesis use with this control algorithm, the movement duration distribution would be the parameter updated with each correctly executed movement. Thus, as the time of use increases, the distribution should become optimal and result in increased movement classification accuracy. The other derived algorithm (HHMM) was ranked overall just behind the SVM ECOC, mostly because of low accuracy when frequency related timebased features were used as inputs.
The results also reveal that reducing number of segments per movement to one does not result in considerable drop in accuracy scores for the derived algorithms (with the exception of VAR feature). This fact has significant impact on the possible implementation of the algorithms as the optimization, even in the case of large time series, for example, tens of minutes, could be executed in a matter of seconds. As the illustration, optimization of algorithm variants with free parameter of segments per movement takes up to 10 minutes per movement (Python 3.6 and Intel i7 processor). Depending on the embedded platform, our estimate is that optimization with one segment per movement could be still done in less than a minute. The other fact that contributes to low processing during optimization and realtime application is that algorithm does not need high order autoregressive models. As presented in Figure 11, using even only one previous sample point could result in the highest accuracy score.
As this study was also aiming at finding out which feature extraction method works the best when coupled with the proposed algorithms and commonly used ones, we carried out a systematic test that included all features and subjects. With the EMG signals that we acquired, using MAV feature as the first step in the classification chain produced the highest accuracy scores among classifiers. The performance of MAV was closely matched by RMS, while WA produced the best results among the frequency related timebased features.
When compared with the reported results of Li et al. [38], the obtained MS and MC times calculated with our data and algorithms show some deviations. This is mostly present in MS values that represent the delay of the first properly classified sample. In the case of Li et al., MS values were in range 0.16–0.33 s, while in our simulations, only HSMM with MAV, RMS, and WA achieved similar performance. This behavior is the result of the difference in the ground truth choice: Li et al. used actual onset of the EMG activity, while we selected visual cue onset. With this in mind, it is expected that the MS times in our case have the superimposed lag of approximately 200–400 ms resulting from the latency between the cue and the muscle activation. That practically means that the HSMM would outperform reported classifiers if used in similar conditions. The same rationale goes for the MC times. Although classifiers tested in this study produced longer MC times on average, if reduced by the visual cue to muscle activity offset, the results would show that HSMM manages to correctly detect 40 samples (which correspond to 10 samples in paper of Li et al.) faster than the classifiers reported by Li et al.
5. Conclusion
The study performed on our set of multichannel surface EMG indicates that using MAV feature coupled with HSMM algorithm leads to movement decoding accuracy higher than other combinations of features and classifiers. This combination also guaranties the shortest MS and MC times, influencing a response of a prosthetic hand to a user intent.
The main advantages inherent to our algorithms compared to the existing methods are the following. (1) Low computation complexity for the execution of the algorithm: following the optimization procedure, which is computationally demanding but executed only once per model, the implementation of our algorithm comes down to the basic algebraic operations that could be implemented even on lowend microcontrollers. This feature also permits true realtime execution with the delay only related to EMG feature extraction and AR depth. Additionally, we tested even faster variations of developed algorithms that have only one segment per movement. With the small drop in performances, these algorithms, especially HSMM with one segment per movement, significantly decreased optimization time and further decreased execution load. (2) Easy expansion: once optimized, models could be stored, and in the case of introducing a new movement class, the optimization could be performed only on the newly added movement class. (3) Noise resilient: in comparison with some weak classifiers and threshold centered rulebased algorithms, the VARHHMM approach implicitly takes into account possible sources of stochastic noise.
The method proposed in this paper is intended for decoding individual finger movements for directly controlling actuated fingers of a prosthetic hand. Although desirable, this approach is not common as the prosthesis use during activities of daily living mostly consists of synergistic movements (grasps). Thus, as the alternative control paradigm, the proposed method could be implemented as the state machine classifiercontroller. In this approach, the separable classes might be used to initiate fixed grasps (pinch/power/lateral/open) as this is what most commercial hands support.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This research is supported by the EUFunded DeTOP Poject (EITICT242015, GA no. 687905) and the Swedish Research Council (6372013444).
References
 M. Zecca, S. Micera, M. C. Carrozza, and P. Dario, “Control of multifunctional prosthetic hands by processing the electromyographic signal,” Critical Reviews in Biomedical Engineering, vol. 30, no. 4–6, pp. 459–485, 2002. View at: Publisher Site  Google Scholar
 N. Jiang, S. Dosen, K.R. Muller, and D. Farina, “Myoelectric control of artificial limbs—is there a need to change focus? [In the Spotlight],” IEEE Signal Processing Magazine, vol. 29, no. 5, pp. 148–152, 2012. View at: Publisher Site  Google Scholar
 A. Fougner, O. Stavdahl, P. J. Kyberd, Y. G. Losier, and P. A. Parker, “Control of upper limb prostheses: Terminology and proportional myoelectric controla review,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 20, no. 5, pp. 663–677, 2012. View at: Publisher Site  Google Scholar
 C. Connolly, “Prosthetic hands from touch bionics,” Industrial Robot, vol. 35, no. 4, pp. 290–293, 2008. View at: Publisher Site  Google Scholar
 S. Došen, C. Cipriani, M. Kostić, M. Controzzi, M. C. Carrozza, and D. B. Popovič, “Cognitive vision system for control of dexterous prosthetic hands: experimental evaluation,” Journal of NeuroEngineering and Rehabilitation, vol. 7, no. 1, article 42, 2010. View at: Publisher Site  Google Scholar
 M. Controzzi, F. Clemente, D. Barone, A. Ghionzoli, and C. Cipriani, “The SSSAMyHand: A dexterous lightweight myoelectric hand prosthesis,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 25, no. 5, pp. 459–468, 2017. View at: Publisher Site  Google Scholar
 Touchbionics, http://www.touchbionics.com/.
 bebionic, http://bebionic.com/.
 Michelangelo, http://www.ottobockus.com/prosthetics/upperlimbprosthetics/solutionoverview/michelangeloprosthetichand/.
 Vincent hand, http://vincentsystems.de/en/.
 M. OrtizCatalan, R. Brånemark, B. Håkansson, and J. Delbeke, “On the viability of implantable electrodes for the natural control of artificial limbs: review and discussion,” Biomedical Engineering Online, vol. 11, article 33, 2012. View at: Publisher Site  Google Scholar
 S. Muceli and D. Farina, “Simultaneous and proportional estimation of hand kinematics from EMG during mirrored movements at multiple degreesoffreedom,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 20, no. 3, pp. 371–378, 2012. View at: Publisher Site  Google Scholar
 A. Phinyomark, C. Limsakul, and P. Phukpattaranont, “A novel feature extraction for robust EMG pattern recognition,” Journal of Computing, vol. 1, no. 1, pp. 71–80, 2009. View at: Google Scholar
 M. A. Oskoei and H. Hu, “Myoelectric control systems—a survey,” Biomedical Signal Processing and Control, vol. 2, no. 4, pp. 275–294, 2007. View at: Publisher Site  Google Scholar
 M. OrtizCatalan, B. Håkansson, and R. Brånemark, “An osseointegrated humanmachine gateway for longterm sensory feedback and motor control of artificial limbs,” Science Translational Medicine, vol. 6, no. 257, Article ID 257re6, 2014. View at: Publisher Site  Google Scholar
 A. Ameri, E. N. Kamavuako, E. J. Scheme, K. B. Englehart, and P. A. Parker, “Support vector regression for improved realtime, simultaneous myoelectric control,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 22, no. 6, pp. 1198–1209, 2014. View at: Publisher Site  Google Scholar
 F. Sebelius, L. Eriksson, C. Balkenius, and T. Laurell, “Myoelectric control of a computer animated hand: a new concept based on the combined use of a treestructured artificial neural network and a data glove,” Journal of Medical Engineering & Technology, vol. 30, no. 1, pp. 2–10, 2006. View at: Publisher Site  Google Scholar
 S. Bitzer and P. Van Der Smagt, “Learning EMG control of a robotic hand: towards active prostheses,” in Proceedings of the IEEE International Conference on Robotics and Automation (ICRA '06), pp. 2819–2823, IEEE, May 2006. View at: Publisher Site  Google Scholar
 G. R. Naik, A. H. AlTimemy, and H. T. Nguyen, “Transradial amputee gesture classification using an optimal number of sEMG sensors: an approach using ICA clustering,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 24, no. 8, pp. 837–846, 2016. View at: Publisher Site  Google Scholar
 A. D. C. Chan and K. B. Englehart, “Continuous myoelectric control for powered prostheses using hidden Markov models,” IEEE Transactions on Biomedical Engineering, vol. 52, no. 1, pp. 121–124, 2005. View at: Publisher Site  Google Scholar
 C. Huang, X. Chen, S. Cao, B. Qiu, and X. Zhang, “An isometric muscle force estimation framework based on a highdensity surface EMG array and an NMF algorithm,” Journal of Neural Engineering, vol. 14, no. 4, Article ID 046005, 2017. View at: Publisher Site  Google Scholar
 R. N. Khushaba, A. H. AlTimemy, A. AlAni, and A. AlJumaily, “A framework of temporalspatial descriptorsbased feature extraction for improved myoelectric pattern recognition,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 25, no. 10, pp. 1821–1831, 2017. View at: Publisher Site  Google Scholar
 A. A. Adewuyi, L. J. Hargrove, and T. A. Kuiken, “An analysis of intrinsic and extrinsic hand muscle EMG for improved pattern recognition control,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 24, no. 4, pp. 485–494, 2016. View at: Publisher Site  Google Scholar
 N. Celadon, S. Došen, I. Binder, P. Ariano, and D. Farina, “Proportional estimation of finger movements from highdensity surface electromyography,” Journal of NeuroEngineering and Rehabilitation, vol. 13, no. 1, article 73, 2016. View at: Publisher Site  Google Scholar
 C. Cipriani, C. Antfolk, M. Controzzi et al., “Online myoelectric control of a dexterous hand prosthesis by transradial amputees,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 19, no. 3, pp. 260–270, 2011. View at: Publisher Site  Google Scholar
 M. Yang, “Some properties of vector autoregressive processes with Markovswitching coefficients,” Econometric Theory, vol. 16, no. 1, pp. 23–43, 2000. View at: Publisher Site  Google Scholar  MathSciNet
 M. Dong and D. He, “Hidden semiMarkov modelbased methodology for multisensor equipment health diagnosis and prognosis,” European Journal of Operational Research, vol. 178, no. 3, pp. 858–878, 2007. View at: Publisher Site  Google Scholar
 S.Z. Yu, “Hidden semiMarkov models,” Artificial Intelligence, vol. 174, no. 2, pp. 215–243, 2010. View at: Publisher Site  Google Scholar  MathSciNet
 N. Malešević, D. Marković, G. Kanitz, M. Controzzi, C. Cipriani, and C. Antfolk, “Decoding of individual finger movements from surface EMG signals using vector autoregressive hierarchical hidden Markov models (VARHHMM),” in Proceedings of the International Conference on Rehabilitation Robotics (ICORR '17), pp. 1518–1523, IEEE, July 2017. View at: Publisher Site  Google Scholar
 H. Huang, T. Li, C. Bruschini et al., “EMG pattern recognition using decomposition techniques for constructing multiclass classifiers,” in Proceedings of the 6th IEEE RAS/EMBS International Conference on Biomedical Robotics and Biomechatronics (BioRob '16), pp. 1296–1301, IEEE, June 2016. View at: Publisher Site  Google Scholar
 C. Antfolk, C. Cipriani, M. Controzzi et al., “Using EMG for realtime prediction of joint angles to control a prosthetic hand equipped with a sensory feedback system,” Journal of Medical and Biological Engineering, vol. 30, no. 6, pp. 399–406, 2010. View at: Publisher Site  Google Scholar
 N. Alves, E. Sejdić, B. Sahota, and T. Chau, “The effect of accelerometer location on the classification of singlesite forearm mechanomyograms,” Biomedical Engineering Online, vol. 9, article 23, 2010. View at: Publisher Site  Google Scholar
 B. Hudgins, P. Parker, and R. N. Scott, “A new strategy for multifunction myoelectric control,” IEEE Transactions on Biomedical Engineering, vol. 40, no. 1, pp. 82–94, 1993. View at: Publisher Site  Google Scholar
 G. D. Forney, “The viterbi algorithm,” Proceedings of the IEEE, vol. 61, no. 3, pp. 268–278, 1973. View at: Publisher Site  Google Scholar
 A. Logothetis and V. Krishnamurthy, “Expectation maximization algorithms for MAP estimation of jump Markov linear systems,” IEEE Transactions on Signal Processing, vol. 47, no. 8, pp. 2139–2156, 1999. View at: Publisher Site  Google Scholar
 M. Johnson, pyhsmmautoregresive, https://github.com/mattjj/pyhsmmautoregressive.
 B. North and A. Blake, “Learning dynamical models using expectationmaximisation,” in Proceedings of the IEEE 6th International Conference on Computer Vision, pp. 384–389, IEEE, 1998. View at: Publisher Site  Google Scholar
 G. Li, A. E. Schultz, and T. A. Kuiken, “Quantifying pattern recognition—based myoelectric control of multifunctional transradial prostheses,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 18, no. 2, pp. 185–192, 2010. View at: Publisher Site  Google Scholar
 M. OrtizCatalan, B. Håkansson, and R. Brånemark, “Realtime and simultaneous control of artificial limbs based on pattern recognition algorithms,” IEEE Transactions on Neural Systems and Rehabilitation Engineering, vol. 22, no. 4, pp. 756–764, 2014. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2018 Nebojša Malešević et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.