Classical and Bayesian inference of the weighted-exponential distribution with an application to insurance data

1.   Introduction

Despite current advances in disease prevention, acute phase treatment and rehabilitation, the global burden of stroke is expected to rise in the future ^[1], and stroke remains the second leading cause of death in the world ^[2]. Restoration of motor function after a stroke is crucial to a patient's recovery. However, due to the large inter-individual variation, it is difficult to accurately predict the outcome of a patient's motor function rehabilitation based on clinical medicine and imaging assessments alone ^[3]. Therefore, there is still a need to identify reliable and inexpensive biomarkers to add predictive information for these patients. Biomarkers can provide clinically useful information when planning a patient's individualized recovery. These biomarkers can also be used for patient selection and stratification, as well as for trials investigating rehabilitation interventions initiated early after a stroke. Ongoing multicenter trials, combined with exercise biomarkers, can help apply them to routine clinical practice. As evidence of neurovascular uncoupling in acute ischemic stroke accumulates, neurophysiological biomarkers appear to be increasingly relevant in predicting the outcomes ^[4].

Since 1976, Cohen et al. ^[5] observed a large increase in low-frequency activity and a decrease in high-frequency activity in EEG signals in stroke patients. Many researchers have begun to use different quantitative electroencephalography (QEEG) methods to explore the biological indicators of stroke patients ^[6]. Finnigan et al. ^[7] studied ischaemic stroke patients by collecting static EEG and calculated QEEG increments for analysis with national institute of health stroke scale (NIHSS) scores, showing that the ratio of δ-band to α-band power and the relative frequency of the α-band were significantly correlated with NIHSS scores. Meanwhile, Finnigan et al. ^[8] studied the consistency with other clinical assessment scales by using the QEEG method in patients with (sub) acute ischaemic stroke. Numerous results on hundreds of patients have shown that the QEEG method can objectively inform clinical management, particularly for the prediction of motor function recovery outcomes. In addition to this, there is a range of studies on stroke that provide strong evidence for the validity of other QEEG indices, such as EEG signal power ^[9], intracortical inhibition ^[10], relative power ^[11], brain symmetry index and transverse coefficient ^[12,13].

As a direct performer of limb movement, especially in stroke patients, muscles are directly affected by motor dysfunction and neglect of the EMG signal will result in a one-sided assessment tool. Since the relationship between EEG signals and EMG signals during exercise was first identified ^[14], communication between brain and muscle has become an important topic of research. To explore the exchange of information between the brain and muscles during exercise, researchers have used various methods of coupled strength calculation and applied these methods to the field of stroke motor dysfunction research. Currently, brain-muscle coupling studies are mainly based on traditional frequency domain coherence methods to explore the exchange of information between the brain drive as well as the muscle response. Hallett et al. in 2001 found that the strength of cerebral muscle coupling in chronic stroke patients was found mainly in the contralateral brain region on the motor side; also, the cerebral muscle coupling on the affected side is lower than that on the healthy side during movement ^[15]. This suggests that the recovery of motor function in stroke may be related to the increased strength of coupling between opposite brain regions and muscles.

However, coherence analysis can only represent the strength of coupling between individual channels and individual channels and is powerless when applied to multi-channel EEG acquisition devices. In order to shift the perspective from individual channels to multiple channels, an analytical approach based on graph theory ^[16] was proposed. When applying the graph theory approach to the EEG signal, what is obtained is a network in which the nodes are the individual channel electrodes and the edges are the strength of the coupling between the individual channels. Graph theory methods provide a powerful approach to the quantitative study of network topologies ^[17], and network characteristics such as shortest path length, clustering coefficients, efficiency and small-world properties can be used to describe networks ^[18]. Coupling strength calculations as a potential generative model of network connectivity can provide insight into the mechanisms by which brain networks transform and process information ^[19]. Fallani et al. ^[20] constructed a functional brain network using imagined EEG signals from unilateral hand movements of stroke patients, and the analysis found that hand movements affected by disease during exercise had significantly reduced small-world properties and local efficiency in the beta band compared to healthy hands. Vecchio et al. ^[21] conducted a study of the correlation between abnormalities in the functional brain network and early clinical outcomes in stroke in order to find potential predictors of functional recovery to address or correct rehabilitation programs, by performing EEG recordings as well as clinical assessments of 139 stroke patients on three scales: the NIHSS, the Barthel index and the action research arm test (ARAT). The results found a significant correlation between the Small World Index in the prediction of recovery of motor function.

Therefore, it is necessary to combine EEG signals and EMG signals to establish a brain-muscle function network through a Graph theory, to study the brain system and muscle system as a whole system and to extract features from the brain-muscle function network for classification by certain feature extraction means, in order to provide a new perspective in quantifying the functional state of the human nervous system during stroke motor rehabilitation, as well as to provide new insights and methods for the assessment of stroke motor function.

In previous study, we concluded that cortical muscle function coupling (FCMC) between motor cortex and muscle can be considered as an assessment mechanism for motor function rehabilitation ^[22], and that FCMC is significantly lower in stroke patients compared to healthy individuals. The reduction in FCMC suggests that damage to the lesioned hemisphere may lead to discontinuities in information transmission in the sensory-motor system. The study also determined changes in paired brain symmetry index (PBSI) and corticomuscular coupling following stroke. In order to improve the understanding of the above-mentioned topics, we have proposed for the first time a symmetry index combining corticomuscular coupling and graph theory to study the recovery process in stroke patients. First, we extracted a number of features based on cortical muscle coupling and graph theory; then we calculated the importance of these features in a classification prediction model; and finally, clinical classification predictions were made based on this importance ranking in combination with multiple classification models to assess the practical significance of these features.

2.   Materials and methods

2.1.   Overview

We analyzed EEG and EMG data during the experiment in 34 subjects. This study consists of the following steps (shown as a schematic shown in Figure 1): first, after preprocessing the EEG and EMG, the coherence and DTW distance are calculated and a weighted adjacency matrix is created using the coherence; then, a threshold value is determined and the weighted matrix is converted into a binary matrix; finally, the brain network degree symmetry index (BNDSI) and corticomuscular coupling symmetry index (CMCSI) are calculated for the final classification analysis.

Figure 1.  Overview of the study workflow. (a) EEG and EMG acquisition, EEG and EMG electrode placement positions; (b) EEG and EMG preprocess; (c) DTW calculation to obtain DTW-EEG and DTW-EMG; (d) Calculation of coherence and creation of adjacency matrix to obtain BNDSI and CMCSI; (e) Classification model performance evaluation.

DownLoad: Full-Size Img PowerPoint

2.2.   Experimental paradigm

We recruited 18 stroke patients (10 females and 8 males, mean age 60.5 years) and 16 healthy volunteers from Hangzhou Mingzhou Naokang Rehabilitation Hospital. Table 1 provides statistical information for each subject (healthy subjects were recorded as controls, patients with grade 6 in Brunnstrom's staging were recorded as group B6, and patients with grade 5 were recorded as group B5). Subjects with other neurological or musculoskeletal disorders, and/or taking illegal drugs and head trauma were excluded. All subjects provided written informed consent in accordance with the Declaration of Helsinki, and the study was approved by the Ethics Committee of Hangzhou Mingzhou Naokang Rehabilitation Hospital. The experimental procedure is shown in Figure 2. The experiment was conducted in a closed room. There was a five-second preparation time before the experiment started, and subjects in groups B5 and B6 were told to place their affected hand on the contralateral waist, and subjects in the healthy group performed two sets of experiments (left hand on the right waist and right hand on the left waist, respectively). After hearing the voice prompt, the subject was required to slowly raise the hand to the chest and hold the position. The whole movement lasted for 5 seconds and the movement was ended after hearing the voice prompt again. After a 30-second rest, the experiment was repeated, with each subject repeating the procedure three times.

Table 1.  Demographics of the participants.

Group	Number of subjects	Age (years)	Gender (M/F)	Affected hand (left/right)
B5	9	59 (± 0.92SE)	(5/4)	(4/5)
B6	9	63 (± 0.75SE)	(3/6)	(5/4)
Control	16	57 (± 0.29SE)	(8/8)	none

---

 | Show Table

DownLoad: CSV

Figure 2.  Experimental setup and paradigm. Flow of the experimental task: each experiment takes 5 seconds, and the interval between each experiment is 30 seconds.

DownLoad: Full-Size Img PowerPoint

2.3.   Data acquisition and pre-processing

The subjects' EEG and EMG data were collected simultaneously during the entire experiment. A 64-channel wireless EEG system (NeuSen.W64, Neuracle, China) was used to collect EEG data with a sampling frequency of 1000 Hz. Referring to the international 10–20 system, we selected 10 channels (F3, F4, FC3, FC4, C3, C4, CP3, CP4, P3 and P4) from the 64 EEG channels for measurement. Simultaneously, a channel EMG electrode (Delsys Inc., Natick, MA, USA) was used to collect EMG data with a sampling frequency of 1000 Hz. Electrodes were placed on the biceps brachii (BB) and front deltoid (FD) muscles of the subject's upper extremity being tested for measurement. The EMG electrode positions are shown in Figure 1(a). The impedance of wireless EEG was kept below 5 kΩ by injecting conductive paste before data acquisition. Furthermore, on the muscle where the EMG electrode was to be applied, each participant shaved their body hair and cleaned the skin with alcohol.

The EELAB toolbox ^[23] in MATLAB was used to pre-process the raw EEG data recorded during the experiment by downscaling the raw signal to 250 Hz. In addition, Cleanline (an EEGLAB plug-in) is used to remove sinusoidal artifacts in the scalp channel that are not effectively removed by trap filtering. The combination of AC power line fluctuations, equipment power and fluorescent lamps may produce sinusoidal artifacts ^[24]. Using the Infomax algorithm in EEGLAB, other artifacts, such as ECG, EMG and EYE, are removed by technical independent component analysis (ICA). Then the artifact-free signal is denoised using wavelet transform. This is an acceptable combination of ICA and wavelet denoising for removing noise from EEG signals ^[25]. The denoised EEG signal and EMG signal was band-pass filtered using the FIR digital filter in EEGLAB to intercept the EEG and EMG signal from 12–45 Hz. Subsequently, the EMG signal was denoised using a combination of empirical mode decomposition (EMD) and wavelet thresholding ^[26]. The db3 wavelet is chosen as the wavelet basis function for denoising. The uniform threshold formula proposed by Donoho and Jonestone ^[27] was used to determine the initial wavelet threshold. We can change the denoising performance by adjusting the two parameters in a modified two-parameter threshold function ^[28] to adjust the threshold value.

During the data interception phase, the length of the collected signals varied due to the different response times of each subject at the beginning and end of the reminder. We used EMG signal power as the interception criterion. First, the resting-state EMG signal was collected from each subject, and then the power of the EMG signal under each window was calculated using a sliding window of variable length from the first 0.5 s to the last 0.5 s of the experimentally collected 5s EMG signal. Finally, the signal range with EMG signal power greater than or equal to 110% of the resting state EMG power and with the longest length is considered as the target range. This range was used as the signal interception range to intercept the corresponding EEG signals and EMG signals as the experimental data for the subsequent study.

3.   Method

3.1.   Dynamic time warping

In this paper, we intercepted data according to different reaction times of each subject, and the signal range with the longest length was used as the interception range when the EMG power was greater than or equal to 110% of the resting-state EMG power. To demonstrate that this data interception method can effectively intercept the EEG signal and EMG signal in the motion state, we select the EEG signal and EMG signal in the resting state of each subject as the reference signal and calculate the distance between the signals using the DTW algorithm. The DTW algorithm is based on dynamic programming and solves the problem of matching data templates of different lengths.

Given these two-time series: ${x}_{l} = {x}_{1}, {x}_{2}, \dots, {x}_{i}, \dots, {x}_{n}$ and ${y}_{k} = {y}_{1}, {y}_{2}, \dots, {y}_{i}, \dots, {y}_{m}$, their lengths are n and m. To align these two series, we need to construct an $n*m$ matrix D. The matrix elements $D(i, j)$ represent the amplitude distance $d({x}_{i}, {y}_{j})$ between the ${x}_{i}$ and ${y}_{j}$ two points, denoted as ${d}_{k}$. The algorithm can be understood as finding the smallest of the paths from the lower left corner to the upper right corner in this matrix. The path is the cumulative addition of the grid points passed. In terms of continuity and monotonicity, each grid point $D(i, j)$, has only three forward directions: $\left(i+1, j\right), \left(i, j+1\right), \left(i+1, j+1\right)$. Our goal is to make the following paths with minimal regularization cost:

where ${d}_{k}$ is the value represented by the matrix lattice on the path, K is the number of matrix lattices through which the path passes, and DTW is the similarity between the control and experimental groups. This path can be obtained by a dynamic programming algorithm.

3.2.   Graph analysis

To quantify the strength of the coupling between the EEG and EMG channels, we calculated the coherence between the channels ^[29].

where, ${P}_{xy}$ represents the cross-spectral density of signals x and y, ${P}_{xx}$ and ${P}_{yy}$ represents the self-spectral density of x and y, i.e., the power spectral density. Then, in order to binarize the weighted adjacency matrix, a cost efficiency (Ce) threshold (th) is used:

Here D is the network density, defined as the ratio of the actual number of edges to the number of all possible edges; Eg represents the global efficiency:

where ${L}_{i, j}$ represents the shortest path length between nodes i and j, and N is the number of nodes in the graph.

Thresholds that maximize cost-effectiveness are used to binarize the weighted matrices. We set separate thresholds for each adjacency matrix, resulting in a directed binary adjacency matrix for further graph analysis.

3.3.   Symmetry index

Van Putten ^[30] was the first to study the paired brain symmetry index (PBSI) in a stroke study and found an association with concomitant acute ischemic stroke (national institutes of health stroke scale, NIHSS). Another study reported that it was sensitive to the brain pathophysiology of subcutaneous stroke compared with health treatment ^[31]. We derive the brain network degree symmetry index (BNDSI) based on graph theory and PBSI. The formula is:

N is the total number of paired channels, and i refer to a pair of channels in N pairs, ${r}_{i}$ and ${l}_{i}$ are the right and left hemispheres of paired channels, respectively.

Similarly, compared to the combined EEG and EMG approach, the above symmetry index has certain drawbacks, ignoring the fact that motor processes are a joint cooperation between brain and muscle. Therefore, based on the PBSI, we proposed the corticomuscular coupled symmetry index (CMCSI):

CMCSI represents the coherence between EEG and EMG channels, N is the total number of pairs of EEG channels, i refers to one of the N pairs of channels, ${R}_{i}$ and ${L}_{i}$ refer to the degree of coupling between the right hemisphere channel and the muscle and the left hemisphere channel and the muscle, respectively, of this pair.

3.4.   Classification model

In this paper, random forest ^[32] is selected to calculate CMCSI, BNDSI, DTW-EEG and DTW-EMG feature importance, and the classification feature data is constructed by combining features in different ways according to the feature importance ranking. Subsequently, these feature data were used as input to the SVM algorithm ^[33] based on 10-fold cross-validation to train the classifier for final recognition. The SVM classifier was implemented by the LIBSVM toolbox in MATLAB, and the parameters were set to default values.

3.5.   Statistical analysis

The statistical analysis software we chose was IBM SPSS statistics 23 (Datanine Software, China). One-way analysis of variance (ANOVA) tests was used to assess statistical differences in subject characteristics. Bonferroni correction was applied to avoid spurious rejections when the variance was chi-squared; Tamhane's T2 method was used to compare between groups when the variance was not chi-squared. The significance level for all statistical analyses was set at P < 0.05. Classification accuracy was used to assess the ability of features to distinguish stroke patients.

4.   Result

4.1.   Dynamic time warping distance

In this paper, the lengths of the intercepted motor state signals were different for different subjects, but the lengths of the EEG signals and EMG signals were the same for the same subjects. To verify that this signal interception method can correctly obtain the motion state signal, we calculated the DTW distance between the motion state signal and the resting state signal for each subject, and the results are shown in Figure 3, with the violin plot on the left of each plot and the corresponding scatter plot for each experimental group on the right.

Figure 3.  DTW in EEG and EMG, respectively. (a) DTW in EEG. There are two representations for each group. For each group, the box plot on the left and the corresponding scatter plot on the right; (b) DTW in EMG. There are two representations for each group. For each group, the box plot on the left and the corresponding scatter plot on the right.

DownLoad: Full-Size Img PowerPoint

Figure 3(a) shows the DTW-EEG of subjects in the healthy, B5 and B6 groups. It can be seen from the figure that the DTW-EEG increases with the degree of stroke, which indicates that stroke does disrupt the exchange of information between brain areas and some extent. We also performed a one-way ANOVA to verify the differences between the groups. The results showed that there was no significant difference between B5 and B6 (p = 0.374), a significant difference between B5 and the healthy group (p = 0.013) and a significant difference between B6 and the healthy group (p = 0.043). Figure 3(b) shows the DTW-EMG of subjects in the healthy, B5 and B6 groups. again, for each experimental group, the violin plot is shown on the left and the corresponding scatter plot on the right. As can be seen in Figure 3(b), for stroke patients, as with the DTW-EEG, the DTW-EMG increased as the stroke deepened. the increase in DTW distance indicates that the coupling between the cerebral cortex and between the muscles of the patient may be reduced and that the coordination function between the muscles of the patient is affected by the stroke. One-way ANOVA results showed significant differences between B5 and B6 (p = 0.017), between B5 and healthy group (p = 0.039) and between B6 and healthy group (p = 0.44).

4.2.   Symmetry index

In this paper, we present the BNDSI and CMCSI based on the brain symmetry index. In this section, the results of the BNDSI and CMCSI calculations are presented and statistically analyzed separately to investigate whether the between-group differences between subjects in different groups are statistically significant.

First, we binarize the weighted networks according to the proposed cost-effective threshold selection method. In the threshold selection process, we selected a threshold for each network independently, which resulted in the unavailability of the arithmetic mean of the representation of the binary network. We randomly selected one subject from each group to draw a schematic representation of the binary network, as shown in Figure 4. Each node is placed on an arc edge of a circle. Different channels correspond to different colors, and the connections between channels are represented by colored edge connections between arc edges. It can be seen from the figure that the coherence between channels in group B5 is lower than that in groups control and B6.

Figure 4.  Binary network representation. Each arc on the circle represents a channel, and the functional coupling between channels is represented by the connections between the arcs.

DownLoad: Full-Size Img PowerPoint

Next, we calculated the BNDSI based on the graph theory proposed, and the results are shown in Figure 5. It is straightforward to see from the figure that the BNDSI of the subjects decreased as the degree of stroke deepened. From the formula of BNDSI, it is known that a high BNDSI represents a large gap in the bilateral channeling degree of the brain region. The results showed that deeper stroke motor dysfunction decreased the degree in the contralateral brain region or increased the degree in the ipsilateral brain region, implying that stroke did change the functional arrangement between the bilateral brain regions. One-way ANOVA results showed significant differences between B5 and B6 (p = 0.017), between B5 and healthy group (p = 0.039) and no significant difference between B6 and healthy group (p = 0.094). Figure 5 also shows the CMCSI in subjects with different stroke levels. In contrast to BNDSI, the CMCSI of the subjects increased with the severity of the stroke. This suggests that stroke may increase the proportion of muscle control in contralateral brain regions or decrease the proportion of muscle control in ipsilateral brain regions. We hypothesize that this may be due to motor dysfunction causing the brain to lose some of its ability to control muscles. The contralateral brain region in stroke patients is required to play a greater role in muscle control compared to normal subjects. One-way ANOVA results showed significant differences between B5 and B6 (p = 0.021), between B5 and healthy group (p = 0.010) and between B6 and healthy group (p = 0.048).

Figure 5.  mean network degree symmetry index and mean cortical muscle coupling symmetry index corresponding to the three groups of subjects.

DownLoad: Full-Size Img PowerPoint

4.3.   Classification

In this section, we calculated the importance of CMCSI, BNDSI, DTW-EEG and DTW-EMG features and performed feature combinations based on feature importance ranking. Finally, we performed a triple classification experiment on the combination of these features to assess their performance in the classification of motor function after stroke.

Table 2 shows the ranking of feature importance calculated by the random forest algorithm. The numbers from number 1 to number 4 represent the order of importance of the four features calculated and selected by the random forest, respectively. Further, according to this ranking order, the features are classified from highest to lowest importance and a set of classification results for the combination of 1 + 2 + 4 are added. Table 3 shows the classification accuracy of all combinations. From the table, it can be seen that the classification accuracy is the best when the combination is 1 + 2 + 3, with an accuracy of 84.4%. Additionally, although the special importance of DTW-EMG is not as good as that of DTW-EEG, the accuracy is higher than that of DTW-EEG when DTW-EMG is used for classification alone.

Table 2.  Feature importance ranking.

Abbreviation	1	2	3	4
Feature	CMCSI	BNDSI	DWT-EEG	DWT-EMG

---

 | Show Table

DownLoad: CSV

Table 3.  Classification accuracy of feature combination classification model.

Feature	ACC(%)	Feature combination	ACC(%)
1	76.7	1 + 2	81.3
2	74.8	1 + 2 + 3	84.4
3	68.4	1 + 2 + 3 + 4	82.8
4	73.1	1 + 2 + 4	82.1
ACC refers to the accuracy of the classification. Numbers 1, 2, 3 and 4 are the same as those in Table 2.

---

 | Show Table

DownLoad: CSV

In order to further analyze the evaluation ability of the proposed method for different motor functions, the class accuracy, recall, F1-score and ROC curve area are calculated for each class under 10-fold cross-validation in the dichotomous test, as shown in Table 4. Figure 6 shows the triclass ROC curve and the dichotomous ROC curve under the best feature combination. The method in this paper showed the highest performance in the evaluation of motor function of B6 subjects, and the accuracy and recall rates reached 95.07% and 91.90%, respectively, which were at a high level, indicating that subjects in B6 group could be correctly identified; for the control group, the area of the ROC curve is 0.94, indicating that its recognition performance for the recognition performance of healthy subjects was also good; while for subjects in group B5, their classification performance was relatively poor, with a ROC curve of only 0.73.

Table 4.  Performance indicators are categorized by category.

Category	ACC(%)	Recall(%)	F1-score	ROC area
control group	85.71	80.00	0.8276	0.94
B5	77.85	82.00	0.7987	0.73
B6	95.07	91.90	0.9346	1

---

 | Show Table

DownLoad: CSV

Figure 6.  ROC curves of the best feature combination classification model. The four curves in the figure represent the overall average ROC curve for the three classifications and the ROC curve for each group.

DownLoad: Full-Size Img PowerPoint

4.4.   The validation result of the feature

In order to verify that the combined features of CMCSI+BNDSI+DTW-EEG can effectively evaluate motor function, this section combined with Brunnstrom motor function staging rating to classify whether there was stroke and three groups of different Brunnstrom motor function staging rating subjects, and sensitivity, specificity and classification accuracy were used to evaluate the classification detection performance.

4.4.1.   Verification results under ideal conditions

Assuming that the EEG and sEMG signals are pure signals without noise after denoising, the pure signals are converted into the combined features of CMCSI+BNDSI+DTW-EEG for feature extraction and feature learning.

Firstly, for the problem of detecting whether stroke is dichotomous, this paper performs dichotomous verification of category data balance and dichotomous verification of category data imbalance. When balancing the data categories, we selected an equal number of healthy participants and stroke participants, in which stroke participants were selected from two Brunnstrom motor function staging rated subjects in the same proportion according to stratified sampling for testing experiments; When the data categories were unbalanced, we selected all participants for the testing experiment. Table 5 summarizes the performance indicators of the detection methods and related stroke patient detection methods in the context of data balance. To make a fair comparison, different training and testing methods were used for validation. For example, compared with other methods using the hold-out method (4:1), the proposed method obtains higher performance, that is, 83.88% sensitivity, 82.07% specificity and 83.91% classification accuracy. For the hold-out method, the method in this article is also tested with different training test data ratios. For example, the hold-out method (19:1) was used to compare the training test data in this paper, and the classification performance of 78.88% sensitivity, 79.00% specificity and 79.08% classification accuracy was obtained. Similarly, using the ratio of training test data by hold-out method (1:2), the proposed method achieves 75.66% sensitivity, 74.03% specificity and 74.56% classification accuracy. More importantly, among the existing stroke dichotomous detection methods, the proposed method achieves the highest classification accuracy of 95.68% under the 10-fold cross-validation method, which is 1.68% higher than the highest accuracy gap reported in the literature. This shows that the method presented in this paper has excellent detection performance for stroke motor dysfunction under the condition of balanced data categories.

Table 5.  Methods in this article and methods in the literature for stroke dichotomous detection results.

Method	Training test methods	Sensitivity (%)	Specificity(%)	ACC(%)
Ella Wahyu Guntari [34]	hold-out method (1:2)	——	——	72.22
1 + 2 + 3	hold-out method (4:1)	83.88	82.07	83.91
1 + 2 + 3	hold-out method (19:1)	78.88	79.00	79.08
1 + 2 + 3	hold-out method (1:2)	75.66	74.03	74.56
Nicolas Vivaldi [35]	10-fold cross-validation	——	——	94.00
Nicolas Vivaldi	10-fold cross-validation	——	——	73.00
1 + 2 + 3	10-fold cross-validation	94.20	95.42	95.68
Nicolas Vivaldi	Independent verification	——	——	76.00
Nicolas Vivaldi	Independent verification	——	——	70.05
Arifah Ummul Fadiyah[36]	Independent verification	——	——	79.69
1 + 2 + 3	Independent verification	87.11	86.16	87.61
Arifah Ummul Fadiyah	Leave-One-Out Cross-Validation	86.10	86.50	86.00
1 + 2 + 3	Leave-One-Out Cross-Validation	84.15	84.38	84.33
Numbers 1, 2, 3 and 4 are the same as those in Table 2.

---

 | Show Table

DownLoad: CSV

For the binary classification of stroke, when the data category is unbalanced, Table 6 shows the performance indicators of the method and other classifiers. To make a fair comparison, each classifier performed detection experiments using different training test methods. For the hold-out method (4:1), 75.87% sensitivity, 76.28% specificity, 76.33% classification accuracy and 79.66% sensitivity, 79.26% specificity and 79.09% classification accuracy were obtained under SVM and KNN classifiers, respectively, while the proposed method was improved on the basis of SVM and KNN classifiers. In addition, when using the 10-fold cross-validation method, the proposed method obtains the best specificity of 92.20%, sensitivity of 92.42% and classification accuracy of 92.68%, which are higher than other methods. This shows that the method proposed in this paper has a slight decrease in detection performance compared with data category imbalance, but still reflects accurate detection performance.

Table 6.  Methods in this article and methods in the literature for stroke dichotomous detection results.

Classifier	Training test methods	Sensitivity (%)	Specificity(%)	ACC(%)
SVM	hold-out method (4:1)	75.87	76.28	76.33
KNN	hold-out method (4:1)	79.66	79.26	79.09
1 + 2 + 3	hold-out method (4:1)	81.88	81.07	81.91
SVM	10-fold cross-validation	84.88	84.26	84.71
KNN	10-fold cross-validation	88.09	88.00	88.59
1 + 2 + 3	10-fold cross-validation	92.20	92.42	92.68

---

 | Show Table

DownLoad: CSV

Finally, in order to verify the performance of the proposed method for the evaluation of different motor functions, the classification problem between three different datasets, health, B6 and B5, is solved. Addressing this problem accurately is more challenging, in fact, few studies have focused on the multi-classification of different motor functions in stroke, but solving this problem is very beneficial for the application of stroke rehabilitation and can provide useful clinical practice methods for medical personnel. This paper applies the proposed methodology, as well as some of the most used classification methods in the field of stroke classification over the past decade, to the collected data, and as in previous comparisons, these methods are compared using different training test methods. The results are shown in Table 7, and the proposed method achieves the highest classification performance under the 10-fold cross-validation test method: sensitivity is 86.7%; specificity is 85.34% and accuracy is 84.4%.

Table 7.  The method in this paper and other classifiers perform tri classification detection results on stroke.

Classifier	Training test methods	Sensitivity (%)	Specificity(%)	ACC(%)
SVM	hold-out method (4:1)	66.37	66.33	67.54
KNN	hold-out method (4:1)	70.96	71.67	71.10
1 + 2 + 3	hold-out method (4:1)	77.67	77.05	77.33
SVM	10-fold cross-validation	75.11	76.50	76.67
KNN	10-fold cross-validation	73.33	74.67	74.51
1 + 2 + 3	10-fold cross-validation	86.70	85.34	84.40

---

 | Show Table

DownLoad: CSV

4.4.2.   Validation results in the presence of noise conditions

To further investigate the robustness of the methods presented in this paper to common noises in evaluating motor function, Matlab is used in this section to attach noises (white noise, muscuation and eye movement) to the pure data in Section 4.4.1. When adding noise, we also adjust the noise amplitude to generate noise data with different SIGNAL-NOISE RATIO (SNR). Figure 7 shows the results of stroke dichotographic classification under a wide range of SNR(−8–24 dB) for data balance of different noise types.

Figure 7.  (a) Accuracy of binary classification under different noise types and levels (data balance); (b) Accuracy of binary classification under different noise types and levels (data imbalance); (c) Accuracy of three-way classification under different noise types and levels.

DownLoad: Full-Size Img PowerPoint

As shown in Figure 7(a), even if the EEG signal and EMG signal are completely immersed in noise, the method proposed in this paper achieves the lowest accuracy of 81.61%. Notably, for data with white noise, the proposed method maintains a higher classification accuracy compared to other noises at almost all SNR levels, achieving a classification accuracy of 87.77% at SNR = −8 dB. The proposed method also retained its high classification performance in the presence of muscle activity noise and eye movement noise: at SNR = −8 dB, the lowest classification accuracy obtained with muscle activity noise and eye movement noise was 81.61% and 85.33%, respectively. These are quite high accuracies for data corrupted by severe noise, and eventually settle in a relatively acceptable range after a small decrease in accuracy as the noise content increases.

Second, for imbalanced binary classification problems, the proposed method is also examined on noisy data contaminated by blinking, muscle artifacts and white noise. For this imbalanced classification problem, the proposed method maintains high classification accuracy even at very low SNR, and Figure 7(b) illustrates the detection results obtained by the method in the presence of each noise type. The results show that when the data are completely immersed in white noise (SNR = −8 dB), the lowest classification accuracy is 79.63%. When SNR > 0 dB, the classification accuracy of the proposed method was higher than 80%. Similar to the data balance case, the proposed method achieved better performance compared to muscle activity in the presence of eye movement artifacts and white noise.

Finally, the performance of the proposed method in the presence of noise for the three-classification problem is evaluated. Figure 7(c) depicts the classification performance of the proposed method under different noise SNR. It is worth noting that for multi-class classification, the proposed method still maintains a high robustness. The detection performance was the worst under muscle activity noise, with the lowest (SNR = −8 dB) of 74.33%. The main reason may be that muscle activity occurs in most of the EEG channels, leading to serious distortion of the EEG shape. In addition, under the eye movement noise, the proposed method maintained a higher accuracy than other noises, with the highest accuracy of 83.88% and the lowest (SNR = −8 dB) of 79.69%.

5.   Discussion

The main objective of this study was to find the biological predictive features between stroke rehabilitation ratings. Therefore, we recruited 34 subjects (9 B5, 9 B6 and 16 healthy individuals, as shown in Table 1) for the trial. Using knowledge such as graph theory, combined with EEG and EMG, four biological characteristics such as DTW and symmetry index were calculated for each subject, focusing on the differences in characteristics between the three groups of participants.

First, due to inter-individual response time bias, we selected the longest signal above 10% baseline as the valid signal based on the power of the EMG signal at resting state. We used the DTW algorithm to avoid the difference in data length. Then, due to the excellent performance of graph theory in brain network studies, we improved on the original PBSI based on signal power and proposed BNDSI based on the degree of brain network. In addition, the cortical muscle function network proposed in our previous study has better representation than brain network during motor control ^[37], therefore, we proposed CMCSI. Further, used random forest-based feature importance to rank these features. Finally, based on this ranking, features are added sequentially to obtain an optimal combination of features and their performance is evaluated.

The DTW scatters plot in Figure 3 shows that there are large individual differences among the three groups of subjects. This is one of the difficulties in measuring the degree of recovery of patients only through simple treatment or physician's experience. However, if the subjects in each experimental group are analyzed as a whole, some significant differences can be found between the different groups, such as the DTW-EEG in group B5 is significantly different from group B6 and the healthy group, respectively; in particular, the DTW-EMG increases with the degree of motor dysfunction (which may be caused by motor dysfunction), and the overall average results of the groups are shown in Figure 3, which is consistent with some previous EMG power spectrum studies ^[38].

From the perspective of BNDSI, the BNDSI gradually decreases with the progressive degree of stroke. Studies have shown that lower PBSI values are associated with better motor function ^[12]. This is contrary to our results. It may be that the PBSI is calculated with different criteria than the network degree symmetry index proposed in this paper. Our binary network built on the basis of cost-effectiveness has chosen an independent threshold for each subject, but the results obtained by the uniform criterion are still valid. Considering the mechanism of motor control by contralateral brain regions ^[39], the decrease in BNDSI indicates a decrease in the dominant role of contralateral brain regions or an increase in the dominant role of ipsilateral brain regions during motor control, suggesting that stroke hemiplegia does affect the dominant role of both brain regions during motor control. For CMCSI outcomes, the CMCSI gradually increases with the degree of stroke. the same findings were reported in the brain symmetry index study by Agius ^[40]. PBSI increased in the subacute phase after stroke compared to 1−2 months after stroke. In addition, the decrease in the mean coupling strength between the channels of each brain region and EMG channels confirms the results of previous studies ^[41].

Finally, for these biomarkers, we calculated their respective importance based on random forest. The results showed that CMCSI had the highest feature importance, followed by BNDSI, then DTW-EEG and finally DTW-EMG. Interestingly, although there was a significant change in DTW-EMG compared with DTW-EEG after stroke, DTW-EMG ranked after DTW-EEG in terms of feature importance and the accuracy rate after adding DTW-EMG did not significantly improved. It is speculated that the difference in DTW-EMG may be due to the influence of post-stroke brain regions. This makes DTW-EMG features and other feature combinations unable to form the largest irrelevant vector group, because the random forest calculates the importance of features only for a certain specific characteristic of adding noise interference ^[42], so the importance of DTW-EMG is the lowest. To test this conjecture, we also added a set of CMCSI+BNDSI+DTW-EMG classification results to the feature combination classification. The improvement in classification accuracy was very small and lower than the CMCSI+BNDSI+DTW-EEG classification accuracy.

In addition, the best classification accuracy was achieved with the combined features of CMCSI+BNDSI+DTW-EEG. There are also some studies ^[43,44] that used different electrophysiological features to obtain better classification accuracy, but most of them used a single feature or a single channel, or only two classifications of stroke, which contributed little to the analysis of the whole process of motor function rehabilitation after stroke. Compared with single features, combined features can aggregate information from multiple sources. In an ideal feature selection, it allows for a more complete representation of the classification model. This is confirmed by our final classification results. In order to verify that the combined features of CMCSI+BNDSI+DTW-EEG can effectively evaluate motor function, combined with Brunnstrom motor function stage ratings, we verified whether the subjects with stroke and three groups of subjects with different Brunnstrom motor function stage ratings were classified. It can be seen from Tables 5–7 that although the method proposed in this paper is more complicated than the methods of other researchers in terms of steps, it has achieved excellent classification results in both two-classification and three-classification. This indicates that the proposed combined features of CMCSI+BNDSI+DTW-EEG have high objectivity and consistency in assessing motor function in different stroke. To further investigate the robustness of the proposed method to common noises when assessing motor function, we added white noise, muscle activity noise and eye movement noise to the processed data. It can be concluded from Figure 7 that the combined features of CMCSI+BNDSI+DTW-EEG proposed in this paper have a reduced ability to evaluate stroke motor function in the presence of noise, but the range of decline is limited. Under ideal conditions and in the presence of noise, excellent accuracy was achieved for the detection of stroke motor dysfunction and for the classification of different motor function levels. This result shows the robustness of our proposed method. Previous studies have also proved that for patients, the fusion of data collected by mobile devices and wearable sensor devices can improve the reliability, robustness and generalization ability of the recognition system ^[45,46,47].

Due to the complexity of motor function damage caused by stroke, the cognition of the activity of human neuromuscular system is still in the preliminary stage of research. There are many deficiencies and improvements to be made in the research work: (1) Improvement of experimental protocol. All the stroke subjects selected in this paper were right-handed, but the affected hand was not uniform, and only 10 EEG channels were selected. In the future work, more EEG leads can be selected, different types of subjects can be recruited, and more standardized and novel experimental paradigms can be designed to study the activity characteristics of the neuromuscular system in stroke from different perspectives. (2) In the process of EEG acquisition experiment, the EEG cap is used in this paper, so there is a volume conduction phenomenon between EEG signals: when the voltage fluctuation on the scalp surface is measured by electrodes, it is actually the result of the joint activity of multiple field potential sources. In future studies, trace-back analysis can be considered to find the real source of EEG signal generation. (3) The Brunnstrom motor function rating scale used in the validation of the brain muscle function network to evaluate the motor function of stroke patients recruited only two levels of subjects, compared with other rating scales, there were fewer categories. In the future study, a more detailed motor function assessment scale can be considered for verification. (4) Our method does not involve the direction characteristics of information flow when measuring the coupling strength, but only reflects the interaction relationship between signals. In future studies, we will increase the exploration of the characteristics of bidirectional information flow between cerebral cortex and muscle on the basis of this paper.

6.   Conclusions

In this paper, by combining corticomuscular coupling and graph theory, we obtained the dynamic time warping (DTW) distance of the signal and two new symmetry indices: the brain network degree symmetry index (BNDSI) and the corticomuscular coupling symmetry index (CMCSI). By analyzing the DTW-EEG and DTW-EMG, we found that stroke disrupts the information exchange between brain regions and the coordination between muscles to some extent; by analyzing the symmetry index, we found that stroke changes the brain-muscle connection during movement to some extent. Subsequently, we used the random forest algorithm to calculate the feature importance of these biological indicators. Finally, the different features were combined and validated by classification based on the results of feature importance. The results showed that the feature importance was from high to low as CMCSI/BNDSI/DTW-EEG/DTW-EMG, and the feature combination with the highest accuracy was CMCSI+BNDSI+DTW-EEG. It can reach 84.4%. Compared with single features, the combination of EEG and EMG features achieved better results in motor function rehabilitation prediction under different stroke levels. Our work also implies that the symmetry index based on graph theory and cortico-muscular coupling can distinguish different motor functions of stroke patients, which is expected to provide a richer evaluation method for the recovery of motor function in stroke patients in the clinical field.

Acknowledgments

This work was supported by the Jinhua Science and Technology Bureau under Grant 2019-3020, Zhejiang Provincial Natural Science Foundation of China(LTGY23H180020), and National Natural Science Foundation of China (62061044).

Conflict of interest

The authors declare there is no conflict of interest.