A muscle synergies-based movements detection approach for recognition of the wrist movements

Abstract

Myoelectric signals are regarded as the control signal for prosthetic limbs. But, the main research challenge is reliable and repeatable movement detection using electromyography. In this study, the analysis of the muscle synergy pattern has been considered as a key idea to cope with this main challenge. The main objective of this research was to provide an analytical tool to recognize six wrist movements through electromyography (EMG) based on analysis of the muscle synergy patterns. In order to design such a system‚ the synergy patterns of the wrist muscles have been extracted and utilized to identify wrist movements. Also, different decision fusion algorithms were used to increase the reliability of the synergy pattern classification. The classification performance was evaluated while no data subject was enrolled. In terms of the achieved performance, using a multi-layer perceptron (MLP) neural network as the fusion algorithm turned out to be the best combination. The classification average accuracy, obtained in an offline manner, was about 99.78 ± 0.45%. While the classification average cross-validation accuracy, obtained in an offline manner, using Bayesian fusion, and Bayesian fuzzy clustering (BFC) fusion algorithm were 99.33 ± 0.80% and 96.43 ± 1.08%, respectively.

1 Introduction

The electromyogram (EMG) signal represents the electrical potentials generated in the muscles during muscle contraction, which shows the important neuromuscular information [1]. The EMG signal recorded from each surface electrode is the total potential of the motor units in the region where the electrode is positioned [2]. Due to the useful application of EMG signal in clinical diagnoses, and biomedical applications as well as rehabilitation, they are considered as one of the best resources of controlling (i.e., prostheses, robots, and human-computer interfaces), recognition of intended limb movements [1, 3]. Extensive research has been done in order to control various functions and increase the efficiency of prostheses (i.e., several degrees-of-freedom (DoF). One class of the known control approaches are based on recognizing the pattern of EMG signals elicited from the residual healthy muscles [3]. Many researchers are engaged in improving the performance of such recognition algorithms in order to improve the efficiency of prostheses. Hence, the estimation algorithms, pattern recognition and regression methods [4, 5], and their combination [6, 7] are among the topics of interest to researchers.

Pattern recognition algorithms are applied to classify the EMG activity patterns in multiple muscles. So, in these researches, the first step is windowing the EMG signal [8, 9] and extracting a set of important features from time-windowed signal in the time domain (TD) [10, 11] and in the frequency domain (FD) [10,11,12,13]. The computational cost of TD features is less than the computational of FD features, yet yield comparable classification accuracy [10].. The classification accuracy in the classification model depends on the type of classification algorithm (i.e., supervised or unsupervised) and the features which are selected [12], [14, 15]. Also, it depends on the type of subject (able-bodied or amputee) [3]. Sang et al. [10] recognized six predetermined tasks of muscle activation patterns using TD features and linear discriminant analysis classifiers in stroke subjects. Similarly, Francesco et al. [16] investigated the recognition of hand and finger movements individually and reported the different accuracies for able-bodied subjects and amputated ones. They used the feed-forward multi-layer perceptron (MLP) and four TD features [16]. In the other study, Chu et al. [13], nine types of wrist movements were detected based on packet wavelet analysis using MLP and self-organizing feature map (SOFM) neural networks with an accuracy of 97%. The research results of the past decades have shown that the pattern recognition-based methodology has promising results [12]. But some impediments such as number of the selected features, the computational cost of the required, and the minimum number of the data required for training the classifier can limit their practical applications. The other challenge is detecting the movements with multi DoF [17, 18], and some research works were focused on detecting the movement with only one DoF [3, 6, 8].

The authors believe that most of these challenges are raised because of specificity movement detection approaches used in the aforementioned researches. In fact, according to a specificity approach, activation dynamics of each involved muscle is analyzed distinctly. In this manner, the muscle coactivation dynamics, characterizing the movement pattern, cannot be revealed. While using the synergy-based analysis, the muscle interaction dynamics can be studied. The synergy-based analysis is a type of movement analysis which addresses the interactive behavior among the involved sub-systems (muscles, joint), instead analyzing the behavior of each sub-system (muscles, joints) separately. Since the real movement emerged from the interaction among the muscles, it is expected that the synergy-based analyses can better reveal and discriminate the different motion dynamics in comparison with the non-synergetic analyses. Therefore, the analyses approaches designed based on muscle synergy patterns can be preferable to the analyses approaches based on the electromyography of the individual muscle distinctly. But most of the related published works have been designed based on electromyography of the individual muscle distinctly [13,14,15,16,17]. Hence, researchers believe that controlling and detecting neural information based on synergistic patterns [3, 8, 18,19,20,21] is preferable than other methods [2, 4,5,6]. Recently, studies have been done on the synergy pattern [3, 8, 18,19,20,21]. For example, Jiang et al. [3] extracted muscle synergy patterns using non-negative matrix analysis in order to two DoF control of wrist. Similarly, Jiaxin et al. [21] used muscle synergy matrixes to control wrist movements (open, close, pronate, and supinate). The main goal of this paper is presenting a new muscle synergy-based wrist movement recognition algorithm. In this study, the extracted muscle synergy patterns were used as the input data for SOFM classifiers instead the TD or FD features. For improving the movement recognition performance, six different classifiers were applied, and different fusion algorithms were used to combine the classification results.

It is worth noting that muscle activation onset can give rise to the presence of transient dynamics which may degrade the performance of movement classifiers. Therefore, in this study, a specific algorithm has been utilized for muscle activation onset detection.

2 Materials and methods

2.1 Experiments and data collection

In this study, eight able-bodied subjects (5 males, 3 females: 26–40 years; all right-handed) without any history of a disease or other skeletal disorders participated in the experiment. Everyone read and signed the consent form before any data was collected and the local ethics committees approved this study. During the experiments, the subject seated comfortably and held a hand horizontally forward and in pronation. For all individuals, six pairs of Ag–AgCl surface bipolar electrodes were placed on the extensor digitorum, flexor pollicis longus, extensor carpi ulnaris, abductor pollicis longus, pronator teres, and supinator forearm muscles. The reference electrode is placed on the ulnar styloid bone. In these subjects, surface EMG signals of wrist movements were recorded by a commercial bio-signal amplifier (16 channels, ME6000, Canada) with a 1000 Hz sampling frequency. Wrist movements included extension, flexion, abduction, adduction, pronation, and supination. Each movement was repeated 10 times; the duration of each repetition and the rest time between the repetitions were recorded 4–5 and 3 s, respectively. Sixty seconds of interruption was considered between movements individually to relieve muscle fatigue. The total duration of each movement was 72 ± 5 s. The EMG data were filtered between 10 and 450 Hz using a fourth-order Butterworth bandpass filter to reduce the effects of noises and movement artifacts. Maximum voluntary contractions were performed before the data were collected to verify the validation of the electrodes and normalize the signal. All processing was done in Python library and MATLAB 64 bit, running on a system with a 1.8 GHz processor and 6 GB of memory. Figure 1 illustrates the experimental setup.

Fig. 1

Experiment setup. a Surface electrode placement on forearm muscles. b Initial wrist position. c Wrist movements. d The muscles activity of wrist extension consist of extensor digitorum (ED), flexor pollicis longus (FPL), extensor carpi ulnaris (ECU), abductor pollicis longus (APL), pronator teres (PT), and supinator (Sup) for 10 times repetition. e Test paradigm which the recording protocols were carried out according to it

The appropriate response time for hand myoelectric control has to be less than 300 ms [22]. Because, in this manner, no delay will be perceived by user [22]. Therefore, the moving window length was set to 256 ms (256 samples at 1000 Hz sampling) in conformity with the previous researches [13, 22]. The more important parameter is window increment length which data processing has to be carried out during this time band. In this study, according to the calculated processing times, the incremental window was set to 10 ms and increased gradually to 150 ms. Each incremental step was 10 ms. In this manner, the influence of the incremental window length on the proposed detection algorithm was evaluated.

2.2 The proposed movement detection algorithm

To address the wrist task discrimination problem, we have proposed a protocol that allows us to use the EMG data for differentiating the tasks from one set. The input of this proposed protocol is EMG data and its outputs are the detected tasks. At the first step, the hierarchical alternating least squares (HALS) algorithm was used for extracting the muscle synergies of the wrist movement from filtered EMG signals. Then, extracted muscle synergies, as the extracted features, were regarded as the input data for classification. Two different classification algorithms were applied and compared. In the first algorithm, a single SOFM network was used as a single classifier. In the second algorithm, a fusion method was utilized alongside six distinct SOFM networks, as six distinct classifiers, for detecting six different wrist movements. The applied SOFM networks were trained using the Kohonen algorithm as an unsupervised algorithm [23]. We evaluated four fusion methods including SOFM, MLP fusion, Bayesian fusion, and Bayesian fuzzy clustering (BFC) fusion algorithm, respectively. The schematic of the proposed scheme is shown in Fig. 2.

Fig. 2

Schematic of the proposed wrist movement detection strategies. In the first algorithm, no fusion algorithms were used and a single SOFM classifier was used. In the second algorithm, the fusion algorithms (MLP based, Bayesian based, BFC based) were used along with six SOFM classifiers

In order to assess the proposed algorithms performance, a K-fold cross-validation (K = 10) method was adopted. For performing this approach, the samples of data recorded during each wrist movement was labeled firstly to determine the type of wrist movement elicited each recorded data. Then, the labeled data was partitioned into 10 subsets. During each round of cross-validation, the training process performed on one subset (called training set), and validating process performed on the other subsets (called the testing set). Multiple rounds of cross-validation were performed using different partitions, and the validation results were combined.

2.3 Identification of the synergy model

EMG signal analysis is a linear combination of time-varying muscle synergies, with each synergy individually having an amplitude and a time shift. The weight coefficient and synergy of each muscle are different for each motion. Each synergy has a weight coefficient and a range of amplitude as well as a time shift, by combining these synergies, the muscle activation pattern is generated [24]. The time-varying synergy model is defined as follows [24]:

$$ M(t)={\sum}_{i=1}^N{c}_i{w}_i\left(t-{t}_i\right)\to M= CW $$

(1)

$$ {\displaystyle \begin{array}{l}W=\left[{w}_1,{w}_2,\dots, {w}_k\right]\in {R}_{+}^{1\times T}\\ {}C=\left[{c}_1,{c}_2,\dots, {c}_k\right]\in {R}_{+}^{1\times T}\end{array}} $$

Here, M is the time-varying muscle activity pattern; N denotes the number of time-varying synergies, c_i is a non-negative coefficient for synergy ith; and t_i is the delay in the start of synergy ith. w is a synergistic matrix and w_i is the muscle activity for ith synergy at time t_i. Here, C and W are non-negative.

The synergy patterns of each movement are extracted by the HALS algorithm [25] from the corresponding EMG signal. In this method, a gradient of a local cost function based on the optimality conditions (KTT) [26] is used to update the stationary points (i.e., the matrix of coefficients and the synergies). The cost function using the Frobenius norm is as follows:

$$ {D}_F^k\left({M}^k\Big\Vert {c}_k{w}_k\right)=\frac{1}{2}\left\Vert {M}^k-\right.{\left.{c}_k{w}_k\right\Vert}_F^2 $$

(2)

Where M^k is the time-varying muscle activity pattern, C_k and W_k are the matrix of synergy matrixes and their corresponding coefficient matrixes, respectively; k denotes the number of current iteration of the HALS algorithm.

Detailed description for the synergy is presented in Table 1.

Table 1 Mathematical description of the HALS algorithm

2.4 Applied classification methods

In this study, two different classification methods have been used. Figure 2 shows the structure of the adopted methods. The first algorithm was designed based on using a single SOFM without using the fusion algorithms. The second algorithm was designed based on using multi SOFM alongside the fusion algorithms. The proposed methods are elaborated in the following subsections.

2.4.1 Classification using single SOFM

In this study, we propose an unsupervised learning and SOFM as a competitive network. During training, the SOFM network, the Euclidean distance between each input training samples vector and the weight vector of each neuron is taken. The Kohonen algorithm was used for training the SOFM network. In performing the Kohonen algorithm, the weight vector of each neuron is updated gradually in a manner that the Euclidean distance between the input training samples vectors and the weight vector of each neuron reduces. In this manner, all the neurons are involved in a competition. At the end of the training process, the neuron which its weight vector has the minimum Euclidean distance to the input training samples vectors is regarded as the winner neuron among all neurons [23]. This winner neuron is called the best matching unit (BMU), and its weight vector is called voronoi vector. The position of the winner neuron within the SOFM network determines the label of the class which the date belongs to it. The position of the winner neuron changes as the class which the input training samples vectors belongs to it, changes. In this manner, the SOFM network can be used for multiclass classification. Accordingly, in this study, the performance of a single SOFM was firstly evaluated for classification of six wrist movements (Fig. 2, the first algorithm). The range of changes is decreased based on the time and distance from the BMU. In this study, a 10 × 10 two-dimensional (2D) lattice SOFM with neighborhood function h(t) with radius of σ(t), a learning-rate parameter η(t), and training of 1000 epochs were used. The learning rate and neighborhood functions are applied as decreasing functions to the SOFM network [13] and defined in (3) and (5), respectively.

$$ h(t)=\exp \left(\frac{-{dist}^2}{2{\sigma}^2(t)}\right) $$

(3)

Here, dist is the neuron distance from BMU in the two-dimensional SOFM network, σ(t) is the width of neighborhood at the instant, and t is defined by an exponential decay function in (4).

$$ \sigma (t)={\sigma}_0{e}^{\left(-\frac{t}{\tau}\right)} $$

(4)

Here, σ₀ is the radius of neighborhood at the instant t₀is 20 [13].

$$ \eta (t)={\eta}_0{e}^{\left(-\frac{t}{\tau}\right)} $$

(5)

Where η₀ is the value η(t) at the instant t₀ is 0.9 [13]. The learning-rate parameter decrease gradually with increasing time t. Constant value τ is 2000 and t is the number of iterations. The t and τ are unit-less constants. The network input is the synergy pattern extracted using HALS algorithm.

In the end, voronoi vectors are applied to the fusion block. Detailed description for the task discrimination with this approach is presented in Algorithm-I

2.4.2 Classification using multi SOFM

In this approach, each SOFM network was used as a two-class classifier to recognize a single movement. Since six different movements must be differentiated, six SOFM networks were applied (Fig. 2, the second algorithm). Then, for improving the performance, different fusion algorithms were used to combine the classification results (Fig. 2, the second algorithm). The applied fusion algorithms are elaborated in the following sub-sections.

MLP-based fusion

In this study, the MLP network has been used to process output data of the SOFM network. The extracted muscle synergy matrixes were regarded as the input data of the MLP network. The MLP network has two hidden layers. Each hidden layer had 10 neurons with the rectified linear activation function (Eq. 6). The number of neurons in the output layer was six, and the softmax was regarded as the activation function of the output neurons (Eq. 7).

$$ F(x)={x}^{+}=\max\;\left(0,x\right) $$

(6)

Where x is the input to a neuron in hidden layer ith.

$$ {O}_k=\frac{\exp \left({z}_k\right)}{\sum_{\mathrm{j}}\exp \left({z}_j\right)} $$

(7)

$$ {z}_k={w}_{0k}+\sum \limits_{j=1}^j{h}_j(x){w}_{kj} $$

In the formula, k represents the index of the output node, and j represents the indexes of all nodes in the group or layer. The Levenberg-Marquardt used as the learning algorithm. The network structure was used as determined by trial and error based on K-fold cross-validation results. Two hidden layers were envisioned which each layer had 10 neurons. Figure 3 shows the structure of the used MLP.

Algorithm-I: Task discrimination (SOFM algorithm)

Step 1: Randomly determine the initial values (t=0) of the network’s weights (from the input i to the output neuron j) (0 ≤ i ≤ n − 1)wij(t) Step 2: Enter the inputs to the network Xt = [x0(t), x1(t), x2(t), …, xn − 1(t)] Step 3: Calculate distances between input vector Xt and each output node j weight vector wij(t) using the Euclidean criterion. ${d}_j={\sum}_{i=0}^{n-1}{\left({x}_i(t)-{w}_{ij}(t)\right)}^2$ Step 4: Select the winner j* neuron with the least distance min(d1, d2, …, dn − 1)=dj* Step 5: update wij(t),σ(t), η(t) The new weight coefficients are calculated using the given equation. wij(t + 1) = wij(t) + η(t)h(t)(x(t)- wij(t)) Step 6: t+1=t, if the network training is not completed, the process is repeated from the step 2.

Fig. 3

Structure of the used MLP artificial neural network with two hidden layers. Where x and y are the input and the output. hⁱ is ith hidden layer and bⁱ is a bias term. Wⁱ are the parameters of the ith layer

Bayesian model-based fusion

The Bayesian classifier is a straight forward based on probabilities inference. There is a probability model for each quantity that can be adopted optimal decisions by assumption about the probability distribution of new data. For estimating the parameters Bayesian model, the maximum likelihood estimation approach is applied. Given the proposed probability model, Bayesian classifier is trained by supervised learning. Let us suppose the X is input data which contains n synergy vectors {x_i : i = 1, 2, …, n}, and Cj shows the jth class which X can belong to it {C_j : j = 1, 2, …, c}. Then, the conditional probability of class C_j to input X is can be calculated as follows [27]:

$$ P\left(X|{C}_j\right)=\frac{P\left({C}_j|X\right).P(X)}{P\left({C}_j\right)} $$

(8)

The conditional probability of class C_j to input X is calculated using minimum error Bayesian classification [27]. The input belongs to a class which has the highest posterior probability.

$$ \underset{c}{\arg \max }P\left({C}_j|X\right)=\underset{c}{\arg \max}\frac{P\left(X|{C}_j\right).P\left({C}_j\right)}{P(X)} $$

(9)

Where

P(C_j): The probability of belonging the input data to the jth class (C_j)

P(C_j| X): The probability of belonging the input data to the jth class (C_j) given that the input data is X

P(X| C_j): The conditional probability of class C_j to input X.

P(X): The probability that the input data is X. P(X) is constant for all classes, only the P(X| C_j)P(C_j) need to be maximized.

Bayesian fuzzy clustering-based fusion

Fuzzy clustering is an unsupervised algorithm which every data point from a data set simultaneously belongs to all clusters with a varying degree of membership [28]. Each data point, X_n for n = 1, 2, …N has a non-negative membership function u_nc in each fuzzy set c = 1.2. …. C. The membership function sum of each data point in all sets is equal to 1, i.e., $\sum \limits_{c=1}^C{u}_{nc}=1$. In the BFC algorithm, first, the clusters center, y_c, is selected as a random variable; then the membership function (u_nc) is defined using (10) [28].

$$ {u}_{nc}=\frac{1/d{\left({X}_n-{y}_n\right)}^2}{\sum_{k=1}^C1/d{\left({X}_n-{y}_n\right)}^2} $$

(10)

Where, d(X_n − y_c)=‖X_n − y_c‖ is Euclidean distance.

The prior distribution for clusters membership functions is defined as follows [28]:

$$ \tilde{p}\left(\mathrm{U}\left|\mathrm{Y}\right.\right)=\left({\prod}_{c=1}^C{u}_{nc}^{\frac{- mD}{2}}\right)\mathrm{Dirichlet}\left({\mathrm{u}}_{\mathrm{n}}\left|\alpha \right.\right) $$

(11)

Where D is the data dimensions. The negative power in $\prod \limits_{c=1}^C{u}_{nc}^{\frac{- mD}{2}}$ produce some values of m and D out of the interval [0, 1]. Therefore, inverse-gamma distributions are used. $\upalpha =\frac{\mathrm{mD}}{2}-1$ have a small value. Dirichlet is defined as follows [28]:

$$ \mathrm{Dirichlet}\left(x\left|\alpha \right.\right)=\frac{\Gamma \left({\sum}_{k=1}^K{a}_k\right)}{\prod_{k=1}^K\Gamma \left({a}_k\right)}{\prod}_{k=1}^K{\chi}_k^{a_{k^{-1}}} $$

(12)

Here, x_k is for all k = 1. …. K x_k ≥ 0 and$\sum \limits_{k=1}^K{x}_k=1$. If α > 1, membership values converge to values (0 or 1) else fuzzy membership values are selected. α = 1_C is a column vector of all ones that clustering output has not been affected. To simplify, conditional distribution of the memberships, p(U|X,Y) is obtained from an uniform symmetric Dirichlet distribution [28].

$$ {u}_n^{+}\sim \mathrm{Dirichlet}\left(\upalpha ={1}_C\right)\kern0.5em $$

(13)

To substitute u_n with the new membership sample ${u}_n^{+}$, first the values u_n and ${u}_n^{+}$ must be obtained by (14), then ${u}_n^{+}$ replaces u_n with the condition (15) [28].

$$ \overset{\sim }{p}\left({x}_n.{u}_n\right|\ Y\left)=p\left({x}_n|{u}_n.Y\right)\overset{\sim }{p}\left({u}_n\right|\ Y\right)\propto \prod \limits_{c=1}^C\exp \left\{\frac{-1}{2}{u}_{nc}^m\left\Vert {X}_n-{\left.{y}_c\right\Vert}^2\right.\right\}{u}_{nc}^{\alpha_c-1} $$

(14)

$$ {a}_u=\min \left\{1.\frac{\overset{\sim }{p}\left({x}_n.{u}_n^{+}\right|\ Y\Big)}{\overset{\sim }{p}\left({x}_n.{u}_n\right|\ Y\Big)}\right\} $$

(15)

Gaussian prior distribution on the clusters center using average and covariance of the clusters center is as follows [28]:

$$ p\left(\mathrm{Y}\right)=\prod \limits_{c=1}^C\mathcal{N}\left({y}_c|{\mu}_y.{\varSigma}_y\right) $$

(16)

$$ {\mu}_y=\frac{1}{N}\sum \limits_{n=1}^N{x}_n\kern0.5em $$

(17)

$$ {\varSigma}_y=\frac{\gamma }{N}\sum \limits_{n=1}^N\left({x}_n-{\mu}_n\right){\left({x}_n-{\mu}_n\right)}^T $$

(18)

Here, the control parameter γ is Gaussian prior variance on the clusters center, in which γ < 1 tends toward average data. For simplicity, the new values of the clusters center are considered as follows [28]:

${y}_n^{+}\sim$ $\mathcal{N}\left({y}_c.\frac{1}{\delta }{\varSigma}_y\right)$ (19)

Here, δ is the sample acceptance rate.

In order to create a new cluster in the samples space first, the cluster values ${y}_n^{+}$ and y_n should be calculated in (20), then the new cluster is replaced by y_n with considering condition (21) [28]:

$$ p\left(\mathrm{X}.{u}_n|\mathrm{U}\right)=p\left(\mathrm{X}|\mathrm{U}.{y}_c\right)p\left({y}_c\right)\propto \exp \left\{\frac{-1}{2}\sum \limits_{n=1}^N{u}_{nc}^m\left\Vert {X}_n-{\left.{y}_c\right\Vert}^2\right.\right\}\times \exp \left\{\frac{-1}{2}{\left({y}_c-{\mu}_y\right)}^T{\varSigma}_y^{-1}\left({y}_c-{\mu}_y\right)\right\} $$

(20)

$$ {a}_y=\min \left\{1.\frac{\overset{\sim }{p}\left({x}_n.{y}_c^{+}\right|\ \mathrm{U}\Big)}{\overset{\sim }{p}\left({x}_n.{y}_c\right|\ \mathrm{U}\Big)}\right\} $$

(21)

Detailed description for the task discrimination with this approach is presented in Algorithm-II.

Algorithm-II: Task discrimination (Bayesian Fuzzy Clustering)

Step 1: Determine the initial values of the clusters (yc) randomly from the samples and calculate the initial membership functions from (10) Step 2: Calculate the mean value and covariance of clusters centers from (17) and (18) Step 3: Create membership functions and new cluster centers for all samples from (13) and (19) Step 4: Calculate the new and initial membership functions using (14) Step 5: If $\overset{\sim }{p}\left({x}_n.{u}_n^{+}\right|\ Y\left)>\overset{\sim }{p}\left({x}_n.{u}_n\right|\ Y\right)$ then un← ${u}_n^{+}$ Step 6: Calculate the new and initial cluster centers using (20) Step 7: If $\overset{\sim }{p}\left({x}_n.{y}_c^{+}\right|\ \mathrm{U}\left)>\overset{\sim }{p}\left({x}_n.{y}_c\right|\ \mathrm{U}\right)$ then yc← ${y}_c^{+}$ Step 8: Step 3 to 7, they are repeated in Niter

3 Results

3.1 Wrist muscle synergies

At the first step, the number of the synergy vectors required to discriminate different tasks was determined. R² criterion was used to determine the number of synergies in (22) [19]. Here, SSE is the sum of the remaining squares, SST is the sum of the remaining squares of the mean activation vector ($\overline{\mathrm{m}}$). Figure 4 shows the number of muscle synergies is used in the HALS algorithm to reach the best reconstruction accuracy. The R² criterion for all individuals and tasks with muscles synergy number was evaluated (1–8) using (22).

$$ {R}^2=1-\frac{\mathrm{SSE}}{\mathrm{SST}}=1-\frac{\sum \limits_s\sum \limits_{k=1}^{k_s}{\left\Vert {m}^s\left({t}_k\left)-{\sum}_i{c}_i^s{w}_i\right({t}_k-{t}_i^s\right)\right\Vert}^2}{\sum \limits_s\sum \limits_{k=1}^{k_s}{\left\Vert {m}^s\left({t}_k\Big)-\overline{m}\right)\right\Vert}^2} $$

(22)

The R² value obtained as the number of synergies were larger than 5 and was nearly stable (R² ≥ 0.96). Figures 4 and 5 show the value of the computed R² changes with respect to the number of the regarded synergy matrixes. As it can be seen, the reconstruction accuracy has significantly increased as the number of synergy matrixes hit 5 and more (p < 0.05). But, no significant improvement was observed when the number of synergy matrixes increased more the 5 (p > 0.05). As shown in Fig. 6a, reconstruction accuracy has been significantly improved due to increasing the number of synergy patterns from one to five (p < 0.05). According to the carried out analyses, the achieved results were repeatable. For example, as shown in Fig. 6b, the same result has been achieved using the data recorded during different recording trials. It can be seen again in Fig. 6 that reconstruction accuracy has been significantly improved due to increasing the number of synergy patterns from one to five.

Fig. 4

R² index variation vs. synergy vectors number plotted for each wrist muscle

Fig. 5

The value of the computed R², as a quantitative measure showing the reconstruction accuracy, with respect to the number of the regarded synergy patterns

Fig. 6

The sample reconstructed muscle activation patterns, extracted as the number of the assigned muscle synergy patterns has been increased, using the collected data during two different trials

3.2 Quantitative comparison of the proposed methods

3.2.1 Evaluation without applying the moving window

At first, the performance of the proposed detection algorithm implemented using different fusion strategies was evaluated while only 70% of the data recorded during performing movement used for processing and detection. In other words, the early and final parts of the data samples, related to muscle activation onset and offset, were removed. In this manner, the data were related to the transient times coincide with the muscle activation onset and offset were removed. As a preliminary evaluation, it can be an acceptable compromise between off-line accuracy and real-time performance [8].

In order to assess algorithms performance, a K-fold cross-validation (K = 10) method was adopted. Then, the quantitative measures computed to compare the efficacy of the applied classification methods using ANOVA test. The statistical analysis showed that there is no significant difference (p > 0.05) among the overall accuracies and sensitivities of the classification methods designed utilizing the MLP, Bayesian, and BFC fusion algorithms along with the classifiers. Overall specificity of BFC-based classifiers significantly differs from two other classifiers designed based on MLP and Bayesian fusion algorithms (p < 0.05). In addition, there is a significant difference (p < 0.05) between the performance of the SOFM-based classification method and used fusion-based classification methods.

Figure 7 shows the overall results and comparison of the proposed algorithms using the data recorded from different subjects. As Fig. 7 shows, using a single SOFM entailed the accuracy value, sensitivity value, and specificity value equal to 87.88 ± 1.54%, 83.38 ± 1.93%, and 74.51 ± 3.03%, respectively. Also, the accuracy values for MLP, Bayesian, and Bayesian fuzzy clustering (BFC) fusion algorithms were 99.78 ± 0.45%, 99.33 ± 0.80%, and 96.43 ± 1.08%, respectively. Also, sensitivity value for MLP, Bayesian, and Bayesian fuzzy clustering (BFC) fusion algorithms were 99.78 ± 0.45%, 99.31 ± 0.79%, and 95.93 ± 1.74%, respectively. In addition, the specificity values for MLP, Bayesian, and Bayesian fuzzy clustering (BFC) fusion algorithms were 98.96 ± 1.06%, 98.56 ± 1.08%, and 90.63 ± 1.56%, respectively. Overall, the average accuracy, average specificity, and average sensitivity achieved through utilizing the fusion algorithms were 97%, 95%, and 96%, respectively. While the average accuracy, average specificity, and average sensitivity achieved through utilizing a single SOFM classifier were less than 88%, less than 75%, and less than 85%, respectively.

Fig. 7

The computed accuracy (a), specificity (b), and sensitivity (c) of each implemented classification method using the recorded data related to 8 subjects without applying the moving window

The presence of learning interference due to using a single SOFM classifier might cause such relatively low efficient classification. In fact, it can be clearly seen that adopting the fusion algorithms has improved the movements separation.

3.3 Evaluation with applying the moving window

After preliminary evaluations, the detection algorithms were evaluated again while a 256 ms moving window was applied. In fact, a sliding window (256 ms) was moved over the data recorded during a movement. The time window was moved by the moving step called window increment length. In this context, an online mechanism for muscle activation onset detection (integrated profile method) has to be used [29]. According to the adopted mechanism, the continuous integration of rectified SEMG was conducted. Then, the difference between the conducted integration and a reference line was computed over the time. The time instance at which the computed difference value was equal to the maximum value and was regarded as the onset time [29]. Table 2 describes the applied detection strategy, briefly.

Table 2 The integrated profile method utilized for muscle activation onset detection [29]

As emphasized previously, the window increment length which data processing has to be carried out during this time band is so important. In fact, the processing time must be less than or at least equal to the envisioned window increment length. Therefore, the processing time needed for performing the proposed algorithm as implemented using different incorporated fusion algorithms. Since the classification performance using a single SOFM was not acceptable, the single SOFM classifier was excluded. Table 3 shows the calculated processing time related to the proposed detection method as implemented using different fusion algorithm.

Table 3 The calculated processing time related to the proposed detection method as implemented using different fusion algorithm

As Table 3 shows, the needed processing times are less than 10 ms. Accordingly, the increment window length must be set to at least 10 ms. In order to determine the optimal increment window length, incremental window length increased gradually from 10 ms to 150 ms. Each incremental step was 10 ms. The statistical analyses (AVOVA) proved that increasing the increment window length did not have a significant influence on the overall performance of the proposed classification methods designed using incorporating the different fusion algorithms (p > 0.05). Therefore, the increment window length was set to the least possible size equals 10 ms. Then, the assessment procedures repeated again using K-fold cross-validation method (K = 10). As Fig. 8 shows, using a single SOFM entailed the accuracy value, sensitivity value, and specificity value equal to 86.54 ± 2.12%, 80.73 ± 1.63%, and 72.98 ± 2.13%, respectively. Also, the accuracy values for MLP, Bayesian, and Bayesian fuzzy clustering (BFC) fusion algorithms were 99.78 ± 0.45%, 99.33 ± 0.80%, and 96.43 ± 1.08%, respectively. Also, sensitivity value for MLP, Bayesian, and Bayesian fuzzy clustering (BFC) fusion algorithms were 98.51 ± 0.83%, 97.52 ± 0.95%, and 95.80 ± 1.51%, respectively. In addition, the specificity values for MLP, Bayesian, and Bayesian fuzzy clustering (BFC) fusion algorithms were 97.77 ± 0.50%, 97.15 ± 0.11%, and 88.24 ± 2.11, respectively. Overall, using the fusion algorithms yielded the average accuracy more than 96%, the average specificity more than 93, and the average sensitivity more than 94%. According to the statistical analysis, it was showed again that there is no significant difference (p > 0.05) among the overall accuracies and sensitivities of three utilized fusion-based classification methods. While overall specificity of BFC-based classifiers significantly differs from two other classifiers designed based on MLP and Bayesian fusion algorithms (p < 0.05).

Fig. 8

The computed accuracy (a), specificity (b), and sensitivity (c) of each implemented classification method using the recorded data related to 8 subjects with applying the moving window (moving window length 256 ms, increment window length 10 ms)

4 Discussion

In this paper, the separation of different wrist movements based on muscle synergy patterns was performed using SOFM and based on fusion (i.e., MLP, Bayesian, and BFC) algorithms. The idea of this approach is to use synergy patterns to cover the coordination between the muscles. Furthermore, the fusion mechanism was used to enhance the separation of movements due to the different synergy patterns. The obtained results can be analyzed from different aspects. We will discuss these aspects as follows.

4.1 Number of the synergy patterns

In this study, the effect of the number of extracted synergy vectors on reconstruction error of the wrist muscle activation pattern was analyzed. The results show that increasing the number of synergy vectors more than 5 did not improve the reconstruction performance. In addition, the performance of the movement detection for all subjects are close to each other using five identified synergy patterns. As a result, it can be considered that muscle activity pattern of each movement is extractable with a limited number of synergy vectors. Such results are in conformity with assumptions related to the muscle synergy theory [18], [30, 31].

4.2 The role of the fusion algorithms

The preliminary results that showed the performance of a single SOFM, as a classifier, were not acceptable. It was attributed to the learning interference could arise from possible mutual dependence among the different EMG signals recorded from different participants. Therefore, adopting the fusion approaches seemed necessary to prevent learning interference [32]. Before applying the fusion algorithm, some pre-processing is recommended to enhance the separability of the input data [32, 33]. In this study, three fusion algorithms with non-linear, probabilistic nature, and fuzzy characteristics have been utilized. The obtained results showed that when a fusion method is used, the results of various movements separation in all participants increased significantly (p < 0.05) compared to a classification using a single SOFM. In fact, instead of using a classification block SOFM for all different movements, a classification block SOFM is considered for each separate movement. Then, the final output of the classification blocks is applied to a fusion block. The statistical analysis showed that there is no significant difference (p > 0.05) among the overall accuracies and sensitivities of three utilized fusion-based classification methods. But, it is worth noting that the average specificity of BFC-based classification method is significantly less than the average specificity of the classification methods designed based on two other fusion algorithms. The BFC fusion algorithm works based on creating the membership functions for each clusters and determining the center of cluster using a probabilistic approach [28]. Therefore, the performance of this fusion algorithm is fully dependent on the clustering accuracy. In other words, mis-clustering can give rise to increase the number of members negatively clustered in the clusters except the desired cluster (false positive (FP)). Thus, such relatively lower specificity can be attributed to relatively high value of FP arises from likely mis-clustering.

4.3 Comparing with the previous works

Movement discrimination based on synergy muscle can potentially be more appropriate than other methods since it contains vital information in relation to the muscle activity pattern generated by neural processes [8]; because synergy vectors represent contractions of low-dimensions muscle activities [34]. Hence, the analysis of synergy vectors with different techniques could be a good strategy to improve the approach of the upper limb prostheses. In the last decade, studies have been conducted for movement separation, especially for controlling prostheses based on muscle synergy analyses [3, 8]. The reported separation accuracies are about 95% [3, 8]. While, in this study, the separation average accuracy was obtained more than 97%. This relatively increased accuracy can be attributed to effect of incorporating the fusion algorithms and using six distinct SOFM neural networks as classifiers. This elucidates again the degrading role of learning interference on the performance of the single SOFM-based classifier.

5 Conclusion

In this study, the separation of different wrist movements based on muscle synergy analysis was investigated. In fact, instead of considering each muscle activity patterns and features of each muscle individually, the synergy pattern of movement was extracted. Since the separation of the activity of each muscle and the elimination of the interactions of other muscles is difficult, the authors believe that synergy vectors can be used directly for controlling the rehabilitation and dependent devices in further innovation. The promising achieved classification accuracy can confirm such believe. Nevertheless, more experimental data has to be recorded and has to be analyzed to prove decisively that the achieved results are generalizable. In addition, though the EMG is a good tool to analyze the muscle electrical activity, but the EMG devices are relatively expensive devices and the EMG preprocessing is time consuming. Thus, using the cheap and reliable muscle sensors is preferable for practical implementation of electromyography-based motion detection system. It is recommended to perform initial evaluations on patients in order to implement it practically on basic microcontroller systems.