Upper Limb Action Identification Based on Physiological Signals and Its Application in Limb Rehabilitation Training

Traditionally, the sEMG signal and acceleration signal from the upper limbs, including shoulder muscles, arm muscles, forward arm muscles and hand muscles, are acquired through segmentation. Despite reflecting the difference between dynamic upper limb actions, the traditional method cannot be directly applied to upper limb action recognizers to classify upper limb actions. Therefore, it is necessary to extract effective upper limb action features from physiological signals of different segments. In this paper, the upper limb actions in limb rehabilitation training are identified by time-domain feature method, frequency-domain feature method, time-frequency domain feature method, and entropy feature method.

2.1 Time-frequency domain feature extraction

The mathematical statistical method was adopted to process the time variation of the waveforms of the collected sEMG signal and acceleration signal. The acquired parameter features of the signals are the time-domain features of upper limb actions. Let a_m={a₁,a₂,...,a_M} be the time series of physiological signal samples collected from a segment, with M being the signal length. The degree of oscillation of the collected physiological signals was characterized by the number N_ZERO of zero points being passed, i.e., the number of times that the acquired signal curve passes through the time axis:

$\begin{aligned}

&N_{\text {ZERO }}=\sum_{m=2}^{M} \operatorname{sgn}\left(-a_{m} a_{m-1}\right) \\

&\operatorname{sgn}(a)=\left\{\begin{array}{l}

1, \text { if } a \geq 0 \\

0, \text { otherwise }

\end{array}\right.

\end{aligned}$        (1)

Then, the mean absolute value of the physiological signal samples from a segment:

${{V}_{ABS}}=\frac{1}{M}\sum\limits_{m=1}^{M}{\left| {{a}_{m}} \right|}$       (2)

To reflect the energy variation of the collected physiological signals between segments, the mean absolute value was added by a weight θ to obtain the weighted mean absolute value:

$\begin{aligned}

&V_{W A B S}=\frac{1}{M} \sum_{m=2}^{M} \theta\left|a_{m}\right| \\

&\theta= \begin{cases}1, \text { if } & 0.25 M \leq m \leq 0.75 M \\

0.5, \text { otherwise } & \end{cases}

\end{aligned}$       (3)

The root mean square (RMS) can be calculated by:

${{V}_{ROOT}}=\sqrt{\frac{1}{M}\sum\limits_{m=1}^{M}{a_{m}^{2}}}$       (4)

In order to characterize the waveform complexity of the collected physiological signals, the absolute values of the difference between the adjacent sample points were accumulated, and the cumulative sum was defined as the waveform length of the collected physiological signals:

${{V}_{WAVE}}=\sum\limits_{m=1}^{M-1}{\left| {{a}_{m+1}}-{{a}_{m}} \right|}$        (5)

To prevent the physiological signal length from affecting the identification result, this paper introduces the mean waveform length feature to reflect the change law of the collected physiological signals:

${{V}_{AW}}=\frac{1}{M}\sum\limits_{m=1}^{M-1}{\left| {{a}_{m=\text{1}}}-{{a}_{m}} \right|}$       (6)

The EMG signal is the superposition between the action potentials of the motor units in various muscle fibers. The level of the action potentials can be characterized by Wilson amplitude, i.e., the number of amplitude differences between two adjacent sample points being greater than the preset threshold FZ. This parameter was introduced to characterize the contraction degree of upper limb muscles:

$\begin{aligned}

&V_{W I L}=\sum_{m=1}^{M-1} g\left(\left|a_{m+1}-a_{m}\right|\right) \\

&g(a)= \begin{cases}1, \text { if } & a \geq F Z \\

0, \text { otherwise } & \end{cases}

\end{aligned}$       (7)

The collected physiological signals were subjected to Fourier analysis, and the power spectrum analysis was conducted to extract frequency domain features. The extracted features can reveal the distribution and variation of the frequency bands of the collected physiological signals. Let γ_j and η_j be the frequency and spectrum intensity of the collected physiological signals, respectively; N be the length of the power spectrum. Then, the median frequency, mean frequency, and peak frequency of the collected physiological signals can be respectively calculated by:

${{F}_{M}}=\sum\limits_{j={{F}_{M}}}^{N}{{{\eta }_{j}}=\frac{1}{2}}\sum\limits_{j=1}^{N}{{{\eta }_{j}}}$       (8)

${{F}_{AV}}=\sum\limits_{j=1}^{N}{{{\gamma }_{j}}{{\eta }_{j}}/}\sum\limits_{j=1}^{N}{{{\eta }_{j}}}$       (9)

${{F}_{PE}}=max\left( {{\eta }_{j}} \right),j=1,....,N$       (10)

The wavelet coefficient set, which contains rich information, can be obtained through wavelet decomposition of different mother wavelets. Let ξ be the decomposition coefficient for the j-th layer component obtained through wavelet decomposition; Q_j,l be the l-th coefficient on the corresponding layer. Then, the wavelet coefficient energy can be taken as the time-frequency domain feature of the collected physiological signals:

$TF{{C}_{j}}=\sqrt{\frac{1}{\xi }\sum\limits_{\xi =1}^{\xi }{Q_{j,\xi }^{2}}}$       (11)

This paper selects db6 wavelet to perform three-layer wavelet decomposition of the collected physiological signals. The TFC_j of each layer was solved, and taken as the time-frequency feature of that layer for analysis.

2.2 Multiscale entropy feature extraction

Approximate entropy, permutation entropy, fuzzy entropy, and sample entropy, as measures of the orderliness of the target system, are widely adopted to analyze physiological signals like EMG signal and electroencephalogram (EEG) signal. The fuzzy entropy FE can measure the fuzziness of the fuzzy domain. Let n be the dimensionality of the collected physiological signals; _jSEM_ij be the similarity between two samples; s be the mean similarity of SEM_ij. Then, FE can be calculated by:

$\begin{aligned}

&\operatorname{FE}(M, n, s)=\lim _{M \rightarrow \infty}\left(\ln \Phi_{n}-\Phi_{n+1}\right) \\

&\Phi_{n}=\frac{1}{(M-n)} \sum_{i=1}^{M-n}\left[\frac{\sum_{j=1, i \neq i}^{M-n} S E M_{i j}^{n}}{(M-n-1)}\right]

\end{aligned}$       (12)

Approximate entropy is often used to quantify the regularity and unpredictability of the fluctuations of the time series of signals. For a signal time series, the approximate entropy is positively correlated with its complexity. Similar to approximate entropy, sample entropy can measure the complexity of signal time series. But the latter boasts two advantages: independence of data length, and high consistency. The sample entropy is under the balanced impacts from the variation in the dimensionality n of the collected physiological signals, and the variation in the mean similarity s of SEM_ij. The sample entropy can be calculated by the following steps:

Step 1. Construct a n-dimensional vector of the collected physiological signals in turn:

${{A}_{n}}\left( i \right)=\left[ a\left( i \right),a\left( i+1 \right),...,a\left( i+n-1 \right) \right]$       (13)

where, i=1,2,...,M-n.

Step 2. Define the distance between vectors A_n(i) and A_n(j) as the largest difference between the corresponding elements:

$u\left[ {{A}_{n}}\left( i \right),{{A}_{n}}\left( j \right) \right]=\underset{l=0,1,...,n-1}{\mathop{max}}\,\left| a\left( i+l \right)-a\left( j+l \right) \right|$       (14)

where, i, j=1,2,...,M-n, and i≠j.

Step 3. Count the number Xⁿ_i(s) of u[A_n(i),A_n(j)] smaller than the mean similarity s, and compute the ratio Yⁿ_i(s) of Xⁿ_i(s) to the total number of distances M-n+1:

$Y_{i}^{n}\left( s \right)=\frac{X_{i}^{n}\left( s \right)}{N-n+1}$        (15)

Step 4. Compute the mean Yⁿ(s) of Y_iⁿ(s):

${{Y}^{n}}\left( s \right)=\frac{1}{M-n}\sum\nolimits_{i=1}^{M-n}{Y_{i}^{n}\left( s \right)}$       (16)

Step 5. Construct a n+1-dimensional vector of the collected physiological signals in turn, and repeat Steps 1-4 to obtain Yⁿ⁺¹(s).

Step 6. Compute the sample entropy of the physiological signal sample series:

$SE\left( n,s \right)=-\underset{M\to \infty }{\mathop{lim}}\,\frac{{{Y}^{n+1}}\left( s \right)}{{{Y}^{n}}\left( s \right)}$       (17)

Step 7. Due to the limited length of the collected physiological signals, denote the natural logarithm as ln, and estimate the sample entropy of the series of M samples, before applying the method to upper limb action identification:

$SE\left( n,s,M \right)=-ln\frac{{{Y}^{n+1}}\left( s \right)}{{{Y}^{n}}\left( s \right)}$       (18)

Any of the above entropy metrics, whether approximate entropy, permutation entropy, fuzzy entropy, or sample entropy, can quantify the regularity and complexity of signal time series. In real-world limb rehabilitation training, however, the physiological signal time series generated by the complex physiological system are extremely complex in time and space dimensions. The upper limb actions may be recognized insufficiently with only one metric. Hence, this paper adopts multiscale entropy to identify the upper limb actions in limb rehabilitation training.

Multiscale entropy is essentially derived from the sample entropies of signal time series on different scales. The derivation process is relatively simple. But multiscale entropy does not apply to the analysis of non-stationary nonlinear signals, and could not reflect the low-frequency components in the collected physiological signals. Therefore, empirical modal decomposition (EMD) was carried out to extract the sample entropy features of the collected physiological signals. Specifically, EMD was implemented to extract the intrinsic mode functions (IMFs) from the signal time series. Then, the sample entropy was calculated from the cumulative sum of the IMFs, and taken as the feature of the signal time series. This strategy fully retains the rich information in the various frequency components of the signal time series, and remains robust to the variation in low-frequency components, providing an efficient feature for identifying upper limb actions. The EMD process is as follows:

Step 1. Extract the maximum and minimum points of the original signal time series a(p).

Step 2. Fit the upper envelope BH_max(p) between maximum points and the lower envelope BH_min(p) between minimum points of the original signal time series, using the third-order spline difference function.

Step 3. Compute the mean n(p)=BH_min(p)+ BH_max(p))/2 of BH_max(p) and BH_min(p).

Step 4. Extract the details of z(p)=a(p)-n(p) the signal time series.

Step 5. Judge whether z(p) satisfies the two conditions of IMFs. If yes, treat z(p) as a component of IMF; otherwise, repeat Steps 1-4 on z(p) until it satisfies the two conditions.

Step 6. Subtract the extracted IMF component from a(p). Repeat Steps 1-5 to obtain the next IMF, until all IMFs are obtained. Define the remaining monotonic function with no extreme as the residual CA. Let C be the number of IMFs Ω. Then, the final decomposition of the EMD can be expressed as:

$a\left( p \right)=\sum\limits_{i=1}^{C}{{{\Omega }_{i}}\left( p \right)}+CA$      (19)

The EMD-based sample entropy feature can be derived from the cumulative sum of the IMFs obtained through the EMD. Starting from the first IMF, the smallest scale, to the end of the complete physiological signals, the total cumulative number QS^k(t) of IMFs can be expressed as:

$QS_{{}}^{c}=\sum\nolimits_{i=1}^{c}{{{\Omega }_{i}}\left( p \right)}$        (20)

Taking residual CA as the C+1-th IMF, QS^C+1 corresponds to the original signal time series a(p). The sample entropy feature based on the EMD can be calculated by:

$SEF\left( c,n,s \right)=SE\left( QS_{IMF}^{c}\left( p \right),n,s \right)$        (21)

The information fusion between sEMG signal and acceleration signal is crucial to the identification of upper limb actions. This paper extracts multiple features from the main feature segments of the upper limb actions of patients. Figure 3 is the scatterplot of the features of the two types of physiological signals. It can be observed that sEMG signal and acceleration signal are not greatly disturbed by external factors, and are highly discriminable. This is because the patient’s action trajectories vary significantly between stages of rehabilitation. Nevertheless, the action features of different patients in different stages may overlap, for different patients have different habits in executing the same actions. It is necessary to assist the action identification with sEMG signal information that accurately mirrors the fine actions of the upper limbs.

The proposed upper limb action identification method, which fuses sEMG signal with acceleration signal, is dominated by sEMG signal. The reason is very simple: the details on hand muscles, wrist muscles, and elbow muscles conveyed by sEMG signal are the key to the successful recognition of upper limb movement.

This paper carries out experiments separately on time-domain features, time-frequency domain features, entropy-domain features, and fused features (which merge the features from different domains). On this basis, the authors analyzed how a single information, i.e., the source of sEMG signal influences the identification of dynamic upper limb actions. Table 1 compares the feature identification accuracy on different patients. The fused features led to the highest accuracy and best effect of upper limb action identification.

Feature extraction was performed on the physiological signals from six patients. The extracted features were divided into a training set and a test set. During the experiments, 140 feature samples were extracted randomly each time. 80 samples were allocated for training, and the remaining 60 for testing. Three different models, namely, KNN, ANN, and SVM, were applied separately to identify upper limb actions. The results in Table 2 suggest that the identification accuracy was significantly improved using the fused signal, resulting in an ideal identification effect. For all three classification models, more than 90% of upper limb actions were identified accurately, using the fused features of sEMG signal and acceleration signal. Therefore, the proposed physiological signal feature extraction approach adapts well to different upper limb action identification models.

Following the proposed flow for upper limb rehabilitation evaluation, the final evaluation results were obtained for the upper limb habilitation of the six patients (Table 3). The proposed evaluation model achieved largely the same results with clinical evaluation, reflecting the feasibility and effectiveness of our model.

3.png

Figure 3. sEMG signal features of different patients

Table 1. Feature identification accuracy on different patients

Table 2. Mean accuracy of upper limb action identification

Note: KNN and ANN are short for k-nearest neighbors and artificial neural network, respectively.

Table 3. Evaluation results on upper limb motor function

Finally, 13 muscles were selected from the shoulder muscles, arm muscles, fore arm muscles, and hand muscles of a patient subjecting to different rehabilitation trainings. The relative activation degree and activation duration of these muscles were tested in each training. The test results are presented in Figures 4 and 5. The relative activation degree of the muscles was characterized by motor function, balance function, and coordination function, while the activation duration of these muscles was characterized by training duration and training intensity. As shown in Figures 4 and 5, after the patient completed rehabilitation trainings (i.e., active training, anti-resistance training, and passive training), both the relative activation degree and activation duration of the muscles were improved. The results are consistent with the evaluation results on the patient’s upper limb motor function. This further confirms the effectiveness of our evaluation model.

4.png

Figure 4. Relative activation degree of upper limb muscles of a patient in different rehabilitation trainings

5.png

Figure 5. Activation duration of upper limb muscles of a patient in different rehabilitation trainings