CN116312951A - 一种基于多模态耦合分析的运动功能评估方法及其系统 - Google Patents
一种基于多模态耦合分析的运动功能评估方法及其系统 Download PDFInfo
- Publication number
- CN116312951A CN116312951A CN202310464893.4A CN202310464893A CN116312951A CN 116312951 A CN116312951 A CN 116312951A CN 202310464893 A CN202310464893 A CN 202310464893A CN 116312951 A CN116312951 A CN 116312951A
- Authority
- CN
- China
- Prior art keywords
- data
- muscle
- brain
- motor
- function
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 230000008878 coupling Effects 0.000 title claims abstract description 160
- 238000010168 coupling process Methods 0.000 title claims abstract description 160
- 238000005859 coupling reaction Methods 0.000 title claims abstract description 160
- 238000004458 analytical method Methods 0.000 title claims abstract description 89
- 238000000034 method Methods 0.000 title claims abstract description 84
- 210000003205 muscle Anatomy 0.000 claims abstract description 240
- 210000004556 brain Anatomy 0.000 claims abstract description 185
- 230000007659 motor function Effects 0.000 claims abstract description 73
- 230000006870 function Effects 0.000 claims description 186
- 230000033001 locomotion Effects 0.000 claims description 162
- 238000012549 training Methods 0.000 claims description 115
- 210000005036 nerve Anatomy 0.000 claims description 98
- QVGXLLKOCUKJST-UHFFFAOYSA-N atomic oxygen Chemical compound [O] QVGXLLKOCUKJST-UHFFFAOYSA-N 0.000 claims description 64
- 229910052760 oxygen Inorganic materials 0.000 claims description 64
- 239000001301 oxygen Substances 0.000 claims description 64
- 238000011156 evaluation Methods 0.000 claims description 51
- 230000000694 effects Effects 0.000 claims description 40
- 238000004364 calculation method Methods 0.000 claims description 36
- 230000011218 segmentation Effects 0.000 claims description 33
- 230000003702 neurovascular coupling effect Effects 0.000 claims description 32
- 238000013210 evaluation model Methods 0.000 claims description 24
- 230000003183 myoelectrical effect Effects 0.000 claims description 22
- 238000007781 pre-processing Methods 0.000 claims description 22
- 238000010801 machine learning Methods 0.000 claims description 21
- 238000012795 verification Methods 0.000 claims description 18
- 238000002329 infrared spectrum Methods 0.000 claims description 14
- 238000005259 measurement Methods 0.000 claims description 14
- 230000001054 cortical effect Effects 0.000 claims description 13
- 238000009826 distribution Methods 0.000 claims description 12
- 238000005315 distribution function Methods 0.000 claims description 12
- 230000005611 electricity Effects 0.000 claims description 11
- 238000001228 spectrum Methods 0.000 claims description 11
- 208000016285 Movement disease Diseases 0.000 claims description 10
- 230000035484 reaction time Effects 0.000 claims description 10
- 230000004044 response Effects 0.000 claims description 10
- 238000013527 convolutional neural network Methods 0.000 claims description 9
- 230000001360 synchronised effect Effects 0.000 claims description 9
- 238000004590 computer program Methods 0.000 claims description 8
- 238000000605 extraction Methods 0.000 claims description 7
- 238000012360 testing method Methods 0.000 claims description 7
- 230000009471 action Effects 0.000 claims description 6
- 238000013145 classification model Methods 0.000 claims description 6
- 230000001427 coherent effect Effects 0.000 claims description 6
- 230000001537 neural effect Effects 0.000 claims description 6
- 238000005070 sampling Methods 0.000 claims description 6
- 238000003062 neural network model Methods 0.000 claims description 4
- 230000002232 neuromuscular Effects 0.000 claims description 4
- 238000005457 optimization Methods 0.000 claims description 3
- 230000008447 perception Effects 0.000 claims description 3
- 238000002360 preparation method Methods 0.000 claims description 3
- 238000007636 ensemble learning method Methods 0.000 claims 1
- 238000012545 processing Methods 0.000 abstract description 6
- 210000003414 extremity Anatomy 0.000 description 91
- 230000006378 damage Effects 0.000 description 16
- 210000003128 head Anatomy 0.000 description 10
- 230000002490 cerebral effect Effects 0.000 description 8
- 206010008190 Cerebrovascular accident Diseases 0.000 description 7
- 208000006011 Stroke Diseases 0.000 description 7
- 230000008569 process Effects 0.000 description 7
- 230000008859 change Effects 0.000 description 5
- 238000006243 chemical reaction Methods 0.000 description 5
- 238000001914 filtration Methods 0.000 description 5
- 230000004220 muscle function Effects 0.000 description 5
- 210000001364 upper extremity Anatomy 0.000 description 5
- 239000008280 blood Substances 0.000 description 4
- 210000004369 blood Anatomy 0.000 description 4
- 238000012937 correction Methods 0.000 description 4
- 210000003141 lower extremity Anatomy 0.000 description 4
- 239000000523 sample Substances 0.000 description 4
- 238000012706 support-vector machine Methods 0.000 description 4
- 208000012661 Dyskinesia Diseases 0.000 description 3
- 108010054147 Hemoglobins Proteins 0.000 description 3
- 102000001554 Hemoglobins Human genes 0.000 description 3
- 208000027418 Wounds and injury Diseases 0.000 description 3
- 208000014674 injury Diseases 0.000 description 3
- 210000003314 quadriceps muscle Anatomy 0.000 description 3
- INGWEZCOABYORO-UHFFFAOYSA-N 2-(furan-2-yl)-7-methyl-1h-1,8-naphthyridin-4-one Chemical compound N=1C2=NC(C)=CC=C2C(O)=CC=1C1=CC=CO1 INGWEZCOABYORO-UHFFFAOYSA-N 0.000 description 2
- 206010061296 Motor dysfunction Diseases 0.000 description 2
- 108010064719 Oxyhemoglobins Proteins 0.000 description 2
- 230000002159 abnormal effect Effects 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 2
- 230000006399 behavior Effects 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 2
- 210000004204 blood vessel Anatomy 0.000 description 2
- 238000004422 calculation algorithm Methods 0.000 description 2
- 230000001364 causal effect Effects 0.000 description 2
- 230000001149 cognitive effect Effects 0.000 description 2
- 230000003111 delayed effect Effects 0.000 description 2
- 108010002255 deoxyhemoglobin Proteins 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000004064 dysfunction Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 230000036284 oxygen consumption Effects 0.000 description 2
- 238000007637 random forest analysis Methods 0.000 description 2
- 238000011084 recovery Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 230000035945 sensitivity Effects 0.000 description 2
- 238000003860 storage Methods 0.000 description 2
- 238000002560 therapeutic procedure Methods 0.000 description 2
- 229920000742 Cotton Polymers 0.000 description 1
- LFQSCWFLJHTTHZ-UHFFFAOYSA-N Ethanol Chemical compound CCO LFQSCWFLJHTTHZ-UHFFFAOYSA-N 0.000 description 1
- 206010061225 Limb injury Diseases 0.000 description 1
- 206010049565 Muscle fatigue Diseases 0.000 description 1
- 206010049816 Muscle tightness Diseases 0.000 description 1
- 208000005392 Spasm Diseases 0.000 description 1
- 208000003443 Unconsciousness Diseases 0.000 description 1
- 230000006978 adaptation Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008344 brain blood flow Effects 0.000 description 1
- 230000003925 brain function Effects 0.000 description 1
- 210000003710 cerebral cortex Anatomy 0.000 description 1
- 230000001684 chronic effect Effects 0.000 description 1
- 238000004140 cleaning Methods 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 238000013500 data storage Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000018109 developmental process Effects 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 210000001061 forehead Anatomy 0.000 description 1
- 230000036541 health Effects 0.000 description 1
- 238000012880 independent component analysis Methods 0.000 description 1
- 230000010354 integration Effects 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 210000000337 motor cortex Anatomy 0.000 description 1
- 230000008035 nerve activity Effects 0.000 description 1
- 230000002035 prolonged effect Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000000284 resting effect Effects 0.000 description 1
- 230000002441 reversible effect Effects 0.000 description 1
- 210000004761 scalp Anatomy 0.000 description 1
- 210000002027 skeletal muscle Anatomy 0.000 description 1
- 210000003625 skull Anatomy 0.000 description 1
- 208000020431 spinal cord injury Diseases 0.000 description 1
- 210000001562 sternum Anatomy 0.000 description 1
- 208000024891 symptom Diseases 0.000 description 1
Images
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/40—Detecting, measuring or recording for evaluating the nervous system
- A61B5/4076—Diagnosing or monitoring particular conditions of the nervous system
- A61B5/4082—Diagnosing or monitoring movement diseases, e.g. Parkinson, Huntington or Tourette
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/103—Detecting, measuring or recording devices for testing the shape, pattern, colour, size or movement of the body or parts thereof, for diagnostic purposes
- A61B5/11—Measuring movement of the entire body or parts thereof, e.g. head or hand tremor, mobility of a limb
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/145—Measuring characteristics of blood in vivo, e.g. gas concentration, pH value; Measuring characteristics of body fluids or tissues, e.g. interstitial fluid, cerebral tissue
- A61B5/14542—Measuring characteristics of blood in vivo, e.g. gas concentration, pH value; Measuring characteristics of body fluids or tissues, e.g. interstitial fluid, cerebral tissue for measuring blood gases
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
- A61B5/369—Electroencephalography [EEG]
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
- A61B5/389—Electromyography [EMG]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N20/00—Machine learning
- G06N20/20—Ensemble learning
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H10/00—ICT specially adapted for the handling or processing of patient-related medical or healthcare data
- G16H10/60—ICT specially adapted for the handling or processing of patient-related medical or healthcare data for patient-specific data, e.g. for electronic patient records
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H20/00—ICT specially adapted for therapies or health-improving plans, e.g. for handling prescriptions, for steering therapy or for monitoring patient compliance
- G16H20/30—ICT specially adapted for therapies or health-improving plans, e.g. for handling prescriptions, for steering therapy or for monitoring patient compliance relating to physical therapies or activities, e.g. physiotherapy, acupressure or exercising
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
- Y02P90/00—Enabling technologies with a potential contribution to greenhouse gas [GHG] emissions mitigation
- Y02P90/30—Computing systems specially adapted for manufacturing
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Public Health (AREA)
- General Health & Medical Sciences (AREA)
- Physics & Mathematics (AREA)
- Biophysics (AREA)
- Veterinary Medicine (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- Pathology (AREA)
- Molecular Biology (AREA)
- Surgery (AREA)
- Animal Behavior & Ethology (AREA)
- Neurology (AREA)
- Physiology (AREA)
- Neurosurgery (AREA)
- Theoretical Computer Science (AREA)
- Software Systems (AREA)
- Primary Health Care (AREA)
- Epidemiology (AREA)
- Optics & Photonics (AREA)
- Data Mining & Analysis (AREA)
- Physical Education & Sports Medicine (AREA)
- Psychiatry (AREA)
- Oral & Maxillofacial Surgery (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Psychology (AREA)
- Evolutionary Computation (AREA)
- Dentistry (AREA)
- Computing Systems (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Developmental Disabilities (AREA)
- Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
Abstract
本发明涉及多模态数据处理技术领域，具体公开了一种基于多模态脑肌耦合分析的运动功能评估方法及其系统，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据，依次进行数据预处理、逐一特征提取和脑肌参数耦合分析，获得单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算肢体运动肌力功能指标，并结合机器学习方法计算运动神经环路功能指标，利用集成学习方法建立运动功能康复评估模型，实现从神经环路和运动能力两个方面对肢体运动障碍进行多模态脑肌耦合分析，以实现康复训练中运动功能状态的快速精准评估。
Description
技术领域
本发明涉及多模态数据处理技术领域，具体涉及一种基于多模态脑肌耦合分析的运动功能评估方法及其系统。
背景技术
脑卒中、脊髓损伤等疾病会导致肢体功能运动降低甚至完全丧失，康复治疗是恢复运动功能的有效途径，而运动功能评定是脑卒中患者肢体运动功能障碍诊断和康复治疗的前提和基础。传统的运动功能评定方法依靠医师的主观经验，其评定结果具有较强主观性，且评定流程复杂且繁琐，对康复医师水平要求居高。运动功能客观评定方法使用传感检测代替医师进行运动功能评定，可直接或间接反映患者运动功能的生理参数用于评定，相较于传统方法更直接客观。其中肌电、肌氧信号可间接反映患者的肌力水平，而脑卒中患者的肢体运动障碍根源是神经环路受损，可检测脑电、脑氧等生理信号体现患者的脑功能状态用于运动功能客观评定。然而目前大多数相关研究采用的生理参数数量少或模态较为单一，利用脑电信号或肌电信号进行单独评价只能间接从单一维度评定运动功能，不能较为准确地反映脑卒中患者的神经环路受损情况。
发明内容
鉴于以上所述现有方法的局限，本发明的目的在于提供一种基于多模态脑肌耦合分析的运动功能评估方法及其系统，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据，对所述多模态生理数据依次进行数据预处理、逐一特征提取和脑肌参数耦合分析，获得单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标，并结合机器学习方法计算得到运动神经环路功能指标，最后利用集成学习方法建立运动功能康复评估模型，实现从神经环路和运动能力两个方面对肢体运动障碍进行多模态脑肌耦合分析，以实现康复训练中运动功能状态的快速精准评估。
为了实现上述目的，根据本发明的一方面，提供一种基于多模态脑肌耦合分析的运动功能评估方法，所述方法包括以下步骤：
S100，采集多模态生理数据，包括神经环路生理数据和运动肌力生理数据；
S200，对所述多模态生理数据进行数据预处理，得到多模态真实数据；对所述多模态真实数据进行逐一特征提取，得到单生理运动参数；
S300，结合所述单生理运动参数，对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数；
S400，根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标；
S500，根据所述单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标；
S600，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型。
进一步地，在S100中，所述多模态生理数据包括神经环路生理数据和运动肌力生理数据，其中所述神经环路生理数据至少包括受试者特定脑区的脑电信号和脑氧数据，所述运动肌力生理数据至少包括受试者多个特定肌肉部位的肌电信号和肌氧数据。
进一步地，在S200中，所述数据预处理至少包括去噪处理、基线校正和去除伪迹；所述单生理运动参数包括：所述脑电信号的功率谱（PSD）、所述肌电信号的最大肌力值（MVC）和均方根（RMS）、所述肌氧数据的氧合血红蛋白浓度（HbO2）中的至少两种（本申请选择脑电、脑氧、肌电和肌氧四种生理参数信息进行多模态耦合分析，其中脑电信号是大脑皮层活动产生的电信号，具有较高的灵敏性，能够快速反映皮层活动状态，脑氧数据能够表征大脑皮层的血氧浓度变化，在进行主动意图和认知活动时会产生显著变化；而肌电和肌氧数据均与肌肉活动状态密切相关，可客观反映肢体运动功能状态）。
进一步地，在S300中，结合所述单生理运动参数，对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数的方法为：所述脑肌耦合运动参数，包括脑电脑氧相关可变量、脑电肌电相干耦合度量和肌电肌氧多模态相关系数，其计算方法包括：
S301，采用滑动窗口方法对所述多模态真实数据进行按段进行分割，按序得到各段分割多模态生理子数据；其中，设定对所述多模态真实数据分割成N个分割多模态生理子数据，N为常数；
S302，对所述各段分割多模态生理子数据进行神经血管耦合分析，包括在各段分割的采集时间段内计算所述脑电脑氧相关可变量记作VarECO，其计算公式为：
其中，x t表示为脑电信号在时刻为t值的采样值，t 为时间变量，各段分割多模态生理子数据对应的时间段范围为t∈[t1,t2]，t1为起始分割时刻，t2为结束分割时刻；y t表示为脑氧数据在t时刻对应的光强度值；F(xt+1, x t, y t)表示为计算脑电信号在t时刻、t+1时刻对应的x t+1、x t值与脑氧数据在t时刻的y t值之间的联合分布函数，F(xt+1 | x t)表示为计算脑电信号在t时刻对应下一时刻出现xt+1值的联合分布概率，F(y t | x t)表示为计算脑电信号在t时刻时对应的脑氧数据出现y t值的联合分布概率，F(xt+1 , y t | x t)表示为计算脑电信号在t+1时刻对应的x t+1值与脑电信号在t时刻对应的脑氧数据出现y t值的概率之间联合分布函数；log2(F(xt+1 , y t | x t) / F(xt+1 | x t)·F(y t | x t))表示为计算以2为底F(xt+1 , y t | x t) 与F(xt+1 | x t)·F(y t | x t)的比值的对数；
S303，对所述各段分割多模态生理子数据进行肌电肌氧耦合分析，包括结合所述单生理运动参数，计算脑电肌电相干耦合度量记作EMCC，其计算公式为：
其中，Pxy(f)表示为脑电信号与肌电数据的互功率谱，Pxx(f)表示为脑电信号的自功率谱，Pyy(f)表示为肌电数据的自功率谱；
S304，对所述各段分割多模态生理子数据进行肌电肌氧耦合分析，包括结合所述单生理运动参数，在时域上各段分割的采集时间段内计算肌电肌氧多模态相关系数记作MOCC，其计算公式为：
其中，x(t)为某一肌肉部位的肌电信号从起始分割时刻到t时刻的均方根值，各段分割多模态生理子数据对应的时间段范围为t∈[t1,t2]，t1为起始分割时刻，t2为结束分割时刻；
为所述肌电信号在t∈[t1,t2]的采集时间段内x(t)的平均值，y(t)为同一肌肉部位的肌氧数据中的HbO2浓度，/>
为所述肌氧数据在t∈[t1,t2]的采集时间段内y(t)的平均值（脑卒中患者患有肢体运动障碍的根本原因是运动神经环路损伤，现有技术方案通常采用单一生理信号数据或运动信号数据，并计算单一信号的特征进行康复评估，无法反映患者的神经环路恢复情况；在步骤S300中，对所述多模态真实数据进行脑肌参数耦合分析，选择脑电、脑氧、肌电和肌氧四种生理参数信息，相较于单一的信号数据能更全面且更完整表征由于脑部神经环路受损导致肢体运动障碍的耦合特征；其中通过计算脑电信号和脑氧数据之间的脑电脑氧相关可变量，能够度量由脑电信号相关源引起的脑氧信号变化的程度，可更准确地表征受试者的运动意图状态；计算脑电信号与肌电数据之间的脑电肌电相干耦合度量，从而表征脑电信号与肌电数据在频域上的关联性，在一定程度上表征大脑皮层神经环路与肢体肌肉肌力之间的功能连接关系；进一步计算肌电信号与肌氧数据之间的动态关联特征如肌电肌氧多模态相关系数，用于度量肌电信号与肌氧数据在时域上的相关性，能够直接体现肢体肌肉的活动状态；利用多模态真实数据耦合反映大脑功能区域与肢体运动肌力之间的病症表达信息，与肢体运动障碍病症的关联度更高，以达到进一步提高对患者运动功能评估准确度的目的）。
进一步地，在S400中，根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标的方法为：其中，所述运动神经环路功能指标包括运动神经环路传导时间、神经肌肉协同程度、主动运动反应时间和主动运动持续时间、以及主动运动参与度，其计算方法具体包括：
S401，所述脑肌耦合运动参数包括所述脑电脑氧相关可变量和脑电肌电相干耦合度量，当所述各段分割多模态生理子数据对应的脑肌耦合运动参数大于等于对应的脑肌耦合特征阈值时，根据该段分割多模态生理子数据对应的起始分割时刻计算得到运动神经环路传导时间记作△T；
S402，根据所述各段分割多模态生理子数据进行运动肌肉活动分析，结合机器学习方法构建运动功能预测子模型，由各段分割多模态生理子数据对应输出的预测值按序构成运动肌肉活动曲线；
S403，根据所述各段分割多模态生理子数据对应的脑肌耦合运动参数，结合机器学习方法构建运动功能预测子模型，由各段分割多模态生理子数据对应输出的预测值按序构成神经血管耦合曲线；
S404，对所述神经血管耦合曲线和运动肌肉活动曲线进行运动耦合分析，计算从神经血管到肌肉活动方向的相关可变量，得到神经肌肉协同程度记作VarNMS，其计算公式为：
其中，ni表示为神经血管耦合曲线在所述各段分割多模态生理子数据的第i个分段对应的预测值，i为分段序号， i∈[1,N]，N为所述分割多模态生理子数据的总个数；mi表示为运动肌肉活动曲线在所述各段分割多模态生理子数据的第i个分段对应的预测值，F(ni+1, ni, mi)表示为计算神经血管耦合曲线在第i个分段、第i+1个分段对应的预测值与运动肌肉活动曲线在第i个分段的预测值之间的联合分布函数，F(ni+1 | ni)表示为计算神经血管耦合曲线在第i个分段对应下一个分段出现ni+1值的联合分布概率，F(mi | ni)表示为计算神经血管耦合曲线在第i个分段为ni值时对应的运动肌肉活动曲线出现mi值的联合分布概率，F(ni+1 , mi | ni)表示为计算神经血管耦合曲线在第i个分段为ni值与其对应的运动肌肉活动曲线在第i个分段出现mi值的概率之间的联合分布函数；log2表示为计算以2为底的对数；
S405，所述脑肌耦合运动参数包括所述肌电肌氧多模态相关系数，当所述各段分割多模态生理子数据对应的肌电肌氧多模态相关系数大于等于对应的肌电肌氧耦合特征阈值时，根据该段分割多模态生理子数据对应的起始分割时刻计算得到主动运动反应时间和主动运动持续时间，分别记作Tr和Td；
S406，计算主动运动持续时间与所述脑肌耦合运动参数对应的采集时间总长的比值，得到主动运动参与度记作Sp，其计算公式为：
其中，T为所述脑肌耦合运动参数对应的采集时间总长（由于脑卒中肢体运动障碍患者易出现无意识运动、痉挛、运动困难等神经-肌肉失调行为，所述运动神经环路传导时间可表征神经环路受损程度，对所述神经血管耦合曲线和运动肌肉活动曲线进行运动耦合分析，得到神经肌肉协同程度可反映神经活动与肌肉行为的脑肌协同耦合程度；以及受试者在进行主动运动康复训练时，神经血管耦合程度和肌电肌氧耦合程度会同步发生变化并保持一段时间，由于脑卒中患者的神经环路损伤和肢体控制能力减弱，计算其主动运动反应时间和主动运动持续时间，能表征皮层肌肉同步变化的持续时间，在一定程度上能反映其运动功能特征；而计算主动运动参与度可表征受试者在整个主动运动过程中，自主发力的时间占比；计算上述特征构成所述运动神经环路功能指标，能更客观准确地反映受试者的运动神经环路损伤程度）。
进一步地，在S401中，所述根据该段分割多模态生理子数据对应的起始分割时刻计算得到运动神经环路传导时间的方法为：
S4011，当所述各段分割多模态生理子数据对应的脑电脑氧相关可变量大于等于预设的神经血管耦合特征阈值时，判定该段出现主动运动意图，记录该段分割多模态生理子数据对应的起始分割时刻为T1；
S4012，当所述各段分割多模态生理子数据对应的脑电肌电相干耦合度量大于等于预设的皮层肌肉耦合特征阈值时，判定该段出现准备运动动作，记录该段分割多模态生理子数据对应的起始分割时刻为T2；
S4013，计算运动神经环路传导时间记作△T，其中△T= | T2-T1 |（脑卒中肢体运动障碍患者由于神经环路损伤，神经环路受损，易导致神经环路传导延迟甚至中断；其主动运动意图产生与实际肌肉运动之间存在会更长的传导时间，因此计算运动神经环路传导时间可作为神经环路受损程度的指标之一）。
进一步地，在S405中，所述根据该段分割多模态生理子数据对应的起始分割时刻计算得到主动运动反应时间和主动运动持续时间的方法为：
S4051，当所述各段分割多模态生理子数据对应的肌电肌氧多模态相关系数大于等于预设的肌电肌氧耦合特征阈值时，判定该段出现实际肌肉活动，记录该段分割多模态生理子数据对应的起始分割时刻为T3；
S4052，根据步骤S4011中判定出现主动运动意图对应的分割多模态生理子数据的起始分割时刻为T1，计算主动运动反应时间记作Tr，Tr = | T3-T1 | ；
S4053，当所述各段分割多模态生理子数据满足对应的肌电肌氧多模态相关系数大于等于对应的肌电肌氧耦合特征阈值，且对应的脑电脑氧相关可变量大于等于预设的神经血管耦合特征阈值时，统计对应的各段分割多模态生理子数据的采集时间总和，得到主动运动持续时间记作Td（受试者在进行主动运动康复训练时，由于存在较强的主观运动意志会导致脑部血流量增加，引起神经血管的耦合程度发生变化；在一段反应时间后，其肌肉发力、肌肉神经和血氧浓度都会发生变化，其中肌氧和肌电信号的耦合程度逐渐增大。因此在主动运动期间，神经血管耦合程度和肌电肌氧耦合程度会同步发生变化并保持一段时间，由于脑卒中患者的神经环路损伤和肢体控制能力减弱，计算其主动运动的反应时间、皮层肌肉同步变化的持续时间在一定程度上能反映其运动功能辅助评估）。
进一步地，在S500中，根据所述单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标的方法为：
S501，根据所述单生理运动参数，利用所述肌电信号的最大肌力值对所述运动肌力生理数据中的肌电信号的幅值进行归一化，得到相对肌力值，并计算在整个采集时间段内所述相对肌力值的总和得到相对肌力值积分；
S502，在整个采集时间段内计算受试者的左侧特定肌肉部位对应的相对肌力值平均值与右侧特定肌肉部位对应的相对肌力值平均值的比值得到发力偏侧性指标；
S503，得到肢体运动肌力功能指标，包括所述相对肌力值积分和发力偏侧性指标（在受试者进行上下肢联动康复训练时，相对肌力值积分可表示多个特定肌肉部位如左右侧上斜方肌、肱二头肌、肱三头肌、股直肌、股二头肌、胫前肌、腓肠肌外侧头肌等肌肉的持续发力能力，而发力偏侧性指标可反映其各特定肌肉部位的肌肉发力是否具有较强偏侧性；计算所述肢体运动肌力功能指标能反映肢体运动障碍的肢体运动能力）。
进一步地，在S600中，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型的方法为：
S601，利用Brunnstrom评定量表对具有不同肢体运动障碍部位及功能程度的多个受试者进行运动功能水平跟踪预测，记录各个受试者在各个康复训练周期对应的Brunnstrom分期水平，并对各个受试者完成步骤S100~S500得到对应的所述运动神经环路功能指标和所述肢体运动肌力功能指标；
S602，根据各个受试者对应的所述运动神经环路功能指标和所述肢体运动肌力功能指标作为数据集，并按照设定比例划分为训练数据集、测试数据集和验证数据集；利用机器学习方法对进行所述训练数据集进行分类训练，构建初级分类模型；
S603，对所述初级分类模型采用十折交叉验证方法或留一法验证进行模型验证，选取每个验证数据集上对应识别准确率最高的模型组成集成初级模型；
S604，利用 BP神经网络模型、多层感知网络(MLP)、卷积神经网络(CNN)作为集成次级模型，将所述集成初级模型的输出端作为所述集成次级模型的输入端，训练得到最终决策模型，即为运动功能康复评估初级模型；
S605，利用所述测试数据集对所述运动功能康复评估初级模型进行优化训练，得到最优的运动功能康复评估模型（现有技术方案用于肢体运动功能障碍的运动功能客观指标评测如最大肌张力、肢体最大屈伸度等，多数是基于力学或运动学信息的，且需要专业人员和设备才能进行评测，使用场景受限于大型医疗机构，操作流程复杂；本申请结合所述脑肌耦合分析得到的运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型，能够实现在康复训练过程中持续对受试者的神经环路损伤程度和肢体运动能力进行客观的运动功能评价，为制定个性化的康复决策提供客观依据，且不需要额外的测试设备，操作流程简单，可为社区医疗机构和居家场景提供高效、便捷的评测服务）。
进一步地，在S600中，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型的方法还包括：
S606，利用验证数据集对所述运动功能康复评估模型进行评估得到各个受试者的运动功能评价结果，所述各个受试者的运动功能评价结果与对应的Brunnstrom分期水平进行比对，进一步根据运动神经环路功能指标和肢体运动肌力功能指标对受试者的康复训练进行自适应调整训练参数。
进一步地，所述一种基于多模态脑肌耦合分析的运动功能评估方法还包括：
S700，将各个受试者在不同康复训练过程中对应的运动神经环路功能指标和肢体运动肌力功能指标输入所述运动功能康复评估模型，得到在不同康复训练过程对应的运动功能评价结果。
根据本发明的另一方面，还提供一种基于多模态脑肌耦合分析的运动功能评估系统，所述系统包括：存储器、处理器以及存储在所述存储器中并在所述处理器上运行的计算机程序；所述处理器执行所述计算机程序运行在以下系统的模块中：
多模态生理数据采集模块，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据；
生理数据预处理模块，用于对所述多模态生理数据进行数据预处理，得到多模态真实数据；
逐一特征提取模块，用于对所述多模态真实数据进行逐一特征提取，得到单生理运动参数；
脑肌参数耦合分析模块，用于对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数；
皮层肌力功能指标计算模块，用于根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标；并结合所述单生理运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标；
运动功能康复评估模块，用于结合所述运动神经环路功能指标和所述肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型；
康复训练自适应调整模块，用于根据所述运动神经环路功能指标和肢体运动肌力功能指标自适应调整康复训练参数，根据所述运动功能康复评估模型自适应调整受试者的康复训练方案。
如上所述，本发明所述的一种基于多模态脑肌耦合分析的运动功能评估方法及系统，具有以下有益效果：（1）多模态生理数据进行耦合分析获得运动神经环路功能指标和肢体运动肌力功能指标，能更全面且更完整表征由于脑部神经环路受损导致肢体运动障碍的耦合特征；（2）相较于单一生理信号特征与肢体运动障碍病症的关联度更高，有效避免单一生理参数测量误差导致的噪声干扰；（3）结合所述脑肌耦合分析得到的运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型，能够实现在康复训练过程中持续对受试者的神经环路损伤程度和肢体运动能力进行客观的运动功能评价；（4）达到进一步提高对患者运动功能评估准确度的目的。
附图说明
通过对结合附图所示出的实施方式进行详细说明，本发明的上述以及其他特征将更加明显，本发明附图中相同的参考标号表示相同或相似的元素，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图，在附图中：
图1所示为一种基于多模态脑肌耦合分析的运动功能评估方法于一实施例中的流程图；
图2所示为一种基于多模态脑肌耦合分析的运动功能评估系统于一实施例中的系统结构图；
图3所示为一种多模态生理数据采集系统于一实施例中受试者在主动运动期间进行同步数据采集获得多模态生理数据的康复训练过程示意图。
具体实施方式
以下将结合实施例和附图对本发明的构思、具体结构及产生的技术效果进行清楚、完整的描述，以充分地理解本发明的目的、方案和效果。需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互组合。
如图1所示为根据本发明的一种基于多模态脑肌耦合分析的运动功能评估方法的流程图，下面结合图1来阐述根据本发明的实施方式的一种基于多模态脑肌耦合分析的运动功能评估方法。本发明提出一种基于多模态脑肌耦合分析的运动功能评估方法，所述方法具体包括以下步骤：
S100，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据，包括神经环路生理数据和运动肌力生理数据；
S200，对所述多模态生理数据进行数据预处理，得到多模态真实数据；对所述多模态真实数据进行逐一特征提取，得到单生理运动参数；
S300，结合所述单生理运动参数，对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数；
S400，根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标；
S500，根据所述单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标；
S600，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型。
进一步地，在S100中，所述多模态生理数据包括神经环路生理数据和运动肌力生理数据，其中所述神经环路生理数据至少包括受试者特定脑区的脑电信号和脑氧数据，所述运动肌力生理数据至少包括受试者多个特定肌肉部位的肌电信号和肌氧数据。
优选地，在本具体实施例中，受试者在主动运动期间利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据；所述多模态生理数据采集系统至少包括脑电信号采集设备、近红外光谱脑氧采集设备、近红外光谱肌氧采集设备和表面肌电采集设备。
进一步地，在S200中，所述数据预处理至少包括去噪处理、基线校正和去除伪迹；所述单生理运动参数包括：所述脑电信号的功率谱（PSD）、所述肌电信号的最大肌力值（MVC）和均方根（RMS）、所述肌氧数据的氧合血红蛋白浓度（HbO2）中的至少两种。
优选地，在本具体实施例中，所述脑电信号是在受试者的主动运动状态下进行全导联采集的，因此相较于静息态采集易出现坏导情况，其中混合肌电、眼电、心电信号的伪迹也较大。因为，对所述脑电信号的数据预处理包括：对全导联脑电信号进行0.4~40Hz的带通滤波，去除直流偏置和非脑电频段（一般脑电信号在uV级别，频率不超过40Hz，因此采集时需要给出较强的直流偏置增大整体幅值来保证信号稳定性，在后续信号处理中需要把直流偏置去除），对幅值异常的数据段进行插值处理，剔除波动异常的通道，并对数据进行重参考，最后使用独立成分分析法去除肌电、眼电、心电伪迹，得到预处理后的脑电信号。
优选地，在本具体实施例中，所述肌电信号相较于脑电信号，其幅值更高，信号更稳定，能量集中在40~150Hz。其中，在所述多模态生理数据采集系统进行数据采集时，保证表面肌电采集设备的采集端（如穿戴式绑带）紧贴皮肤，提高数据传输质量。对所述肌电信号进行数据预处理包括：基线校正、50Hz工频陷波、51~150Hz的带通滤波，之后按打标点对预处理后的肌电信号进行分段，得到预处理后的肌电信号。具体地，若肌电信号中发现坏导，则考虑坏导是否出于患肢及其对侧肢的肌电数据，若是，则本次数据不可使用，若不是，则去除坏导后可以继续进行分析。
优选地，在本具体实施例中，所述脑氧数据由高精度的近红外光谱采集设备进行采集，使用连续波照射模式，仅测量通过头部皮层散射的光强度值，然后通过校正的朗伯-比尔定律计算出氧合血红蛋白和脱氧血红蛋白的浓度变化得到脑氧数据。所述脑氧数据的数据预处理包括：使用小波滤波去除运动伪影，使用0.5Hz低通滤波和0.01Hz高通滤波去除生理噪声，得到预处理后的脑氧数据。
优选地，在本具体实施例中，所述肌氧数据可使用相对较低精度的近红外光谱采集设备进行采集，一般不影响数据分析结果。其中，所述肌氧数据的数据预处理包括：对所述肌氧数据进行基线校正，采用自适应滤波器去除运动干扰，利用样条插值对伪迹进行矫正，得到预处理后的肌氧数据。
具体地，所述多模态真实数据包括预处理后的脑电信号、肌电信号、脑氧数据和肌氧数据。
进一步地，在S300中，结合所述单生理运动参数，对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数的方法为：所述脑肌耦合运动参数，包括脑电脑氧相关可变量、脑电肌电相干耦合度量和肌电肌氧多模态相关系数，其计算方法包括：
S301，采用滑动窗口方法对所述多模态真实数据进行按段进行分割，按序得到各段分割多模态生理子数据；其中，设定对所述多模态真实数据分割成N个分割多模态生理子数据，N为常数；
S302，对所述各段分割多模态生理子数据进行神经血管耦合分析，包括在各段分割的采集时间段内计算所述脑电脑氧相关可变量记作VarECO，其计算公式为：
其中，x t表示为脑电信号在时刻为t值的采样值，t 为时间变量，各段分割多模态生理子数据对应的时间段范围为t∈[t1,t2]，t1为起始分割时刻，t2为结束分割时刻；y t表示为脑氧数据在t时刻对应的光强度值，F(xt+1, x t, y t)表示为计算脑电信号在t时刻、t+1时刻对应的x t+1、x t值与脑氧数据在t时刻的y t值之间的联合分布函数，F(xt+1 | x t)表示为计算脑电信号在t时刻对应下一时刻出现xt+1值的联合分布概率，F(y t | x t)表示为计算脑电信号在t时刻时对应的脑氧数据出现y t值的联合分布概率，F(xt+1 , y t | x t)表示为计算脑电信号在t+1时刻对应的x t+1值与脑电信号在t时刻对应的脑氧数据出现y t值的概率之间联合分布函数；
S303，对所述各段分割多模态生理子数据进行肌电肌氧耦合分析，包括结合所述单生理运动参数，计算脑电肌电相干耦合度量记作EMCC，其计算公式为：
其中，Pxy(f)表示为脑电信号与肌电数据的互功率谱，Pxx(f)表示为脑电信号的自功率谱，Pyy(f)表示为肌电数据的自功率谱；
S304，对所述各段分割多模态生理子数据进行肌电肌氧耦合分析，包括结合所述单生理运动参数，在时域上各段分割的采集时间段内计算肌电肌氧多模态相关系数记作MOCC，其计算公式为：
其中，x(t)为某一肌肉部位的肌电信号从起始分割时刻到t时刻的均方根值，各段分割多模态生理子数据对应的时间段范围为t∈[t1,t2]，t1为起始分割时刻，t2为结束分割时刻；
为所述肌电信号在t∈[t1,t2]的采集时间段内x(t)的平均值，y(t)为同一肌肉部位的肌氧数据中的HbO2浓度，/>
为所述肌氧数据在t∈[t1,t2]的采集时间段内y(t)的平均值。
进一步地，在S400中，根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标的方法为：其中，所述运动神经环路功能指标包括运动神经环路传导时间、神经肌肉协同程度、主动运动反应时间和主动运动持续时间、以及主动运动参与度，其计算方法具体包括：
S401，所述脑肌耦合运动参数包括所述脑电脑氧相关可变量和脑电肌电相干耦合度量，当所述各段分割多模态生理子数据对应的脑肌耦合运动参数大于等于对应的脑肌耦合特征阈值时，根据该段分割多模态生理子数据对应的起始分割时刻计算得到运动神经环路传导时间记作△T；
S402，根据所述各段分割多模态生理子数据进行运动肌肉活动分析，结合机器学习方法构建运动功能预测子模型，由各段分割多模态生理子数据对应输出的预测值按序构成运动肌肉活动曲线；
S403，根据所述各段分割多模态生理子数据对应的脑肌耦合运动参数，结合机器学习方法构建运动功能预测子模型，由各段分割多模态生理子数据对应输出的预测值按序构成神经血管耦合曲线；
S404，对所述神经血管耦合曲线和运动肌肉活动曲线进行运动耦合分析，计算从神经血管到肌肉活动方向的相关可变量，得到神经肌肉协同程度记作VarNMS，其计算公式为：
其中，ni表示为神经血管耦合曲线在所述各段分割多模态生理子数据的第i个分段对应的预测值，i为分段序号， i∈[1,N]，N为所述分割多模态生理子数据的总个数；mi表示为运动肌肉活动曲线在所述各段分割多模态生理子数据的第i个分段对应的预测值，F(ni+1, ni, mi)表示为计算神经血管耦合曲线在第i个分段、第i+1个分段对应的预测值与运动肌肉活动曲线在第i个分段的预测值之间的联合分布函数，F(ni+1 | ni)表示为计算神经血管耦合曲线在第i个分段对应下一个分段出现ni+1值的联合分布概率，F(mi | ni)表示为计算神经血管耦合曲线在第i个分段为ni值时对应的运动肌肉活动曲线出现mi值的联合分布概率，F(ni+1 , mi | ni)表示为计算神经血管耦合曲线在第i个分段为ni值与其对应的运动肌肉活动曲线在第i个分段出现mi值的概率之间的联合分布函数；log2表示为计算以2为底的对数；
S405，所述脑肌耦合运动参数包括所述肌电肌氧多模态相关系数，当所述各段分割多模态生理子数据对应的肌电肌氧多模态相关系数大于等于对应的肌电肌氧耦合特征阈值时，根据该段分割多模态生理子数据对应的起始分割时刻计算得到主动运动反应时间和主动运动持续时间，分别记作Tr和Td；
S406，计算主动运动持续时间与所述脑肌耦合运动参数对应的采集时间总长的比值，得到主动运动参与度记作Sp，其计算公式为：
其中，T为所述脑肌耦合运动参数对应的采集时间总长。
优选地，所述机器学习方法包括但不限于支持向量机(SVM)、随机森林、K近邻算法（KNN）、多层感知网络(MLP)、卷积神经网络(CNN)等算法，用于构建运动功能预测子模型。
进一步地，在S401中，所述根据该段分割多模态生理子数据对应的起始分割时刻计算得到运动神经环路传导时间的方法为：
S4011，当所述各段分割多模态生理子数据对应的脑电脑氧相关可变量大于等于预设的神经血管耦合特征阈值时，判定该段出现主动运动意图，记录该段分割多模态生理子数据对应的起始分割时刻为T1；
S4012，当所述各段分割多模态生理子数据对应的脑电肌电相干耦合度量大于等于预设的皮层肌肉耦合特征阈值时，判定该段出现准备运动动作，记录该段分割多模态生理子数据对应的起始分割时刻为T2；
S4013，计算运动神经环路传导时间记作△T，其中△T= | T2-T1 |。
进一步地，在S405中，所述根据该段分割多模态生理子数据对应的起始分割时刻计算得到主动运动反应时间和主动运动持续时间的方法为：
S4051，当所述各段分割多模态生理子数据对应的肌电肌氧多模态相关系数大于等于预设的肌电肌氧耦合特征阈值时，判定该段出现实际肌肉活动，记录该段分割多模态生理子数据对应的起始分割时刻为T3；
S4052，根据步骤S4011中判定出现主动运动意图对应的分割多模态生理子数据的起始分割时刻为T1，计算主动运动反应时间记作Tr，Tr = | T2-T1 |；
S4053，当所述各段分割多模态生理子数据满足对应的肌电肌氧多模态相关系数大于等于对应的肌电肌氧耦合特征阈值，且对应的脑电脑氧相关可变量大于等于预设的神经血管耦合特征阈值时，统计对应的各段分割多模态生理子数据的采集时间总和，得到主动运动持续时间记作Td。
优选地，所述神经血管耦合特征阈值、皮层肌肉耦合特征阈值和肌电肌氧耦合特征阈值可通过置信度计算得到，其计算公式为：
其中，K表示为阈值大小，N所述分割多模态生理子数据的总个数，α为置信度，α的取值范围在0~1之间，不同受训者的取值不尽相同，具体数值需预先进行主、被动运动实验数据对比确定。
进一步地，在S500中，根据所述单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标的方法为：
S501，根据所述单生理运动参数，利用所述肌电信号的最大肌力值对所述运动肌力生理数据中的肌电信号的幅值进行归一化，得到相对肌力值，并计算在整个采集时间段内所述相对肌力值的总和得到相对肌力值积分；
S502，在整个采集时间段内计算受试者的左侧特定肌肉部位对应的相对肌力值平均值与右侧特定肌肉部位对应的相对肌力值平均值的比值得到发力偏侧性指标；
S503，得到肢体运动肌力功能指标，包括所述相对肌力值积分和发力偏侧性指标。
进一步地，在S600中，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型的方法为：
S601，利用Brunnstrom评定量表对具有不同肢体运动障碍部位及功能程度的多个受试者进行运动功能水平跟踪预测，记录各个受试者在各个康复训练周期对应的Brunnstrom分期水平，并对各个受试者完成步骤S100~S500得到对应的所述运动神经环路功能指标和所述肢体运动肌力功能指标；
S602，根据各个受试者对应的所述运动神经环路功能指标和所述肢体运动肌力功能指标作为数据集，并按照设定比例划分为训练数据集、测试数据集和验证数据集；利用机器学习方法对进行所述训练数据集进行分类训练，构建初级分类模型；优选地，所述机器学习方法包括但不限于支持向量机(SVM)、随机森林、K近邻算法（KNN）、多层感知网络(MLP)、卷积神经网络(CNN)等分类算法；
S603，对所述初级分类模型采用十折交叉验证方法或留一法验证进行模型验证，选取每个验证数据集上对应识别准确率最高的模型组成集成初级模型；
S604，利用 BP神经网络模型、多层感知网络(MLP)、卷积神经网络(CNN)作为集成次级模型，将所述集成初级模型的输出端作为所述集成次级模型的输入端，训练得到最终决策模型，即为运动功能康复评估初级模型；
S605，利用所述测试数据集对所述运动功能康复评估初级模型进行优化训练，得到最优的运动功能康复评估模型。
进一步地，在S600中，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型的方法还包括：
S606，利用验证数据集对所述运动功能康复评估模型进行评估得到各个受试者的运动功能评价结果，与所述各个受试者对应的Brunnstrom分期水平进行比对，进一步根据运动神经环路功能指标和肢体运动肌力功能指标对受试者的康复训练进行自适应调整训练参数，提高受试者主动运动康复训练的舒适度和效率。
优选地，在本具体实施例中，在步骤S606中，进一步根据运动神经环路功能指标和肢体运动肌力功能指标对受试者的康复训练进行自适应调整训练参数的具体步骤包括：
S6061，调整康复训练处方：比对各个受试者在上一次康复训练对应的Brunnstrom分期水平与各个受试者本次康复训练对应的运动功能评价结果，当Brunnstrom分期水平与运动功能评价结果不一致时，上传本次康复训练对应的多模态真实数据和运动功能评价结果至上位机模块，由康复医师进行进一步测评，若康复医师的测评结果与所述运动功能评价结果一致，则依据受试者在本次康复训练对应的Brunnstrom分期水平选择预先制定好的康复训练处方；若测评结果与所述运动功能评价结果不一致，则以测评结果为准，并将本次康复训练对应的多模态真实数据保存为勘误数据集，用于对运动功能康复评估模型进行优化（现有的康复训练处方制定方法是由康复医师通过运动功能测评，了解患者的肢体运动功能损伤程度和损伤部位后进行制定，并按康复周期进行调整；但各个患者的康复进度存在较大的个体差异性，需要对个体定制康复训练处方；而现有方案中缺少对患者运动功能恢复的跟踪评测方法，康复训练处方难以及时作出个体化调整。本申请根据受试者在康复训练中的生理状态表现，预测其目前的运动功能水平，实现个体的运动功能跟踪测评，在运动功能水平发生变化时能及时调整训练处方提高康复效率，在运动功能水平未变化时，针对受试者的神经反应能力和薄弱肢体部位，制定个性化专项训练任务，加快康复进度）；
S6062，调整主动模式训练时长：根据所述运动神经环路功能指标中的主动运动参与度在允许范围内适当调整主动模式训练时长，若检测到一段时间内主动运动参与度过低，则延长主动模式训练时长；具体地，同步进行语音提示，提醒受训者保持专注训练；
S6063，调整主动运动阻抗：根据所述肢体运动肌力功能指标中各个特定肌肉部位对应的相对肌力值积分，适当调整主动运动阻抗，若一段时间内的整体相对肌力值积分过高，则适当降低阻抗；若一段时间内的整体相对肌力值积分过低，则适当加大阻抗；具体地，在进行阻抗调整前，提前30秒使用语音进行预告；
S6064，调整特定肢体专项训练：根据所述肢体运动肌力功能指标中各个特定肌肉部位对应的发力偏侧性指标，在完成本次康复训练任务后，对发力较少的肢体进行专项训练，逐步提升受训者对偏弱肢体的控制能力；
S6065，调整反应能力专项训练：设计康复训练游戏对受试者进行反应能力训练，逐步提高受训者对薄弱肢体的控制能力，实现运动神经环路的恢复与重构；具体地，受试者需在康复训练游戏在屏幕出现特定画面停留的一段时间内，完成一次运动动作，成功完成则加分，延迟完成或未完成则不加分，根据所述运动神经环路功能指标中的运动神经环路传导时间和主动运动反应时间来确定特定画面的停留时间；其中运动神经环路传导时间和主动运动反应时间越短，则停留时间越短，反之则越长（现有的康复训练方案采用固定的训练程序进行训练，其各个训练参数无法根据受试者的实际训练状态进行实时调整，容易出现受试者过疲劳和训练量不够的情况，且在受试者注意力不集中时无法及时提醒，训练效率低；本申请提出了根据运动功能评价结果结合运动神经环路功能指标和肢体运动肌力功能指标对康复训练参数进行自适应调整的方法，根据实时测量得到的相对肌力值积分和主动运动参与度，自适应调整主动运动阻抗，防止受训者过疲劳或训练量不达标；并在主动运动参与度较低的时候，发出语音提醒，并适当延长主动运动训练持续时间，提高受训者的训练舒适度，训练效率高）。
如图2所示为本发明的一种基于多模态脑肌耦合分析的运动功能评估系统结构图，在本具体实施例中，提出了一种基于多模态脑肌耦合分析的运动功能评估系统，利用多模态生理数据采集模块、生理数据预处理模块、逐一特征提取模块和脑肌参数耦合分析模块进行多模态脑肌耦合分析，可根据受训者不同肢体部位、障碍程度以及训练表现，利用脑肌参数耦合分析模块和皮层肌力功能指标计算模块对肢体运动功能水平进行评测，能够实现康复训练处方和训练参数的自适应调整，以及个性化专项训练内容推荐，进行个性化的康复训练服务，提高受训者康复训练的舒适度和康复效率，适用于康复医师不足的社区医疗机构、老人院及居家场景。
其中，所述系统包括：存储器、处理器以及存储在所述存储器中并在所述处理器上运行的计算机程序；所述处理器执行所述计算机程序运行在以下系统的模块中：
多模态生理数据采集模块，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据；
生理数据预处理模块，用于对所述多模态生理数据进行数据预处理，得到多模态真实数据；
逐一特征提取模块，用于对所述多模态真实数据进行逐一特征提取，得到单生理运动参数；
脑肌参数耦合分析模块，用于对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数；
皮层肌力功能指标计算模块，用于根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标；并结合所述单生理运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标；
运动功能康复评估模块，用于结合所述运动神经环路功能指标和所述肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型；
康复训练自适应调整模块，用于根据所述运动神经环路功能指标和肢体运动肌力功能指标自适应调整康复训练参数，根据所述运动功能康复评估模型自适应调整受试者的康复训练方案。
优选地，在本具体实施例中，所述多模态生理数据采集模块，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据；所述多模态生理数据采集系统至少包括脑电信号采集设备、近红外光谱脑氧采集设备、近红外光谱肌氧采集设备和表面肌电采集设备。其中，所述多模态生理数据采集系统通过局域网连接上位机建立数据传输，由上位机进行无线同步控制和数据存储，所述多模态生理数据采集模块获取的多模态生理数据传输至上位机的数据预处理模块进行信号预处理，得到多模态真实数据。
具体地，在本具体实施例中，所述多模态生理数据采集系统将脑电、脑氧、肌电和肌氧的采集设备集成一体化，均满足非侵入、空间小、采集成本低的特点，可方便集成入现有的康复训练设备中，实现采集、评估、训练一体化，操作简单更适合社区康复和居家康复，能够满足普通社区康复医疗工作人员经过简单培训便可操作使用的需求。
具体地，在本具体实施例中，所述数据预处理模块用于对采集到的脑电、脑氧、肌电和肌氧信号数据分别进行预处理，根据每种生理信号特点进行不同的预处理方法，能够去除信号数据的各种伪迹，保留有效的信号数据用于进行逐一特征提取和脑肌参数耦合分析。
优选地，在另一种可能的实施例中，在进行获取受试者的多模态生理数据之前，应对受试者的基本身体信息和既往病史等进行登记，包括年龄、性别、身高、体重、运动功能量表评分、既往慢性病史、近期不适症状、是否有肢体运动不协调现象等，用于受试者的长期随时日记和肢体功能障碍病症发展的评估监护使用。
优选地，在另一种可能的实施例中，所述脑电信号采集设备的电极装置的脑电采集电极头与所述近红外光谱脑氧采集设备中采集脑氧信号的近红外探头分别交叉安放，集成在一个头戴式的EEG-fNIRS联合采集帽上，以覆盖待测者的大脑全头皮区域，也可优化安置区域集中在大脑左右半球的前额叶运动皮层区域，因为该特定脑区主要控制自主运动，与主动运动意图的产生高度相关。其中所述脑电采集电极头包括分为湿式电极头和干式电极头，优选使用湿式电极头。
具体地，脑电信号是大脑皮层活动产生的电信号，具有较高的灵敏性，能够快速反映皮层活动状态，目前研究发现脑电信号的gamma波段与持续主动发力有较强关联。脑氧数据用于表征大脑皮层的血氧浓度变化，在受试者产生主动运动意图和认知活动时，大脑活动消耗氧气增加，导致氧合血红蛋白浓度HbO2下降，因此脑氧数据更加稳定，且与脑电信号具有因果关系，对脑电信号和脑氧信号进行因果耦合分析，可更准确地表征受试者的运动意图状态。
具体地，在本具体实施例中，所述脑电信号采集设备采用32导联脑电帽进行数据采集，采用3厘米的源-探测器距离，设置采样率为1000Hz。所述表面肌电采集设备选择16导联肌电采集系统，设置采样率为1000Hz，将肌电电极可选择地放置于人体上斜方肌、肱二头肌、肱三头肌、股直肌、股二头肌、胫前肌、腓肠肌外侧头肌，用于采集肌电信号，将参考电极放置在胸骨部位。所述近红外光谱脑氧采集设备和近红外光谱肌氧采集设备选择采样率为10Hz的24通道，使用连续波照射模式，近红外探头仅测量通过头部的光衰减，采集受试者头部的散射光强度，然后通过校正的朗伯-比尔定律计算出氧合血红蛋白和脱氧血红蛋白的浓度变化，得到脑氧数据和肌氧数据；其中，采集脑氧数据时将近红外探头对受试者头部的颅骨顶叶部位进行散射光强探测；为了保持生理数据的关联性，采集脑氧信号时将近红外探头对受试者的肌肉采集部位与所述表面肌电采集设备的肌电电极的采集位置相同。
优选地，在另一种可能的实施例中，受试者在主动运动康复训练期间利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据如图3所示。具体地，所述多模态生理数据采集系统与具有主动运动模式的康复训练设备搭配使用，包括但不限于上肢康复训练设备、下肢康复训练设备、上下肢联动康复训练设备。其中，根据不同的训练设备提供对应的肢体训练部位，相应调整肌电电极的采集位置，如使用上肢康复训练设备时，只需采集上肢肌电肌氧信号，不需要采集下肢肌电肌氧信号。
具体地，在本具体实施例中，由所述EEG-fNIRS联合采集帽对受试者头部的脑电信号EEG和脑氧信号fNIRS进行同步采集；所述表面肌电采集设备和所述近红外光谱肌氧采集设备集成在穿戴式的肌电肌氧采集绑带上，如图3所示。在未确认受试者的肢体损伤部位时，可选择地将所述肌电肌氧采集绑带戴在受训者的左右侧上斜方肌、肱二头肌、肱三头肌、股直肌、股二头肌、胫前肌、腓肠肌外侧头肌部位。在确认受试者的肢体损伤部位后，可选择只采集左右侧健、患肢体肌肉数据。具体地，为了保证采集肌电信号质量，佩戴前需使用酒精棉签对皮肤表面进行清洁。由于所述近红外光谱脑氧采集设备、近红外光谱肌氧采集设备在实际使用时易受外源环境光干扰，优选采集的理想环境在暗室内进行，或应当保持较暗的灯光环境，尽量避免光线直射。
优选地，在本具体实施例中，在所述逐一特征提取模块中，计算单生理运动参数包括最大肌力值，需要在康复训练前需采集受试者在康复训练状态下各部位肌肉的最大肌力，其具体方式为：在多模态生理数据采集系统进行肌电信号采集前，将康复训练设备的电机制动，让受试者采用与康复训练一样的动作进行发力，采集此时发力的肌电信号，并调整不同旋转角度，改变其主要发力肌肉，继续采集肌电信号。其中，最大肌力值的测量无需每次训练都进行，可依据受试者状态进行调整，优选测量的训练间隔在5~15次范围内，避免最大肌力值长时间不更新可能导致评估存在偏差。
具体地，肌电信号和肌氧数据均与肌肉活动状态密切相关，但肌电信号倾向于表征肌肉发力程度，由于每个人肌肉力量不同，其肌电信号幅值具有一定差异性，在一种可能的实施例中，先测量最大肌力值时的肌电幅值MVC，后续分析计算其相对肌力值。而肌氧数据更倾向于表征肌肉负荷程度，尤其在不同阻抗运动中，肌氧数据能够出现显著差异。另外，在受试者进行主动运动康复训练时，肌肉发力出现肌电信号幅值增大，同时消耗血氧导致肌氧信号的HbO2下降，在康复训练前期时肌电和肌氧的信号变化具有较强的关联性，而在训练后期由于肌肉疲劳导致关联性略微下降，每个人运动功能不同也会导致这种关联性出现不同。因此，对肌电肌氧数据进行耦合分析能够直接体现目前肌肉活动状态，也可表征受训者的运动能力。而在康复训练时，受试者产生运动意图再出现实际肌肉活动，对脑电肌电信号进行耦合分析，可一定程度上反映神经环路受阻情况，所述皮层肌力功能指标计算模块计算得到的运动神经环路功能指标和肢体运动肌力功能指标，相较于单一生理信号特征与肢体运动障碍程度的关联度更高，更能表征受训者的运动功能。
优选地，在另一种可能的实施例中，在所述皮层肌力功能指标计算模块中，计算运动神经环路功能指标需要提前确定受试者出现主动运动意图、准备运动动作以及实际肌肉活动对应的神经血管耦合特征阈值、皮层肌肉耦合特征阈值和肌电肌氧耦合特征阈值，因此需要提前采集一次短时间被动训练期间和一次短时间主动训练期间的多模态生理数据用于计算对应的耦合特征阈值；具体地，设定计算耦合特征阈值的频率与计算最大肌力值的频率保持一致。
优选地，在本具体实施例中，所述皮层肌力功能指标计算模块通过计算受试者主动运动期间的运动神经环路功能指标，包括运动神经环路传导时间、神经肌肉协同程度、主动运动反应时间和主动运动参与度，用于间接反映反映运动神经环路损伤程度；还计算受试者的肢体运动肌力功能指标，包括发力偏侧性指标、相对肌力值积分，用于间接反映肢体运动能力损伤程度；具体地，所述运动神经环路功能指标能更客观全面地从生理层面反映受试者真实的运动能力，相较于单一信号提取的信号特征量化评估肢体功能障碍患者的运动功能更加符合个人的个体化特征和提高评估准确性。
优选地，在本具体实施例中，所述运动功能康复评估模块，用于结合所述运动神经环路功能指标和所述肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型；优选地，使用决策矩阵选择最优特征分类路径，在同时再嵌入一个反向学习神经网络模型，进一步提升模型的识别正确率。
优选地，在本具体实施例中，所述康复训练自适应调整模块中，根据所述皮层肌肉功能活动指标自适应调整受试者的康复训练参数，包括主动运动时长、主动运动阻抗；具体地，根据所述运动神经环路传导时间、主动运动反应时间和发力偏侧性，定制肢体运动和神经反应专项训练，根据所述运动功能康复评估模块得到的Brunnstrom分期水平用于自适应调整受试者的康复训练处方。
具体地，在另一种可能的实施例中，所述康复训练自适应调整模块中，根据所述运动功能康复评估模型自适应调整受试者的康复训练方案包括设定个性化专项训练任务，如特定肢体专项训练和反应能力专项训练。其中特定肢体专项训练是根据所述指标测评模块得到的各部位肌肉的发力偏侧性指标，对发力较少的肢体进行专项训练，逐步提升受训者对偏弱肢体的控制能力；反应能力专项训练主要是通过康复训练游戏对受训者进行反应能力训练，受训者需在屏幕出现特定画面停留的一段时间内，完成一次运动动作，成功完成则加分，延迟完成或未完成则不加分，根据所述指标测评模块得到的运动神经环路传导时间指标来确定特定画面的停留时间。传导时间和反应时间越短，则停留时间越短，反之则越长，以此逐步提高受训者的反应能力，实现运动神经环路的恢复和重构。
所述一种基于多模态脑肌耦合分析的运动功能评估系统可以运行于桌上型计算机、笔记本、掌上电脑及云端数据中心等计算设备中。所述一种基于多模态脑肌耦合分析的运动功能评估系统包括，但不仅限于，处理器、存储器。本领域技术人员可以理解，所述例子仅仅是一种基于多模态脑肌耦合分析的运动功能评估方法及系统的示例，并不构成对一种基于多模态脑肌耦合分析的运动功能评估方法及系统的限定，可以包括比例子更多或更少的部件，或者组合某些部件，或者不同的部件，例如所述一种基于多模态脑肌耦合分析的运动功能评估系统还可以包括输入输出设备、网络接入设备、总线等。
所称处理器可以是中央处理单元(Central Processing Unit，CPU)，还可以是其他通用处理器、数字信号处理器 (Digital Signal Processor，DSP)、专用集成电路(Application Specific Integrated Circuit，ASIC)、现场可编程门阵列 (Field-ProgrammHCO3-le Gate Array，FPGA) 或者其他可编程逻辑器件、分立元器件门电路或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等，所述处理器是所述一种基于多模态脑肌耦合分析的运动功能评估系统的控制中心，利用各种接口和线路连接整个一种基于多模态脑肌耦合分析的运动功能评估系统的各个分区域。
所述存储器可用于存储所述计算机程序和/或模块，所述处理器通过运行或执行存储在所述存储器内的计算机程序和/或模块，以及调用存储在存储器内的数据，实现所述一种基于多模态脑肌耦合分析的运动功能评估方法及系统的各种功能。所述存储器可主要包括存储程序区和存储数据区，其中，存储器可以包括高速随机存取存储器，还可以包括非易失性存储器，例如硬盘、内存、插接式硬盘，智能存储卡（Smart Media Card, SMC），安全数字（Secure Digital, SD）卡，闪存卡（Flash Card）、至少一个磁盘存储器件、闪存器件、或其他易失性固态存储器件。
尽管本发明的描述已经相当详尽且特别对几个所述实施例进行了描述，但其并非旨在局限于任何这些细节或实施例或任何特殊实施例，从而有效地涵盖本发明的预定范围。此外，上文以发明人可预见的实施例对本发明进行描述，其目的是为了提供有用的描述，而那些目前尚未预见的对本发明的非实质性改动仍可代表本发明的等效改动。
Claims (10)
1.一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，所述方法包括以下步骤：
S100，采集多模态生理数据，包括神经环路生理数据和运动肌力生理数据；
S200，对所述多模态生理数据进行数据预处理，得到多模态真实数据；对所述多模态真实数据进行逐一特征提取，得到单生理运动参数；
S300，结合所述单生理运动参数，对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数；
S400，根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标；
S500，根据所述单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标；
S600，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型。
2.根据权利要求1所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S100中，所述多模态生理数据包括神经环路生理数据和运动肌力生理数据，其中所述神经环路生理数据至少包括受试者特定脑区的脑电信号和脑氧数据，所述运动肌力生理数据至少包括受试者多个特定肌肉部位的肌电信号和肌氧数据。
3.根据权利要求1所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S300中，结合所述单生理运动参数，对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数的方法为：所述脑肌耦合运动参数，包括脑电脑氧相关可变量、脑电肌电相干耦合度量和肌电肌氧多模态相关系数，其计算方法包括：
S301，采用滑动窗口方法对所述多模态真实数据进行按段进行分割，按序得到各段分割多模态生理子数据；其中，设定对所述多模态真实数据分割成N个分割多模态生理子数据，N为常数；
S302，对所述各段分割多模态生理子数据进行神经血管耦合分析，包括在各段分割的采集时间段内计算所述脑电脑氧相关可变量记作VarECO，其计算公式为：
其中，x t表示为脑电信号在时刻为t值的采样值，t 为时间变量，各段分割多模态生理子数据对应的时间段范围为t∈[t1,t2]，t1为起始分割时刻，t2为结束分割时刻；y t表示为脑氧数据在t时刻对应的光强度值，F(xt+1, x t, y t)表示为计算脑电信号在t时刻、t+1时刻对应的x t+1、x t值与脑氧数据在t时刻的y t值之间的联合分布函数，F(xt+1 | x t)表示为计算脑电信号在t时刻对应下一时刻出现xt+1值的联合分布概率，F(y t | x t)表示为计算脑电信号在t时刻时对应的脑氧数据出现y t值的联合分布概率，F(xt+1 , y t | x t)表示为计算脑电信号在t+1时刻对应的x t+1值与脑电信号在t时刻对应的脑氧数据出现y t值的概率之间联合分布函数；log2表示为计算以2为底的对数；
S303，对所述各段分割多模态生理子数据进行肌电肌氧耦合分析，包括结合所述单生理运动参数，计算脑电肌电相干耦合度量记作EMCC，其计算公式为：
其中，Pxy(f)表示为脑电信号与肌电数据的互功率谱，Pxx(f)表示为脑电信号的自功率谱，Pyy(f)表示为肌电数据的自功率谱；
S304，对所述各段分割多模态生理子数据进行肌电肌氧耦合分析，包括结合所述单生理运动参数，在时域上各段分割的采集时间段内计算肌电肌氧多模态相关系数记作MOCC，其计算公式为：
其中，x(t)为某一肌肉部位的肌电信号从起始分割时刻到t时刻的均方根值，各段分割多模态生理子数据对应的时间段范围为t∈[t1,t2]，t1为起始分割时刻，t2为结束分割时刻；
为所述肌电信号在t∈[t1,t2]的采集时间段内x(t)的平均值，y(t)为同一肌肉部位的肌氧数据中的HbO2浓度，/>
为所述肌氧数据在t∈[t1,t2]的采集时间段内y(t)的平均值。
4.根据权利要求3所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S400中，根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标的方法为：其中，所述运动神经环路功能指标包括运动神经环路传导时间、神经肌肉协同程度、主动运动反应时间和主动运动持续时间、以及主动运动参与度，其计算方法具体包括：
S401，所述脑肌耦合运动参数包括所述脑电脑氧相关可变量和脑电肌电相干耦合度量，当所述各段分割多模态生理子数据对应的脑肌耦合运动参数大于等于对应的脑肌耦合特征阈值时，根据该段分割多模态生理子数据对应的起始分割时刻计算得到运动神经环路传导时间记作△T；
S402，根据所述各段分割多模态生理子数据进行运动肌肉活动分析，结合机器学习方法构建运动功能预测子模型，由各段分割多模态生理子数据对应输出的预测值按序构成运动肌肉活动曲线；
S403，根据所述各段分割多模态生理子数据对应的脑肌耦合运动参数，结合机器学习方法构建运动功能预测子模型，由各段分割多模态生理子数据对应输出的预测值按序构成神经血管耦合曲线；
S404，对所述神经血管耦合曲线和运动肌肉活动曲线进行运动耦合分析，计算从神经血管到肌肉活动方向的相关可变量，得到神经肌肉协同程度记作VarNMS，其计算公式为：
其中，ni表示为神经血管耦合曲线在所述各段分割多模态生理子数据的第i个分段对应的预测值，i为分段序号， i∈[1,N]，N为所述分割多模态生理子数据的总个数；mi表示为运动肌肉活动曲线在所述各段分割多模态生理子数据的第i个分段对应的预测值，F(ni+1,ni, mi)表示为计算神经血管耦合曲线在第i个分段、第i+1个分段对应的预测值与运动肌肉活动曲线在第i个分段的预测值之间的联合分布函数，F(ni+1 | ni)表示为计算神经血管耦合曲线在第i个分段对应下一个分段出现ni+1值的联合分布概率，F(mi | ni)表示为计算神经血管耦合曲线在第i个分段为ni值时对应的运动肌肉活动曲线出现mi值的联合分布概率，F(ni+1 , mi | ni)表示为计算神经血管耦合曲线在第i个分段为ni值与其对应的运动肌肉活动曲线在第i个分段出现mi值的概率之间的联合分布函数；log2表示为计算以2为底的对数；
S405，所述脑肌耦合运动参数包括所述肌电肌氧多模态相关系数，当所述各段分割多模态生理子数据对应的肌电肌氧多模态相关系数大于等于对应的肌电肌氧耦合特征阈值时，根据该段分割多模态生理子数据对应的起始分割时刻计算得到主动运动反应时间和主动运动持续时间，分别记作Tr和Td；
S406，计算主动运动持续时间与所述脑肌耦合运动参数对应的采集时间总长的比值，得到主动运动参与度记作Sp，其计算公式为：
其中，T为所述脑肌耦合运动参数对应的采集时间总长。
5.根据权利要求4所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S401中，所述根据该段分割多模态生理子数据对应的起始分割时刻计算得到运动神经环路传导时间的方法为：
S4011，当所述各段分割多模态生理子数据对应的脑电脑氧相关可变量大于等于预设的神经血管耦合特征阈值时，判定该段出现主动运动意图，记录该段分割多模态生理子数据对应的起始分割时刻为T1；
S4012，当所述各段分割多模态生理子数据对应的脑电肌电相干耦合度量大于等于预设的皮层肌肉耦合特征阈值时，判定该段出现准备运动动作，记录该段分割多模态生理子数据对应的起始分割时刻为T2；
S4013，计算运动神经环路传导时间记作△T，其中△T= | T2-T1 |。
6.根据权利要求5所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S405中，所述根据该段分割多模态生理子数据对应的起始分割时刻计算得到主动运动反应时间和主动运动持续时间的方法为：
S4051，当所述各段分割多模态生理子数据对应的肌电肌氧多模态相关系数大于等于预设的肌电肌氧耦合特征阈值时，判定该段出现实际肌肉活动，记录该段分割多模态生理子数据对应的起始分割时刻为T3；
S4052，根据步骤S4011中判定出现主动运动意图对应的分割多模态生理子数据的起始分割时刻为T1，计算主动运动反应时间记作Tr，Tr = | T2-T1 |；
S4053，当所述各段分割多模态生理子数据满足对应的肌电肌氧多模态相关系数大于等于对应的肌电肌氧耦合特征阈值，且对应的脑电脑氧相关可变量大于等于预设的神经血管耦合特征阈值时，统计对应的各段分割多模态生理子数据的采集时间总和，得到主动运动持续时间记作Td。
7.根据权利要求3所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S500中，根据所述单生理运动参数和脑肌耦合运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标的方法为：
S501，根据所述单生理运动参数，利用所述肌电信号的最大肌力值对所述运动肌力生理数据中的肌电信号的幅值进行归一化，得到相对肌力值，并计算在整个采集时间段内所述相对肌力值的总和得到相对肌力值积分；
S502，在整个采集时间段内计算受试者的左侧特定肌肉部位对应的相对肌力值平均值与右侧特定肌肉部位对应的相对肌力值平均值的比值得到发力偏侧性指标；
S503，得到肢体运动肌力功能指标，包括所述相对肌力值积分和发力偏侧性指标。
8.根据权利要求1~7中任一所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S600中，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型的方法为：
S601，利用Brunnstrom评定量表对具有不同肢体运动障碍部位及功能程度的多个受试者进行运动功能水平跟踪预测，记录各个受试者在各个康复训练周期对应的Brunnstrom分期水平，并对各个受试者完成步骤S100~S500得到对应的所述运动神经环路功能指标和所述肢体运动肌力功能指标；
S602，根据各个受试者对应的所述运动神经环路功能指标和所述肢体运动肌力功能指标作为数据集，并按照设定比例划分为训练数据集、测试数据集和验证数据集；利用机器学习方法对进行所述训练数据集进行分类训练，构建初级分类模型；
S603，对所述初级分类模型采用十折交叉验证方法或留一法验证进行模型验证，选取每个验证数据集上对应识别准确率最高的模型组成集成初级模型；
S604，利用 BP神经网络模型、多层感知网络、卷积神经网络作为集成次级模型，将所述集成初级模型的输出端作为所述集成次级模型的输入端，训练得到最终决策模型，即为运动功能康复评估初级模型；
S605，利用所述测试数据集对所述运动功能康复评估初级模型进行优化训练，得到最优的运动功能康复评估模型。
9.根据权利要求8所述的一种基于多模态脑肌耦合分析的运动功能评估方法，其特征在于，在S606中，输入所述运动神经环路功能指标和肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型的方法还包括：
S606，利用验证数据集对所述运动功能康复评估模型进行评估得到各个受试者的运动功能评价结果，所述各个受试者的运动功能评价结果与对应的Brunnstrom分期水平进行比对，进一步根据运动神经环路功能指标和肢体运动肌力功能指标对受试者的康复训练进行自适应调整训练参数。
10.一种基于多模态脑肌耦合分析的运动功能评估系统，其特征在于，所述系统包括：存储器、处理器以及存储在所述存储器中并在所述处理器上运行的计算机程序；所述处理器执行所述计算机程序是实现权利要求1至9中的任意一种所述基于多模态脑肌耦合分析的运动功能评估方法中的步骤，具体包括：
多模态生理数据采集模块，利用多模态生理数据采集系统进行同步数据采集获得多模态生理数据；
生理数据预处理模块，用于对所述多模态生理数据进行数据预处理，得到多模态真实数据；
逐一特征提取模块，用于对所述多模态真实数据进行逐一特征提取，得到单生理运动参数；
脑肌参数耦合分析模块，用于对所述多模态真实数据进行脑肌参数耦合分析，得到脑肌耦合运动参数；
皮层肌力功能指标计算模块，用于根据所述脑肌耦合运动参数结合机器学习方法构建多个运动功能预测子模型，计算得到运动神经环路功能指标；并结合所述单生理运动参数进行肢体运动功能分析，计算得到肢体运动肌力功能指标；
运动功能康复评估模块，用于结合所述运动神经环路功能指标和所述肢体运动肌力功能指标，利用集成学习方法建立运动功能康复评估模型；
康复训练自适应调整模块，用于根据所述运动神经环路功能指标和肢体运动肌力功能指标自适应调整康复训练参数，根据所述运动功能康复评估模型自适应调整受试者的康复训练方案；
其中，所述多模态生理数据采集系统至少包括脑电信号采集设备、近红外光谱脑氧采集设备、近红外光谱肌氧采集设备和表面肌电采集设备；
所述一种基于多模态脑肌耦合分析的运动功能评估系统运行于桌上型计算机、掌上电脑中。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310464893.4A CN116312951B (zh) | 2023-04-27 | 2023-04-27 | 一种基于多模态耦合分析的运动功能评估方法及其系统 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310464893.4A CN116312951B (zh) | 2023-04-27 | 2023-04-27 | 一种基于多模态耦合分析的运动功能评估方法及其系统 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN116312951A true CN116312951A (zh) | 2023-06-23 |
| CN116312951B CN116312951B (zh) | 2023-08-04 |
Family
ID=86797995
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310464893.4A Active CN116312951B (zh) | 2023-04-27 | 2023-04-27 | 一种基于多模态耦合分析的运动功能评估方法及其系统 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116312951B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117009876A (zh) * | 2023-10-07 | 2023-11-07 | 长春光华学院 | 基于人工智能的运动状态量评估方法 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2012060874A2 (en) * | 2010-11-01 | 2012-05-10 | Udaya Sankar Devanaboyina | Methods for diagnosis and treatment of disorders of the gastrointestinal tract. and apparatus for use therewith |
| US20130150687A1 (en) * | 2010-05-31 | 2013-06-13 | Toshinori Kato | Apparatus and program for evaluating biological function |
| US20170296048A1 (en) * | 2014-12-11 | 2017-10-19 | Indian Institute of Technology, Gandhinagar | Smart eye system for Visuomotor dysfunction diagnosis and its operant conditioning |
| CN112353407A (zh) * | 2020-10-27 | 2021-02-12 | 燕山大学 | 一种基于神经功能康复主动训练的评估系统及方法 |
| CN113499084A (zh) * | 2021-07-09 | 2021-10-15 | 燕山大学 | 脑卒中运动功能评价模型建立方法、评价方法及系统 |
-
2023
- 2023-04-27 CN CN202310464893.4A patent/CN116312951B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130150687A1 (en) * | 2010-05-31 | 2013-06-13 | Toshinori Kato | Apparatus and program for evaluating biological function |
| WO2012060874A2 (en) * | 2010-11-01 | 2012-05-10 | Udaya Sankar Devanaboyina | Methods for diagnosis and treatment of disorders of the gastrointestinal tract. and apparatus for use therewith |
| US20170296048A1 (en) * | 2014-12-11 | 2017-10-19 | Indian Institute of Technology, Gandhinagar | Smart eye system for Visuomotor dysfunction diagnosis and its operant conditioning |
| CN112353407A (zh) * | 2020-10-27 | 2021-02-12 | 燕山大学 | 一种基于神经功能康复主动训练的评估系统及方法 |
| CN113499084A (zh) * | 2021-07-09 | 2021-10-15 | 燕山大学 | 脑卒中运动功能评价模型建立方法、评价方法及系统 |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117009876A (zh) * | 2023-10-07 | 2023-11-07 | 长春光华学院 | 基于人工智能的运动状态量评估方法 |
| CN117009876B (zh) * | 2023-10-07 | 2024-01-09 | 长春光华学院 | 基于人工智能的运动状态量评估方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN116312951B (zh) | 2023-08-04 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Pincus et al. | Physiological time-series analysis: what does regularity quantify? | |
| Roy et al. | A combined sEMG and accelerometer system for monitoring functional activity in stroke | |
| US11071494B2 (en) | Wearable technologies for joint health assessment | |
| Jordanic et al. | Spatial distribution of HD-EMG improves identification of task and force in patients with incomplete spinal cord injury | |
| Prashant et al. | Brain computer interface: A review | |
| CN112955066A (zh) | 治疗空间评估 | |
| CN111386068A (zh) | 基于摄像机的压力测量系统与方法 | |
| Maura et al. | Literature review of stroke assessment for upper-extremity physical function via EEG, EMG, kinematic, and kinetic measurements and their reliability | |
| CN116312951B (zh) | 一种基于多模态耦合分析的运动功能评估方法及其系统 | |
| Baghdadi et al. | Dasps: a database for anxious states based on a psychological stimulation | |
| Bonato | Advances in wearable technology for rehabilitation | |
| CN113855053A (zh) | 基于肌电的佩戴式肌肉阈值监测系统 | |
| Senanayake et al. | A knowledge-based intelligent framework for anterior cruciate ligament rehabilitation monitoring | |
| CN116746885A (zh) | 麻醉深度监测系统及其方法 | |
| CN117133404B (zh) | 一种针对胸廓出口综合征的智能化康复护理装置 | |
| Örücü et al. | Design and Validation of Multichannel Wireless Wearable SEMG System for Real‐Time Training Performance Monitoring | |
| Patel et al. | Using wearable sensors to predict the severity of symptoms and motor complications in late stage Parkinson's Disease | |
| Li et al. | Detection of muscle fatigue by fusion of agonist and synergistic muscle semg signals | |
| Chan et al. | Estimating heart rate from seismocardiogram signal using a novel deep dominant frequency regressor and domain adversarial training | |
| CN116186502B (zh) | 一种多模态视觉神经功能检测方法及其系统 | |
| Chee et al. | Gait acceleration-based diabetes detection using hybrid deep learning | |
| Zhou et al. | Upper-limb functional assessment after stroke using mirror contraction: A pilot study | |
| CN112656430A (zh) | 基于站立位失衡诱发脑电的卒中平衡康复评估方法 | |
| Moghadam et al. | Quantification of Parkinson tremor intensity based on EMG signal analysis using fast orthogonal search algorithm | |
| CN106983511A (zh) | 肌力与肌张力状态突变点的识别方法及识别装置 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |