CN117017323B - 基于盲源分离的高密度表面膈肌肌电采集与预处理方法 - Google Patents
基于盲源分离的高密度表面膈肌肌电采集与预处理方法 Download PDFInfo
- Publication number
- CN117017323B CN117017323B CN202311186323.XA CN202311186323A CN117017323B CN 117017323 B CN117017323 B CN 117017323B CN 202311186323 A CN202311186323 A CN 202311186323A CN 117017323 B CN117017323 B CN 117017323B
- Authority
- CN
- China
- Prior art keywords
- ith
- signal
- unmixed
- iteration
- diaphragmatic
- 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.)
- Active
Links
- 238000000926 separation method Methods 0.000 title claims abstract description 13
- 238000002203 pretreatment Methods 0.000 title claims abstract description 7
- 230000000241 respiratory effect Effects 0.000 claims abstract description 43
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 14
- 238000001914 filtration Methods 0.000 claims abstract description 10
- 239000013598 vector Substances 0.000 claims description 50
- 230000029058 respiratory gaseous exchange Effects 0.000 claims description 28
- 238000000034 method Methods 0.000 claims description 19
- 238000005314 correlation function Methods 0.000 claims description 11
- 230000003183 myoelectrical effect Effects 0.000 claims description 10
- 238000004364 calculation method Methods 0.000 claims description 9
- 239000011159 matrix material Substances 0.000 claims description 9
- 230000006870 function Effects 0.000 claims description 8
- 238000004458 analytical method Methods 0.000 claims description 5
- 238000004590 computer program Methods 0.000 claims description 5
- 238000007781 pre-processing Methods 0.000 claims description 5
- 238000012217 deletion Methods 0.000 claims description 3
- 230000037430 deletion Effects 0.000 claims description 3
- 210000003205 muscle Anatomy 0.000 claims description 3
- 238000010606 normalization Methods 0.000 claims description 3
- 230000017105 transposition Effects 0.000 claims description 3
- 239000000470 constituent Substances 0.000 claims 1
- 230000000694 effects Effects 0.000 abstract description 2
- 238000002567 electromyography Methods 0.000 abstract 1
- 239000000284 extract Substances 0.000 abstract 1
- 238000012545 processing Methods 0.000 description 3
- QVLMCRFQGHWOPM-ZKWNWVNESA-N N-arachidonoyl vanillylamine Chemical compound CCCCC\C=C/C\C=C/C\C=C/C\C=C/CCCC(=O)NCC1=CC=C(O)C(OC)=C1 QVLMCRFQGHWOPM-ZKWNWVNESA-N 0.000 description 2
- 230000008602 contraction Effects 0.000 description 2
- 210000003238 esophagus Anatomy 0.000 description 2
- 238000000605 extraction Methods 0.000 description 2
- 230000003434 inspiratory effect Effects 0.000 description 2
- 230000033764 rhythmic process Effects 0.000 description 2
- 208000019693 Lung disease Diseases 0.000 description 1
- 206010047700 Vomiting Diseases 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000000740 bleeding effect Effects 0.000 description 1
- 238000009530 blood pressure measurement Methods 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000012880 independent component analysis Methods 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000005399 mechanical ventilation Methods 0.000 description 1
- 230000004007 neuromodulation Effects 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000002685 pulmonary effect Effects 0.000 description 1
- 238000011002 quantification Methods 0.000 description 1
- 210000003019 respiratory muscle Anatomy 0.000 description 1
- 208000023504 respiratory system disease Diseases 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 210000001562 sternum Anatomy 0.000 description 1
- 238000009423 ventilation Methods 0.000 description 1
- 230000003519 ventilatory effect Effects 0.000 description 1
- 230000002747 voluntary effect Effects 0.000 description 1
- 230000008673 vomiting Effects 0.000 description 1
Classifications
-
- 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]
-
- 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]
- A61B5/397—Analysis of electromyograms
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7203—Signal processing specially adapted for physiological signals or for diagnostic purposes for noise prevention, reduction or removal
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7235—Details of waveform analysis
- A61B5/725—Details of waveform analysis using specific filters therefor, e.g. Kalman or adaptive filters
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/72—Signal processing specially adapted for physiological signals or for diagnostic purposes
- A61B5/7235—Details of waveform analysis
- A61B5/7253—Details of waveform analysis characterised by using transforms
- A61B5/7257—Details of waveform analysis characterised by using transforms using Fourier transforms
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61M—DEVICES FOR INTRODUCING MEDIA INTO, OR ONTO, THE BODY; DEVICES FOR TRANSDUCING BODY MEDIA OR FOR TAKING MEDIA FROM THE BODY; DEVICES FOR PRODUCING OR ENDING SLEEP OR STUPOR
- A61M16/00—Devices for influencing the respiratory system of patients by gas treatment, e.g. mouth-to-mouth respiration; Tracheal tubes
- A61M16/021—Devices for influencing the respiratory system of patients by gas treatment, e.g. mouth-to-mouth respiration; Tracheal tubes operated by electrical means
- A61M16/022—Control means therefor
- A61M16/024—Control means therefor including calculation means, e.g. using a processor
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61M—DEVICES FOR INTRODUCING MEDIA INTO, OR ONTO, THE BODY; DEVICES FOR TRANSDUCING BODY MEDIA OR FOR TAKING MEDIA FROM THE BODY; DEVICES FOR PRODUCING OR ENDING SLEEP OR STUPOR
- A61M16/00—Devices for influencing the respiratory system of patients by gas treatment, e.g. mouth-to-mouth respiration; Tracheal tubes
- A61M16/021—Devices for influencing the respiratory system of patients by gas treatment, e.g. mouth-to-mouth respiration; Tracheal tubes operated by electrical means
- A61M16/022—Control means therefor
- A61M16/024—Control means therefor including calculation means, e.g. using a processor
- A61M16/026—Control means therefor including calculation means, e.g. using a processor specially adapted for predicting, e.g. for determining an information representative of a flow limitation during a ventilation cycle by using a root square technique or a regression analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/10—Pre-processing; Data cleansing
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61M—DEVICES FOR INTRODUCING MEDIA INTO, OR ONTO, THE BODY; DEVICES FOR TRANSDUCING BODY MEDIA OR FOR TAKING MEDIA FROM THE BODY; DEVICES FOR PRODUCING OR ENDING SLEEP OR STUPOR
- A61M2230/00—Measuring parameters of the user
- A61M2230/08—Other bio-electrical signals
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2218/00—Aspects of pattern recognition specially adapted for signal processing
- G06F2218/02—Preprocessing
- G06F2218/04—Denoising
Landscapes
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Animal Behavior & Ethology (AREA)
- Veterinary Medicine (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- Public Health (AREA)
- General Health & Medical Sciences (AREA)
- Biophysics (AREA)
- Surgery (AREA)
- Molecular Biology (AREA)
- Medical Informatics (AREA)
- Pathology (AREA)
- Mathematical Physics (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Signal Processing (AREA)
- Data Mining & Analysis (AREA)
- Psychiatry (AREA)
- Physiology (AREA)
- Emergency Medicine (AREA)
- Pulmonology (AREA)
- Anesthesiology (AREA)
- Hematology (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- General Engineering & Computer Science (AREA)
- Evolutionary Biology (AREA)
- Computing Systems (AREA)
- Evolutionary Computation (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
Abstract
本发明公开了一种基于盲源分离的膈肌高密度表面肌电采集与预处理方法，其步骤包括：1、使用高密度表面阵列式肌电采集系统从胸廓皮肤表面采集含有表面膈肌肌电图的电生理信号，并使用基本频域滤波器进行降噪；2、执行盲源分离算法滤除心电等干扰；3、使用带时域约束性FastICA算法检测去噪和提取结果的可靠性；4、根据滤除干扰后的表面膈肌肌电图sEMGdi信号提取膈肌电活动EAdi曲线。本发明能以方便、无创的方式获取膈肌高密度表面肌电信号，对其进行预处理去除心电干扰和其他噪声，提取能够表现膈肌呼吸驱动强度的EAdi曲线，以展示受试者的呼吸信息。
Description
技术领域
本发明属于肌电信号处理技术领域，具体涉及一种膈肌肌电的无创采集和预处理方法，实现受试者呼吸节律的识别与呼吸状态评估。
背景技术
在临床操作中，机械通气的一个关键方面是准确识别和量化患者的自主呼吸力。传统方法采用流量和压力测量来实现此目的，但吸气肌收缩产生的压力必须对抗气道阻力和肺泡内源性呼气末正压(PEEPi)后，才能在气道产生触发呼吸机的流量或压力变化，因此，从吸气肌收缩到气道产生流量或压力变化存在一定的延迟，导致呼吸机送气的时刻比患者实际吸气时刻会相对滞后。相比之下，人体呼吸时膈肌产生的肌电信号(EMGdi)比流量和压力信号更早出现，并受到较少来自肺部和气道疾病等因素的干扰。作为主要的呼吸肌群之一，膈肌在呼吸过程中负责60％-80％的呼吸运动，因此通过测量膈肌肌电信号并进行处理可用于触发呼吸机送气，也可作为膈肌肌电强弱的评判标准之一。
近年来，一种神经调节通气辅助((Neurally Adjusted Ventilatory Assist，NAVA)的呼吸机利用食道膈肌肌电信号实现了这一目标，然而，通过食道采集EMGdi信号实施上比较困难，对采集条件和使用环境要求较高，对患者不方便，同时也具有一定的侵入性，可能引发不适感(如出血、呕吐等),且由单通道采集到的数据易受电极接触问题导致的伪影影响，信号整体鲁棒性和稳定性较差。还有一种通过表面电极获得的EMGdi信号的手段，但由于表面电极常常受到各种干扰和噪声的影响，如心电干扰(ECG)，电极移动伪影、食道蠕动、工频干扰等噪声，对于降噪有更高的要求，目前流行的降噪手段无法达到理想效果。所以即使食道采集存在缺陷，目前市面上的NAVA呼吸机仍以食道采集EMGdi为主。
同时，在采集到的混合信号中，ECG干扰信号的能量主要集中在0-70Hz的频率范围内，而EMGdi信号的能量主要分布在20-250Hz的频带范围内，两者存在频谱重叠现象。为了方便对EMGdi信号进行后续处理，必须采取措施去除ECG干扰以及其他高频噪声。目前常用的的处理方法，如频域分析法，使用快速傅里叶变换获取信号的频谱特征，然后对EMGdi信号进行高通滤波处理。尽管这种方法简单易行，但只能对EMGdi信号做出粗略估计，因为它只能去除ECG信号的一小部分，若进一步滤波将导致EMGdi信息丢失较大。
发明内容
本发明为了克服现有技术的不足之处，提出了一种基于盲源分离的膈肌高密度表面肌电采集与预处理方法，以期能够获取去噪后纯净的多通道膈肌生理信号，从而能够准确提取受试者的呼吸节律并判断其呼吸状态。
本发明为达到上述发明目的，采用如下技术方案：
本发明一种基于盲源分离的膈肌高密度表面肌电采集与预处理方法的特点在于，是按如下步骤进行：
步骤一：采集表面膈肌肌电信号，并使用频域分析法进行去噪；
步骤1.1：采用a行×b列的阵列式电极片采集多通道含有噪声和干扰的表面膈肌肌电信号C＝[c1,c2,…,ci,…,cn]T，其中，ci表示第i个通道的含有噪声和干扰的表面膈肌肌电信号，m表示通道总数，且m＝a×b；
步骤1.2：对含有噪声和干扰的表面膈肌肌电信号C＝[c1,c2,…,ci,…,cn]T进行离散傅里叶变换，得到滤除低频和高频噪声后的表面膈肌肌电信号X＝[x1,x2,…,xi,…,xn]T，其中，xi表示滤除低频和高频噪声后的第i个通道的表面膈肌肌电信号，T表示转置；
步骤二：执行盲源分离算法；
步骤2.1：设置算法执行的总次数为n，利用式(1)对滤除低频和高频噪声后的表面膈肌肌电信号X＝[x1,x2,…,xi,…,xn]T进行扩展，得到扩展信号
式(1)中，K是延迟因子；xm(t-K+1)表示第m通道的延时信号；
定义并初始化膈肌呼吸信号源集合γ及其对应的解混向量集合ε为空集；
定义并初始化残差信号初始化i＝1；定义/>为求解每个解混向量的迭代次数；
步骤2.2：初始化j＝0；
步骤2.3.：定义并初始化残差信号中第i个解混向量在第j次的计算结果wi,j为一个均值为0，方差为1的高斯信号；利用式(2)计算第i个解混向量在第j+1次的计算结果wi,j+1：
式(2)中，G'表示非多项式函数G的一阶导数；G”表示G的二阶导数；E是数学期望；代表第i个解混向量在第j次的计算结果wi,j的转置；
步骤2.4：如果则令j+1赋值给j，返回至步骤2.3；若i＜n且/>则表示得到第i个解混向量wi＝wi,j，令i+1赋值给i后，返回至步骤2.2；若i≥n且/>则表示得到n个解混向量及其组成的解混矩阵W＝[w1,w2,…,wi,…,wn]T，并利用式(3)得到正交化后的解混矩阵W'：
W'＝(WWT)-1/2W (3)
步骤2.5：计算滤除干扰后的呼吸源信号且Y＝[y1,y2,…,yi,…,yn]T，其中，yi代表第i个解混向量wi所对应的呼吸源信号；
从Y中得到信噪比最高的一条呼吸源信号yh；
步骤三：使用带时域约束性FastICA算法检测去噪和提取结果的可靠性；
步骤3.1：定义当前迭代次数为d，并初始化d＝1；定义表示总迭代次数；
定义并初始化第d次迭代的第i个拉格朗日乘子μid＝0、第d次迭代的第i个惩罚因子δid为随机小数；
初始化第d次迭代的第i个呼吸源信号yi,d＝yi、第d次迭代的第i个解混向量wi,d＝wi，其中，解混向量wi为解混矩阵W中的第i个向量；
定义yh和待测的第i个呼吸源信号yi之间的相关性函数为其中，ξ是一个预设的相关性下界，E是数学期望；
初始化第i通道，第d次迭代的第i个相关性下界ξi,d＝1、第d次迭代的第i个相关性函数
步骤3.2：通过式(4)得到第d+1次迭代的第i个解混向量wi,d+1：
式(9)中，g'i,d为第d次迭代的第i个相关性函数gi,d的一阶导数，G'表示非多项式函数G的一阶导数；
步骤3.3：将wi,d+1/||wi,d+1||2赋予wi,d+1，从而对第d+1次迭代的第i个解混向量wi,d+1进行标准化，其中，||wi,d+1||2是wi,d+1的2-范数；
步骤3.4：根据式(10)、式(11)和式(12)分别得到第d+1次迭代的第i个拉格朗日乘子μi,d+1、第d+1次迭代的第i个惩罚因子δi,d+1、第d+1次迭代的第i个相关性下界ξi,d+1：
μi,d+1＝max{0,μi,d+δi,dgi,d} (5)
δi,d+1＝α×δi,d (6)
ξi,d+1＝β×ξi,d (7)
式(5)、式(6)和式(7)中，α、β都为常数；
步骤3.5：更新第d+1次迭代的第i个呼吸源信号更新第d+1次迭代的第i个相关性函数/>
步骤3.6：若则令d+1赋值给d，返回步骤4.2；否则，表示得到第i个呼吸源信号yi＝yi,d+1，并步骤4.7；
步骤3.7：采用步骤一到步骤三的过程得到另一个更新后的呼吸源信号y'i，y'i表示另一个更新后的第i个通道的呼吸源信号；
计算yi和y'i的互相关系数ci；
若ci>C，则表示得到的第i个呼吸源信号yi可靠，并将新的呼吸源信号y'i和其对应的解混向量wi,d+1分别放入膈肌呼吸信号源集合γ和解混向量集合ε中；否则，将新的呼吸源信号y'i和其对应的解混向量wi,d+1直接删除；并第i+1将赋值给第i后，返回步骤3.1，从而对第i+1个呼吸源信号yi+1继续检测更新，直至i>m为止；最终获得更新后检验无误的呼吸源信号γ＝[y'1,y'2,…,y'i,…,y'q]T，其中，y'i表示其中的第i个呼吸源信号；q表示呼吸源信号的数量；C表示信号间的相关度界限；
步骤四：提取EAdi曲线；
从滤除干扰后的膈肌肌电信号γ＝[y'1,y'2,…,y'i,…,y'q]T中得到信噪比最高的膈肌肌电信号y'h；提取y'h的方均根包络线sh并作为最终得到的EAdi曲线。
本发明一种电子设备，包括存储器以及处理器的特点在于，所述存储器用于存储支持处理器执行所述高密度表面膈肌肌电采集与预处理方法的程序，所述处理器被配置为用于执行所述存储器中存储的程序。
本发明一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序的特点在于，所述计算机程序被处理器运行时执行所述高密度表面膈肌肌电采集与预处理方法的步骤。
与已有技术相比，本发明的有益效果体现在：
1、本发明相对于传统的呼吸机，减少了人机对抗，直接反应的受试者的呼吸意图，而相对于食道测量肌电的NAVA呼吸机，使用了方便、无创的肌电采集方式，获取了高密度、时空信息丰富的sEMGdi信号，简化了信号采集步骤，减少了受试者在采集过程中的痛苦。
2、本发明将盲源分离应用于高密度表面膈肌肌电分析，去除了采集信号中的干扰与噪声，并充分保留了其中的sEMGdi信号中的信息，获取EAdi曲线，以准确判断受试者的呼吸状态。
3、本发明使用带时域约束性FastICA算法检测并优化去噪和提取结果，提高了预处理过程的可靠性。
附图说明
图1为本发明的流程图；
图2为本发明采集并基础滤波后一个通道的原始信号图；
图3为本发明去噪后得到的sEMGdi信号图；
图4为本发明提取出的EAdi信号图。
具体实施方式
本实施例中，一种基于盲源分离的膈肌高密度表面肌电(surface DiaphragmElectromyogram,sEMGdi)采集与预处理方法，如图1所示，是按如下步骤进行：
步骤一：采集表面膈肌肌电信号，并使用频域分析法进行去噪；
步骤1.1：采用a行×b列的阵列式电极片采集多通道含有噪声和干扰的表面膈肌肌电信号C＝[c1,c2,…,ci,…,cn]T，其中，ci表示第i个通道的含有噪声和干扰的表面膈肌肌电信号，m表示通道总数，且m＝a×b；
本实施例中，数据采集设备采用32通道的高密度表面肌电采集设备，阵列式电极片贴在人体右侧腋前线与第7-9根肋骨相交的位置，此外还有一个公共的参考电极，将其贴在胸骨上缘的位置。系统采样频率为1000Hz，电极排列为4行×8列，即a＝4，b＝8，m＝32,采集时长为20秒。
步骤1.2：对含有噪声和干扰的表面膈肌肌电信号C＝[c1,c2,…,ci,…,cn]T进行离散傅里叶变换，得到滤除低频噪声后的表面膈肌肌电信号X＝[x1,x2,…,xi,…,xn]T，其中，xi表示滤除低频噪声后的第i个通道的表面膈肌肌电信号，T表示转置；
本实施例中，由于sEMGdi信号的能量主要分布在20-250Hz的频带范围内，因此采用20-250Hz的六阶巴特沃斯带通滤波器去除可能存在的低频和高频噪声，并采用一个50Hz的陷波滤波器消除工频干扰。进行滤波处理后的信号波形如图2所示。
步骤二：执行盲源分离算法，本实施例中，利用盲源分离算法中的FastICA算法进行独立成分分析；
步骤2.1：设置算法执行的总次数为n；为了尽可能迭代得到较多的解混向量，一般n会取得较大，示例性的，可以设置n＝200；
利用式(1)对滤除低频噪声后的表面膈肌肌电信号X＝[x1,x2,…,xi,…,xn]T进行扩展，得到扩展信号
式(1)中，K是延迟因子；本实施例中，K设置为5；xm(t-K+1)表示第m通道的延时信号。
定义并初始化膈肌呼吸信号源集合γ及其对应的解混向量集合ε为空集；
定义并初始化残差信号初始化i＝1；
步骤2.2：初始化j＝0；
步骤2.3：定义并初始化残差信号中第i个解混向量在第j次的计算结果wi,j为一个均值为0，方差为1的高斯信号；利用式(2)计算第i个解混向量在第j+1次的计算结果wi,j+1：
式(2)中，G'表示非多项式函数G的一阶导数；G”表示G的二阶导数；E是数学期望；代表第i个解混向量在第j次的计算结果wi,j的转置；本实施例中的非多项式函数通常设置为G(x)＝log(cosh(x))。
步骤2.4：如果其中/>为求解每个解混向量的迭代次数，/>通常设置为40；则令j+1赋值给j，返回至步骤2.3；若i＜n且/>则得到第i个解混向量wi＝wi,j，则令i+1赋值给i，返回至步骤2.2；若i≥n且/>则得到n个解混向量及其组成的解混矩阵W＝[w1,w2,…,wi,…,wn]T，并利用式(3)得到正交化后的解混矩阵W'：
W'＝(WWT)-1/2W (3)
步骤2.5：计算滤除干扰后的呼吸源信号且Y＝[y1,y2,…,yi,…,yn]T，其中，yi代表第i个解混向量wi所对应的呼吸源信号；以膈肌肌电信号的幅值作为生理特征以提取出信噪比最高的膈肌肌电信号，从多通道的呼吸源信号Y中得到信噪比最高的一条信号yh，如图3所示。
步骤三：使用带时域约束性FastICA算法检测去噪和提取结果的可靠性，带时域约束的FastICA算法是为了检测分解的呼吸源信号的可靠性，通常认为，一个可靠的信号源收敛至尽可能相似于待测信号源的输出，这样判定了待测源信号的可靠性，同时通过更新解混向量更正了可能的错误或遗漏的信号。
步骤3.1：定义当前迭代次数为d，并初始化d＝1；定义表示总迭代次数；
定义并初始化第d次迭代的第i个拉格朗日乘子μid＝0、第d次迭代的第i个惩罚因子δid为随机小数；本实施例的惩罚因子通常初始设置为一个接近于0的小数，示例性的，可以设置为10-5；
初始化第d次迭代的第i个呼吸源信号yi,d＝yi、第d次迭代的第i个解混向量wi,d＝wi，其中，解混向量wi为解混矩阵W中的第i个向量；
定义yh和待测的第i个呼吸源信号yi之间的相关性函数为其中，ξ是一个预设的相关性下界，E是数学期望；
初始化第i通道，第d次迭代的第i个相关性下界ξi,d＝1、第d次迭代的第i个相关性函数
步骤3.2：通过式(4)得到第d+1次迭代的第i个解混向量wi,d+1：
式(9)中，g'i,d为第d次迭代的第i个相关性函数gi,d的一阶导数，G'表示非多项式函数G的一阶导数；
步骤3.3：将wi,d+1/||wi,d+1||2赋予wi,d+1，从而对第d+1次迭代的第i个解混向量wi,d+1进行标准化，其中，||wi,d+1||2是wi,d+1的2-范数。
步骤3.4：根据式(10)、式(11)和式(12)分别得到第d+1次迭代的第i个拉格朗日乘子μi,d+1、第d+1次迭代的第i个惩罚因子δi,d+1、第d+1次迭代的第i个相关性下界ξi,d+1：
μi,d+1＝max{0,μi,d+δi,dgi,d} (5)
δi,d+1＝α×δi,d (6)
ξi,d+1＝β×ξi,d (7)
式(5)、式(6)和式(7)中，α、β都为常数；为保证收敛，α通常设置为靠近1但比1大的数值，β通常设置为靠近1但比1小的数值，示例性的，可以设置为α＝1.02，β＝0.97。
步骤3.5：更新第d+1次迭代的第i个呼吸源信号更新第d+1次迭代的第i个相关性函数/>
步骤3.6：若则令d+1赋值给d，返回步骤4.2；否则，表示得到第i个呼吸源信号yi＝yi,d+1，并步骤4.7；本实施例中，迭代次数/>通常设置为20。
步骤3.7：采用步骤一到步骤三的过程得到另一个更新后的呼吸源信号y'i，y'i表示另一个更新后的第i个通道的呼吸源信号；
计算yi和y'i的互相关系数ci；
若ci>C，则表示得到的第i个呼吸源信号yi可靠，并将新的呼吸源信号y'i和其对应的解混向量wi,d+1分别放入膈肌呼吸信号源集合γ和解混向量集合ε中；否则，将新的呼吸源信号y'i和其对应的解混向量wi,d+1直接删除；并第i+1将赋值给第i后，返回步骤3.1，从而对第i+1个呼吸源信号yi+1继续检测更新，直至i>m为止；其中，C表示信号间的相关度界限；最终全部更新后获得检验无误的呼吸源信号γ＝[y'1,y'2,…,y'i,…,y'q]T，也就是纯净的膈肌肌电信号，其中，q表示呼吸源信号的数量，y'i表示第i个呼吸源信号。
本实施例中的C是两段信号之间的相关度界限，通常设置为0.4，若两段信号的互相关系数超过了0.4，则认为这两段信号来自于同一呼吸源，基于此原理，步骤3.7用来判断更新前后呼吸信号是否发生了质的变化，最终得到q＝10。
步骤四：提取膈肌电活动EAdi(electrical activityofdiaphragm)曲线；
从滤除干扰后的呼吸源信号γ＝[y'1,y'2,…,y'i,…,y'q]T得到信噪比最高的膈肌肌电信号y'h；对于y'h提取其方均根包络线sh，以作为最终得到的EAdi曲线。
本实施例中，以膈肌肌电信号，也就是滤除干扰后的呼吸源信号的幅值作为生理特征以得到信噪比最高的膈肌肌电信号；最终提取出作为EAdi曲线的sh如图4所示。
本实施例中，一种电子设备，包括存储器以及处理器，该存储器用于存储支持处理器执行上述方法的程序，该处理器被配置为用于执行该存储器中存储的程序。
本实施例中，一种计算机可读存储介质，是在计算机可读存储介质上存储有计算机程序，该计算机程序被处理器运行时执行上述方法的步骤。
Claims (3)
1.一种基于盲源分离的膈肌高密度表面肌电采集与预处理方法，其特征在于，是按如下步骤进行：
步骤一：采集表面膈肌肌电信号，并使用频域分析法进行去噪；
步骤1.1：采用a行×b列的阵列式电极片采集多通道含有噪声和干扰的表面膈肌肌电信号C＝[c1，c2，...，ci，...，cn]T，其中，ci表示第i个通道的含有噪声和干扰的表面膈肌肌电信号，m表示通道总数，且m＝a×b；
步骤1.2：对含有噪声和干扰的表面膈肌肌电信号C＝[c1，c2，...，ci，...，cn]T进行离散傅里叶变换，得到滤除低频和高频噪声后的表面膈肌肌电信号X＝[x1，x2，...，xi，...，xn]T，其中，xi表示滤除低频和高频噪声后的第i个通道的表面膈肌肌电信号，T表示转置；
步骤二：执行盲源分离算法；
步骤2.1：设置算法执行的总次数为n，利用式(1)对滤除低频和高频噪声后的表面膈肌肌电信号X＝[x1，x2，...，xi，...，xn]T进行扩展，得到扩展信号
式(1)中，K是延迟因子；xm(t-K+1)表示第m通道的延时信号；
定义并初始化膈肌呼吸信号源集合γ及其对应的解混向量集合ε为空集；
定义并初始化残差信号初始化i＝1；定义/>为求解每个解混向量的迭代次数；
步骤2.2：初始化j＝0；
步骤2.3：定义并初始化残差信号中第i个解混向量在第j次的计算结果wi，j为一个均值为0，方差为1的高斯信号；利用式(2)计算第i个解混向量在第j+1次的计算结果wi，j+1：
式(2)中，G′表示非多项式函数G的一阶导数；G″表示G的二阶导数；E是数学期望；代表第i个解混向量在第j次的计算结果wi，j的转置；
步骤2.4：如果则令j+1赋值给j，返回至步骤2.3；若i＜n且/>则表示得到第i个解混向量wi＝wi，j，令i+1赋值给i后，返回至步骤2.2；若i≥n且/>则表示得到n个解混向量及其组成的解混矩阵W＝[w1，w2，...，wi，...，wn]T，并利用式(3)得到正交化后的解混矩阵W′：
W′＝(WWT)-1/2W (3)
步骤2.5：计算滤除干扰后的呼吸源信号且Y＝[y1，y2，...，yi，...，yn]T，其中，yi代表第i个解混向量wi所对应的呼吸源信号；以膈肌肌电信号的幅值作为生理特征以提取出信噪比最高的膈肌肌电信号，从Y中得到信噪比最高的一条呼吸源信号yh；
步骤三：使用带时域约束性FastICA算法检测去噪和提取结果的可靠性；
步骤3.1：定义当前迭代次数为d，并初始化d＝1；定义表示总迭代次数；
定义并初始化第d次迭代的第i个拉格朗日乘子μid＝0、第d次迭代的第i个惩罚因子δid为随机小数；
初始化第d次迭代的第i个呼吸源信号yi，d＝yi、第d次迭代的第i个解混向量wi，d＝wi，其中，解混向量wi为解混矩阵W中的第i个向量；
定义yh和待测的第i个呼吸源信号yi之间的相关性函数为其中，ξ是一个预设的相关性下界，E是数学期望；
初始化第i通道，第d次迭代的第i个相关性下界ξi，d＝1、第d次迭代的第i个相关性函数
步骤3.2：通过式(4)得到第d+1次迭代的第i个解混向量wi，d+1：
式(9)中，g′i，d为第d次迭代的第i个相关性函数gi，d的一阶导数，G′表示非多项式函数G的一阶导数；
步骤3.3：将wi，d+1/||wi，d+1||2赋予wi，d+1，从而对第d+1次迭代的第i个解混向量wi，d+1进行标准化，其中，||wi，d+1||2是wi，d+1的2-范数；
步骤3.4：根据式(10)、式(11)和式(12)分别得到第d+1次迭代的第i个拉格朗日乘子μi，d+1、第d+1次迭代的第i个惩罚因子δi，d+1、第d+1次迭代的第i个相关性下界ξi，d+1：
μi，d+1＝max{0，μi，d+δi，dgi，d} (5)
δi，d+1＝α×δi，d (6)
ξi，d+1＝β×ξi，d (7)
式(5)、式(6)和式(7)中，α、β都为常数；
步骤3.5：更新第d+1次迭代的第i个呼吸源信号更新第d+1次迭代的第i个相关性函数/>
步骤3.6：若则令d+1赋值给d，返回步骤4.2；否则，表示得到第i个呼吸源信号yi＝yi，d+1，并步骤4.7；
步骤3.7：采用步骤一到步骤三的过程得到另一个更新后的呼吸源信号y′i，y′i表示另一个更新后的第i个通道的呼吸源信号；
计算yi和y′i的互相关系数ci；
若ci>C，则表示得到的第i个呼吸源信号yi可靠，并将新的呼吸源信号y′i和其对应的解混向量wi，d+1分别放入膈肌呼吸信号源集合γ和解混向量集合ε中；否则，将新的呼吸源信号y′i和其对应的解混向量wi，d+1直接删除；并将第i+1将赋值给第i后，返回步骤3.1，从而对第i+1个呼吸源信号yi+1继续检测更新，直至i>m为止；最终获得更新后检验无误的呼吸源信号γ＝[y′1，y′2，...，y′i，...，y′q]T，其中，y′i表示其中的第i个呼吸源信号；q表示呼吸源信号的数量；C表示信号间的相关度界限；
步骤四：提取EAdi曲线；
从滤除干扰后的膈肌肌电信号γ＝[y′1，y′2，...，y′i，...，y′q]T中得到信噪比最高的膈肌肌电信号y′h；提取y′h的方均根包络线sh并作为最终得到的EAdi曲线。
2.一种电子设备，包括存储器以及处理器，其特征在于，所述存储器用于存储支持处理器执行权利要求1所述高密度表面膈肌肌电采集与预处理方法的程序，所述处理器被配置为用于执行所述存储器中存储的程序。
3.一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序，其特征在于，所述计算机程序被处理器运行时执行权利要求1所述高密度表面膈肌肌电采集与预处理方法的步骤。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311186323.XA CN117017323B (zh) | 2023-09-14 | 2023-09-14 | 基于盲源分离的高密度表面膈肌肌电采集与预处理方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311186323.XA CN117017323B (zh) | 2023-09-14 | 2023-09-14 | 基于盲源分离的高密度表面膈肌肌电采集与预处理方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117017323A CN117017323A (zh) | 2023-11-10 |
| CN117017323B true CN117017323B (zh) | 2024-03-29 |
Family
ID=88635599
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311186323.XA Active CN117017323B (zh) | 2023-09-14 | 2023-09-14 | 基于盲源分离的高密度表面膈肌肌电采集与预处理方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117017323B (zh) |
Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1733330A (zh) * | 2004-08-10 | 2006-02-15 | 郑则广 | 利用食道电极膈肌肌电图触发呼吸机送气的方法 |
| CN102631195A (zh) * | 2012-04-18 | 2012-08-15 | 太原科技大学 | 人体表面肌电信号单通道盲源分离法 |
| WO2013043157A2 (en) * | 2011-09-20 | 2013-03-28 | Vivaquant Llc | Physiological signal denoising |
| CN110680317A (zh) * | 2019-10-15 | 2020-01-14 | 中国科学技术大学 | 基于独立向量分析的高密度表面肌电信号消噪方法 |
| CN112494054A (zh) * | 2020-11-26 | 2021-03-16 | 天津大学 | 基于多导肌电与脑电相干的卒中下肢运动康复评估方法 |
| CN114028671A (zh) * | 2021-11-04 | 2022-02-11 | 湖南万脉医疗科技有限公司 | 一种呼吸机控制触发方法 |
| CN114767132A (zh) * | 2022-04-20 | 2022-07-22 | 中国科学技术大学 | 一种基于表面肌电分解的运动单位识别方法 |
| CN114969637A (zh) * | 2022-04-18 | 2022-08-30 | 华南理工大学 | 基于线性预测的膈肌肌电实时处理方法及存储介质 |
| CN116089851A (zh) * | 2022-11-07 | 2023-05-09 | 中国科学技术大学 | 一种自适应更新的在线运动单位识别方法 |
| CN116601721A (zh) * | 2020-06-05 | 2023-08-15 | 瑞思迈传感器技术有限公司 | 用于识别用户接口的系统和方法 |
| CN116584960A (zh) * | 2023-04-04 | 2023-08-15 | 北京雅果科技有限公司 | 一种膈肌肌电信号降噪方法 |
Family Cites Families (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE102011016804B4 (de) * | 2011-04-12 | 2016-01-28 | Drägerwerk AG & Co. KGaA | Vorrichtung und Verfahren zur Datenverarbeitung physiologischer Signale |
| US20140207017A1 (en) * | 2013-01-23 | 2014-07-24 | Altec, Inc. | Signal quality monitor for electromyographic sensors |
| US20180014784A1 (en) * | 2015-01-30 | 2018-01-18 | New York University | System and method for electrophysiological monitoring |
| WO2017201538A1 (en) * | 2016-05-20 | 2017-11-23 | The Regents Of The University Of California | Devices, systems and methods for monitoring physiological functions from surface electrophysiological sensors |
| EP4125597A4 (en) * | 2020-03-30 | 2024-05-22 | X-Trodes Ltd | SYSTEM AND METHOD FOR MAPPING MUSCLE ACTIVATION |
-
2023
- 2023-09-14 CN CN202311186323.XA patent/CN117017323B/zh active Active
Patent Citations (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1733330A (zh) * | 2004-08-10 | 2006-02-15 | 郑则广 | 利用食道电极膈肌肌电图触发呼吸机送气的方法 |
| WO2013043157A2 (en) * | 2011-09-20 | 2013-03-28 | Vivaquant Llc | Physiological signal denoising |
| CN102631195A (zh) * | 2012-04-18 | 2012-08-15 | 太原科技大学 | 人体表面肌电信号单通道盲源分离法 |
| CN110680317A (zh) * | 2019-10-15 | 2020-01-14 | 中国科学技术大学 | 基于独立向量分析的高密度表面肌电信号消噪方法 |
| CN116601721A (zh) * | 2020-06-05 | 2023-08-15 | 瑞思迈传感器技术有限公司 | 用于识别用户接口的系统和方法 |
| CN112494054A (zh) * | 2020-11-26 | 2021-03-16 | 天津大学 | 基于多导肌电与脑电相干的卒中下肢运动康复评估方法 |
| CN114028671A (zh) * | 2021-11-04 | 2022-02-11 | 湖南万脉医疗科技有限公司 | 一种呼吸机控制触发方法 |
| CN114969637A (zh) * | 2022-04-18 | 2022-08-30 | 华南理工大学 | 基于线性预测的膈肌肌电实时处理方法及存储介质 |
| CN114767132A (zh) * | 2022-04-20 | 2022-07-22 | 中国科学技术大学 | 一种基于表面肌电分解的运动单位识别方法 |
| CN116089851A (zh) * | 2022-11-07 | 2023-05-09 | 中国科学技术大学 | 一种自适应更新的在线运动单位识别方法 |
| CN116584960A (zh) * | 2023-04-04 | 2023-08-15 | 北京雅果科技有限公司 | 一种膈肌肌电信号降噪方法 |
Non-Patent Citations (6)
| Title |
|---|
| 《下颌表面肌电提取和评价呼吸相关时相型成分》;张达等;《仪器仪表学报》;第35卷(第8期);1875-1882 * |
| 《基于卷积混合盲源分离技术的表面肌电信号分解研究》;李强等;《中国生物医学工程学报》;第25卷(第4期);404-410 * |
| 《基于模式识别和肌力估计的高密度肌电控制方法研究》;李斌；陈香;博士论文 * |
| 《基于盲源分离技术的高密度表面肌电分解》;陈茂启;《医药卫生科技辑》;E080-2 * |
| 《肌电检测中消除工频干扰的方法》;钱晓进等;《中国医疗器械杂志》;第27卷(第4期);260-273 * |
| 《融合盲源分离与小波变换的EMGdi信号降噪研究》;杨智;硕士论文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117017323A (zh) | 2023-11-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US7025729B2 (en) | Apparatus for detecting sleep apnea using electrocardiogram signals | |
| US10945628B2 (en) | Apparatus and method for processing electromyography signals related to respiratory activity | |
| US6415174B1 (en) | ECG derived respiratory rhythms for improved diagnosis of sleep apnea | |
| CN106108889A (zh) | 基于深度学习算法的心电图分类方法 | |
| US20160256063A1 (en) | Analyte assessment and arrhythmia risk prediction using physiological electrical data | |
| CN106214145A (zh) | 一种基于深度学习算法的心电图分类方法 | |
| Sun et al. | An improved morphological approach to background normalization of ECG signals | |
| EP3270774A2 (en) | Continuous non-invasive monitoring of a pregnant human subject | |
| US20210169359A1 (en) | Intra-cardiac pattern matching | |
| CN103961089A (zh) | 基于分段直线拟合的窦性心率震荡趋势检测方法 | |
| CN114449947B (zh) | 用于确定呼吸信号和/或心源性信号的方法和设备 | |
| CN109157220B (zh) | 一种基于多通道动态监测的呼吸指数提取系统及工作方法 | |
| Azad et al. | Seismocardiographic signal variability during regular breathing and breath hold in healthy adults | |
| CN106955101A (zh) | 从心电信号中提取呼吸信号的方法及装置 | |
| CN117017323B (zh) | 基于盲源分离的高密度表面膈肌肌电采集与预处理方法 | |
| CN105943040A (zh) | 一种治疗中枢性睡眠呼吸暂停症的方法及装置 | |
| CN113116328A (zh) | 一种基于心冲击图的心率检测方法 | |
| CN106725487A (zh) | 一种呼吸率提取方法及装置 | |
| CN108784683B (zh) | 包括用于特征检测的滤波器的ecg机器 | |
| Stantic et al. | Detecting abnormal ecg signals utilising wavelet transform and standard deviation | |
| Tan et al. | Analysis of signal-averaged electrocardiogram performance for body surface recordings | |
| CN118044829A (zh) | 基于盲源分离与剥离策略的表面膈肌肌电去噪方法 | |
| Swamy et al. | Pulsation artifact removal from intra-operatively recorded local field potentials using sparse signal processing and data-specific dictionary | |
| Bear et al. | Impact of baseline drift removal on ECG beat classification and alignment | |
| Li et al. | Extraction of Weak Surface Diaphragmatic Electromyogram Using Modified Progressive FastICA Peel-Off |
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 |