# Decoding of Motor Coordination Imagery Involving the Lower Limbs by the EEG-Based Brain Network

Authors: Yunfa Fu; Zhouzhou Zhou; Anmin Gong; Qian Qian; Lei Su; Lei Zhao
Journal: Computational Intelligence and Neuroscience
 (2021)

---

## Abstract

Compared with the efficacy of traditional physical therapy, a new therapy utilizing motor imagery can induce brain plasticity and allows partial recovery of motor ability in patients with hemiplegia after stroke. Here, we proposed an updated paradigm utilizing motor coordination imagery involving the lower limbs (normal gait imagery and hemiplegic gait imagery after stroke) and decoded such imagery via an electroencephalogram- (EEG-) based brain network. Thirty subjects were recruited to collect EEGs during motor coordination imagery involving the lower limbs. Time-domain analysis, power spectrum analysis, time-frequency analysis, brain network analysis, and statistical analysis were used to explore the neural mechanisms of motor coordination imagery involving the lower limbs. Then, EEG-based brain network features were extracted, and a support vector machine was used for decoding. The results showed that the two employed motor coordination imageries mainly activated sensorimotor areas; the frequency band power was mainly concentrated within theta and alpha bands, and brain functional connections mainly occurred in the right forehead. The combination of the network attributes of the EEG-based brain network and the spatial features of the adjacency matrix had good separability for the two kinds of gait imagery (p < 0.05), and the average classification accuracy of the combination feature was 92.96% ± 7.54%. Taken together, our findings suggest that brain network features can be used to identify normal gait imagery and hemiplegic gait imagery after stroke.


---

## Body

## 1. Introduction

Approximately 45% of poststroke hemiplegic patients depend on wheelchairs to move [1–4]. Unfortunately, traditional physical therapy has only a limited effect on rehabilitation of the lower limbs of poststroke hemiplegic patients [5]. In contrast, motor imagery therapy can induce plasticity in the brain [6–10], which can partially restore lower limb movement in poststroke hemiplegic patients [11–14]. At present, most motor imagery therapies have focused on improving the performance of the upper limbs [15, 16], whereas only a few studies have aimed at improving performance of the lower limbs [17]. In addition, current and past studies on motor coordination have mainly focused on control of redundant degrees of freedom in motor coordination [18–23], quantification of motor coordination [24–29], external environments [30], and evaluation methods of motor coordination [31–34]. However, few studies have investigated the neural mechanisms of motor coordination imagery, and there has not been much focus on decoding motor coordination imagery involving the lower limbs.For the decoding of motor coordination imagery, conventional features have mainly adopted time-domain features [35], frequency-domain features [36], time-frequency features [37], spatial-domain features [38], and fusion of time-frequency-spatial multidomain features [39]. These traditional features mainly reflect the activation level of the corresponding brain area. However, motor coordination imagery of the lower limbs involves not only the activation of related brain regions but also the coordinated control among related brain regions [36]. Therefore, it is necessary to use brain network features to more comprehensively analyze the neural correlates of motor coordination imagery.In the present study, we proposed a new paradigm (normal gait imagery and hemiplegic gait imagery after stroke) of motor coordination imagery involving the lower limbs. In our paradigm, normal gait imagery is abbreviated as motor coordination imagery (MCI), and poststroke hemiplegic gait imagery is abbreviated as imaging exceptional gait (IEG). In addition to traditional time-domain analysis, power spectrum analysis, and time-frequency analysis, a method for brain network analysis was also used to investigate the neural mechanisms and decoding of motor coordination imagery involving the lower limbs.The rest of the paper is organized as follows. Details of the method are illuminated in Section2 including experiments configuration and data analysis. Section 3 displays the brain topographic maps, functional brain networks during MCI, IEG, and REST, and comparison results of the average classification accuracy of MCI and IEG. In Section 4, a brief discussion of the results is given. Finally, Section 5 concludes the paper.


## 2. Methods

### 2.1. Subjects

Previous studies have shown that brain activity patterns evoked by motor imagery and actual exercise are similar [33]. Hence, our present study used the Kinesthetic and Visual Imagery Questionnaire (KVIQ) [38–50] to recruit 30 healthy subjects (numbered s1–s30) with strong motor imagery abilities (questionnaire score ≥ 70 points) to serve as participants for data collection. All recruited subjects were male, right-handed, without psychological disease and were 25 ± 1 years old. This study was approved by the Medical Ethics Committee of Kunming University of Science and Technology. Each subject signed an informed written consent before the start of the experiment.


### 2.2. Experimental Paradigm

Coordinated movement is a type of exercise to restore and strengthen motor coordination. It usually involves the coordinated movement of multiple joints and limbs. Specifically, it includes motor coordination of the following: the upper and lower limbs, the limb trunk, and both symmetrical coordination and asymmetrical coordination of the left and right limbs. Among these components, coordinated movement of the upper limbs and hands can train the accuracies, reaction speeds, and rhythms of actions, while coordinated movement of the lower limbs can train correct gaits and motor coordination of the upper and lower limbs. In order to establish a method for recovering motor abilities of poststroke hemiplegic patients, the present study focused on motor coordination imagery of the lower limbs. We designed MCI and IEG. The timing of each trial is shown in Figure1. The screen at the beginning of the experiment showed a fixed cross for 5 s, which prompted subjects to prepare for the experiment. When the fixed cross disappeared, the screen showed a prompt for normal gait imagery for 2 s, which prompted subjects to prepare for normal gait imagery. Once the prompt disappeared, the screen became blank and required the subjects to imagine a normal gait for 3 s. When the imagery period ended, the screen showed a star prompt for 5 s, indicating the subjects to rest. After the star prompt disappeared, the screen showed a prompt for poststroke hemiplegic gait imagery for 2 s, which prompted subjects to prepare for poststroke hemiplegic gait imagery. After the prompt disappeared, the screen became blank and required subjects to imagine a poststroke hemiplegic gait for 3 s. At the end of the imagery period, the screen showed a star prompt for 5 s, prompting the subjects to rest. This entire sequence comprised one experimental block, and each subject completed a total of 50 blocks.Figure 1
The timing of a single trial.


### 2.3. Experimental Setup

During the experiment, subjects were approximately 70 cm away from the prompt screen, and their hands were placed on the desktop. A TCL 24-inch LCD monitor was used for task prompts. A Lenovo ThinkPad computer running Windows 10 was used for data acquisition and processing via MATLAB software. The electroencephalogram (EEG) amplifier consisted of an NT9200 series (Beijing Zhongke Xintuo Instrument Co., Ltd.) with a sampling rate of 1000 Hz and 45 Hz low-pass filtering. At the beginning of the experiment, each subject completed the task according to the prompts in Figure1. During the experiment, the subjects were required to avoid movements and blinking as much as possible during the visual imagery task. Figure 2 shows the schematic diagram of the electrode layout. The collected EEG data consisted of the following channels: Fp1, Fp2, F7, F3, Fz, F4, F8, FT7, FC3, FCz, FC4, FT8, T7, C3, Cz, C4, T8, TP7, TP8, CP3, CP4, CPz, P7, P3, Pz, P4, P8, PO7, PO8, O1, O2, and Oz. The electrode layout included spatial positioning over the frontal lobe, parietal lobe, occipital lobe, and temporal lobe. The reference electrodes are denoted as A1 and A2.Figure 2
Schematic diagram of the electrode layout.


### 2.4. Data Analysis

#### 2.4.1. Data Preparation

Visual inspection was carried out on the collected EEG data to remove any data segments with substantial noise in the EEG signal. Then, the EEG data of the 30 subjects under different task conditions were extracted.


#### 2.4.2. Pretreatment

First, baseline drift correction was performed on the extracted EEG signals to eliminate deviations of the EEG signal from the baseline [51]. An elliptic filter was used for 8–30 Hz band-pass filtering. Then, independent component analysis (ICA) was used to remove electrooculogram artifacts and electromyographical artifacts. Additionally, alpha, beta, gamma, zeta, delta, and theta bands were extracted.


#### 2.4.3. Time-Domain Analysis

After preprocessing, the data were divided into three categories: MCI, IEG, and process of rest (REST). Then, these three categories of data were stacked and averaged, and we compared the movement-related cortical potentials (MRCPs) among them in the time domain and compared the MRCP average across channels and the grand average across subjects among these three categories of data.


#### 2.4.4. Analysis of Brain Topography

In order to determine the activation location of the brain region under the three tasks of MCI, IEG, and REST, the three categories of data were preprocessed and then averaged, after which the brain topographies during the three tasks were determined at 1, 100, 200, 300, 400, and 500 ms.


#### 2.4.5. Power Spectrum Analysis

In order to analyze the EEG power associated with the task and reduce the impact of absolute power between subjects and electrodes, the absolute powers,ABSP, of MCI and IEG were obtained by using the power of MCI minus the power of REST, as well as by using the power of IEG minus the power of REST. The quantitative method for determining the power spectral density was via a fast Fourier transform algorithm. ABSP was calculated as follows [52]:(1)ABSPMCI=POWERMCI−POWERREST,ABSPIEG=POWERIEG−POWERREST.


#### 2.4.6. Brain Network Analysis

In order to determine the connectivity of the brain function network under the three tasks of MCI, IEG, and REST, a multivariable autoregressive (MVAR) model was established to obtain a time-varying direct directed transfer function (dDTF). Because dDTF can eliminate the pseudocausal relationship derived from other leads [53], dDTF was used as a measure of causality, and Granger causality analysis was performed on the three types of preprocessed data to generate a dDTF weight matrix. The threshold value used in the present study was set to 0.2. Then, we generated an adjacency matrix under the selected threshold. Finally, we obtained the brain function network under the three types of tasks.Specifically, the original EEG data,Xt, of the N lead can be expressed as follows:(2)Xt=X1t,X2t,…,XNtT.Among them,Xit is the EEG time series of the ith lead. MVAR can be constructed according to Xt [54]:(3)Xt=∑k=1pAkXt−k+Et.Among them,Ak is a matrix of model coefficients of N×N, which represents the dependence of the time delay, k;p is the model order, and Et is the random noise. Our present study used the Bayesian information criterion to determine the model order, p.By Fourier transform ofAk, we obtain(4)Af=−∑k=0pAke−i2πfk.Furthermore,Xf in the frequency domain can be obtained as follows:(5)Xf=Af−1Ef=HfEf.Among them,Hf is the transfer matrix of the system.Furthermore, the dDTF can be obtained as [55](6)χij2f=Fij2fCij2f.Among them,Fij2fis the full frequency DTF, and the calculation method is as follows:(7)Fij2f=Hijf2∑f∑k=1NHikf2,where Cijf is the partial coherence. The calculation method is as follows:(8)Cijf=MijfMiifMjjf,where Mij is a minor of spectral matrix (matrix of spectra and cross-spectra) with the i-th row and j-th column removed.


#### 2.4.7. Measurement Analysis of the Brain Function Network

In order to analyze the connectivity of the brain function network in detail, four kinds of network measures (characteristic path length, global efficiency, clustering coefficient, and local efficiency) were selected for comparative analysis of the three types of tasks. The characteristic path length, PL, and global efficiency,Eglob, measure the global transmission capability of the network as follows [56]:(9)PL=1NN−1∑i,j∈N,i≠jdij,(10)Eglob=1NN−1∑i,j∈N,i≠jdij.In equations (9) and (10), the distance between two nodes, i and j, in the network is represented by dij, which represents the edge number information of the shortest path, and N is the total number of nodes in the network. It can be determined from formulas (9) and (10) that the calculation method of global efficiency is similar to the characteristic path length. The difference is that the global efficiency, Eglob, does not need to consider the isolated points in the network, which is suitable for a regular network, small world network, or a random network.The local efficiencyEloc measures the local transmission capacity of the network. The local efficiency Eloc of the network is the average of the global efficiency of the subgraphs corresponding to all nodes in the network:(11)Eloc=1N∑i∈NEglobGi.In equation (11), Gi represents the subnetwork formed by all nodes directly connected with node i.The cluster coefficientCLa measures the degree of grouping of the network. The cluster coefficient CLa of a network is the mean value of the cluster coefficient cli of all nodes in the network, where the cluster coefficient of node i is defined as follows:(12)cli=2eikiki−1=∑j,maijaimajmkiki−1.In equation (12), ∑j,maijaimajm is the number of edges included in the subnetwork composed of node i and is directly connected to nodes j and m. kiki−1 is the maximum number of edges of the subnetwork, and ki is the degree of node i. Therefore, the clustering coefficient CLa of a network is as follows:(13)CLa=1N∑i=1Ncli.In equation (13), N represents the total number of nodes in the network.


### 2.5. Classification Verification

In order to verify the difference between normal motor coordination imagery involving the lower limbs and that of poststroke hemiplegia gait imagery, the time-domain MRCP feature, power spectrum feature, time-frequency feature, network attribute feature, combined features of network attributes, and the adjacency matrix space were extracted. Since the SVM classifier has a good classification effect on small sample data [57], SVM was used for classification verification.Specifically, the MRCP feature was used to extract the peak and phase of the 32 electrodes and to form a 32-dimensional feature matrix. The feature of the power spectrum was used to for the AR model to estimate the power spectrum and to obtain a 32-dimensional feature matrix. The time-frequency feature was to extract the time-frequency energy of the 32-lead EEGs by HHT and to form a 32-dimensional feature matrix. The network attribute feature was used to calculate the characteristic path length of the 32 electrodes to form a 32-dimensional characteristic matrix. The spatial feature of the adjacency matrix was used to extract the 16-dimensional feature matrix of the adjacency matrix using a spatial filter.


### 2.6. Statistical Analysis

In order to test the significance of the task differences, the characteristic path lengths of 15 channels for the three types of tasks were compared in pairs by repeated-measures analysis of variance. The statistical significance level was set asp < 0.05.


## 2.1. Subjects

Previous studies have shown that brain activity patterns evoked by motor imagery and actual exercise are similar [33]. Hence, our present study used the Kinesthetic and Visual Imagery Questionnaire (KVIQ) [38–50] to recruit 30 healthy subjects (numbered s1–s30) with strong motor imagery abilities (questionnaire score ≥ 70 points) to serve as participants for data collection. All recruited subjects were male, right-handed, without psychological disease and were 25 ± 1 years old. This study was approved by the Medical Ethics Committee of Kunming University of Science and Technology. Each subject signed an informed written consent before the start of the experiment.


## 2.2. Experimental Paradigm

Coordinated movement is a type of exercise to restore and strengthen motor coordination. It usually involves the coordinated movement of multiple joints and limbs. Specifically, it includes motor coordination of the following: the upper and lower limbs, the limb trunk, and both symmetrical coordination and asymmetrical coordination of the left and right limbs. Among these components, coordinated movement of the upper limbs and hands can train the accuracies, reaction speeds, and rhythms of actions, while coordinated movement of the lower limbs can train correct gaits and motor coordination of the upper and lower limbs. In order to establish a method for recovering motor abilities of poststroke hemiplegic patients, the present study focused on motor coordination imagery of the lower limbs. We designed MCI and IEG. The timing of each trial is shown in Figure1. The screen at the beginning of the experiment showed a fixed cross for 5 s, which prompted subjects to prepare for the experiment. When the fixed cross disappeared, the screen showed a prompt for normal gait imagery for 2 s, which prompted subjects to prepare for normal gait imagery. Once the prompt disappeared, the screen became blank and required the subjects to imagine a normal gait for 3 s. When the imagery period ended, the screen showed a star prompt for 5 s, indicating the subjects to rest. After the star prompt disappeared, the screen showed a prompt for poststroke hemiplegic gait imagery for 2 s, which prompted subjects to prepare for poststroke hemiplegic gait imagery. After the prompt disappeared, the screen became blank and required subjects to imagine a poststroke hemiplegic gait for 3 s. At the end of the imagery period, the screen showed a star prompt for 5 s, prompting the subjects to rest. This entire sequence comprised one experimental block, and each subject completed a total of 50 blocks.Figure 1
The timing of a single trial.


## 2.3. Experimental Setup

During the experiment, subjects were approximately 70 cm away from the prompt screen, and their hands were placed on the desktop. A TCL 24-inch LCD monitor was used for task prompts. A Lenovo ThinkPad computer running Windows 10 was used for data acquisition and processing via MATLAB software. The electroencephalogram (EEG) amplifier consisted of an NT9200 series (Beijing Zhongke Xintuo Instrument Co., Ltd.) with a sampling rate of 1000 Hz and 45 Hz low-pass filtering. At the beginning of the experiment, each subject completed the task according to the prompts in Figure1. During the experiment, the subjects were required to avoid movements and blinking as much as possible during the visual imagery task. Figure 2 shows the schematic diagram of the electrode layout. The collected EEG data consisted of the following channels: Fp1, Fp2, F7, F3, Fz, F4, F8, FT7, FC3, FCz, FC4, FT8, T7, C3, Cz, C4, T8, TP7, TP8, CP3, CP4, CPz, P7, P3, Pz, P4, P8, PO7, PO8, O1, O2, and Oz. The electrode layout included spatial positioning over the frontal lobe, parietal lobe, occipital lobe, and temporal lobe. The reference electrodes are denoted as A1 and A2.Figure 2
Schematic diagram of the electrode layout.


## 2.4. Data Analysis

### 2.4.1. Data Preparation

Visual inspection was carried out on the collected EEG data to remove any data segments with substantial noise in the EEG signal. Then, the EEG data of the 30 subjects under different task conditions were extracted.


### 2.4.2. Pretreatment

First, baseline drift correction was performed on the extracted EEG signals to eliminate deviations of the EEG signal from the baseline [51]. An elliptic filter was used for 8–30 Hz band-pass filtering. Then, independent component analysis (ICA) was used to remove electrooculogram artifacts and electromyographical artifacts. Additionally, alpha, beta, gamma, zeta, delta, and theta bands were extracted.


### 2.4.3. Time-Domain Analysis

After preprocessing, the data were divided into three categories: MCI, IEG, and process of rest (REST). Then, these three categories of data were stacked and averaged, and we compared the movement-related cortical potentials (MRCPs) among them in the time domain and compared the MRCP average across channels and the grand average across subjects among these three categories of data.


### 2.4.4. Analysis of Brain Topography

In order to determine the activation location of the brain region under the three tasks of MCI, IEG, and REST, the three categories of data were preprocessed and then averaged, after which the brain topographies during the three tasks were determined at 1, 100, 200, 300, 400, and 500 ms.


### 2.4.5. Power Spectrum Analysis

In order to analyze the EEG power associated with the task and reduce the impact of absolute power between subjects and electrodes, the absolute powers,ABSP, of MCI and IEG were obtained by using the power of MCI minus the power of REST, as well as by using the power of IEG minus the power of REST. The quantitative method for determining the power spectral density was via a fast Fourier transform algorithm. ABSP was calculated as follows [52]:(1)ABSPMCI=POWERMCI−POWERREST,ABSPIEG=POWERIEG−POWERREST.


### 2.4.6. Brain Network Analysis

In order to determine the connectivity of the brain function network under the three tasks of MCI, IEG, and REST, a multivariable autoregressive (MVAR) model was established to obtain a time-varying direct directed transfer function (dDTF). Because dDTF can eliminate the pseudocausal relationship derived from other leads [53], dDTF was used as a measure of causality, and Granger causality analysis was performed on the three types of preprocessed data to generate a dDTF weight matrix. The threshold value used in the present study was set to 0.2. Then, we generated an adjacency matrix under the selected threshold. Finally, we obtained the brain function network under the three types of tasks.Specifically, the original EEG data,Xt, of the N lead can be expressed as follows:(2)Xt=X1t,X2t,…,XNtT.Among them,Xit is the EEG time series of the ith lead. MVAR can be constructed according to Xt [54]:(3)Xt=∑k=1pAkXt−k+Et.Among them,Ak is a matrix of model coefficients of N×N, which represents the dependence of the time delay, k;p is the model order, and Et is the random noise. Our present study used the Bayesian information criterion to determine the model order, p.By Fourier transform ofAk, we obtain(4)Af=−∑k=0pAke−i2πfk.Furthermore,Xf in the frequency domain can be obtained as follows:(5)Xf=Af−1Ef=HfEf.Among them,Hf is the transfer matrix of the system.Furthermore, the dDTF can be obtained as [55](6)χij2f=Fij2fCij2f.Among them,Fij2fis the full frequency DTF, and the calculation method is as follows:(7)Fij2f=Hijf2∑f∑k=1NHikf2,where Cijf is the partial coherence. The calculation method is as follows:(8)Cijf=MijfMiifMjjf,where Mij is a minor of spectral matrix (matrix of spectra and cross-spectra) with the i-th row and j-th column removed.


### 2.4.7. Measurement Analysis of the Brain Function Network

In order to analyze the connectivity of the brain function network in detail, four kinds of network measures (characteristic path length, global efficiency, clustering coefficient, and local efficiency) were selected for comparative analysis of the three types of tasks. The characteristic path length, PL, and global efficiency,Eglob, measure the global transmission capability of the network as follows [56]:(9)PL=1NN−1∑i,j∈N,i≠jdij,(10)Eglob=1NN−1∑i,j∈N,i≠jdij.In equations (9) and (10), the distance between two nodes, i and j, in the network is represented by dij, which represents the edge number information of the shortest path, and N is the total number of nodes in the network. It can be determined from formulas (9) and (10) that the calculation method of global efficiency is similar to the characteristic path length. The difference is that the global efficiency, Eglob, does not need to consider the isolated points in the network, which is suitable for a regular network, small world network, or a random network.The local efficiencyEloc measures the local transmission capacity of the network. The local efficiency Eloc of the network is the average of the global efficiency of the subgraphs corresponding to all nodes in the network:(11)Eloc=1N∑i∈NEglobGi.In equation (11), Gi represents the subnetwork formed by all nodes directly connected with node i.The cluster coefficientCLa measures the degree of grouping of the network. The cluster coefficient CLa of a network is the mean value of the cluster coefficient cli of all nodes in the network, where the cluster coefficient of node i is defined as follows:(12)cli=2eikiki−1=∑j,maijaimajmkiki−1.In equation (12), ∑j,maijaimajm is the number of edges included in the subnetwork composed of node i and is directly connected to nodes j and m. kiki−1 is the maximum number of edges of the subnetwork, and ki is the degree of node i. Therefore, the clustering coefficient CLa of a network is as follows:(13)CLa=1N∑i=1Ncli.In equation (13), N represents the total number of nodes in the network.


## 2.4.1. Data Preparation

Visual inspection was carried out on the collected EEG data to remove any data segments with substantial noise in the EEG signal. Then, the EEG data of the 30 subjects under different task conditions were extracted.


## 2.4.2. Pretreatment

First, baseline drift correction was performed on the extracted EEG signals to eliminate deviations of the EEG signal from the baseline [51]. An elliptic filter was used for 8–30 Hz band-pass filtering. Then, independent component analysis (ICA) was used to remove electrooculogram artifacts and electromyographical artifacts. Additionally, alpha, beta, gamma, zeta, delta, and theta bands were extracted.


## 2.4.3. Time-Domain Analysis

After preprocessing, the data were divided into three categories: MCI, IEG, and process of rest (REST). Then, these three categories of data were stacked and averaged, and we compared the movement-related cortical potentials (MRCPs) among them in the time domain and compared the MRCP average across channels and the grand average across subjects among these three categories of data.


## 2.4.4. Analysis of Brain Topography

In order to determine the activation location of the brain region under the three tasks of MCI, IEG, and REST, the three categories of data were preprocessed and then averaged, after which the brain topographies during the three tasks were determined at 1, 100, 200, 300, 400, and 500 ms.


## 2.4.5. Power Spectrum Analysis

In order to analyze the EEG power associated with the task and reduce the impact of absolute power between subjects and electrodes, the absolute powers,ABSP, of MCI and IEG were obtained by using the power of MCI minus the power of REST, as well as by using the power of IEG minus the power of REST. The quantitative method for determining the power spectral density was via a fast Fourier transform algorithm. ABSP was calculated as follows [52]:(1)ABSPMCI=POWERMCI−POWERREST,ABSPIEG=POWERIEG−POWERREST.


## 2.4.6. Brain Network Analysis

In order to determine the connectivity of the brain function network under the three tasks of MCI, IEG, and REST, a multivariable autoregressive (MVAR) model was established to obtain a time-varying direct directed transfer function (dDTF). Because dDTF can eliminate the pseudocausal relationship derived from other leads [53], dDTF was used as a measure of causality, and Granger causality analysis was performed on the three types of preprocessed data to generate a dDTF weight matrix. The threshold value used in the present study was set to 0.2. Then, we generated an adjacency matrix under the selected threshold. Finally, we obtained the brain function network under the three types of tasks.Specifically, the original EEG data,Xt, of the N lead can be expressed as follows:(2)Xt=X1t,X2t,…,XNtT.Among them,Xit is the EEG time series of the ith lead. MVAR can be constructed according to Xt [54]:(3)Xt=∑k=1pAkXt−k+Et.Among them,Ak is a matrix of model coefficients of N×N, which represents the dependence of the time delay, k;p is the model order, and Et is the random noise. Our present study used the Bayesian information criterion to determine the model order, p.By Fourier transform ofAk, we obtain(4)Af=−∑k=0pAke−i2πfk.Furthermore,Xf in the frequency domain can be obtained as follows:(5)Xf=Af−1Ef=HfEf.Among them,Hf is the transfer matrix of the system.Furthermore, the dDTF can be obtained as [55](6)χij2f=Fij2fCij2f.Among them,Fij2fis the full frequency DTF, and the calculation method is as follows:(7)Fij2f=Hijf2∑f∑k=1NHikf2,where Cijf is the partial coherence. The calculation method is as follows:(8)Cijf=MijfMiifMjjf,where Mij is a minor of spectral matrix (matrix of spectra and cross-spectra) with the i-th row and j-th column removed.


## 2.4.7. Measurement Analysis of the Brain Function Network

In order to analyze the connectivity of the brain function network in detail, four kinds of network measures (characteristic path length, global efficiency, clustering coefficient, and local efficiency) were selected for comparative analysis of the three types of tasks. The characteristic path length, PL, and global efficiency,Eglob, measure the global transmission capability of the network as follows [56]:(9)PL=1NN−1∑i,j∈N,i≠jdij,(10)Eglob=1NN−1∑i,j∈N,i≠jdij.In equations (9) and (10), the distance between two nodes, i and j, in the network is represented by dij, which represents the edge number information of the shortest path, and N is the total number of nodes in the network. It can be determined from formulas (9) and (10) that the calculation method of global efficiency is similar to the characteristic path length. The difference is that the global efficiency, Eglob, does not need to consider the isolated points in the network, which is suitable for a regular network, small world network, or a random network.The local efficiencyEloc measures the local transmission capacity of the network. The local efficiency Eloc of the network is the average of the global efficiency of the subgraphs corresponding to all nodes in the network:(11)Eloc=1N∑i∈NEglobGi.In equation (11), Gi represents the subnetwork formed by all nodes directly connected with node i.The cluster coefficientCLa measures the degree of grouping of the network. The cluster coefficient CLa of a network is the mean value of the cluster coefficient cli of all nodes in the network, where the cluster coefficient of node i is defined as follows:(12)cli=2eikiki−1=∑j,maijaimajmkiki−1.In equation (12), ∑j,maijaimajm is the number of edges included in the subnetwork composed of node i and is directly connected to nodes j and m. kiki−1 is the maximum number of edges of the subnetwork, and ki is the degree of node i. Therefore, the clustering coefficient CLa of a network is as follows:(13)CLa=1N∑i=1Ncli.In equation (13), N represents the total number of nodes in the network.


## 2.5. Classification Verification

In order to verify the difference between normal motor coordination imagery involving the lower limbs and that of poststroke hemiplegia gait imagery, the time-domain MRCP feature, power spectrum feature, time-frequency feature, network attribute feature, combined features of network attributes, and the adjacency matrix space were extracted. Since the SVM classifier has a good classification effect on small sample data [57], SVM was used for classification verification.Specifically, the MRCP feature was used to extract the peak and phase of the 32 electrodes and to form a 32-dimensional feature matrix. The feature of the power spectrum was used to for the AR model to estimate the power spectrum and to obtain a 32-dimensional feature matrix. The time-frequency feature was to extract the time-frequency energy of the 32-lead EEGs by HHT and to form a 32-dimensional feature matrix. The network attribute feature was used to calculate the characteristic path length of the 32 electrodes to form a 32-dimensional characteristic matrix. The spatial feature of the adjacency matrix was used to extract the 16-dimensional feature matrix of the adjacency matrix using a spatial filter.


## 2.6. Statistical Analysis

In order to test the significance of the task differences, the characteristic path lengths of 15 channels for the three types of tasks were compared in pairs by repeated-measures analysis of variance. The statistical significance level was set asp < 0.05.


## 3. Results

Figure3(a) shows the average EEG waveforms across channels and the grand average across subjects during MCI. The colored line represents the average EEG data during MCI across all subjects, and the black line represents the grand average EEG data of MCI across the 32 channels. It can be seen from Figure 3(a) that a negative ERP was generated at about 200 ms. Figure 3(b) shows the average and grand average waveforms during IEG for all subjects. Figure 3(c) shows a comparison of the stack average of EEG data after preprocessing under MCI, IEG, and REST tasks. Figure 3(d) shows the phase difference of EEG data between MCI and REST, IEG and REST, and MCI and IEG. The EEG waveform trends of IEG and REST states were relatively consistent, while the EEG waveform trends of MCI and REST states were almost opposite to one another; additionally, the EEG waveform trends of MCI and IEG states were also opposite to one another (Figure 3(d)).Figure 3
The average EEG waveforms across subjects, the average EEG across channels, and the phase difference of EEG data during MCI, IEG, and REST.
(a)(b)(c)(d)Figure4(a) shows the MRCP difference between the MCI and IEG on the C3 and C4 channels, and Figure 4(b) shows the MRCP phase difference between the MCI and IEG on the C3 and C4 channels. Figure 4 shows that, for the MCI MRCPs, the waveforms of the C3 and C4 channels were similar, but there was a significant difference at about 200 ms. For the IEG MRCPs, the waveforms of the C3 and C4 channels were nearly opposite to one another, and the waveforms of the C3 channel were half a cycle ahead of those of the C4 channel, while they reached a consistency at about 200 ms.Figure 4
Differences in MRCPs and MRCP phases between MCI and IEG.
(a)(b)Figure5 shows the brain topographic maps at 1, 100, 200, 300, 400, and 500 ms during MCI, IEG, and REST. At about 200 ms, during MCI, the brain topographic map was mainly activated in the right brain sensorimotor area (C4 channel); during IEG, the brain topographic map was mainly activated in the left brain sensorimotor area (C3 channel) and the parietal lobe (Pz channel), while during REST the brain topographic map was mainly activated in the left and right sensorimotor areas (C3 and C4).Figure 5
The brain topographic maps at 1, 100, 200, 300, 400, and 500 ms during MCI, IEG, and REST. The left, middle, and right columns show the brain topographic maps during MCI, IEG, and REST tasks, respectively. The brain topographic maps, from top to bottom, occurred at 1, 100, 200, 300, 400, and 500 ms during the corresponding task.In terms of timing, MCI first activated the right frontal lobe (F4 channel), left sensorimotor area (C3 channel), and parietal lobe (Pz channel) at 100 ms; the right sensorimotor area (C4 channel) was activated at 200 ms, whereas the left sensorimotor area (C3 channel) and parietal lobe (Pz channel) were activated at 300 ms. The central frontal area (FCz channel), central sensorimotor area (Cz channel), and right parietal lobe (P4 channel) were activated at 400 ms; the central areas of right frontal lobe (FC4 channel) and left parietal lobe (P3 channel) were activated at 500 ms. IEG first activated the left and right parietal central areas (CP3 and CP4 channels) at 100 ms; the left sensorimotor area (C3 channel) and the parietal lobe (Pz channel) were activated at 200 ms; the frontal central area (FCz channel) was activated at 300 ms. The left frontal central areas (FC3 channel), middle frontal lobe (FCz channel), and left parietal lobes (P3 channel) were activated at 400 ms; the left frontal lobe (F3 channel) was activated at 500 ms. REST first activated the central area of the left and right parietal lobes (CP3 and CP4 channels) at 100 ms; the left and right sensorimotor regions (C3 and C4) were activated at 200 ms. The left frontal lobes (F3 channel), left sensorimotor region (C3 channel), and the central area of the parietal lobes (Pz channel) were activated at 400 ms.Figure6 shows the power spectrum during MCI, IEG, and REST tasks, as well as the absolute power spectrum of MCI and IEG. According to changes in the absolute power spectrum of MCI and IEG in Figures 6(d) and 6(e), the power distributions of MCI and IEG were in the θ wave (4–7 Hz), α wave (8–13 Hz), and β wave (14–16 Hz). The frequency band corresponding to the absolute power showed that the power related to MCI was mainly concentrated in the α wave, while the power related to IEG was mainly concentrated in the θ wave.Figure 6
Comparison of absolute power changes during MCI and IEG.
(a)(b)(c)(d)(e)Figure7 shows the brain function networks of MCI, IEG, and REST at 1, 100, 200, and 300 ms. The main core nodes of the brain network were F7, FT7, T8, CP3, Pz, and PO8 at 1 ms for MCI, and the network connection took place from the left frontal lobe to the right occipital lobe (Figure 7). The main core nodes of the brain network were FP1, Fz, Oz, and PO8 at 100 ms, and the network connection took place from the forehead to the occipital lobe. The main core nodes of the brain network were F8, FT8, FCz, CP4, TP8, PO7, and O1 at 200 ms. The main core nodes of the brain network were FP2, F3, FT7, T8, TP8, PO8, and Oz at 300 ms, and the network connection occurred from the right forehead to the right occipital lobe. At 1 ms, the main core nodes during IEG were F7, FP2, F4, FCz, C4, P3, and O2, and the network connection took place in the left and right hemispheres; at 100 ms, the main core nodes of the brain network were F7, Fz, P4, and PO8, and the network connection took place from the left forehead to the right occipital lobe. At 200 ms, the main core nodes of the brain network were F4, FCZ, CP4, TP7, and Pz, and the network connection took place in the right hemisphere; at 300 ms, the main core nodes of the brain network were FCz, FT8, P7, Pz, and O1, and the network connection occurred from right forehead to left occipital lobe. When REST was 1 ms, the main core nodes of brain network were FP1, FP2, F7, FC3, T8, O1, and O2, and the network connection took place in the left prefrontal lobe and the right occipital lobe; when REST was 100 ms, the main core nodes of the brain network were PO7, CPz, CP4, P4, and TP8, and the network connection took place in the left and right occipital lobes. At 200 ms during REST, the main core nodes of the brain network were F7, T8, and O2, and the network connection took place in the left prefrontal lobe; when REST was 300 ms, the main core nodes of the brain network were FP2, CPZ, and PO8, and the network connection occurred from the right forehead to the occipital lobe.Figure 7
Functional brain networks at 1, 100, 200, and 300 ms during MCI, IEG, and REST. The upper, middle, and lower rows show the brain function networks during MCI, IEG, and REST, respectively. From left to right, the brain function networks occur at 1, 100, 200, and 300 ms during the corresponding task, respectively.Figure8 shows the comparison of brain network measures among MCI, IEG, and REST tasks. The four selected brain network measures were as follows: characteristic path length, global efficiency, clustering coefficient, and local efficiency. Among these four measures, the difference in the characteristic path lengths was the largest across MCI, IEG, and REST, while the other three measures exhibited only small differences across groups (Figure 8). REST had the largest characteristic path length, followed by MCI and finally IEG. The characteristic path length showed that the speed of information transmission between nodes in IEG was the fastest, followed by MCI and REST.Figure 8
Comparison of network measures during MCI, IEG, and REST.In order to statistically analyze differences among MCI, IEG, and REST, the characteristic path lengths of the 15 channels during MCI, IEG, and REST were selected for repeated-measures single-factor analysis of variance. Table1 shows the statistical results, where significant differences are marked in red. The difference in the characteristic path length between MCI and REST was significant, and only the O2 channel of MCI was not significantly different between MCI and REST (Table 1). There was no significant difference between IEG and REST. Only the p values of P3, C3, CP4, and Pz channel were small during IEG, and there was a significant difference between MCI and IEG.Table 1
Results of repeated-measures one-way analysis of variance of characteristic path lengths of the 15 channels during MCI, IEG, and REST.
p valuesRESTF3F4FC3FC4C3CzC4CP3CP4P3PzP4O1O2OzMCIF30.0000.0000.0010.0000.0000.0020.0010.0000.0000.0000.0000.0000.0000.0000.000F40.0000.0020.0070.0000.0020.0140.0090.0000.0000.0010.0000.0010.0020.0000.000FC30.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000FC40.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000C30.0000.0000.0000.0000.0000.0010.0000.0000.0000.0000.0000.0000.0000.0000.000Cz0.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000C40.0000.0000.0000.0000.0000.0010.0000.0000.0000.0000.0000.0000.0000.0000.000CP30.0000.0000.0010.0000.0000.0020.0010.0000.0000.0000.0000.0000.0000.0000.000CP40.0000.0030.0070.0000.0030.0120.0090.0000.0000.0010.0000.0010.0020.0000.000P30.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000Pz0.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000P40.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000O10.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.0000.000O20.0310.1600.2770.0000.1740.3590.3120.0630.0440.1130.0450.1540.1990.0110.022Oz0.0010.0080.0170.0000.0080.0250.0210.0020.0010.0040.0010.0060.0090.0000.001IEGF30.0000.0050.0120.0000.0050.0190.0160.0010.0000.0020.0000.0020.0050.0000.000F40.0020.0390.1040.0000.0410.1620.1270.0050.0030.0170.0030.0230.0480.0010.003FC30.0500.4830.9490.0000.5270.7960.9360.1330.0790.3150.0820.4880.6860.0180.039FC40.0340.2260.4220.0000.2440.5490.4780.0740.0480.1480.0510.2130.2990.0140.025C30.0000.0030.0090.0000.0030.0170.0120.0000.0000.0010.0000.0010.0030.0000.000Cz0.1180.8520.5880.0000.9250.3650.4800.3100.1840.6210.1950.9330.8270.0400.080C40.0420.4260.8610.0000.4650.8870.9730.1110.0650.2700.0670.4190.6040.0160.033CP30.0210.2140.4540.0000.2320.6220.5310.0480.0310.1210.0330.1810.2910.0090.017CP40.0000.0050.0140.0000.0050.0230.0190.0010.0000.0020.0000.0020.0060.0000.000P30.0490.3720.7020.0000.4040.8970.7870.1160.0730.2440.0770.3630.5070.0190.037Pz0.0000.0060.0170.0000.0070.0290.0230.0010.0000.0020.0000.0030.0070.0000.000P40.0190.2570.5950.0000.2850.8330.6990.0460.0300.1340.0310.2130.3630.0080.017O10.0740.6440.8210.0000.7030.5680.7010.1980.1170.4400.1220.6790.9110.0260.054O20.0280.1560.2730.0000.1660.3480.3090.0590.0400.1080.0420.1480.1990.0110.019Oz0.0220.1470.2750.0000.1590.3560.0150.0150.0150.0150.0150.0150.0150.0150.01Figure9 shows the average classification results obtained by classifying MCI and IEG. The classifier adopted features of MRCPs, time-frequency analysis (Time-freq), power spectrum analysis (Power spec), network attribute analysis (Net-attr), a combination of the network attribute features, and the spatial features of the adjacency matrix (Net-attr + Adj-mat-space). The average classification accuracies of MRCPs, time-frequency analysis, power spectrum analysis, network attribute analysis, a combination of the network attribute features, and the spatial features of the adjacency matrix were 67.89% ± 9.14%, 78.12% ± 8.91%, 81.08% ± 8.99%, 89.12% ± 9.92%, and 92.96% ± 7.54%, respectively. According to the average classification results, the combination of the network attribute feature and the spatial feature of the adjacency matrix had the highest classification accuracies.Figure 9
Comparison of the average classification results of MCI and IEG classified by MRCPs, time-frequency analysis, power spectrum analysis, network attribute analysis, a combination of the network attribute features, and the spatial features of the adjacency matrix.


## 4. Discussion

Compared with the efficacy of traditional physical therapy, motor imagery therapy can induce brain plasticity, which can partially restore motor abilities of poststroke hemiplegic patients. At present, motor imagery therapy is mostly used for the recovery of upper-limb motor abilities but rarely for the recovery of lower-limb motor abilities. In addition, although there have been studies on the actual motor coordination of the lower limbs, few studies have investigated the neural mechanisms and decoding of imagined motor coordination of the lower limbs. Therefore, the focus of our present study was distinct from that of traditional studies on motor coordination, such that we focused on motor coordination imagery involving the lower limbs.First, our present study was different from the simple motor imagery used traditionally for a unilateral limb, such as studies that focused on the upper limbs [58–64], as well as lower limbs [65–67]. Prior to our present study, no related studies on EEG-based lower-limb motor coordination have been published. In order to determine a rehabilitation method for lower-limb motor coordination in poststroke hemiplegic patients, we used a paradigm of motor coordination imagery involving the lower limbs that consisted of both normal gait imagery and poststroke hemiplegic gait imagery. Simple motor imagery of a unilateral limb mainly activates the contralateral brain area [68], while motor coordination imagery involving both lower limbs involves coordination of the left and right lower limbs, which requires the coordinated participation of relevant brain regions on the left and right sides of the brain. This kind of motor coordination imagery training may be beneficial for inducing plasticity of brain network structure/function or affected patients. In future research, we will design and improve our experimental research paradigm for poststroke hemiplegic patients with problems in motor coordination of the lower limbs.The EEG pattern induced by motor coordination imagery involving the lower limbs is closely related to the performance of motor imagery psychological activities of subjects [68]. In order to ensure the effective completion of the motor coordination imagery involving the lower limbs and the reliability of the corresponding data, before the experiment, subjects were selected according to their scores on a motor-vision imagery questionnaire; subjects with scores higher than 70 were selected for participation in our present study. After the experiment, for the designed imagery task, each subject was required to fill in a questionnaire about the implementation of psychological activities and the degree of difficulty of imagination during completion of the experiments in our present study. The subjects reported that they could complete MCI and IEG in a controllable way, but that for MCI was easier to imagine compared to that for IEG. In future research, we intend to use online neurofeedback to further improve subjects’ motor coordination imagery performance of the lower limbs.In terms of analysis content, previous studies have mainly focused on the control of redundant degrees of freedom, quantification of motor coordination, external environments, and evaluation of motor coordination, while our present study focused on the neural mechanisms and decoding of motor coordination imagery involving the lower limb. In the EEG analysis method, in addition to traditional sensorimotor rhythm analysis [58–67], we also used MRCPs and brain networks to study motor coordination imagery involving the lower limbs. Previous studies have shown that MRCPs can represent neural activities related to planning before each exercise, the execution of each exercise, and the end of each exercise [64]. Previous studies have also shown that features of brain networks can represent the coordination among brain regions during activities such as sensation, perception, representation, and cognition [53]. The motor coordination imagery in our present study involved planning, execution, and termination of motor imagery, as well as coordination between brain regions.In the time domain, the EEG patterns of IEG and REST were more consistent, while the EEG patterns of MCI and REST were almost opposite to one another, and the EEG patterns of MCI and IEG were also nearly opposite to one another. These findings may be due to the fact that the subjects were more familiar with resting states and walking states but were not as familiar with poststroke hemiplegic gaits. Unusual imagery can trigger different EEG patterns [36]. In our present study, the MRCPs of MCI and IEG were significantly different in terms of their EEG signals in the C3 and C4 channels. The EEG waveforms of MCI in the C3 and C4 channels were similar, but there was a significant difference at about 200 ms. The EEG waveforms of IEG in the C3 and C4 channels were nearly opposite to one another. The waveform of the C3 channel was a half cycle ahead of the waveform of the C4 channel, whereas they reached a consistent value at about 200 ms. This result may be due to MCI mainly involving bilateral lower limbs, while IEG mainly involves unilateral lower limbs. MRCPs will trigger ERP waveforms around 200 ms [35], so the ERP waveform is formed in the contralateral brain.As demonstrated by brain topography, motor coordination imagery involving the lower limbs is mainly activated in sensorimotor areas, which is due to motor coordination imagery of the lower limbs representing a form of motor imagery, the latter of which has been shown to activate sensorimotor areas [10]. In addition to activating sensorimotor areas, the visual area located in the occipital lobe also has corresponding activation, which may be because the cues provided to the subjects are processed through the visual pathway before being transmitted to the sensorimotor area [68].In the present study, power spectrum analysis showed that the power distribution related to MCI was distributed inα waves, while the power related to IEG was distributed in θ waves. Simultaneously, the absolute power of the IEG was much higher than the absolute power of the MCI. This result may be due to MCI being easier to complete than IEG, resulting in MCI inducing lower energy.Under the condition of motor coordination or imagery, external stimulation acts on sensory organs and causes nerve excitation. The activation of brain-related areas and the coordinated control between brain areas transmit such nerve excitation to synergetic and antagonistic muscles, such that each muscle’s activation is timed appropriately and ensures smooth completion of motor actions. Not only the activation of brain area but also the features of the brain network to represent the coordinated motion or imagination need to be calculated. In the present study, brain network analysis showed that brain functional connections during MCI, IEG, and REST mainly occurred in the right forehead, which may have been due to these psychological activities being related to the subconscious since human subconscious behavior is mainly related to the right forehead network [37, 53, 65, 69]. Our present brain network analysis also showed that the characteristic path length was a good distinguishing feature for MCI and IEG and that the other three measures were not. In addition, because the characteristic path length and the global efficiency determined the efficiency of the entire information transmission of the network [53], the characteristic path length and the global efficiency of the three types of tasks were quite different. IEG had the smallest characteristic path length and the highest global efficiency, so its network-information transfer rate was the fastest, followed by MCI; the slowest was REST, which may have been due to IEG needing to recruit more neural resources.In the present study, the features of MRCPs, time-frequency analysis, power spectrum analysis, network attribute analysis, a combination of the network attribute features, and the spatial features of the adjacency matrix in MCI and IEG were classified by SVM, and the average classification accuracies were 67.89 ± 9.14%, 78.12 ± 8.91%, 81.08 ± 8.99%, 89.12 ± 9.92%, and 92.96 ± 7.54%, respectively. The results showed that the average classification accuracy of the combination of the network attribute features and the spatial features of the adjacency matrix was 92.96%, confirming that these were the two best classification features. This result may be due to features representing collaborative information among related brain regions; also, the feature dimensionality was high, so the relevant features were extracted and retained comprehensively.Statistical analysis showed that there was a significant difference in the characteristic path lengths between MCI and REST, for which only the O2 channel was not statistically significant. There was no significant difference between the characteristic path lengths between IEG and REST, whereas thep values of P3, C3, CP4, and Pz channels of IEG were small, but the difference between MCI and IEG was significant. This result may be due to the fact that imagining MCI and REST is a relatively easy psychological activity, while imagining IEG is a relatively difficult psychological activity.In Figure5, there are some transition areas which are not clear by uncertainty. This is often more frequent in reality than one might think. We may segment each image obtained in order to delimit each activated area using a fuzzy procedure [70, 71].


## 5. Conclusion

In order to determine a rehabilitation method for lower-limb motor coordination in poststroke hemiplegic patients, a novel paradigm using lower-limb motor coordination imagery (normal gait imagery and poststroke hemiplegic gait imagery) was employed in the present study. Time-domain analysis, power spectrum analysis, time-frequency analysis, brain network analysis, and statistical analysis were used to investigate the neural mechanisms and decoding methods of motor coordination imaging of the lower limbs. The results showed that motor coordination imagery involving the lower limbs was mainly activated in sensorimotor areas and that the brain functional connection mainly occurred in the right forehead. The results also showed that the average classification accuracies of the combination feature of the network attributes and the spatial features of the adjacency matrix were 92.96% for normal gait imagery and poststroke hemiplegic gait imagery. Taken together, our findings may provide ideas for the rehabilitation of movement in poststroke hemiplegic patients based on motor coordination imagery involving the lower-limb brain-computer interface (BCI). Finally, our future directions are as follows: (1) EEG combined with fMRI to study the neural mechanisms of motor coordination imagery involving the lower limbs and (2) online verification of motor coordination imagery via BCI.



---

*Source: 5565824-2021-06-23.xml*