Motion unit discharge sequence long-time online decoding method based on surface myoelectricity

Abstract

The invention discloses a surface myoelectricity-based motion unit discharge sequence long-time online decoding method, which relates to the field of biomedical signal processing and comprises the following steps of: identifying an initial recruitment motion unit in an offline training stage; establishing an initial recruitment motion unit decoding model in an offline training stage; real-time decoding of the motion unit discharge sequence in the online decoding stage is carried out; post-processing the motion unit discharge sequence in the online decoding stage; in the online decoding stage, motion unit recognition is collected; establishing a new recruitment motion unit decoding model in an online decoding stage; and adjusting the characteristic change motion unit decoding model in the online decoding stage. According to the method, a decoding model updating framework is constructed, a data enhancement strategy is designed, and the discharge sequence of the motion unit can be accurately and fully decoded.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of biomedical signal processing, and particularly relates to a long-time online decoding method for motor unit discharge sequence based on surface electromyography. BACKGROUND

[0002] The movement command in the natural movement process of human body originates from the cerebral cortex, is conducted to the spinal cord to cause motor neuron discharge, and then drives the innervated muscle fiber to produce contraction to realize the movement intention. Therefore, a single motor neuron and its innervated muscle fiber can be regarded as a functional whole, which is often referred to as a motor unit. The discharge information of the motor unit represents the neural drive signal transmitted from the spinal cord to the muscle, and the accurate online decoding of such signal plays a key role in the scenes of movement rehabilitation, prosthesis control, human-computer interaction and other people's livelihood concerns.

[0003] The surface electromyography signal is obtained by attaching the electrode to the skin surface, and can non-invasively record the discharge information of most motor units in the covered area. Therefore, real-time reverse decoding of the surface electromyography signal is an effective means to obtain rich neural drive signals. The motor unit discharge sequence decoding model is mostly established in the offline training stage, and then directly acts on the online electromyography signal to obtain real-time decoding results. However, in the long-time surface electromyography signal acquisition process, many non-steady-state factors make the online process and offline process have differences beyond the generalization ability of the decoding model, and then cause challenges to the online decoding of the motor unit discharge sequence.

[0004] These non-steady-state factors can be divided into noise-related and motor unit-related two categories. The first category, long-time electrode attachment will cause sweating on the skin surface, reduce the attachment quality of the electrode, and increase the noise level in the collected signal. The second category, long-time electromyography acquisition will cause muscle fatigue, on the one hand, the previously recruited motor units will change in characteristics due to the body protection mechanism, on the other hand, the recruitment of new motor units is triggered to drive the musculoskeletal system to maintain the expected movement and force. In summary, the long-time online decoding of motor unit discharge sequence based on surface electromyography needs to adapt to the changes caused by the three non-steady-state factors of increasing noise level, characteristic changes of the originally recruited motor units, and recruitment of new motor units.

[0005] The existing methods for online decoding of motor unit discharge sequences based on surface electromyography can be divided into two categories. The first category of methods is based on blind source separation technology (Y. Zheng and P. Zhou, "Adaptive Real-Time Decomposition of Electromyogram During Sustained Muscle Activation: A Simulation Study," in IEEE Transactions on Biomedical Engineering, vol. 69, no. 2, pp. 645-653, 2022, doi: 10.1109 / TBME.2021.3102947), which relies on whitening operation to ensure decoding accuracy. However, whitening operation is sensitive to noise, which makes the blind source separation-based method less robust to increased noise levels. The second category of methods is based on deep learning technology (A. K. Clarke et al., "Deep Learning for Robust Decomposition of High-Density Surface EMG Signals," in IEEE Transactions on Biomedical Engineering, vol. 68, no. 2, pp. 526-534, 2021, doi: 10.1109 / TBME.2020.3006508), which improves the robustness of the online decoding process to increased noise levels through noise-oriented data augmentation strategies. However, the architecture of the deep learning-based decoding model used is complex, and the training time is in the order of hours, so the establishment of the decoding model is limited to the offline stage, and the two types of changes related to motor units that occur in the online decoding stage cannot be trained and learned. Based on the above analysis, the existing methods cannot achieve long-time online decoding of motor unit discharge sequences.

[0006] Therefore, the skilled person in the art is committed to developing a long-time online decoding method for motor unit discharge sequences based on surface electromyography. The method can perform online decoding of long-time electromyography signals under the coexistence of three non-steady-state factors: increased noise, changes in motor unit characteristics, and recruitment of new motor units, and can accurately obtain neural discharge information for a long time. SUMMARY

[0007] In view of the above defects of the prior art, the technical problem to be solved by the present application is how to accurately and sufficiently decode the discharge sequence of motor units during long-time electromyography acquisition under the coexistence of multiple non-steady-state factors.

[0008] To achieve the above object, the application provides a surface electromyography-based motor unit discharge sequence long-time online decoding method, which comprises the following steps: Step 1, identification of initial recruitment motor units in the offline training stage, comprising: Step 1.1, 25-second electromyography signals collected offline are processed by a fast independent component analysis algorithm to obtain a candidate independent component set { E Offline} and a corresponding separation vector set { I Candidate}; T Candidate Step 1.2, a candidate discharge sequence set { I Candidate} is extracted from the candidate independent component set { S Candidate}; Step 1.3, the candidate discharge sequences are screened based on the design index of motor unit discharge characteristics to obtain an effective discharge sequence set { S Valid}} and a corresponding separation vector set { T Valid}; Step 1.4, repeated items in the effective discharge sequence set { S Valid} are removed to obtain a discharge sequence set { S Initial} of the initial recruitment motor units, and the corresponding separation vector T Valid constitutes a set { T Initial}; Step 1.5, the action potential waveform of each initial recruitment motor unit is extracted by using the pulse trigger superposition average method to constitute an action potential waveform set { Q Initial}; Step 1.6, the total separation vector set { T} is initialized as { T Initial}, and the total motor unit action potential waveform set { Q} is { Q Initial}; Step 2, establishment of a decoding model of the initial recruitment motor units in the offline training stage, comprising the following steps: Step 2.1, generation of enhanced training data based on noise injection; Step 2.2, establishment of a decoding model based on deep learning; ​Step 2.3, initialize the total decoding model set M} as M Initial}; Step 3, real-time decoding of motor unit discharge sequence in online decoding stage, including: every 0.125 seconds, using all models in decoding model set M} to decode 0.125 seconds of electromyographic signals collected online to obtain the prediction result of real-time discharge sequence; Step 4, post-processing of motor unit discharge sequence in online decoding stage: post-processing each predicted discharge sequence obtained in step 3 to remove redundant discharge pulses; Step 5, identification of recruited motor units in online decoding stage; Step 6, establishment of decoding model for newly recruited motor units in online decoding stage; Step 7, adjustment of decoding model for characteristic change motor units in online decoding stage.

[0009] Further, the step 1.1 further comprises: Step 1.1.1, performing band-pass filtering and 50 Hz comb filtering on E Offline to obtain filtered electromyographic signals E Filtering ; Step 1.1.2, performing 4 times of time domain delay shift with data points of 1, 2, 3 and 4 respectively on E Filtering ; E Filtering and 4 time domain shift signals to obtain extended electromyographic data E Extend ; Step 1.1.3, performing zero-phase component analysis whitening on the extended electromyographic data to obtain whitened electromyographic data E Whiten , denoted as E Whiten , the first dimension number of n Whiten ; Step 1.1.4, performing 50 times of iterative analysis on E Whiten to obtain 50 candidate independent components I Candidate and corresponding separation vectors T Candidate .

[0010] Further, the step 1.2 further comprises: Step 1.2.1, screening ICandidate all local extrema in the sequence; Step 1.2.2, K-means clustering on all local extrema, the class corresponding to the maximum average value is the discharge class; Step 1.2.3, initialization S Candidate a full zero vector of the same length, let I Candidate a full zero vector of the same length, let S Candidate all discharge class element indexes in the sequence are 1, that is, the candidate discharge sequence in the form of binary corresponding to discharge state 1 and rest state 0 is obtained.

[0011] Further, the step 1.5 further comprises: Taking the discharge pulse as the trigger point, each discharge pulse is intercepted for a period of time before and after the discharge pulse n MUAP The electromyographic signals of the sampling points are sampled, and then all the electromyographic signals obtained by the interception are averaged to obtain the action potential waveform.

[0012] Further, the step 2.1 further comprises: Step 2.1.1, zero-mean Gaussian white noise addition; add a zero-mean Gaussian white noise signal to the offline surface electromyographic signal to obtain a noisy electromyographic signal, and the level of added noise is 5 decibels; Step 2.1.2, motor unit main energy channel determination; analyze the 10 channels with the largest amplitude of the current motor unit action potential waveform, select the smallest rectangle that can envelope the 10 channels as the main energy region, and all channels in the main energy region are motor unit main energy channels; Step 2.1.3, point-by-point sliding window generation of training data; the offline surface electromyographic signal and its noisy signal at the motor unit main energy channel are subjected to sliding window processing with a window length of 0.02 seconds and a step length of 0.0005 seconds to obtain training input data, and the label of each training data is determined by the discharge state of the discharge sequence obtained in the step 1.4 at the middle time of the window, the discharge state corresponds to the label 1, and the rest state corresponds to the label 0; Step 2.1.4, inter-class imbalance mitigation; the data of a single round of training process is composed of all training data corresponding to the discharge state and randomly extracted training data corresponding to the rest state, and the proportion of the two is controlled to be about 1:8.

[0013] Further, the step 2.2 further comprises: Step 2.2.1, initialization of the decoding model architecture, the core of which is a three-dimensional convolution layer for feature extraction and a fully connected layer for feature mapping; Step 2.2.2, training the decoding model using the training data obtained in the step 2.1, and using the adaptive moment estimation optimizer to update the model weights; Step 2.2.3: Validate the decoding models obtained in each round of training, and retain the model with the highest decoding performance index, including: The 25-second offline surface electromyography (EMG) signal at the main energy channel of the motor unit determined in step 2.1.2 is used to obtain sequential EMG window data by sliding a window with a window length of 0.02 seconds and a step size of 0.0005 seconds. The sequential electromyography window data is input into the decoding model to obtain the predicted discharge sequence; The consistency between the predicted discharge sequence and the reference discharge sequence obtained in step 1.3 is calculated, which corresponds to the decoding performance of the model. The consistency is evaluated by the F1-score, and the calculation formula is as follows:

[0014] Wherein, TP represents the number of correctly predicted discharges in the predicted discharge sequence, FP represents the number of incorrectly predicted discharges in the predicted discharge sequence, and FN represents the number of missed predicted discharges in the predicted discharge sequence.

[0015] Furthermore, step 4 also includes: Step 4.1: Group all discharge pulses in the original discharge sequence from front to back according to their proximity, with the interval between the first and last discharge pulses in each group not exceeding 0.003 seconds; Step 4.2: For a single discharge pulse group, calculate the confidence level of each discharge pulse and retain the discharge pulse with the highest confidence level as the unique discharge pulse of the current discharge pulse group. Step 4.3: Repeat step 4.2 until the discharge pulse is unique for all groups in step 4.1.

[0016] Furthermore, step 5 also includes: Step 5.1: For a 25-second online surface electromyography signal, obtain the set of discharge sequences corresponding to the currently recruited motor units through the same steps as in Step 1. S Current}, Set of action potential waveforms for motor units { Q Current}; Step 5.2: Analyze the set of motion unit waveforms { Q}and{ Q Current Calculate the similarity between the elements in} Q and Q Current Overall similarity MUSI; Step 5.3, for single elements Q Current If it is with all Qthe similarity of which is not more than 0.9, it is considered that the Q Current corresponding to the newly recruited motor unit, the element in { Q Current and its corresponding discharge sequence S Current is added to { Q} and { S New}, respectively; Step 5.4, for single-element Q , if there is an element with a similarity exceeding 0.9 to it Q Current , it is considered that the Q corresponding to the motor unit characteristics has changed, and the element in { Q} is replaced by the element corresponding to the maximum similarity Q Current , and the Q Current corresponding discharge sequence S Current is added to { S Vary}.

[0017] Further, the step 6 further comprises: Step 6.1, taking the 25-second online surface electromyography signal and the discharge sequence of the current motor unit identified in the step 5.3 as source data, generating noise-enhanced training data through the same step as the step 2.1; Step 6.2, completing the initialization, training and verification of the decoding model of the current motor unit by the same step as the step 2.2; Step 6.3, adding the established decoding model to the total decoding model set { M}.

[0018] Further, the step 7 further comprises: Step 7.1, taking the 25-second online surface electromyography signal and the discharge sequence of the current motor unit identified in the step 5.4 as source data, generating noise-enhanced training data through the same step as the step 2.1; Step 7.2, adjusting the decoding model based on the training data: fixing the weight parameters of the three-dimensional convolution layer in the decoding model, and optimizing the weight parameters of the remaining layers; verifying the weight parameter combination of each round through the same step as the step 2.2.3, and retaining the optimal combination as the adjustment result; Step 7.3, replacing the original decoding model corresponding to the total decoding model set { M} with the adjusted decoding model.

[0019] Compared with the prior art, the present application has at least the following beneficial technical effects: 1、 The present application constructs a movement unit discharge sequence decoding model update framework that fuses transfer learning and full initialization learning, designs a data enhancement strategy based on noise injection, realizes the adjustment of a characteristic change movement unit decoding model under a noise robustness constraint and the establishment of a newly recruited movement unit decoding model, and makes the online decoding accuracy increase by 20% and the effective decoding quantity increase by 40%; 2、 The present application designs a discharge sequence lightweight decoding model focusing on the main energy channel of the movement unit, reduces the online training time consumption of the decoding model, meets the real-time requirement of the online process, and makes the required floating point operation number decrease by 2.6 orders of magnitude; 3、 The present application designs a discharge reliability index based on multi-channel electromyographic energy, post-processes the decoded discharge sequence, realizes the removal of redundant discharges under high time resolution, can remove 33% of the redundant discharges, and makes the decoding accuracy increase by 26%.

[0020] The concept, specific structure and generated technical effects of the present application will be further described below in combination with the drawings, so as to fully understand the purpose, features and effects of the present application. BRIEF DESCRIPTION OF DRAWINGS

[0021] Figure 1 is a method flowchart of a preferred embodiment of the present application; Figure 2 is a decoding model enhancement training data generation flowchart of a preferred embodiment of the present application based on noise injection; Figure 3 is a movement unit discharge sequence decoding model architecture based on deep learning of a preferred embodiment of the present application; Figure 4 is a movement unit discharge sequence post-processing flowchart of a preferred embodiment of the present application. DETAILED DESCRIPTION

[0022] The following reference description of the drawings introduces a plurality of preferred embodiments of the present application, so that the technical content thereof is clearer and easier to understand. The present application can be embodied in many different forms of embodiments, and the protection scope of the present application is not limited to the embodiments mentioned herein.

[0023] In the drawings, components of the same structure are denoted by the same reference numerals, and components having similar structures or functions are denoted by similar reference numerals. The size and thickness of each component shown in the drawings are arbitrarily shown, and the present application does not limit the size and thickness of each component. In order to make the drawing clearer, the thickness of some components is appropriately exaggerated in some places in the drawing.

[0024] The embodiment provides a surface electromyography-based motor unit discharge sequence long-time online decoding method, which is used for realizing reliable online decoding of motor unit discharge information in a long-time electromyography collection process in which multiple disturbance factors coexist. The method comprises an offline training stage and an online decoding stage. As shown in the figure, the method comprises the following steps: Figure 1 Step 1, identification of the initial recruitment motor unit in the offline training stage.

[0025] The process obtains a discharge sequence set of the initial recruitment motor unit, S Initial a separation vector set, T Initial a motor unit action potential waveform set, Q Initial and forms a total separation vector set T and a total motor unit action potential waveform set Q correspondingly. The specific steps are as follows: Step 1.1, performing fast independent component analysis algorithm processing on the 25-second electromyography signal collected offline to obtain a candidate independent component set E Offline and a corresponding separation vector set I Candidate . The specific steps are as follows: T Candidate Step 1.1.1, performing 20-500 Hz band-pass filtering and 50 Hz comb filtering on Offline . E Filtering E

[0026] Step 1.1.2, performing 4 times of time domain delay shift with data points of 1, 2, 3 and 4 lengths on E Filtering . E Filtering E Extend .

[0027] Step 1.1.3, performing zero-phase component analysis whitening on the extended electromyography data to obtain whitened electromyography data E Whiten , and recording the first dimension number of E Whiten . n Whiten .

[0028] Step 1.1.4, performing independent component analysis on E ​​​​Whiten Fifty candidate independent components were obtained through 50 iterations of analysis. I Candidate and the corresponding separation vector T Candidate The specific steps for a single analysis are as follows: Step 1.1.4.1: Randomly initialize the length as... n Whiten vector w 1 and w 2.

[0029] Step 1.1.4.2, Calculation w i+1 :

[0030] Step 1.1.4.3 w i+1 Orthogonalization:

[0031] In the formula, B It is a matrix formed by concatenating the separation vectors of all the obtained candidate independent components through channels.

[0032] Step 1.1.4.4 w i+1 Normalization:

[0033] In the formula, norm(·) represents the calculation of the Euclidean norm.

[0034] Step 1.1.4.5, Calculation w i and w i+1 If the dot product is less than 1.0001, the current iteration converges and proceeds to step 1.1.4.6; otherwise, proceed to step 1.1.4.2.

[0035] Step 1.1.4.6: Convergence obtained in this iteration w That is, the separating vector T Candidate Then, the independent components are obtained through the following formula. I Candidate :

[0036] Step 1.2: From the set of candidate independent components { I Candidate The candidate discharge sequence set is extracted from} S Candidate}。Specific steps are as follows: Step 1.2.1, screening I Candidate all local extreme values in.

[0037] Step 1.2.2, K-means clustering is performed on all local extreme values, and the class corresponding to the maximum average value is the discharge class.

[0038] Step 1.2.3, initialization S Candidate is a zero vector with the same length as, and let I Candidate is a zero vector with the same length as, and let S Candidate all discharge class element indexes in are 1, that is, the candidate discharge sequence in binary form corresponding to the discharge state 1 and the rest state 0 is obtained.

[0039] Step 1.3, based on the discharge characteristics of the movement unit, the index is designed to screen the candidate discharge sequence, and the effective discharge sequence set { S Valid} and the corresponding separation vector set { T Valid} are obtained. The screening process uses the following index: Step 1.3.1, the discharge frequency needs to be lower than 35 hertz.

[0040] Step 1.3.2, the discharge interval coefficient of variation c ISI needs to be lower than 0.45. The calculation formula of the discharge interval coefficient of variation is as follows:

[0041] , wherein is the standard deviation of adjacent discharge intervals of the candidate discharge sequence, is the mean of adjacent discharge intervals of the candidate discharge sequence.

[0042] Step 1.3.3, the pulse noise ratio p PNR needs to be higher than 20 decibels. The calculation formula of the pulse noise ratio is as follows:

[0043] , wherein i is the data point index.

[0044] Step 1.4, the repeated items in { S Valid} are removed, and the obtained non-repeated item set is the initial recruitment movement unit discharge sequence set { S Initial}, and the corresponding separation vector TValid constitute a set { T Initial}. The specific steps are as follows: Step 1.4.1, for a single discharge sequence S Valid , judge its consistency with all the rest of the discharge sequences, if the number of synchronous discharges does not exceed 80%, then the two discharge sequences are inconsistent, otherwise they are consistent.

[0045] Step 1.4.2, if there is no discharge sequence consistent with S Valid , add S Valid to the set { S Initial}, add the corresponding T Valid to the set { T Initial}, and remove it from { S Valid}.

[0046] Step 1.4.3, if there is a discharge sequence consistent with S Valid , calculate S Valid and the pulse noise ratio of these discharge sequences according to step 1.3.3. Add the discharge sequence corresponding to the highest pulse noise ratio to the set { S Initial}, add the corresponding T Valid to the set { T Initial}, and remove S Valid and its consistent discharge sequence from { S Valid}.

[0047] Step 1.4.4, loop the above arrangement until { S Valid} is empty.

[0048] Step 1.5, extract the action potential waveform of each initial recruitment motor unit using the pulse trigger superposition averaging method, and constitute an action potential waveform set { Q Initial}.

[0049] This method takes the discharge pulse as the trigger point, and intercepts the electromyographic signal of each discharge pulse at the front and back n MUAP data points, and then averages all the intercepted electromyographic signal segments to obtain the action potential waveform, the calculation formula is as follows:

[0050] wherein, E Filtering is the filtered surface electromyography signal, n MUAP is set to 62, i is the data point index, j is the row channel index, k is the column channel index.

[0051] Step 1.6, initialization of the total separation vector set{ T} is{ T Initial} the total motor unit action potential waveform set{ Q} is{ Q Initial}.

[0052] Step 2, establishment of the initial recruitment motor unit decoding model in the offline training phase.

[0053] A decoding model is established for each initial recruitment motor unit to form a decoding model set{ M Initial} and a total decoding model set{ M} is formed accordingly. The steps for establishing a single motor unit decoding model are as follows: Step 2.1, enhanced training data generation based on noise injection, by adding zero-mean Gaussian white noise which has a greater impact on decoding accuracy to the training data, the diversity of the mapping between electromyography signal and discharge sequence is increased, and the robustness of the decoding model obtained by training to noise increase in the forward inference process is improved.

[0054] As shown in Figure 2 , it is a flow chart of the training data generation process, including the following steps: Step 2.1.1, zero-mean Gaussian white noise addition. Add zero-mean Gaussian white noise signal to the offline surface electromyography signal to obtain noisy electromyography signal, the level of added noise is 5 decibels.

[0055] Step 2.1.2, determination of motor unit main energy channel. Analyze the 10 channels with the largest amplitude of the motor unit action potential waveform, select the smallest rectangle that can envelope these 10 channels as the main energy region, and all channels within the main energy region are the main energy channels of the motor unit. From the perspective of feature description, the discharge behavior of the motor unit is mainly reflected in the channels within its main energy region, therefore, the decoding process selects the electromyography signal on these channels as the input to achieve discharge feature preservation and input space compression at the same time.

[0056] Step 2.1.3, point-by-point sliding window generation of training data. The offline surface electromyography signals at the main energy channel of motor units and their noisy signals are slid by a window length of 0.02 seconds and a step length of 0.0005 seconds to obtain training input data. The label of each training data is determined by the discharge state of the discharge sequence obtained in step 1.4 at the middle time of the window. The discharge state corresponds to a label of 1, and the rest state corresponds to a label of 0.

[0057] Step 2.1.4, inter-class imbalance mitigation. Due to the discharge characteristics of motor units, the number of rest states in the discharge sequence is much larger than that of discharge states, resulting in inter-class imbalance in the training data, which affects the convergence direction of the decoding model training process. Therefore, the data in each round of the training process is composed of training data corresponding to discharge states and randomly selected training data corresponding to rest states, and the proportion of the two is controlled to be about 1:8.

[0058] Step 2.2, establishment of a deep learning-based decoding model, including the following steps: Step 2.2.1, decoding model initialization. As shown in Figure 3 , the decoding model architecture is based on deep learning. The core is a three-dimensional convolution layer for feature extraction and a fully connected layer for feature mapping. Considering that the main feature type contained in the input data is spatio-temporal feature, three-dimensional convolution blocks and fully connected blocks are selected for feature extraction and feature mapping, respectively, so that the decoding model architecture has the characteristics of few parameters and high parallelization, thereby ensuring the feasibility of online updating of the decoding model in the online decoding stage from the aspects of feature fitting capability and training efficiency.

[0059] Step 2.2.2, decoding model training. The training data obtained in step 2.1 are used to train the decoding model. The number of rounds of the training process is 20, the batch size is 300, the loss function is the mean square error function, and the adaptive moment estimation optimizer with an initial learning rate of 0.001 is used for model weight update.

[0060] Step 2.2.3, decoding model verification. The decoding models obtained in each round of step 2.2.2 are verified, and the model corresponding to the highest decoding performance index is retained for use in the online decoding stage. The verification of a single model includes the following steps: Step 2.2.3.1, sliding window of 0.02 seconds and step length of 0.0005 seconds on the 25-second offline surface electromyography signal at the main energy channel of the motor unit determined in step 2.1.2 to obtain sequence electromyography window data.

[0061] Step 2.2.3.2, sequentially inputting the sequence electromyography window data into the decoding model to obtain the predicted discharge sequence.

[0062] Step 2.2.3.3, calculate the consistency between the predicted discharge sequence and the reference discharge sequence obtained in step 1.3, which corresponds to the decoding performance of the model. The consistency is evaluated by the F1-score, and the calculation formula is as follows:

[0063] Where TP is the number of correctly predicted discharges in the predicted discharge sequence, FP is the number of incorrectly predicted discharges in the predicted discharge sequence, and FN is the number of missed predicted discharges in the predicted discharge sequence. When calculating TP and FP, multiple discharges in the predicted discharge sequence with a span of no more than 0.005 seconds are considered as a whole, and only one is counted.

[0064] Step 2.3, initialize the total decoding model set M} = { M Initial}.

[0065] Step 3, real-time decoding of the motor unit discharge sequence in the online decoding stage.

[0066] Every 0.125 seconds, all models in the decoding model set M} are used to decode the 0.125 seconds of online collected electromyographic signals to obtain the predicted results of the real-time discharge sequence. The specific steps are as follows: Step 3.1, the 0.125 seconds of online electromyographic signals at the motor unit main energy channel determined in step 2.1.2 are windowed with a 0.02 second window length and a 0.0005 second step length to obtain a sequence of electromyographic window data.

[0067] Step 3.2, sequentially input the sequence of electromyographic window data into each model in the total decoding model set M} to obtain a set of predicted discharge sequences S online}.

[0068] Step 4, post-processing of the motor unit discharge sequence in the online decoding stage.

[0069] The small step size setting of the decoding process is beneficial to improve the time resolution of the obtained discharge sequence, but at the same time it will cause adjacent items in the sequence input to have high similarity, and then the predicted discharge sequence will have redundant discharges near the expected discharge pulse. According to the energy aggregation characteristics of the motor unit action potential, the surface electromyographic signal at the true discharge time will show a local extremum at the main energy channel of the action potential waveform, therefore, the total energy of the instantaneous electromyographic signal at the main energy channel is used as a discharge reliability indicator to post-process each predicted discharge sequence obtained in step 3.2 to eliminate redundant discharge pulses.

[0070] The post-processing process is shown in Figure 4 , and the specific steps are as follows: Step 4.1, Group all discharge pulses in the original discharge sequence in a proximity from front to back, the interval between the first and last discharge pulse in each group is no more than 0.003 seconds.

[0071] Step 4.2, For a single discharge pulse group, calculate the credibility of each discharge pulse, and keep the discharge pulse with the highest credibility as the only discharge pulse of the current discharge pulse group. The calculation of the credibility is as follows: Step 4.2.1, According to the motor unit action potential waveform corresponding to the discharge sequence being post-processed, select the 3 channels with the highest amplitude as the calculation channels.

[0072] Step 4.2.2, Calculate the average energy of the electromyographic signal at all calculation channels at the discharge pulse moment, which is the credibility of the discharge pulse being processed.

[0073] Step 4.3, Cycle step 4.2 until the discharge pulse is unique for all groups in step 4.1.

[0074] Step 5, Recruit motor units in the online decoding phase.

[0075] Every 175 seconds, analyze the currently recruited motor units based on 25 seconds of online surface electromyography to obtain a new set of motor unit discharge sequences{ S New} and a set of characteristic changing motor unit discharge sequences{ S Vary}, and update the total motor unit action potential set{ Q}. The specific steps are as follows: Step 5.1, For 25 seconds of online surface electromyography, obtain the discharge sequence set corresponding to the currently recruited motor units{ S Current}, the motor unit action potential waveform set{ Q Current} through the same steps as step 1.

[0076] Step 5.2, Analyze the similarity of each element in the motor unit waveform set{ Q} and{ Q Current}. Calculate the overall similarity MUSI of Q and Q Current by the following formula:

[0077] In the formula, j is the row channel index, k is the column channel index.

[0078] Step 5.3, for single element Q Current , if its similarity to all Q is no more than 0.9, it is considered that the Q Current corresponds to a newly recruited motor unit. Then, the Q Current and its corresponding discharge sequence S Current are added to Q and S New respectively.

[0079] Step 5.4, for single element Q , if there exists an element with similarity more than 0.9 to it Q Current , it is considered that the Q corresponds to a motor unit with changed properties. Then, the element in Q is replaced by the element Q Current with the largest similarity, and the Q Current and its corresponding discharge sequence S Current are added to S Vary .

[0080] Step 6, establishment of decoding model for newly recruited motor unit in online decoding phase.

[0081] For all newly recruited motor units in step 5.3, decoding models are established for them through full initialization learning to recognize their real-time discharge sequences in the subsequent online decoding process. The steps to establish decoding model for a single newly recruited motor unit are as follows: Step 6.1, generate noise-enhanced training data through the same steps as step 2.1, with 25 seconds of online sEMG signal and the discharge sequence of the current motor unit identified in step 5.3 as source data.

[0082] Step 6.2, complete the initialization, training and validation of the decoding model for the current motor unit through the same steps as step 2.2.

[0083] Step 6.3, add the established decoding model to the total decoding model set M .

[0084] Step 7, adjustment of decoding model for motor unit with changed properties in online decoding phase.

[0085] For all the characteristic-changing motor units in step 5.4, the weights of their decoding models are adjusted through transfer learning, which improves the performance of the decoding model in the subsequent decoding process while improving the efficiency of the adjustment process by reducing the number of training parameters. The adjustment steps of the decoding model of a single characteristic-changing motor unit are as follows: Step 7.1, using 25 seconds of online surface electromyography signal and the discharge sequence of the current motor unit identified in step 5.4 as source data, generate noise-enhanced training data through the same steps as step 2.1.

[0086] Step 7.2, adjust the decoding model based on the training data, the specific steps are as follows: Step 7.2.1, fix the weight parameters of the three-dimensional convolution layer in the decoding model, and optimize the weight parameters of the remaining layers. The number of rounds, batch size, loss function, and optimizer settings are the same as in step 2.2.2.

[0087] Step 7.2.2, verify the weight parameter combination of each round through the same steps as step 2.2.3, and keep the optimal combination as the adjustment result.

[0088] Step 7.3, replace the original decoding model in the total decoding model set M} with the adjusted decoding model.

[0089] Step 8, evaluate the data.

[0090] The evaluation data is obtained by simulation, with a duration of 750 seconds. The non-steady-state factors considered in the simulation process include motor unit characteristic changes, new motor unit recruitment, and increased noise. The specific steps are as follows: Step 8.1, motor unit action potential waveform generation, including initial waveform extraction and time-dependent waveform change, as follows: Step 8.1.1, the initial action potential waveform of the motor unit is obtained by independent component analysis method and pulse superposition average method from the experimental electromyography signal. The number of motor units is 100, the waveform channel dimension is 8x8, and the channel spacing is 10 millimeters.

[0091] Step 8.1.2, the action potential waveform of the motor unit is changed in amplitude and time course based on the initial waveform. The change coefficient is linearly related to the recruitment duration, and the coefficient extreme value follows a uniform distribution U (1.2, 1.4).

[0092] Step 8.2, motor unit discharge sequence generation, including excitation curve definition, initial discharge sequence generation, and rotation calculation, as follows: Step 8.2.1, the simulation data contains 60 segments of stimulation curves, and each segment of stimulation curve follows a trapezoidal pattern with a duration of 12.5 seconds, including a 1-second pre-rest, a 2-second linear rise, a 6.5-second constant stimulation at a 40% level, a 2-second linear decline, and a 1-second post-rest.

[0093] Step 8.2.2, the recruitment order of motor units is related to the amplitude of their action potentials, and motor units with small amplitudes are recruited first. The initial firing frequency of motor units is 8 Hz, and then increases at a rate of 0.3 Hz per 1% stimulation level until reaching a maximum of 35 Hz. The interval between adjacent discharges follows a Gaussian distribution with a coefficient of variation of 20%.

[0094] Step 8.2.3, group motor units according to the similarity of action potential amplitude, each group contains 2-3 motor units. In the simulation process, only one motor unit in a single group is recruited at any time. The probability of switching the recruitment state within the group increases exponentially with the cumulative recruitment duration.

[0095] Step 8.3, surface electromyography signal generation.

[0096] Convolve the time-varying action potential waveform of all motor units and the discharge sequence considering rotation to obtain noise-free raw electromyography data. Add zero-mean Gaussian white noise to the raw electromyography data to obtain the final simulation electromyography data. The method provided in this embodiment is first run on electromyography data with a signal-to-noise ratio of 30 decibels, the first 25 seconds correspond to the offline training phase, and the last 725 seconds correspond to the online decoding phase, to obtain all decoding models. Then test the decoding performance of the method on electromyography data with a signal-to-noise ratio of 5 decibels to verify the robustness to increased noise.

[0097] Step 9, evaluation results.

[0098] The effectiveness of the method provided in this embodiment is evaluated from three aspects: the online decoding accuracy of the long-time electromyography acquisition process, the computational complexity of the decoding model, and the necessity of post-processing the decoded discharge sequence. The specific steps are as follows: Step 9.1, the online decoding accuracy of the long-time electromyography acquisition process is evaluated by the F1-score between the decoded discharge sequence and the simulation true value, as shown in Table 1. (Note: bold indicates that the decoding result is valid, and — indicates that it cannot be decoded.)

[0099] The first comparison method is a blind source separation-based method (Y. Zheng and P. Zhou, "Adaptive Real-Time Decomposition of Electromyogram During Sustained Muscle Activation: A Simulation Study," in IEEE Transactions on Biomedical Engineering, vol. 69, no. 2, pp. 645-653, 2022, doi: 10.1109 / TBME.2021.3102947), and the second comparison method is a deep learning-based method (A. K. Clarke et al., "Deep Learning for Robust Decomposition of High-Density Surface EMG Signals," in IEEE Transactions on Biomedical Engineering, vol. 68, no. 2, pp. 526-534, 2021, doi: 10.1109 / TBME.2020.3006508).

[0100] Taking whether the F1-score reaches 0.7 as the criterion for judging whether the decoding result is valid, the number of valid decoding results of the embodiment is 21 in total, which is better than 15 of the first comparison method and 10 of the second comparison method. For the initial recruitment motor units, the average decoding accuracy of the embodiment is 0.925 ± 0.062, which is better than 0.839 ± 0.182 of the first comparison method and 0.897 ± 0.050 of the second comparison method. For the newly recruited motor units in the online decoding stage, the average decoding accuracy of the embodiment is 0.820 ± 0.080, which is better than 0.631 ± 0.284 of the first comparison method, and the second comparison method cannot decode the discharge sequences of these motor units. Therefore, for the online decoding of long-time electromyography signals coexisting with three non-steady-state factors of noise increase, motor unit characteristic change, and new motor unit recruitment, the embodiment has higher decoding accuracy and decoding quantity.

[0101] Table 1 Online decoding accuracy of motor unit discharge sequences in long-time electromyography acquisition process

[0102] Step 9.2, the computational complexity of the decoding model is evaluated by the number of floating point operations required for forward inference, and the comparison method is based on deep learning (A. K. Clarke et al., "Deep Learning for Robust Decomposition of High-Density Surface EMG Signals," in IEEE Transactions on Biomedical Engineering, vol. 68, no. 2, pp. 526-534, 2021, doi: 10.1109 / TBME.2020.3006508). Analysis shows that the average number of floating point operations of the decoding model established for all motor units in this embodiment is 170650, and the average number of floating point operations of the comparison method is 75393199. It can be seen that the lightweight characteristics of the decoding model provided in this embodiment are obvious.

[0103] Step 9.3, the effectiveness of the decoding discharge sequence post-processing is evaluated by comparing the online decoding accuracy of this embodiment with and without the post-processing step. As shown in Table 2, the comparison results of decoding accuracy. Analysis shows that the proportion of redundant discharges in the direct decoding results of the decoding model is 33%, the average decoding accuracy of the present method when performing post-processing is 0.868 ± 0.089, and the average decoding accuracy when not performing post-processing is 0.689 ± 0.145. It can be seen that the execution of the post-processing step is obviously improved for decoding accuracy, which at the same time shows the necessity and effectiveness of the decoding discharge sequence post-processing.

[0104] Table 2 Effect of decoding discharge sequence post-processing on decoding accuracy

[0105] The preferred embodiments of the present application are described in detail above. It should be understood that those skilled in the art can make many modifications and changes without creative labor based on the concept of the present application. Therefore, any technical solutions obtained by logical analysis, reasoning or limited experiments based on the prior art according to the concept of the present application shall be within the protection scope determined by the claims.

Claims

1. A surface myoelectric-based motor unit discharge sequence length online decoding method, characterized in that, The method comprises the following steps: Step 1, identification of initial recruited motor units in offline training phase, comprising: Step 1.1, on-line acquisition of 25 seconds of EMG signal E Offline performing a fast independent component analysis algorithm to obtain a candidate independent component set { I Candidate} and a corresponding separating vector set { T Candidate}; Step 1.2, extracting a candidate discharge sequence set { I Candidate} from a candidate independent component set { S Candidate}; Step 1.3, screening the candidate discharge sequences based on the discharge characteristics of motor units to obtain an effective discharge sequence set S Valid} and a corresponding set of separated vectors T Valid} Step 1.

4. Remove the duplicates in the set of effective discharge sequences { S Valid} to obtain the set of discharge sequences of initial recruitment motor units { S Initial}, the corresponding separation vectors T Valid constitute the set { T Initial}. Step 1.5, the action potential waveform of each initial recruitment motor unit is extracted by using the pulse trigger superimposed average method to form an action potential waveform set Q Initial} Step 1.6, initialize total separation vector set{ T} to{ T Initial}, total motion unit action potential waveform set{ Q} to{ Q Initial}; Step 2, establishment of initial recruited motor unit decoding model in offline training phase, comprising the following steps: Step 2.1, noise injection-based enhanced training data generation; Step 2.2, establishment of deep learning-based decoding model; Step 2.3, initialize the total set of decoding models M} for { M Initial}; Step 3, real-time decoding of the online decoding stage motor unit discharge sequence, comprising: every 0.125 seconds, decoding the 0.125 second electromyogram signal obtained online using all models in the decoding model set { M} to obtain the prediction result of the real-time discharge sequence; Step 4, post-processing of motor unit firing sequence in online decoding phase: post-processing each predicted firing sequence obtained in the step 3 to remove redundant firing pulses; Step 5, identification of recruited motor units in online decoding phase; Step 6, establishment of new recruited motor unit decoding model in online decoding phase; Step 7, adjustment of characteristic change motor unit decoding model in online decoding phase.

2. The sEMG-based motor unit discharge sequence length online decoding method of claim 1, wherein, The step 1.1 further comprises: Step 1.1.1, to E Offline Band-pass filtering with a bandwidth of 20-500 Hz and 50-Hz comb filtering are performed on the filtered electromyography signal E Filtering ; Step 1.1.2, to E Filtering perform 4 time domain delay shifts with 1, 2, 3, 4 data points in length, respectively, to E Filtering and channel concatenation with the 4 time domain shifted signals to obtain extended electromyography data E Extend ; Step 1.1.3, performing zero-phase component analysis whitening on the extended electromyography data to obtain whitened electromyography data E Whiten , record E Whiten The first dimension number of n Whiten ; Step 1.1.4, on E Whiten 50 iterations of the analysis resulted in 50 candidate independent components I Candidate and corresponding separation vectors T Candidate .

3. The sEMG-based motor unit discharge sequence length online decoding method of claim 1, wherein, The step 1.2 further comprises: Step 1.2.1, screening I Candidate all local extrema in Step 1.2.2, K-means bimodality clustering of all local extrema, and a class corresponding to the maximum average value is a firing class; Step 1.2.3, Initialization S Candidate For each discharge class element index i, let I Candidate For each discharge class element index i, let S Candidate The candidate discharge sequence is obtained by setting the value of all discharge class element indices to 1, i.e., the binary form of the candidate discharge sequence is 1 for discharge state and 0 for rest state.

4. The sEMG-based motor unit discharge sequence length online decoding method of claim 1, wherein, The step 1.5 further comprises: With the discharge pulse as the trigger point, intercept each discharge pulse before and after each n MUAP The data points of the electromyogram, and then the average of all electromyogram segments obtained by intercepting the action potential waveform.

5. The sEMG-based motor unit discharge sequence length online decoding method of claim 1, wherein, The step 2.1 further comprises: Step 2.1.1, zero-mean Gaussian white noise addition; adding a zero-mean Gaussian white noise signal to the offline surface electromyogram to obtain a noisy electromyogram, and the level of added noise is 5 decibels; Step 2.1.2, determination of motor unit main energy channel; analyzing the 10 channels with the maximum amplitude of the motor unit action potential waveform, selecting the smallest rectangle that can envelope the 10 channels as the main energy region, and all channels within the main energy region are motor unit main energy channels; Step 2.1.3, point-by-point sliding window generation of training data; the offline surface electromyogram and its noisy signal at the motor unit main energy channel are slid with a window length of 0.02 seconds and a step length of 0.0005 seconds to obtain training input data, and the label of each training data is determined by the firing state of the firing sequence obtained in the step 1.4 at the middle time of the window, with a firing state corresponding to a label of 1 and a rest state corresponding to a label of 0; Step 2.1.4, inter-class imbalance mitigation; the data of a single round of training process is composed of all training data corresponding to the firing state and randomly extracted training data corresponding to the rest state, and the proportion of the two is controlled to be about 1:

8.

6. The sEMG-based motor unit discharge sequence length online decoding method of claim 1, wherein, The step 2.2 further comprises: Step 2.2.1, initialization of the decoding model architecture, the core of which is a three-dimensional convolution layer for feature extraction and a fully connected layer for feature mapping; Step 2.2.2, training the decoding model using the training data obtained in the step 2.1, and using the adaptive moment estimation optimizer to update the model weights; Step 2.2.3, validating each decoding model obtained in each round of the training process, and retaining the model corresponding to the highest decoding performance index, comprising: Sliding the 25-second offline surface electromyogram at the motor unit main energy channel determined in the step 2.1.2 with a window length of 0.02 seconds and a step length of 0.0005 seconds to obtain sequence electromyogram window data; Inputting the sequence electromyogram window data into the decoding model in sequence to obtain a predicted firing sequence; Calculating the consistency between the predicted firing sequence and the reference firing sequence obtained in the step 1.3, which corresponds to the decoding performance of the model, and the consistency is evaluated by F1-score, and the calculation formula is as follows: Wherein, TP is the number of correctly predicted discharges in the predicted discharge sequence, FP is the number of incorrectly predicted discharges in the predicted discharge sequence, and FN is the number of missed predicted discharges in the predicted discharge sequence.

7. The sEMG-based motor unit discharge sequence length online decoding method of claim 1, wherein, The step 4 further comprises: Step 4.1, grouping all discharge pulses in the original discharge sequence in a front-to-back order according to their proximity, wherein the interval between the first and last discharge pulse in each group is no more than 0.003 seconds; Step 4.2, for each single discharge pulse group, calculating the confidence of each discharge pulse, and retaining the discharge pulse with the highest confidence as the only discharge pulse of the current discharge pulse group; Step 4.3, repeating the step 4.2 until the discharge pulse uniqueness is achieved for all groups in the step 4.

1.

8. The sEMG-based motor unit discharge sequence length online decoding method of claim 6, wherein, The step 5 further comprises: Step 5.1, for 25 seconds of online surface electromyography signals, obtain a set of discharge sequences corresponding to the currently recruited motor units by the same procedure as Step 1 S Current} a set of motor unit action potential waveforms Q Current} Step 5.2, analyzing the set of motion unit waveforms { Q} and { Q Current} the similarity of the elements of Q with Q Current the overall similarity MUSI of Step 5.3, for single elements Q Current If it is with all Q If the similarity is no more than 0.9, then it is considered that... Q Current Corresponding to the newly recruited movement units, this Q Current and its corresponding discharge sequence S Current Add to { Q }and{ S New }middle; Step 5.4, for single elements Q , if there is an element with similarity over 0.9 Q Current , it is considered that the Q corresponding motion unit characteristics have changed, then replace the element in Q with the element corresponding to the maximum similarity Q Current , and add the Q Current corresponding discharge sequence S Current to S Vary .

9. The sEMG-based motor unit discharge sequence length online decoding method of claim 8, wherein, The step 6 further comprises: Step 6.1, taking the 25-second online surface electromyography signal and the discharge sequence of the current motor unit identified in the step 5.3 as source data, and generating noise-enhanced training data through the same steps as the step 2.1; Step 6.2, completing the initialization, training and verification of the decoding model of the current motor unit through the same steps as the step 2.2; Step 6.

3. Add the established decoding model to the total decoding model set M}.

10. The sEMG-based motor unit discharge sequence length online decoding method of claim 8, wherein, The step 7 further comprises: Step 7.1, taking the 25-second online surface electromyography signal and the discharge sequence of the current motor unit identified in the step 5.4 as source data, and generating noise-enhanced training data through the same steps as the step 2.1; Step 7.2, adjusting the decoding model based on the training data: fixing the weight parameters of the three-dimensional convolution layer in the decoding model, and optimizing the weight parameters of the remaining layers; verifying the weight parameter combination of each round through the same steps as the step 2.2.3, and retaining the optimal combination as the adjustment result; Step 7.

3. Replace the corresponding original decoding model in the total decoding model set with the adjusted decoding model. M}

Images