CN116049735A - 一种基于多元经验模态分解的子宫肌电信号分类方法 - Google Patents

一种基于多元经验模态分解的子宫肌电信号分类方法 Download PDF

Info

Publication number
CN116049735A
CN116049735A CN202310052408.2A CN202310052408A CN116049735A CN 116049735 A CN116049735 A CN 116049735A CN 202310052408 A CN202310052408 A CN 202310052408A CN 116049735 A CN116049735 A CN 116049735A
Authority
CN
China
Prior art keywords
uterine
feature
channel
time
signals
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.)
Pending
Application number
CN202310052408.2A
Other languages
English (en)
Inventor
张旭
崔嘉雯
骆轩玉
李心慧
陈香
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Science and Technology of China USTC
Original Assignee
University of Science and Technology of China USTC
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by University of Science and Technology of China USTC filed Critical University of Science and Technology of China USTC
Priority to CN202310052408.2A priority Critical patent/CN116049735A/zh
Publication of CN116049735A publication Critical patent/CN116049735A/zh
Pending legal-status Critical Current

Links

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/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/43—Detecting, measuring or recording for evaluating the reproductive systems
    • A61B5/4306—Detecting, measuring or recording for evaluating the reproductive systems for evaluating the female reproductive systems, e.g. gynaecological evaluations
    • A61B5/4318—Evaluation of the lower reproductive system
    • A61B5/4325—Evaluation of the lower reproductive system of the uterine cavities, e.g. uterus, fallopian tubes, ovaries
    • 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/7264—Classification of physiological signals or data, e.g. using neural networks, statistical classifiers, expert systems or fuzzy systems
    • 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/7264—Classification of physiological signals or data, e.g. using neural networks, statistical classifiers, expert systems or fuzzy systems
    • A61B5/7267—Classification of physiological signals or data, e.g. using neural networks, statistical classifiers, expert systems or fuzzy systems involving training the classification device

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Artificial Intelligence (AREA)
  • General Health & Medical Sciences (AREA)
  • Biophysics (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Medical Informatics (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Pathology (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Biomedical Technology (AREA)
  • Signal Processing (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Physiology (AREA)
  • Psychiatry (AREA)
  • Evolutionary Computation (AREA)
  • Fuzzy Systems (AREA)
  • Mathematical Physics (AREA)
  • Reproductive Health (AREA)
  • Gynecology & Obstetrics (AREA)
  • Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)

Abstract

本发明公开了一种基于多元经验模态分解的子宫肌电信号分类方法，包括：1、使用多通道表面电极采集一段宫缩时期内的子宫肌电信号并进行降噪预处理；2、采用多元经验模态分解将降噪预处理后的子宫肌电信号扩展到时频域，获取包含子宫肌电信号不同时间尺度的本征模态分量，根据方差贡献度和相关系数选择其中的有效分量；3、对筛选后的本征模态分量进行滑动分窗，并计算六个特征值；4、采用卡方检验预筛选和序列双向搜索组合的特征选择算法，挑选最优的特征子集；5、基于支持向量机模型，分类子宫肌电信号中的正常和异常活动。本发明能利用子宫肌电信号的时频域局部信息，滤除无关的噪声和冗余特征，从而完成对子宫肌电信号中的异常活动识别。

Description

一种基于多元经验模态分解的子宫肌电信号分类方法
技术领域
本发明属于肌电信号处理及分析领域，具体的说是一种基于多元经验模态分解的子宫肌电信号分类方法。
背景技术
子宫肌电(Eelctrohysterogram,EHG)信号作为一种伴随子宫收缩产生的电生理信号，是无数个子宫平滑肌细胞单个电活动的综合表现。临近分娩时，子宫平滑肌细胞缝隙连接的增加为动作电位的传播提供了快速途径，体现了细胞兴奋性的增强，当动作电位沿肌层间歇性传播时，会导致子宫收缩，而规律和同步宫缩可作为判断孕妇是否分娩的重要指标。随着肌电采集技术的进步，高密度表面电极可以从孕妇腹部以无创的电生理监测方式记录子宫收缩时产生的电活动，获取到包含丰富时空信息的EHG信号，客观地描述子宫收缩的持续时间、频率和强度，进而识别出宫缩期间子宫肌电信号的异常活动，为孕妇分娩提供预警。
由于EHG信号的产生和传播是由无数子宫平滑肌细胞复杂的电活动引起的，具有很强的非平稳性和非线性。因此，对EHG信号整个频率范围的分析可能无法充分描述妊娠的演变，并使后续分类变得困难。近年来，许多不同的信号时频处理技术已被用于在信号不同尺度上进行分析来提升特征表征能力。其中，相比于短时傅里叶变换在全局范围内的时间分辨率和频率分辨率是固定的，以及小波变换需要人为地选择基函数和分解层数，经验模态分解作为一种完全数据驱动的方法，可以自适应地将非平稳时间序列分解为不同时间尺度的本征模态分量。然而，对于多通道信号，经验模态分解只能在单个通道上逐一执行，导致了本征模态分量在不同通道间的数量和频率分布不一致以及频谱混叠问题，这对后续的信号分析及异常的子宫肌电活动识别是一个巨大的挑战。
发明内容
本发明是为了解决上述现有技术存在的不足之处，提出一种基于多元经验模态分解的子宫肌电信号分类方法，以期能够自适应将多通道的EHG信号分解成不同时间尺度的本征模态分量，同时保证不同通道间的数量和频谱一致性，并滤除无关的噪声和冗余特征，从而能解决模态失调和模态混叠问题，以完成对子宫肌电信号中的异常活动识别。
本发明为达到上述发明目的，采用如下技术方案：
本发明一种基于多元经验模态分解的子宫肌电信号分类方法的特点在于，是按如下步骤进行：
步骤1、信号采集：
步骤1.1、利用子宫肌电测量设备和Q通道的高密度表面电极采集某一待测对象在t时刻子宫收缩时的表面子宫肌电信号数据，记为：x(t)＝[x1(t)，x2(t)，...，xq(t)，...，xQ(t)]T，1≤q≤Q，1≤t≤T1，其中，xq(t)表示第q通道上第t时刻的子宫肌电信号数据，Q为高密度电极的总通道数，T1为记录的总时长；
步骤1.2、在采样频率fs下，得到子宫肌电信号x(t)的离散样本点数据，记为：x(n)＝[x1(n)，x2(n)，...，xq(n)，...，xQ(n)]T，1≤n≤N，其中，xq(n)表示第q通道上第n个样本点的子宫肌电信号数据，N为记录的总样本点数；从而经过T1时长的测量后，得到样本长度为N的子宫肌电时间序列x＝[x1，x2，...，xq，...，xQ]T，其中，xq表示第q通道的子宫肌电时间序列，且xq＝{xq(n)|1≤n≤N}，令子宫肌电时间序列x的类别标签为y∈{1，0}，其中，y＝1表示正常，y＝0代表异常；
步骤2、信号预处理：
使用带通滤波器对子宫肌电时间序列x进行降噪，滤除干扰的生理信号和无关噪声后，得到降噪处理后的Q通道的子宫肌电时间序列x*＝[x* 1，x* 2...，x* q，...，x* Q]T，其中，x* q表示降噪处理后的第q通道的宫肌电时间序列；
步骤3、信号时频扩展：
步骤3.1、采用式(1)所示的多元经验模态分解方法将子宫肌电信号时间序列x*在Q通道上联合分解为I个本征模态分量和残差项：
Figure BDA0004058790450000021
式(1)中，
Figure BDA0004058790450000022
表示第q通道的第i个本征模态分量，rq表示第q通道的残差项；
步骤3.2、根据式(2)和式(3)分别计算
Figure BDA0004058790450000023
在第i个本征模态分量上的方差贡献度Varci和相关系数Cori：
Figure BDA0004058790450000024
Figure BDA0004058790450000025
式(2)和式(3)中，
Figure BDA0004058790450000026
和var(x* q)分别表示
Figure BDA0004058790450000027
和x* q的方差，
Figure BDA0004058790450000028
表示
Figure BDA0004058790450000029
和x* q的互方差；
根据{Varci|1≤i≤I}和{Cori|1≤i≤I}，选取I0个对子宫收缩有贡献且相关的本征模态分量，记为
Figure BDA00040587904500000210
其中，
Figure BDA00040587904500000211
表示第q通道的第i0个本征模态分量；
步骤4、特征提取：
设置窗口长度为Nw和窗口重叠为NL，对每一个长度为N的
Figure BDA0004058790450000031
进行滑动窗分帧处理后，共得到((N-Nw)/NL+1)个滑动窗口时间序列，并将其中一个滑动窗口时间序列记为
Figure BDA0004058790450000032
且
Figure BDA0004058790450000033
其中，
Figure BDA0004058790450000034
表示长度为Nw的时间序列上的第n0个样本点数据；
分别计算每个滑动窗口时间序列的六个特征值，包括：时域上的均方根幅度和自相关过零数，频域上的峰值频率和中值频率，非线性域上的样本熵和模糊熵；并对所有滑动窗口时间序列求得的六个特征值分别取平均，从而得到三维特征矩阵f＝[f(i0，j)，f2(i0，j)，...，fq(i0，j)，...，fQ(i0，j)]T，1≤i0≤I0，1≤j≤6，其中，fq(i0，j)表示第q个通道上的第i0个本征模态分量的第j个特征值；
步骤5、特征选择：
步骤5.1、将三维特征矩阵f转换为一维的特征集合F＝{F(k)|k＝1，2，…，Q×10×6}，并通过卡方检验计算第k个特征F(k)的卡方分数，并对所有卡方分数进行降序排序后，选取大于预设阈值th的K1个特征，从而得到初筛选后的特征集合F*＝{F*(k′)|k″＝1，2，…，K1}；其中，F*(k′)表示初筛选后的第k′个特征，K1表示初筛选特征集合中的特征个数；
步骤5.2、基于序列双向搜索算法对初筛选后的特征集合F*进行特征优化，得到最优特征子集F**＝{F**(k″)|k″＝1，2，…，K2}，其中，F**(k″)表示双向搜索算法选择的第k″个特征，K2表示最优特征子集中的特征个数；
步骤6、按照步骤1-步骤5的过程得到M个待测对象的最优特征子集并形成数据集R＝{F** 1(k″)，F** 2(k″)，...，F** M(k″)|k″＝1，2，…，K2}，其中，F** M(k″)表示第M个待测对象的最优特征子集；
将数据集R分为Ma个已知标签的训练数据
Figure BDA0004058790450000035
及其对应类别标签
Figure BDA0004058790450000036
Mb个未知标签的测试数据
Figure BDA0004058790450000037
使用支持向量机对训练数据
Figure BDA0004058790450000038
建立与类别标签
Figure BDA0004058790450000039
之间的映射关系，找到能够分不同类别的超平面γ，从而基于所述支持向量机模型的超平面γ对测试数据
Figure BDA00040587904500000310
进行分类，并得到预测的分类结果
Figure BDA00040587904500000311
本发明所述的基于多元经验模态分解的子宫肌电信号分类方法的特点也在于，所述步骤5.2包括：
步骤5.2.1、设置算法的执行总次数为T，从特征集合F*中随机生成一个特征子集F*(k0)；
步骤5.2.2、从特征集合F*中依次选择一个特征F*(k+)加入特征子集F*(k0)中，并判断F*(k0)的F1分数是否达到最优，若是，则执行步骤5.2.3，否则，继续执行步骤5.2.2；
步骤5.2.3、依次从特征子集F*(k0)中删除一个特征F(k-)，并判断删除后的特征子集F*(k0)的F1分数是否降低，若是，则在特征子集F*(k0)中保留特征F(k-)，否则继续执行步骤5.2.3；从而得到最优特征子集F**＝{F**(k″)|k″＝1，2，…，K2}，其中，F**(k″)表示双向搜索算法选择的第k″个特征，K2表示最优特征子集中的特征个数。
本发明一种电子设备，包括存储器以及处理器的特点在于，所述存储器用于存储支持处理器执行所述子宫肌电信号分类方法的程序，所述处理器被配置为用于执行所述存储器中存储的程序。
本发明一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序的特点在于，所述计算机程序被处理器运行时执行所述子宫肌电信号分类方法的步骤。
与现有技术相比，本发明的有益效果在于：
1、本发明步骤三对子宫肌电信号进行了自适应地时频分解，获取了分布在不同频率区间上的本征模态分量，步骤四中通过提取能够反映子宫肌电信号的强度、频率和复杂度特征，从而刻画子宫收缩时电活动的传播趋势，步骤五中通过滤除与子宫收缩不相关的特征分量和冗余噪声，提升对子宫肌电信号中异常活动的识别性能。
2、本发明步骤三中通过多元经验模态分解，提高了子宫肌电信号不同尺度的时频分辨率，并保证不同通道间得到的本征模态分量的数量和频谱一致性，解决了模态失调和模态混叠问题。
3、本发明获取了各个分量上的时域、频域、非线性域特征，并通过步骤五中的卡方检验和序列双向搜索组合的特征选择算法，滤除了不相关的特征分量，消除冗余信息和无关噪声，更好的定位反映子宫肌电收缩活动的信息。
附图说明
图1为本发明的一种基于多元经验模态分解的子宫肌电信号分类方法流程图；
图2为本发明的多通道EHG信号经MEMD分解后的本征模态分量波形图和对应的功率谱密度分布图
图3为本发明的各级IMF对原始EHG信号的方差贡献度的箱线图
图4为本发明的各级IMF对原始EHG信号的相关系数的箱线图；
图5为本发明的特征选择算法中预设阈值th和迭代次数T对F1分数的影响趋势图；
图6为本发明的基于不同早产预测方法下的性能对比柱状图。
具体实施方式
本实施例中，一种基于多元经验模态分解的子宫肌电信号分类方法，如图1所示，包括如下步骤：
步骤1、信号采集：
步骤1.1、利用子宫肌电测量设备和Q通道的高密度表面电极采集某一待测对象在t时刻子宫收缩时的表面子宫肌电信号数据，记为：x(t)＝[x1(t)，x2(t)，...，xq(t)，...，xQ(t)]T，1≤q≤Q，1≤t≤T1，其中，xq(t)表示第q通道上第t时刻的子宫肌电信号数据，Q为高密度电极的总通道数，T1为记录的总时长；
步骤1.2、在采样频率fs下，得到子宫肌电信号x(t)的离散样本点数据，记为：x(n)＝[x1(n)，x2(n)，...，xq(n)，...，xQ(n)]T，1≤n≤N，其中，xq(n)表示第q通道上第n个样本点的子宫肌电信号数据，N为记录的总样本点数；从而经过T1时长的测量后，得到样本长度为N的子宫肌电时间序列x＝[x1，x2，...，xq，...，xQ]T，其中，xq表示第q通道的子宫肌电时间序列，且xq＝{xq(n)|1≤n≤N}，令子宫肌电时间序列x的类别标签为y∈{1，0}，其中，y＝1表示正常，y＝0代表异常；
具体实施中包括：
(1)招募M位胎龄在30周以上的受试者，引导每位受试者选择舒适的姿势平躺，将四个电极对称贴附于受试者肚脐下方的腹部皮肤表面上，参考电极和接地电极分别放置于胯部两侧，设备以采样频率fs采集子宫平滑肌的表面子宫肌电信号，为了降低由外部设备和身体移动带来的噪声干扰，对四个电极采用两两差分操作形成Q个双极通道，电极阵列中相邻电极中心间距为d。示例性的，可以设置：M＝300，fs＝20Hz，Q＝6，d＝7cm；
(2)逐一采集受试者经历子宫收缩期间的连续时序信号，手动去除明显噪声干扰和异常片段，得到M段时长为T1的子宫肌电信号，经过采样离散化后对应M组子宫肌电信号的离散样本点数据，并对每一个受试者的子宫肌电数据添加类别标签y∈{1，0}，从而将其存入子宫肌电数据集。示例性的，可以设置：T1＝30min，类别标签y＝1表示正常，y＝0代表异常；
步骤2、信号预处理：
使用带通滤波器对子宫肌电时间序列x进行降噪，滤除干扰的生理信号和无关噪声后，得到降噪处理后的Q通道的子宫肌电时间序列x*＝[x* 1，x* 2...，x* q，...，x* Q]T，其中，x* q表示降噪处理后的第q通道的宫肌电时间序列；
本实施例中，采用带通截止频率为0.3～3Hz的四级巴特沃斯因果IIR滤波器，滤除EHG信号中常常混有母体心电、胎儿心电、母体呼吸等生理信号，以及采集过程中易引入50Hz工频干扰和基线漂移等噪声。进一步地，为了避免由于滤波带来的瞬态效应，去除整段时长为T1的信号的前后T0分钟片段，中间的T1-2T0时长的子宫肌电信号用作后续分析，对应的离散子宫肌电数据的总样本点数为N。示例性的，可以设置：T0＝5min，N＝24000。
步骤3、对降噪后的子宫信号进行时频扩展，将其分解为不同频率区间的单独的分量进行分析，能够更好的定位有用信息。其中，多元经验模态分解作为完全数据驱动的时频分解方法，能够自适应的将复杂的非平稳信号分解为包含不同时间尺度的本征模态分量，同时保证各通道间的能量分布一致性。本发明通过执行步骤3.1至步骤3.2得到分解后的I个本征模态分量并选择其中对子宫收缩频率相关的I0个作后续分析，排除无关频率分量的干扰。
步骤3.1、采用式(1)所示的多元经验模态分解方法将子宫肌电信号时间序列x*在Q通道上联合分解为I个本征模态分量和残差项：
Figure BDA0004058790450000061
式(1)中，
Figure BDA0004058790450000062
表示第q通道的第i个本征模态分量，rq表示第q通道的残差项；
多元经验模态分解的具体实现过程如下：对于Q通道的离散输入的子宫肌电信号x*＝[x* 1，x* 2...，x* q，...，x* Q]T，在Q-1维超球面上沿着角度
Figure BDA0004058790450000063
所给出的方向进行投影，得到一组k维的投影向量，定义为
Figure BDA0004058790450000064
接着采用三次样条插值，将投影向量在每个方向上的极值点连接起来，形成k个包络并对其取平均得到多元信号的局部均值，从多元输入信号x*中迭代减去局部均值，直到它不包含足够多的极值点，且满足成为IMF的两个条件就可被定义为一个IMF；重复上述步骤，直到残差小于预设阈值，IMF的筛选过程结束。需要说明的是，定义IMF的两个条件分别为：局部极值点与零点交叉点的数目必须相等或最多相差一个；局部最大值(上包络)和局部最小值(下包络)的平均值必须为零。
本实施例中，如图2所示，一段子宫肌电信号经过MEMD分解后在六个通道上获得了相同数量的I个本征模态分量IMF，并在频率尺度上逐一对齐。IMF各分量的波动趋势按分解顺序由剧烈到平缓，对应由高频到低频的变化，为了观察本征模态分量的局部细节，截取整段信号的其中五分钟时长作展示。通常I的取值范围在10～15内，示例性的，可以设置为I＝12。
步骤3.2、根据式(2)和式(3)分别计算
Figure BDA0004058790450000065
在第i个本征模态分量上的方差贡献度Varci和相关系数Cori：
Figure BDA0004058790450000066
Figure BDA0004058790450000071
式(2)和式(3)中，
Figure BDA0004058790450000072
和var(x* q)分别表示
Figure BDA0004058790450000073
和x* q的方差，
Figure BDA0004058790450000074
表示
Figure BDA0004058790450000075
和x* q的互方差；
根据{Varci|1≤i≤I}和{Cori|1≤i≤I}，选取I0个对子宫收缩有贡献且相关的本征模态分量，记为
Figure BDA0004058790450000076
其中，
Figure BDA0004058790450000077
表示第q通道的第i0个本征模态分量；
本实施例中，对步骤3.1所获得的12个IMF在每一个分解水平上分别计算方差贡献度和相关系数，并分别绘制箱线图如图3和图4所示，结果表明，所有IMF的前5个水平反映了较高的Varc和Cor值，这表明对EHG信号的重要贡献和相关性很强。我们假设，由于其余6～12级的IMFs两个指标的值较低，这包含了与早产分类无关的频率分量。因此，选择前五个IMF进行后续特征提取。
步骤4、特征提取：
设置窗口长度为Nw和窗口重叠为NL，对每一个长度为N的
Figure BDA0004058790450000078
进行滑动窗分帧处理后，共得到((N-Nw)/NL+1)个滑动窗口时间序列，并将其中一个滑动窗口时间序列记为
Figure BDA0004058790450000079
且
Figure BDA00040587904500000710
其中，
Figure BDA00040587904500000711
表示长度为Nw的时间序列上的第n0个样本点数据；
分别计算每个滑动窗口时间序列的六个特征值，包括：时域上的均方根幅度和自相关过零数，频域上的峰值频率和中值频率，非线性域上的样本熵和模糊熵；并对所有滑动窗口时间序列求得的六个特征值分别取平均，从而得到三维特征矩阵f＝[f(i0，j)，f2(i0，j)，...，fq(i0，j)，...，fQ(i0，j)]T，1≤i0≤I0，1≤j≤6，其中，fq(i0，j)表示第q个通道上的第i0个本征模态分量的第j个特征值；
本实施例中，采用窗口长度为Nw和窗口重叠为NL的滑动窗口分析，提取本征模态分量细节部分的特征信息。窗口长度和重叠尺寸的选择是计算复杂度和特征表示连续性之间的综合考虑。示例性的，可以设置：Nw＝1200，NL＝300。对于每一个通道的特定本征模态分量，从每个窗口提取六个特征，并将所有窗口的平均值作为最终特征结果。因此对于每一个受试者，都会得到涵盖通道、本征模态分量、特征三个维度的矩阵f。其中，给定输入信号为
Figure BDA00040587904500000712
六个特征的具体计算过程如下：
1)均方根幅度(RMS)：通过计算EHG信号的功率来反映子宫收缩的强度，定义为所有样本点的值求平方除以样本点个数再开根号：
Figure BDA0004058790450000081
2)自相关过零点(AuZC)：自相关函数R(τ)测量信号在不同时间的相似性，可用于区分EHG信号中的随机活动和周期活动，定义为：
Figure BDA0004058790450000082
因此，AuZC被定义为从信号
Figure BDA0004058790450000083
的自相关函数R(τ)的峰值开始的第一个过零点：
AuZc＝min{τ，τ∈τxx}，R(τxx)＝0   (6)
3)峰值频率(PF)：通过快速傅里叶变换得到输入信号的功率谱密度P(n0)。PF定义为P(n0)取最大幅值时所对应的频率值：
Figure BDA0004058790450000084
4)中值频率(MDF)：定义为使信号功率谱P(n)上下部分之和相同的频率：
Figure BDA0004058790450000085
5)样本熵(SampEn)：描述了一个有限长度的时间序列的不规则性，时间序列越复杂，熵
值就越高。可定义为两个序列在容忍值r内相似的条件概率的负对数：
Figure BDA0004058790450000086
其中，Am+1(r)和Bm(r)分别表示两个序列匹配m+1和m个点的概率。示例性的，可以设置入维数m＝3，容忍值r＝0.15SD(SD为时间序列的标准差)。
6)模糊熵FuzzyEn：测量了当时间序列的维度发生变化时产生新模式的概率，采用指数函数
Figure BDA0004058790450000087
对相似度的测量进行模糊化，使熵值随着参数的变化连续而平稳地变化，适用于每类之间存在边缘模糊的情况。FuzzyEn被计算为：
Figure BDA0004058790450000088
其中，φm+1(r)和ψm(r)分别表示两个序列在权重相似度s下匹配m+1和m个点的概率。示例性的，可以设置嵌入维数m＝2，容忍值r＝0.2SD，权重相似度s＝3。
步骤5、特征选择：
步骤5.1、将三维特征矩阵f转换为一维的特征集合F＝{F(k)|k＝1，2，…，Q×I0×6}，并通过卡方检验计算第k个特征F(k)的卡方分数，并对所有卡方分数进行降序排序后，选取大于预设阈值th的K1个特征，从而得到初筛选后的特征集合F*＝{F*(k′)|k′＝1，2，…，K1}；其中，F*(k′)表示初筛选后的第k′个特征，K1表示初筛选特征集合中的特征个数；
步骤5.2、基于序列双向搜索算法的特征优化，并设置算法的执行总次数为T：
步骤5.2.1、从特征集合F*中随机生成一个特征子集F*(k0)；
步骤5.2.2、从特征集合F*中依次选择一个特征F*(k+)加入特征子集F*(k0)中，并判断F*(k0)的F1分数是否达到最优，若是，则执行步骤5.2.3，否则，继续执行步骤5.2.2；
步骤5.2.3、依次从特征子集F*(k0)中删除一个特征F(k-)，并判断删除后的特征子集F*(k0)的F1分数是否降低，若是，则在特征子集F*(k0)中保留特征F(k-)，否则继续执行步骤5.2.3；从而得到最优特征子集F**＝{F**(k″)|k″＝1，2，…，K2}，其中，F**(k″)表示双向搜索算法选择的第k″个特征，K2表示最优特征子集中的特征个数；
本实施例中，对高维的特征特征集合f实施特征选择算法去除冗余项和无关分量，提高特征质量，设计了两步式的特征选择算法：5.1特征预筛选和5.2特征优化。步骤5.1中，假设特征和类别之间是相互独立的，卡方值越大，特征和类别之间的相关性就越强，步骤5.2中，使用包装器方法，为分类器定制基于机器学习的最优特征子集，基于启发式搜索思想，设计序列双向搜索特征选择算法(SBS)，通过F1分数的性能反馈进行特征优化。其中，F1分数是分类器模型精确率和召回率的一种加权平均，顾全了查准率和查全率，计算如下：
Figure BDA0004058790450000091
上述特征选择过程中的卡方检验阈值th和序列双向搜索算法的迭代次数T两个参数是可控的，这意味着可以通过最大化验证数据的F1分数来选择全局最优解。设置th的变化范围为0到0.4，增量为0.1，T的变换范围为1到5，增量为1，根据不同的参数设定得到的F1分数变化曲线如图5所示，示例性的，可以设置th＝0.2，T＝4，此时F1分数达到最优性能。
步骤6、按照步骤1-步骤5的过程得到M个待测对象的最优特征子集并形成数据集R＝{F** 1(k″)，F** 2(k″)，...，F** M(k″)|k″＝1，2，…，K2}，其中，F** M(k″)表示第M个待测对象的最优特征子集。
将数据集R分为Ma个已知标签的训练数据
Figure BDA0004058790450000092
及其对应类别标签
Figure BDA0004058790450000093
Mb个未知标签的测试数据
Figure BDA0004058790450000094
使用支持向量机对训练数据
Figure BDA0004058790450000095
建立与类别标签
Figure BDA0004058790450000096
之间的映射关系，找到能够分不同类别的超平面γ，从而基于支持向量机模型的超平面γ对测试数据
Figure BDA0004058790450000097
进行分类，并得到预测的分类结果
Figure BDA0004058790450000098
本实施例中，对由M个待测对象得到的最优特征集合组成的数据集R，按8：2划分训练集和测试集，且保证两个集合中不同类别的比例是一致的，得到Ma个已知标签的训练数据和Mb个未知标签的测试数据。由于EHG信号的特征参数呈现线性不可分割的状态，支持向量机模型中选择RBF核将输入特征映射到高维空间，找到能分离不同类别的超平面，进而对测试数据进行分类预测。此外，模型中惩罚系数C和核尺度gamma等超参数的设置对分类器的性能影响较大，可通过网格参数搜索进行优化。示例性的，可以设置Ma＝240，Mb＝60，C＝{2-3，2-2，2-1，21，22}，gamma＝{5-4，5-3，5-2，5-1}。
本实施例中，一种电子设备，包括存储器以及处理器，该存储器用于存储支持处理器执行上述子宫肌电信号分类方法的程序，该处理器被配置为用于执行该存储器中存储的程序。
本实施例中，一种计算机可读存储介质，是在计算机可读存储介质上存储有计算机程序，该计算机程序被处理器运行时执行上述子宫肌电信号分类方法的步骤。
为了量化评估本发明的分类性能，使用常用的评估指标包括敏感性、特异性、F1分数、准确度、曲线下面积(AUC)对早产预测性能进行评价，通过300位不同受试者的实验数据将本发明方法与传统的子宫肌电信号分类方法相比较，所有的数据处理过程在软件MATLAB2020a和上运行。
四种对比方法如下：
基于短时傅里叶变换的子宫肌电信号分类方法，过程如下：对采集到的表面子宫肌电信号进行短时傅里叶变换，扩展到不同频率区间。为了满足对高频和时间分辨率的需求，窗口的长度设置为n＝200，步长设置为s＝100。当EHG信号的采样频率为20Hz时，短时傅里叶变换后的时间和频率分辨率分别为5s和0.1Hz；接着计算不同频率成分的近似熵和样本熵特征，并应用主成分分析选择前20个主成分作为优化后的特征子集；最后选择合适的机器学习模型，完成子宫肌电信号分类。
基于遗传算法的子宫肌电信号分类方法，过程如下：对采集到的表面子宫肌电信号计算四个频带(0.1～4Hz，0.2～34Hz，0.34～4Hz，0.34～1Hz)上的时域、频域、非线性域特征，加上临床信息共222个特征；接着使用遗传算法优化特征子集，降低特征维度，它是一种基于种群的启发式搜索方法，通过遗传功能模拟自然进化过程，完成迭代式的种群优化；最后基于机器学习模型，完成子宫肌电信号分类。
基于EMD和小波包分解的子宫肌电信号分类方法，过程如下：对采集到的表面子宫肌电信号进行EMD分解，得到11个本征模态函数，再通过小波包分解在每个IMF上得到12个低频和高频分量；进一步地，提取包含分形维数、平均能量、样本熵、标准偏差等十个特征；接着粒子群优化算法被用于特征选择，使用巴塔查里亚距离对所选择的重要特征进行排序；最后利用分类器对所有排序的特征进行处理，完成子宫肌电信号分类。
基于希尔伯特黄变换的子宫肌电信号分类方法，过程如下：对采集到的表面子宫肌电信号进行希尔伯特黄变换，这包括EMD和希尔伯特变换两个步骤，首先通过EMD获得一组IMFs，分别计算前十个IMF的瞬时振幅和瞬时频率的熵值，将两者的比值作为最终特征；最后基于机器学习算法建立模型，对输入样本进行训练，完成子宫肌电信号分类。
本发明提供的方法以及对比方法在子宫肌电信号分类任务下后得到的性能表现的均值标准差柱状图如图6所示。实验结果证明，在五项指标下，基于多元经验模态分解的子宫肌电信号分类方法都能够提供较高的异常子宫活动预测精度和稳定性。
综上所述，本实施例中的一种基于多元经验模态分解的子宫肌电信号分类方法，包括：获取待测对象子宫收缩时的表面子宫肌电信号；根据子宫收缩频率范围对信号预处理降噪；通过多元经验模态分解将原始的子宫肌电时序信号扩展到不同尺度的频率分量上，并保证各通道间的频谱一致性，更好的观察宫缩时子宫电活动的细节信息；从时域、频域、时频域和非线性域提取反映肌电信号强度、复杂度的六个指标，形成对子宫肌肉活动的宏观描述；融合卡方检验预筛选和序列双向搜索的特征选择算法得到优化后的特征子集，能够降低特征维度，避免了模型过拟合，去除冗余项和无关噪声分量；基于机器学习算法找到分离正常和异常子宫肌电信号类别的边界平面，挖掘了EHG特征描述与子宫电生理活动之间的潜在关系，实现了客观、有效的异常子宫活动评估。本发明的提出对于临床宫缩监测、早产筛查和产后大出血评估具有重要意义。

Claims (4)

1.一种基于多元经验模态分解的子宫肌电信号分类方法，其特征在于，是按如下步骤进行：
步骤1、信号采集：
步骤1.1、利用子宫肌电测量设备和Q通道的高密度表面电极采集某一待测对象在t时刻子宫收缩时的表面子宫肌电信号数据，记为：x(t)＝[x1(t),x2(t),…,xq(t),…,xQ(t)]T，1≤q≤Q，1≤t≤T1，其中，xq(t)表示第q通道上第t时刻的子宫肌电信号数据，Q为高密度电极的总通道数，T1为记录的总时长；
步骤1.2、在采样频率fs下，得到子宫肌电信号x(t)的离散样本点数据，记为：x(n)＝[x1(n),x2(n),…,xq(n),…,xQ(n)]T，1≤n≤N，其中，xq(n)表示第q通道上第n个样本点的子宫肌电信号数据，N为记录的总样本点数；从而经过T1时长的测量后，得到样本长度为N的子宫肌电时间序列x＝[x1,x2,…,xq,…,xQ]T，其中，xq表示第q通道的子宫肌电时间序列，且xq＝{xq(n)|1≤n≤N}，令子宫肌电时间序列x的类别标签为y∈{1,0}，其中，y＝1表示正常，y＝0代表异常；
步骤2、信号预处理：
使用带通滤波器对子宫肌电时间序列x进行降噪，滤除干扰的生理信号和无关噪声后，得到降噪处理后的Q通道的子宫肌电时间序列x*＝[x* 1,x* 2…,x* q,…,x* Q]T，其中，x* q表示降噪处理后的第q通道的宫肌电时间序列；
步骤3、信号时频扩展：
步骤3.1、采用式(1)所示的多元经验模态分解方法将子宫肌电信号时间序列x*在Q通道上联合分解为I个本征模态分量和残差项：
Figure FDA0004058790440000011
式(1)中，
Figure FDA0004058790440000012
表示第q通道的第i个本征模态分量，rq表示第q通道的残差项；
步骤3.2、根据式(2)和式(3)分别计算
Figure FDA0004058790440000013
在第i个本征模态分量上的方差贡献度Varci和相关系数Cori：
Figure FDA0004058790440000014
Figure FDA0004058790440000015
式(2)和式(3)中，
Figure FDA0004058790440000016
和var(x* q)分别表示
Figure FDA0004058790440000017
和x* q的方差，
Figure FDA0004058790440000018
表示
Figure FDA0004058790440000021
和x* q的互方差；
根据{Varci|1≤i≤I}和{Cori|1≤i≤I}，选取I0个对子宫收缩有贡献且相关的本征模态分量，记为
Figure FDA0004058790440000022
其中，
Figure FDA0004058790440000023
表示第q通道的第i0个本征模态分量；
步骤4、特征提取：
设置窗口长度为Nw和窗口重叠为NL，对每一个长度为N的
Figure FDA0004058790440000024
进行滑动窗分帧处理后，共得到((N-Nw)/NL+1)个滑动窗口时间序列，并将其中一个滑动窗口时间序列记为
Figure FDA0004058790440000025
且
Figure FDA0004058790440000026
其中，
Figure FDA0004058790440000027
表示长度为Nw的时间序列上的第n0个样本点数据；
分别计算每个滑动窗口时间序列的六个特征值，包括：时域上的均方根幅度和自相关过零数，频域上的峰值频率和中值频率，非线性域上的样本熵和模糊熵；并对所有滑动窗口时间序列求得的六个特征值分别取平均，从而得到三维特征矩阵f＝[f(i0,j),f2(i0,j),…,fq(i0,j),…,fQ(i0,j)]T,1≤i0≤i0，1≤j≤6，其中，fq(i0,j)表示第q个通道上的第i0个本征模态分量的第j个特征值；
步骤5、特征选择：
步骤5.1、将三维特征矩阵f转换为一维的特征集合F＝{F(k)|k＝1,2,…,Q×I0×6}，并通过卡方检验计算第k个特征F(k)的卡方分数，并对所有卡方分数进行降序排序后，选取大于预设阈值th的K1个特征，从而得到初筛选后的特征集合F*＝{F*(k')|k'＝1,2,…,K1}；其中，F*(k')表示初筛选后的第k'个特征，K1表示初筛选特征集合中的特征个数；
步骤5.2、基于序列双向搜索算法对初筛选后的特征集合F*进行特征优化，得到最优特征子集F**＝{F**(k”)|k”＝1,2,…,K2}，其中，F**(k”)表示双向搜索算法选择的第k”个特征，K2表示最优特征子集中的特征个数；
步骤6、按照步骤1-步骤5的过程得到M个待测对象的最优特征子集并形成数据集R＝{F** 1(k”),F** 2(k”),…,F** M(k”)|k”＝1,2,…,K2}，其中，F** M(k”)表示第M个待测对象的最优特征子集；
将数据集R分为Ma个已知标签的训练数据
Figure FDA0004058790440000028
及其对应类别标签
Figure FDA0004058790440000029
Mb个未知标签的测试数据
Figure FDA00040587904400000210
使用支持向量机对训练数据
Figure FDA00040587904400000211
建立与类别标签
Figure FDA00040587904400000212
之间的映射关系，找到能够分不同类别的超平面γ，从而基于所述支持向量机模型的超平面γ对测试数据
Figure FDA00040587904400000213
进行分类，并得到预测的分类结果
Figure FDA00040587904400000214
2.根据权利要求1所述的基于多元经验模态分解的子宫肌电信号分类方法，其特征在于，所述步骤5.2包括：
步骤5.2.1、设置算法的执行总次数为T，从特征集合F*中随机生成一个特征子集F*(k0)；
步骤5.2.2、从特征集合F*中依次选择一个特征F*(k+)加入特征子集F*(k0)中，并判断F*(k0)的F1分数是否达到最优，若是，则执行步骤5.2.3，否则，继续执行步骤5.2.2；
步骤5.2.3、依次从特征子集F*(k0)中删除一个特征F(k-)，并判断删除后的特征子集F*(k0)的F1分数是否降低，若是，则在特征子集F*(k0)中保留特征F(k-)，否则继续执行步骤5.2.3；从而得到最优特征子集F**＝{F**(k”)|k”＝1,2,…,K2}，其中，F**(k”)表示双向搜索算法选择的第k”个特征，K2表示最优特征子集中的特征个数。
3.一种电子设备，包括存储器以及处理器，其特征在于，所述存储器用于存储支持处理器执行权利要求1或2所述子宫肌电信号分类方法的程序，所述处理器被配置为用于执行所述存储器中存储的程序。
4.一种计算机可读存储介质，计算机可读存储介质上存储有计算机程序，其特征在于，所述计算机程序被处理器运行时执行权利要求1或2所述子宫肌电信号分类方法的步骤。
CN202310052408.2A 2023-02-02 2023-02-02 一种基于多元经验模态分解的子宫肌电信号分类方法 Pending CN116049735A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310052408.2A CN116049735A (zh) 2023-02-02 2023-02-02 一种基于多元经验模态分解的子宫肌电信号分类方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310052408.2A CN116049735A (zh) 2023-02-02 2023-02-02 一种基于多元经验模态分解的子宫肌电信号分类方法

Publications (1)

Publication Number Publication Date
CN116049735A true CN116049735A (zh) 2023-05-02

Family

ID=86123498

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310052408.2A Pending CN116049735A (zh) 2023-02-02 2023-02-02 一种基于多元经验模态分解的子宫肌电信号分类方法

Country Status (1)

Country Link
CN (1) CN116049735A (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117033965A (zh) * 2023-08-11 2023-11-10 湖北工业大学 生物疫苗数据特征选择方法、装置、设备及介质
CN117874434A (zh) * 2024-03-11 2024-04-12 湖南工程学院 一种环境温湿度智能检测过程中数据清洗方法
CN117982155A (zh) * 2024-04-07 2024-05-07 天津市肿瘤医院（天津医科大学肿瘤医院） 一种胸外科检查数据智能处理方法

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117033965A (zh) * 2023-08-11 2023-11-10 湖北工业大学 生物疫苗数据特征选择方法、装置、设备及介质
CN117874434A (zh) * 2024-03-11 2024-04-12 湖南工程学院 一种环境温湿度智能检测过程中数据清洗方法
CN117874434B (zh) * 2024-03-11 2024-05-24 湖南工程学院 一种环境温湿度智能检测过程中数据清洗方法
CN117982155A (zh) * 2024-04-07 2024-05-07 天津市肿瘤医院（天津医科大学肿瘤医院） 一种胸外科检查数据智能处理方法
CN117982155B (zh) * 2024-04-07 2024-06-11 天津市肿瘤医院（天津医科大学肿瘤医院） 一种胸外科检查数据智能处理方法

Similar Documents

Publication Publication Date Title
CN108714026B (zh) 基于深度卷积神经网络和在线决策融合的细粒度心电信号分类方法
Jager et al. Characterization and automatic classification of preterm and term uterine records
CN116049735A (zh) 一种基于多元经验模态分解的子宫肌电信号分类方法
Sangaiah et al. An intelligent learning approach for improving ECG signal classification and arrhythmia analysis
Rubin et al. Recognizing abnormal heart sounds using deep learning
Ren et al. Machine learning driven approach towards the quality assessment of fresh fruits using non-invasive sensing
CN109998525B (zh) 一种基于判别式深度置信网络的心律失常自动分类方法
Zheng et al. A novel hybrid energy fraction and entropy-based approach for systolic heart murmurs identification
CN104523266B (zh) 一种心电信号自动分类方法
CN104887224A (zh) 面向癫痫脑电信号的特征提取与自动识别方法
CN111449644A (zh) 一种基于时频变换和数据增强技术的生物电信号分类方法
Zeng et al. Classification of myocardial infarction based on hybrid feature extraction and artificial intelligence tools by adopting tunable-Q wavelet transform (TQWT), variational mode decomposition (VMD) and neural networks
CN113158964B (zh) 一种基于残差学习和多粒度特征融合的睡眠分期方法
CN109948396B (zh) 一种心拍分类方法、心拍分类装置及电子设备
CN112932498A (zh) 一种基于深度学习的强泛化能力的t波形态分类系统
Diery et al. Automated ECG diagnostic P-wave analysis using wavelets
Lo et al. Deep learning for detection of fetal ECG from multi-channel abdominal leads
Feng et al. Recovery of missing samples in fetal heart rate recordings with Gaussian processes
CN116361688A (zh) 用于心电节律自动分类的多模态特征融合模型构建方法
Suhail et al. Cardiac disease detection from ECG signal using discrete wavelet transform with machine learning method
Al-Dujaili et al. Novel approach for reinforcement the extraction of ECG signal for twin fetuses based on modified BSS
Ziani Enhancing fetal electrocardiogram classification: A hybrid approach incorporating multimodal data fusion and advanced deep learning models
CN113925495B (zh) 结合统计学习和时频分析的动静脉瘘异常震颤信号识别系统及方法
CN114947850A (zh) 基于脉搏波高斯模型特征的脑力负荷等级客观检测方法
CN115429284A (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