Movement disorder treatment effect evaluation system based on HD-sEMG and neuronal analysis

Abstract

The present invention belongs to the technical field of medical diagnosis and treatment evaluation, specifically a movement disorder treatment effect evaluation system based on HD-sEMG and neuron analysis. The present invention induces pathological tremor symptoms in the hands of patients with movement disorders, collects HD-sEMG data of the forearm, and analyzes and evaluates the improvement effect of pathological tremor symptoms in patients with movement disorders before and after DBS treatment; the system includes: a data acquisition platform, which records surface electromyography data of forearm flexors and extensors; a data preprocessing module, which removes noise and baseline drift by filtering to improve signal quality; a motor unit decomposition and identification module, which decomposes and identifies the motor units of forearm flexors and extensors; a micro-neuron characteristic analysis module, which analyzes neuronal discharge characteristics and reflects changes in neural control; by comparing the changes in neuronal characteristics before and after deep brain stimulation treatment, the improvement of tremor symptoms is evaluated, and doctors are assisted in optimizing treatment plans, thereby improving treatment effects and improving patients' quality of life.

Description

Technical Field

[0001] The present invention belongs to the technical field of medical diagnosis and treatment evaluation, and specifically relates to a system for evaluating the therapeutic effect of deep brain electrical stimulation on pathological tremor symptoms in the hands of patients with movement disorders. Background Art

[0002] Movement disorders are a class of neurological conditions that affect motor control, typically resulting in abnormalities in both voluntary and involuntary movements. Tremor is a common movement disorder characterized by involuntary, rhythmic tremors in a particular part of the body. Parkinson's disease is one of the most common pathological tremor disorders. It manifests as a resting tremor due to impaired motor control caused by degeneration of the dopamine system. Dystonia can also cause pathological tremor, which is often exacerbated by specific postures or movements. Pathological tremor significantly impacts patients' daily functioning, social interactions, psychological well-being, and overall quality of life. Therefore, early diagnosis and comprehensive assessment strategies are crucial to mitigate the negative impact of pathological tremor on patients' lives. Deep brain stimulation (DBS) is a treatment method that delivers electrical stimulation to specific areas of the brain using electrodes implanted in these areas, thereby modulating abnormal neural activity and improving patients' pathological tremor. DBS has shown significant efficacy in improving movement disorder symptoms, particularly for pathological tremor that is refractory to medication.

[0003] Currently, the diagnosis and assessment of pathological tremor mostly rely on neurological examination results and tremor assessment scales. However, these traditional assessment techniques and methods may have some limitations and disadvantages: (1) Invasiveness and discomfort: Neurological examinations may require patients to undergo relatively invasive examinations. For example, detailed neurological examinations may require some patients to endure discomfort. For example, acupuncture electromyography (EMG) examinations. (2) Limited to specific environments: Some patients with movement disorders need to undergo magnetic resonance imaging (MRI) examinations to further investigate the patient's intracranial conditions. However, MRI examinations require patients to remain still. Patients with head or limb tremors are prone to imaging artifacts, which can affect the assessment of intracranial tremor-related conditions. Some auxiliary examinations, such as magnetic resonance imaging (MRI), require specific facilities and environments, which limits their application in outpatient or mobile diagnostic environments. (3) Subjectivity: Tremor assessment scales may be affected by the patient's subjective report and the doctor's subjective evaluation. There are certain issues of subjectivity and measurement inaccuracy, and there is a lack of detailed analysis of the microscopic level of neuromuscular activity.

[0004] To overcome these shortcomings, assessment methods based on surface electromyography (sEMG) technology have emerged, aiming to provide more objective and accurate data. sEMG records the electrical signals generated by neurons during muscle activity, which are transmitted through muscle fibers and subcutaneous tissue to the skin surface. Pathological tremor is a typical manifestation of the peripheral nervous system in patients with movement disorders, manifesting as sustained muscle excitation and contraction. Therefore, direct recording of muscle electrical activity can be used to assess the severity of pathological tremor symptoms in patients with movement disorders. sEMG technology, which records muscle electrical activity by placing electrodes on the skin surface, is non-invasive, more comfortable, and safer for patients. Tremor measurements can be repeated at different time points and under different conditions, regardless of the environment. Furthermore, objective physiological signals avoid the limitations of subjective assessments and, by capturing subtle changes in neuromuscular activity, provide detailed analysis of motor unit activity at the microscopic level.

[0005] Currently, researchers have proposed various methods for evaluating changes in muscle activity characteristics before and after DBS from an sEMG perspective, focusing primarily on macroscopic EMG analysis. Commonly used methods include spectrum analysis, signal morphology analysis, and nonlinear dynamics analysis. These methods often use monopolar or bipolar electrodes to evaluate the effects of different DBS parameter settings on pathological tremor symptoms. Existing solutions mainly include:

[0006] 1. A method for analyzing the differences in sEMG signal characteristics of patients with Parkinson's disease pathological tremor under different DBS conditions. This method uses bipolar electrodes and analyzes the electromyographic characteristics of 13 patients with Parkinson's disease pathological tremor based on spectrum, signal morphology, and nonlinear dynamics. The therapeutic effect of DBS is evaluated by quantifying the electromyographic signal characteristics under different DBS parameter settings (pulse amplitude, frequency, and width) [1].

[0007] 2. A method for evaluating the effects of DBS on Parkinson's disease tremor based on principal component tracking. This method uses monopolar or bipolar electrodes (depending on the patient) and analyzes sEMG data from 13 patients with Parkinson's disease tremor in the DBS-on and DBS-off states and 13 healthy controls based on correlation dimension and recurrence rate characteristics. This method evaluates the therapeutic effect by quantifying the effect of DBS on neuromuscular function in patients with Parkinson's disease tremor[2].

[0008] These existing technical solutions mainly analyze the tremor improvement effect before and after DBS by using the amplitude, frequency, and nonlinear characteristics of sEMG signals. However, they have the following shortcomings:

[0009] Limitation 1: Limitations of high-density surface EMG acquisition. Existing methods often use monopolar or bipolar electrodes, which have limited signal coverage and the number of channels, often insufficient to fully capture the EMG activity of the target muscle group. Due to the limited number of electrodes, ensuring the comprehensiveness and accuracy of the data is difficult, and some key EMG signal features may be missed, affecting the accuracy of the analysis results.

[0010] Flaw 2: Limitations of macroscopic EMG signal characteristics. Existing methods primarily focus on macroscopic EMG signal characteristics, analyzing parameters such as amplitude and frequency based on spectrum, signal morphology, and nonlinear dynamics. However, they lack in-depth understanding of the microscopic characteristics and patterns of neuronal firing. This limitation limits our comprehensive understanding of changes in neuromuscular activity before and after DBS stimulation. Summary of the Invention

[0011] The purpose of the present invention is to provide a movement disorder treatment effect evaluation system based on HD-sEMG and neuronal analysis to evaluate the improvement of tremor symptoms before and after DBS treatment, and to help doctors formulate and optimize treatment plans.

[0012] The movement disorder treatment efficacy evaluation system provided by this invention, based on HD-sEMG and neuron analysis, induces pathological tremor symptoms in the hands of movement disorder patients, collects HD surface electromyography data from the forearm, and analyzes its micro-neuron characteristics, thereby accurately evaluating the improvement of pathological tremor symptoms in movement disorder patients before and after DBS treatment. Specifically, it includes: a data acquisition module for pathological tremor data in movement disorder patients' hands, a data preprocessing module, a motor unit decomposition and identification module, and a micro-neuron characteristic analysis module; wherein:

[0013] (1) The data acquisition platform is a high-definition sEMG acquisition device used to record 256 channels of surface electromyography data from the forearm flexor and extensor muscles of patients with movement disorders. This high-density data acquisition capability allows for more precise analysis of muscle activity patterns and motor unit characteristics.

[0014] Specifically, four electrode arrays were used to collect HD-sEMG signals from a single forearm. Each array consisted of 64 gel-coated circular electrodes, aligned roughly with the muscle fibers of the patient with movement disorders. Two electrode arrays were placed on the flexor muscles, one above the radial extensor carpi and the other above the digitorum extensor muscles. Two electrode arrays were placed on the extensor muscles, one above the radial extensor carpi and the other above the digitorum extensor muscles. The collected EMG signals were stored on an internal SD card in each acquisition device to ensure data integrity and stability.

[0015] (2) The data preprocessing module preprocesses the raw EMG signals through preliminary filtering, LMS filtering, baseline drift correction, and signal normalization. Preliminary filtering removes power frequency interference and preserves the effective frequency range. The LMS filter adaptively adjusts the coefficient to eliminate dynamic noise. Baseline drift correction removes chronic drift through polynomial fitting. Finally, normalization ensures signal amplitude consistency, thereby improving signal clarity and accuracy.

[0016] Specifically:

[0017] (1) Preliminary filtering: Apply a low-pass filter and a band-pass filter to the raw EMG signal (10 Hz to 500 Hz) to remove power frequency interference and retain the effective signal frequency range. The low-pass filter is used to remove high-frequency noise above the signal frequency range. The band-pass filter is used to retain the main components of the EMG signal and remove irrelevant low-frequency and high-frequency noise.

[0018] (2) LMS filter: After the initial filtering, the least mean square (LMS) adaptive filter is used to adjust the filter coefficients in real time to minimize the error between the output signal and the desired signal, thereby eliminating dynamic noise and improving signal quality. The error signal is determined by the difference between the desired signal and the filter output. The core of LMS is to gradually adjust the filter coefficients ω based on the gradient descent method. i (n) to minimize the mean square value of the error signal. By adjusting the objective function J(n) for the filter coefficient ω i (n) Calculate the partial derivative, obtain the gradient, and use the gradient descent method to update the filter coefficients.

[0019] (3) Baseline drift correction: Baseline drift is a low-frequency component caused by changes in the electrode-skin interface impedance or other chronic interference. Polynomial fitting techniques can be used to correct for nonlinear trends in the signal and remove baseline drift. The given signal is fitted using the least squares method to obtain the polynomial coefficients, which are then subtracted from the original signal.

[0020] (4) Signal normalization: The signal amplitude is normalized to the minimum-maximum range and the signal amplitude is normalized to the range of [0, 1] to ensure the scale consistency of the signal and reduce the impact of the impedance difference of the electrode-skin interface.

[0021] The above-mentioned preprocessing module can effectively cope with complex signal environments, improve signal clarity and accuracy while removing noise and drift, and provide reliable input for subsequent signal analysis and processing.

[0022] (3) The motor unit decomposition and identification module is used to perform motor unit decomposition on the flexor muscle groups and extensor muscle groups of patients with movement disorders, including:

[0023] (1) Demeaning, radial basis function transformation, and principal component analysis submodules; by demeaning, the baseline offset of the signal is eliminated to ensure the accuracy of the analysis; using radial basis function transformation, the signal is mapped to a nonlinear high-dimensional space to better capture complex signal characteristics; applying principal component analysis to reduce the dimensionality of high-dimensional data, extract the main features, reduce the influence of noise, and improve the efficiency and accuracy of subsequent analysis;

[0024] (2) The motor unit extraction submodule uses the FastICA algorithm combined with the q-non-Gaussian metric to decompose the EMG signal after dimensionality reduction. The q-non-Gaussian metric provides better robustness by introducing Tsallis entropy and adjusting parameters, and can effectively handle noise and outliers. The weight vector is optimized through fixed-point iteration and updated to maximize the q-non-Gaussian metric. In each iteration, Gram-Schmidt orthogonalization is performed to ensure the orthogonality between independent components and maintain the validity of the decomposition results.

[0025] (3) An effective motor unit identification submodule uses Morlet wavelet to decompose the electromyographic signal based on the decomposed independent components. By detecting the significant local extreme values ​​in the wavelet coefficients and binarizing them, a 01 value time series of the motor unit discharge moment is generated.

[0026] (4) Screening and post-processing submodule, i.e., screening and post-processing the obtained discharge time sequence; including extracting the main frequency components in the frequency range of 4 Hz to 35 Hz by Fourier transform, and retaining the motor unit sources that meet this frequency range; identifying the motor unit source pairs with synchronous discharge by calculating the cross-correlation function and constructing the synchronization matrix; for the identified synchronization sources, comprehensively scoring them through signal energy and stability evaluation, and finally retaining the synchronization source with the highest score, thereby ensuring the quality and stability of the selected motor unit source.

[0027] Further:

[0028] (1) Demeaning, radial basis function transformation and principal component analysis submodules; among which:

[0029] The main purpose of deaveraging is to eliminate baseline drift in the signal and make the signal data centralized at zero. Specifically, the mean of each channel is calculated and subtracted from each data point of each channel to obtain the deaveraged sEMG signal.

[0030] The radial basis function (RBF) transform maps the de-averaged sEMG signal to a high-dimensional feature space, thereby capturing the complex nonlinear characteristics of the signal, making the nonlinear relationship between independent components more significant, and better reversely solving the discharge pulse signal of motor neurons, so that subsequent processing (such as principal component analysis or motor unit decomposition) can more effectively identify and separate different motor units.

[0031] Principal component analysis involves normalizing the RBF-transformed data so that each feature has zero mean and unit variance to ensure that different features have the same scale. The covariance matrix of the normalized data matrix is ​​calculated and regularization terms are added to the covariance matrix to improve its stability and reduce the effects of noise.

[0032] Perform eigenvalue decomposition on the regularized covariance matrix to obtain eigenvalues ​​and eigenvectors. Select the first k eigenvectors with the largest eigenvalues ​​as principal components, which form the new feature space. Use the selected principal components to project the standardized data to achieve dimensionality reduction.

[0033] (2) The motor unit extraction submodule uses the FastICA algorithm combined with the q-non-Gaussian metric to perform reverse calculation on the reduced dimensionality electromyographic signals to separate the individual motor units.

[0034] Traditional FastICA algorithms typically use negentropy as a measure of non-Gaussianity. Negentropy uses the difference in entropy to measure the non-Gaussianity of a signal. However, negentropy is sensitive to outliers and noise, which may lead to a decrease in the accuracy of independent component extraction in the presence of noise. The q-non-Gaussianity metric provides a more generalized non-Gaussianity measure by introducing Tsallis entropy. Compared with traditional non-Gaussianity measures (such as those based on negentropy), the q-non-Gaussianity metric introduces an adjustable parameter q, which can be adjusted to control robustness to outliers. Using a fixed-point iterative algorithm, the weight vector is gradually updated to maximize the q-non-Gaussian metric. In each iteration, to ensure orthogonality between the independent components, the updated weight vector is orthogonalized using the Gram-Schmidt orthogonalization process.

[0035] (3) An effective motor unit identification submodule uses the Morlet wavelet function to decompose the electromyographic signal into wavelet coefficients of different scales based on the decomposed independent components, so as to extract significant features for the detection of discharge events.

[0036] Specifically, significant local maxima are detected in the wavelet coefficients as potential discharge events. Local maxima can be determined by comparing the values ​​within the wavelet coefficient neighborhood. Local extreme values ​​in the wavelet coefficients are binarized, with the portion exceeding the threshold marked as 1 (discharge moment) and the remainder marked as 0 (background noise). The binarized results of all discharge events are integrated to generate a time series of 0-1 values ​​that mark the discharge moments of the motor unit.

[0037] (4) Screening and post-processing submodule, i.e., screening and post-processing the obtained discharge time sequence. This includes extracting the main frequency components between 4 Hz and 35 Hz by Fourier transform and retaining the motor unit sources that meet this frequency range; identifying the motor unit source pairs with synchronous discharge by calculating the cross-correlation function and constructing the synchronization matrix; for the identified synchronization sources, comprehensively scoring them through signal energy and stability evaluation, and finally retaining the synchronization source with the highest score, thereby ensuring the quality and stability of the selected motor unit sources.

[0038] Specifically, a Fourier transform is applied to the discharge time sequence of each motor unit to obtain the spectrum and calculate the power spectral density of the signal. The main frequency components in the frequency range of 4Hz to 35Hz are selected, and the motor unit sources within this range are retained. For each pair of motor unit sources A and B, the cross-correlation function between their discharge sequences is calculated to evaluate the degree of synchronization between the two signal sequences. The synchronization matrix is ​​constructed using the synchronization calculation results of all motor unit sources. Statistical tests are performed on the synchronization matrix to identify the synchronous discharge phenomenon between different motor unit sources. Motor unit source pairs with significant synchronization are screened out and marked as sources with significant synchronous discharge phenomenon. For each synchronous source, the quality of each synchronous source is evaluated by calculating the signal energy and signal stability (standard deviation). A comprehensive score is calculated for each synchronous source, and the synchronous source with the highest score is selected based on the comprehensive score. The synchronous source with the highest score is selected as the motor unit source that is finally retained, which effectively selects the most representative source from multiple synchronous sources. This method ensures that the retained sources have the best quality and stability, providing a reliable data basis for subsequent analysis.

[0039] (IV) The micro-neuron characteristic analysis module includes analyzing the discharge characteristics and patterns of neurons at the micro level by calculating the instantaneous discharge rate of motor units and statistically analyzing the discharge patterns of motor units. By analyzing the differential changes in the micro-level motion characteristics of electromyographic signals before and after DBS, the improvement effect of tremor symptoms in patients with movement disorders before and after deep brain stimulation treatment can be more accurately assessed. Specifically, it includes:

[0040] (1) Motion unit instantaneous discharge rate analysis submodule. The calculation formula of the instantaneous discharge rate between the discharge pulses of the motion unit is:

[0041]

[0042] Where Δn is the number of sampling points between two consecutive discharges, and fs is the sampling rate. Based on this data, the motor unit decomposition results and instantaneous discharge rates of different muscle groups in patients with movement disorders were analyzed before surgery, under DBS-on, and DBS-off conditions. This was used to assess changes in neuronal activity in patients with movement disorders under different conditions. These changes directly reflect the modulatory effects of DBS on neuronal activity. Understanding and adjusting motor unit discharge patterns through microscopic analysis can achieve the goal of improving symptoms through DBS regulation.

[0043] (2) Motor unit discharge pattern analysis submodule. By counting the instantaneous discharge rates of all motor units and plotting them into a histogram, the intervention effect of DBS on the discharge pattern is demonstrated. The statistics and plotting of the histogram can help researchers and clinicians observe the changes in neural activity after DBS interference and adjust pathological neural activity to a pattern that is more consistent with normal physiological rhythms. By optimizing these discharge patterns, tremor symptoms can be effectively improved.

[0044] Compared with the prior art, the present invention has the following advantages:

[0045] (1) Comprehensiveness and accuracy of high-density surface electromyography acquisition

[0046] Compared with existing movement disorder assessment technologies such as invasive lumbar puncture to obtain cerebrospinal fluid samples, it is non-invasive and more comfortable and safe for patients. Compared with the simple application of tremor assessment scales for manual empirical naked eye observation scoring, the present invention significantly improves the comprehensiveness, accuracy and objectivity of data collection by using 256-channel HD-sEMG acquisition equipment. The forearm flexor and extensor muscle groups closely related to hand movements are selected as target muscle groups, and high-density surface electromyography data are collected by inducing tremor symptoms, ensuring the breadth and detail of data coverage. This enables the present invention to obtain richer and more reliable electromyography data, which has obvious advantages over traditional monopolar or bipolar electrode technology.

[0047] (2) Accuracy of microscopic neuron characteristic analysis

[0048] By calculating the instantaneous discharge rate and discharge pattern of motor units, the discharge characteristics and laws of neurons are deeply analyzed, and the condition and treatment effect of pathological tremors in patients with movement disorders are evaluated at the microscopic level. This precise microscopic analysis method enables the present invention to not only accurately capture the dynamic changes in neuronal activity, but also systematically count and analyze the discharge patterns of motor units. The results show that before and after DBS treatment, neuronal activity has changed significantly at the microscopic level. DBS can reduce the rhythmic discharge of motor units, reduce the synchronization of discharges between motor units, and inhibit the activity of pathological motor units. This precise microscopic analysis not only improves the accuracy of diagnosis and treatment, but also provides a scientific basis for the formulation of personalized treatment plans, helping doctors to more effectively formulate and adjust treatment plans, thereby improving treatment effects. BRIEF DESCRIPTION OF THE DRAWINGS

[0049] Figure 1 Illustration of the HD-sEMG data acquisition platform.

[0050] Figure 2 Figure 1 is an illustration of the data acquisition software.

[0051] Figure 3 is the motor unit discharge pulse and its instantaneous discharge rate.

[0052] Figure 4 The distribution of motor unit discharge patterns. DETAILED DESCRIPTION

[0053] The present invention is further described below with reference to the accompanying drawings. A system for evaluating the therapeutic effect of deep brain stimulation on hand movement disorder symptoms in patients with movement disorders based on HD-sEMG specifically comprises:

[0054] 1: HD-sEMG data acquisition system

[0055] HD-sEMG signals were acquired using a portable EMG collector, the Sessantaquattro system (OT Bioelectronica, Torino, Italy). This device transmits signals wirelessly, eliminating the need for connection to fixed equipment. This allows the experimenter to move the device to coordinate with the subject's movements, facilitating experimental execution. The EMG sampling rate was 2000 Hz. This EMG collector can acquire up to 64 channels of EMG data. Four EMG collectors were used simultaneously, each powered by an internal battery. To reduce the impedance between the skin and the electrodes, the skin of the forearm on the acquisition side was cleaned with a scrub and then rinsed with alcohol before data acquisition. High-density EMG signals were acquired from the unilateral forearm using four 5×13 HD electrode arrays (Adhesive Matrix ELSCH064NM2 Pin Out, OT Bioelectronica, Torino, Italy). Each electrode array consisted of 64 gel-coated circular electrodes, one of which had an electrode channel missing from one corner. The center-to-center distance between adjacent electrode channels is 8 mm, and the electrode diameter is 3 mm. The electrode arrays are roughly aligned with the muscle fibers of the test subject. Two electrode arrays are placed on the flexor muscles, and two electrode arrays are placed on the extensor muscles. A moistened wristband is wrapped around the forearm as a reference electrode.

[0056] Data collection platforms such as Figure 1 As shown in the figure, the directed lines in the figure represent the data signal transmission flow. The platform uses a laptop computer for host computer control. A portable router establishes a network hotspot and connects the host computer and four EMG data collectors via the TCP protocol. Within the local area network, the laptop computer and EMG data collectors can exchange signals, enabling functions such as parameter reading and configuration, starting and ending data acquisition, and online data reading and offline storage. Due to the time delay between devices, and the high latency requirements for EMG signal analysis (on the order of tens of milliseconds), this platform uses a hardware synchronization module to synchronize the signals. A mechanical button serves as the input of the synchronization circuit, and the output of the hardware circuit is connected to the auxiliary synchronization signal input ports of the four HD-sEMG devices. Pressing the button twice in succession starts acquisition, and pressing it three times in succession ends it. Through these steps, the HD-sEMG data acquisition system can comprehensively and accurately record high-density EMG signals from the forearm flexor and extensor muscles of patients with movement disorders, providing high-quality data for subsequent analysis of micro-neuronal characteristics.

[0057] Data acquisition is controlled by software. The software control block diagram is as follows: Figure 2As shown, the system includes functions such as sensor module connection (displaying and managing the connection status of the HD-sEMG sensor modules), start and stop acquisition (providing control buttons to start and end acquisition), acquisition timing (displaying the duration of data readout and storage), single-module sEMG data display (displaying real-time data from a specific HD-sEMG module), and data saving (managing file save settings for acquired data). Below the timing display, the connection status of each of the four HD-sEMG sensor modules is displayed. When the host computer initiates an acquisition command, four indicator lights turn green and red, respectively, indicating successful or failed connection to the corresponding module. These indicators allow you to check the status of the corresponding module in the event of a connection failure, ensuring smooth data acquisition. The two input boxes below, from left to right, are used to enter the file name prefix for saved data files and the number of the sEMG module displaying data, respectively. The "Start" and "End" buttons at the bottom are used to initiate acquisition start and end commands. Due to the large volume of EMG data, EMG signals are only sent to the host computer when the signal is displayed. During synchronous acquisition, EMG signals are stored on the internal SD card of each collector to ensure data integrity and stability. Through the above functions, the data acquisition software can effectively manage and control the HD-sEMG data acquisition process, monitor the connection status of the sensor module in real time, and provide precise time reference and data storage functions to ensure the accuracy and integrity of the collected data.

[0058] 2. Experimental Paradigm

[0059] This study recruited eight patients with movement disorders (4 males and 4 females; age range 50-66 years, mean age 61.6±5.3 years). Data were collected in a dedicated EMI / EMC-shielded room in the Functional Neurosurgery Department of a Shanghai hospital. The experiment included three conditions: preoperative, postoperative DBS-on (DBS-on), and postoperative DBS-off (DBS-off). Participants sat in a comfortable chair with armrests, resting their arm on the side with the most severe tremor symptoms, maintaining a relaxed and stable posture. Before data collection, the experimenter instructed the participants to perform a countdown task with intervals of two to distract them and increase their mental workload, thereby inducing or exacerbating tremor symptoms. Data collection began when the patient's tremor became visibly noticeable or after the countdown task lasted for one minute. In each experimental condition, participants continued the countdown task for at least 45 seconds while HD-sEMG signals were collected. Multiple trials were collected for each experimental condition. The number of trials varied among participants due to factors such as discomfort caused by discontinuation of medication, electrode contact issues, and environmental interference. The average number of trials under the three experimental conditions was 10.6 ± 4.3, 7.7 ± 3.4, and 7.3 ± 3.2, respectively. This method allows for comprehensive collection of HD-sEMG data from patients with Parkinson's disease (PD) under different experimental conditions, enabling analysis of the effects of DBS on hand tremor symptoms and providing a reliable data foundation for subsequent data processing and analysis.

[0060] 3: Data preprocessing module

[0061] (1) Preliminary filtering: Apply a low-pass filter and a band-pass filter to the raw EMG signal (10 Hz to 500 Hz) to remove power frequency interference and retain the effective signal frequency range. The low-pass filter is used to remove high-frequency noise above the signal frequency range, with a cutoff frequency of 500 Hz. The band-pass filter is used to retain the main components of the EMG signal, with a frequency range of 10 Hz to 500 Hz, and remove irrelevant low-frequency and high-frequency noise.

[0062] (2) LMS filter: After preliminary filtering, a least mean square (LMS) adaptive filter is used to adjust the filter coefficients in real time to minimize the error between the output signal and the desired signal, thereby eliminating dynamic noise and improving signal quality. The LMS filter is based on the finite impulse response (FIR) filter structure, and its output y(n) is expressed as:

[0063]

[0064] Where y(n) is the output of the filter, M is the filter order, ω i (n) is the filter coefficient, and x(n) is the input signal. The error signal is determined by the difference between the desired signal and the filter output:

[0065] e(n)=d(n)-y(n)

[0066] Where: e(n) is the error signal, d(n) is the error signal. The core of the LMS algorithm is to gradually adjust the filter coefficient ω based on the gradient descent method. i (n) to minimize the mean square value of the error signal:

[0067] J(n)=E[e 2 (n)] = E[(d(n) - y(n)) 2 ]

[0068] Where E[·] represents the expected value. By using the objective function J(n) to calculate the filter coefficient ω i (n) Take the partial derivative and get the gradient:

[0069]

[0070] And update the filter coefficients using gradient descent:

[0071] ω i (n+1)=ω i (n)+2μe(n)·x(ni)

[0072] Where μ is the step size factor (learning rate), which determines the magnitude of the coefficient update. e(n) x(ni) is the product of the error signal and the input signal, which is used to adjust the coefficients. Choosing a step size factor μ that is too large may cause the algorithm to diverge, while choosing a step size factor that is too small may slow the convergence speed. Therefore, the step size factor should be selected to maximize the convergence speed while ensuring stability. To ensure algorithm stability, the step size factor must meet the following conditions:

[0073]

[0074] λ max is the largest eigenvalue of the autocorrelation matrix of the input signal.

[0075] (3) Baseline drift correction: Baseline drift is a low-frequency component caused by changes in the electrode-skin interface impedance or other chronic interference. Polynomial fitting technology can be used to correct the nonlinear trend of the signal and remove the baseline drift. The baseline drift of the signal can be expressed by a low-order polynomial:

[0076] b(t)=a0+a1t+a2t 2 +…+a n t n

[0077] Where b(t) represents the drift signal, a0, a1,…, a nis the fitting coefficient. Use the least squares method to fit the given signal x(t) and get the polynomial coefficient a i . This drift signal is then subtracted from the original signal:

[0078] x corrected =x(t)-b(t)

[0079] (4) Signal normalization: The signal amplitude is normalized to the minimum-maximum range and the signal amplitude is normalized to the range of [0, 1] to ensure the scale consistency of the signal and reduce the impact of the impedance difference of the electrode-skin interface.

[0080] The above-mentioned preprocessing module can effectively cope with complex signal environments, improve signal clarity and accuracy while removing noise and drift, and provide reliable input for subsequent signal analysis and processing.

[0081] 4: Motion unit decomposition module

[0082] The patient's flexor and extensor muscle groups are decomposed into motor units. Specifically,

[0083] (1) Demeaning, radial basis function transformation, and principal component analysis submodules. The main purpose of demeaning is to eliminate baseline drift in the signal so that the signal data is concentrated at zero. The mean of the signal of each channel is calculated and subtracted from each data point. The original sEMG signal X is an N×T matrix, where N is the number of signal channels and T is the number of time points. Each row represents the time series signal of a channel:

[0084]

[0085] For each channel (each row of data), calculate its mean:

[0086]

[0087] in, is the mean of the jth channel, x jt is the signal value of the jth channel at time point t. Subtract the mean of each channel from each data point of each channel:

[0088]

[0089] Where x′ jtis the signal value after de-averaging. The de-averaged signals of all channels are combined into a new matrix X′. sEMG is a linear combination of the motor neuron input pulse signal and its delay superimposed with noise. In order to better reversely solve the discharge pulse signal of the motor neuron, it is necessary to apply the radial basis function (RBF) transformation to the de-averaged sEMG signal and map it to a high-dimensional feature space, thereby capturing the complex nonlinear characteristics of the signal and making the nonlinear relationship between independent components more significant. Using Gaussian RBF, its form is:

[0090]

[0091] μ is the center of the RBF, σ is the width (standard deviation) of the Gaussian function, and controls the width and shape of the function. Perform RBF transformation on the de-averaged signal to obtain the RBF feature matrix:

[0092]

[0093] Where, M is the number of RBFs, φ m (x t ) is the value of the mth RBF at time point t. The RBF transform can capture the nonlinear characteristics of the original signal and map the signal to a higher-dimensional space, so that subsequent processing (such as principal component analysis or motion unit decomposition) can more effectively identify and separate different motion units.

[0094] The data after RBF transformation is standardized so that each feature has zero mean and unit variance to ensure that different features have the same scale:

[0095] X″=X′ RBF -mean(X′ RBF )

[0096] Among them, X″ is the standardized data matrix, mean(X′ RBF ) is the mean of each column feature. Calculate the covariance matrix C of the standardized data matrix:

[0097]

[0098] Where N is the number of samples. Add a regularization term to the covariance matrix C to improve its stability and reduce the impact of noise:

[0099] C reg =C+λI

[0100] Where λ is the regularization parameter and I is the identity matrix. The regularization term λI helps to avoid problems that occur during eigenvalue decomposition, such as eigenvalues ​​that are too small or the matrix is ​​close to singularity.reg Perform eigenvalue decomposition to obtain eigenvalues ​​and eigenvectors:

[0101] C reg =VΛ

[0102] Where V is the eigenvector matrix, Λ is a diagonal matrix, and the elements on the diagonal are eigenvalues. The first k eigenvectors with the largest eigenvalues ​​are selected as principal components, and these eigenvectors form the new feature space:

[0103] W=[v1v2…v k ]

[0104] Among them, v i is the eigenvector corresponding to the first k eigenvalues. Project the normalized data using the selected principal components to complete the dimensionality reduction:

[0105] X PCA =X″W

[0106] Among them, X PCA is the data matrix after dimensionality reduction, and W is the eigenvector matrix containing the first k principal components.

[0107] (2) Using the FastICA algorithm submodule combined with the q-non-Gaussianity metric, the dimensionality-reduced electromyographic signal is reversely solved to separate the independent motor units. The traditional FastICA algorithm usually uses negative entropy as a measure of non-Gaussianity. Negative entropy uses the difference in entropy to measure the non-Gaussianity of a signal. However, negative entropy is more sensitive to outliers and noise, which may lead to a decrease in the accuracy of independent component extraction in the presence of noise. The q-non-Gaussianity metric provides a more generalized non-Gaussianity metric by introducing Tsallis entropy. Compared with traditional non-Gaussianity metrics (such as metrics based on negative entropy), the q-non-Gaussianity metric introduces an adjustable parameter q. By adjusting this parameter, the robustness to outliers can be controlled. The q-non-Gaussianity metric formula for the signal is:

[0108]

[0109] w is the weight vector of the independent components to be estimated, x is the observation data matrix, and q is the adjustment parameter. The weight vector w is gradually updated through the fixed point iteration algorithm to maximize the q-non-Gaussianity measure. The iteration formula is:

[0110] w new =Ε[xg(w T x)]-Ε[g′(w T x)]w

[0111] Where g(·) is a nonlinear function and w is the weight vector of the current iteration. In each iteration, to ensure the orthogonality between the independent components, the updated weight vector is orthogonalized through the Gram-Schmidt orthogonalization process.

[0112] (3) The effective motor unit identification submodule decomposes the electromyographic signal into wavelet coefficients of different scales based on the Morlet wavelet function to extract significant features for the detection of discharge events:

[0113]

[0114] Among them, W(a,b) is the wavelet coefficient at scale a and position b, is a wavelet basis function. Significant local maxima are detected in the wavelet coefficients as potential discharge events. Local maxima can be determined by comparing values ​​within the wavelet coefficient neighborhood. Local extreme values ​​in the wavelet coefficients are binarized, with those exceeding the threshold marked as 1 (discharge event) and the remainder as 0 (background noise). The binarized results of all discharge events are integrated to generate a time series of 0-1 values ​​that mark the discharge moments of motor units.

[0115] (4) Screening and post-processing submodule, which screens and post-processes the obtained discharge time sequence. Fourier transform is applied to the discharge time sequence of each motor unit to obtain the spectrum and calculate the power spectrum density of the signal. The main frequency components in the frequency range of 4Hz to 35Hz are selected, and the motor unit sources within this range are retained. This frequency range is concentrated in the main frequency band of the electromyographic signal, which can effectively extract the most representative discharge frequency while reducing the interference of high-frequency noise. For each pair of motor unit sources A and B, the cross-correlation function between their discharge sequences is calculated to evaluate the degree of synchronization of the two signal sequences. The synchronization matrix is ​​constructed based on the synchronization calculation results of all motor unit sources. The synchronization matrix is ​​statistically tested to identify the synchronous discharge phenomenon between different motor unit sources. The motor unit source pairs with significant synchronization are screened out, and these sources are marked as sources with significant synchronous discharge phenomenon. For each synchronization source, the quality of each synchronization source is evaluated by calculating the signal energy and signal stability (standard deviation). The signal energy formula is:

[0116]

[0117] Among them, A i (t) is the signal value of the ith synchronization source at time t, and T is the total number of time points. The signal stability formula is:

[0118]

[0119] Among them A iis the mean value of the ith synchronization source signal. Calculate the comprehensive score of each synchronization source and select the synchronization source with the highest score based on the comprehensive score. The comprehensive score S of each synchronization source is:

[0120]

[0121] The highest-scoring synchronization source is selected as the final retained motor unit source, which effectively selects the most representative source from multiple synchronization sources. This method ensures that the retained sources have the best quality and stability, providing a reliable data foundation for subsequent analysis.

[0122] 5. Micro-neuron Analysis Module

[0123] (1) Motion unit instantaneous discharge rate analysis submodule. The calculation formula of the instantaneous discharge rate between the discharge pulses of the motion unit is:

[0124]

[0125] Among them, Δn is the number of sampling points between two consecutive discharge moments, and fs is the sampling rate. Figure 3 The results of motor unit decomposition and transient discharge rates of different muscle groups in representative subjects under three experimental conditions: preoperative, DBS-on, and DBS-off. The results of motor unit discharge pulse trains showed that under preoperative and DBS-off conditions, motor units were highly rhythmic, with double or multiple continuous spike discharges. The transient discharge rate before surgery was concentrated around a stable low frequency of approximately 5 Hz, corresponding to the typical PD tremor frequency. Some motor units showed higher transient discharge rates, which corresponded to double or triple peaks in the tremor action potential. Under DBS-on conditions, the transient discharge rates of motor units were distributed around 10 Hz, corresponding to the typical physiological tremor frequency. Under DBS-off conditions, the transient discharge rates of motor units fluctuated within the range of 3-12 Hz. Compared with preoperative conditions, more motor units spread to a higher discharge rate range, and more motor units showed more active short-term continuous discharges.

[0126] (2) Motor unit discharge pattern analysis submodule: The analysis of the motor unit discharge pattern is described by counting the instantaneous discharge rates of all motor units and drawing them into a histogram. Figure 4The distribution of motor unit discharge patterns in eight patients is presented. Preoperatively, the discharge pattern clustered in the pathological tremor band (approximately 5 Hz) with relatively high, sharp peaks, indicating that a large number of motor units were firing synchronously within this frequency range. This high-peak discharge pattern corresponds to the typical tremor symptoms of dyskinesia and reflects pathological neural activity. With DBS on, the discharge pattern exhibited two diffuse peaks at higher frequencies (approximately 10 Hz and 15 Hz, respectively), with flatter and broader peak profiles, indicating that motor unit discharges were more dispersed and less concentrated at these frequencies. These frequencies correspond to the normal physiological tremor band, demonstrating that DBS effectively modulates the pathological motor unit discharge pattern to a more normal physiological state. With DBS off, some motor units re-clusted in the symptomatic tremor band (approximately 5 Hz), with a sharp peak, but the peak value was lower than that of the preoperative state. Other motor units clustered in the normal physiological tremor band (approximately 12 Hz), indicating that some motor unit discharge patterns remained within the higher frequency range. This suggests that after DBS was turned off, some motor units resumed their pathological firing rhythms, but the overall condition was still better than before the surgery. Specifically, DBS improves tremor symptoms by interfering with pathological neural activity and adjusting the firing patterns of motor units to patterns more consistent with normal physiological rhythms.

[0127] References:

[0128] [1]. Rissanen SM, Ruonala V, Pekkonen E, et al. Signal features of surface electromyography in advanced Parkinson's disease during different settings of deep brain stimulation [J]. Clinical Neurophysiology, 2015, 126(12): 2290-2298.

[0129] [2].Rissanen SM, M,Tarvainen MP,et al.Analysis of EMG andacceleration signals for quantifying the effects of deep brain stimulation inParkinson'sdisease[J].IEEE transactions on biomedical engineering,2011,58(9):2545-2553.

Claims

1. A movement disorder treatment effect evaluation system based on HD-sEMG and neuronal analysis, characterized in that: By inducing pathological tremor symptoms in the hands of patients with movement disorders, collecting HD-sEMG data from the forearm, and analyzing its micro-neuron characteristics, the improvement effect of pathological tremor symptoms in patients with movement disorders before and after DBS treatment can be accurately evaluated. The system specifically includes a pathological tremor data acquisition platform for patients with movement disorders, a data preprocessing module, a motor unit decomposition and identification module, and a micro-neuron characteristic analysis module. Among them: (1) The data acquisition platform is an HD-sEMG acquisition device used to record 256-channel surface electromyography data of the flexor and extensor muscles of the forearm of patients with movement disorders; (2) The data preprocessing module preprocesses the raw electromyographic signal through preliminary filtering, LMS filter, baseline drift correction and signal normalization; the preliminary filtering is used to remove power frequency interference and retain the effective frequency range, the LMS filter adaptively adjusts the coefficient to eliminate dynamic noise, and the baseline drift correction removes chronic drift through polynomial fitting; normalization ensures the consistency of signal amplitude, thereby improving the clarity and accuracy of the signal; (3) The motor unit decomposition and identification module is used to perform motor unit decomposition on the patient's flexor muscle group and extensor muscle group respectively; specifically including: (1) Demeaning, radial basis function transformation, and principal component analysis submodules; by demeaning, the baseline offset of the signal is eliminated to ensure the accuracy of the analysis; using radial basis function transformation, the signal is mapped to a nonlinear high-dimensional space to better capture complex signal characteristics; applying principal component analysis to reduce the dimensionality of high-dimensional data, extract the main features, reduce the influence of noise, and improve the efficiency and accuracy of subsequent analysis; (2) The motor unit extraction submodule uses the FastICA algorithm combined with the q-non-Gaussian metric to decompose the EMG signal after dimensionality reduction. The q-non-Gaussian metric provides better robustness by introducing Tsallis entropy and adjusting parameters, and can effectively handle noise and outliers. The weight vector is optimized through fixed-point iteration and updated to maximize the q-non-Gaussian metric. In each iteration, Gram-Schmidt orthogonalization is performed to ensure the orthogonality between independent components and maintain the validity of the decomposition results. (3) An effective motor unit identification submodule uses Morlet wavelet to decompose the electromyographic signal based on the decomposed independent components. By detecting the significant local extreme values ​​in the wavelet coefficients and binarizing them, a 01 value time series of the motor unit discharge moment is generated. (4) Screening and post-processing submodule, i.e., screening and post-processing the obtained discharge time sequence; including extracting the main frequency components in the frequency range of 4 Hz to 35 Hz by Fourier transform, retaining the motor unit sources that meet this frequency range; identifying the motor unit source pairs with synchronous discharge by calculating the cross-correlation function and constructing the synchronization matrix; for the identified synchronization sources, comprehensively scoring them through signal energy and stability evaluation, and finally retaining the synchronization source with the highest score, thereby ensuring the quality and stability of the selected motor unit sources; (4) The micro-neuron characteristic analysis module includes analyzing the discharge characteristics and patterns of neurons at the micro level by calculating the instantaneous discharge rate of motor units and counting the discharge patterns of motor units; and accurately evaluating the improvement effect of tremor symptoms in patients with movement disorders before and after deep brain stimulation treatment by analyzing the changes in micro-neuron characteristics of electromyographic signals before and after deep brain stimulation treatment.

2. The evaluation system according to claim 1, wherein: The data acquisition platform uses a four-electrode array to collect high-density electromyographic signals from a single forearm. Each electrode array consists of 64 gel-coated circular electrodes. The electrode arrays are aligned with the patient's muscle fibers. Two electrode arrays are placed on the flexor muscles, one above the flexor carpi and the other above the superficial flexor digitorum. Two electrode arrays are placed on the extensor muscles, one above the radial extensor carpi and the extensor digitorum. The collected electromyographic signals are stored on SD cards built into each acquisition device to ensure data integrity and stability.

3. The evaluation system according to claim 1, wherein: In the data preprocessing module: (1) Preliminary filtering: Apply low-pass filter and band-pass filter to perform preliminary filtering on the original electromyographic signal to remove power frequency interference and retain the effective signal frequency range; low-pass filter is used to remove high-frequency noise higher than the signal frequency range; band-pass filter is used to retain the main components of the electromyographic signal and remove irrelevant low-frequency and high-frequency noise; (2) LMS filter: After preliminary filtering, the least mean square (LMS) adaptive filter is used to adjust the filter coefficients in real time to minimize the error between the output signal and the desired signal, thereby eliminating dynamic noise and improving signal quality. The error signal is determined by the difference between the desired signal and the filter output. LMS is based on the gradient descent method to gradually adjust the filter coefficients ω. i (n) to minimize the mean square value of the error signal; by adjusting the filter coefficient ω for the objective function J(n) i (n) Calculate the partial derivative, obtain the gradient, and use the gradient descent method to update the filter coefficients; (3) Baseline drift correction: Baseline drift is a low-frequency component caused by changes in the electrode-skin interface impedance or other chronic interference. The nonlinear trend of the signal is corrected and the baseline drift is removed by polynomial fitting technology. The given signal is fitted using the least squares method to obtain the polynomial coefficients, and then the drift signal is subtracted from the original signal. (4) Signal normalization: The signal amplitude is normalized to the minimum-maximum range and the signal amplitude is normalized to the range of [0, 1] to ensure the scale consistency of the signal and reduce the impact of the impedance difference of the electrode-skin interface.

4. The evaluation system according to claim 3, wherein: In the motor unit decomposition module, the specific process of the screening and post-processing submodule is as follows: applying Fourier transform to the discharge time sequence of each motor unit to obtain the spectrum and calculate the power spectral density of the signal; selecting the main frequency components between 4 Hz and 35 Hz and retaining the motor unit sources within this range; for each pair of motor unit sources A and B, calculating the cross-correlation function between their discharge sequences to evaluate the degree of synchronization between the two signal sequences; The synchronization matrix is ​​constructed by calculating the synchronization degree of all motion unit sources; Perform statistical tests on the synchronization matrix to identify the synchronous discharge phenomenon between different motor unit sources; screen out pairs of motor unit sources with significant synchronization, and mark these sources as sources with significant synchronous discharge phenomenon; for each synchronization source, evaluate the quality of each synchronization source by calculating signal energy and signal stability; calculate the comprehensive score of each synchronization source, and select the synchronization source with the highest score based on the comprehensive score; select the synchronization source with the highest score as the final retained motor unit source to provide a reliable data basis for subsequent analysis.

5. The evaluation system according to claim 4, wherein: The micro-neuron characteristic analysis module specifically includes: (1) Motor unit instantaneous discharge rate analysis submodule; the calculation formula of the instantaneous discharge rate between motor unit discharge pulses is: Where Δn is the number of sampling points between two consecutive discharge moments, and fs is the sampling rate. Based on this, the motor unit decomposition results and instantaneous discharge rates of different muscle groups in patients with movement disorders under three conditions (pre-surgery, DBS-on, and DBS-off) were analyzed to assess changes in neuronal activity in patients with movement disorders under different conditions. Through microscopic analysis, the discharge patterns of motor units can be understood and adjusted to achieve the goal of improving symptoms through DBS regulation. (2) Motor unit discharge pattern analysis submodule: By counting the instantaneous discharge rates of all motor units and plotting them into a histogram, the intervention effect of DBS on the discharge pattern is demonstrated; the statistics and plotting of the histogram can help researchers and clinicians observe the changes in neural activity after DBS interference and adjust pathological neural activity to a pattern that is more consistent with normal physiological rhythms; by optimizing these discharge patterns, tremor symptoms can be effectively improved.

Images