CN108288058B - 一种改进的小波阈值膝关节摆动信号去噪算法 - Google Patents

Abstract

本发明涉及生物医学信号领域，由于采集的膝关节摆动信号存在一些噪声的干扰，为了得到高质量的VAG信号，我们需要对信号做一些预处理，有利于后续的特征分析，这对计算机的辅助诊断和治疗是必要的。本发明根据已有的经验模态分解方法，提出一种改进的小波阈值膝关节摆动信号去噪算法，并对小波阈值去噪函数做了改进。实验仿真结果证明，本发明提出的算法和其他已有的算法相比，得到了较高的信噪和较小的均方根误差，使膝关节摆动信号的质量得到改善。

Description

一种改进的小波阈值膝关节摆动信号去噪算法

技术领域

本发明涉及膝关节摆动信号领域，涉及到集成经验模态分解算法和小波阈值去噪算法，特别涉及到一种改进的小波阈值膝关节摆动信号去噪算法。

背景技术

膝关节损伤和疾病的早期诊断和治疗，对于保障人口健康具有重要意义。目前临床上主要采用有创或微创方法来检测和评估膝关节的损伤与疾病情况，其中，MRI和CT等无创检测方法能够较为有效地实现对膝关节损伤疾病的检测，但不能提供早期诊断，且这些大型仪器设备成本高，不便于日常检查。作为微创手段的关节镜检查能提供诊断的详细信息，但其对膝关节的反复切口，会对患者造成痛苦，且不利于损伤疾病的恢复。另一方面，这些现行的检测方法在展现膝关节动态特征方面均存在一定的缺陷。因此，基于膝关节摆动(VAG)信号分析方法对膝关节损伤病变进行检测和辅助诊断是近年来兴起的一种新型无创检测方法。通过对膝关节摆动信号的特征进行分析，采用机器学习等分类方法，能够有效区分与病理学有关的正常与异常VAG信号，从而可使医生清楚地了解膝关节内部的情况，实现对膝关节损伤和病变的诊断与监测。

由加速度传感器记录的膝关节摆动信号会受一些干扰，主要包括由低频移动产生的噪声，由医疗器械的热效应和半导体引起的随机噪声和膝关节运动时的肌肉收缩干扰。为了得到高质量的VAG信号，我们需要对信号做一些预处理，有利于后续的特征分析，这对计算节的辅助诊断和治疗是必要的。本发明根据已有的经验模态分解减少VAG信号随机噪声的方法，提出了一种改进的集成经验模态分解和小波阈值联合的去噪方法，并对小波阈值去噪函数做了改进，得到了较高的信噪比并减小了均方误差，使VAG信号的质量得到改善。

发明内容

由于采集的膝关节摆动信号存在一些噪声的干扰，为了得到高质量的VAG信号，需要对信号做一些预处理，有利于后续的特征分析，这对计算节的辅助诊断和治疗是必要的。本发明的主要目的就是为了得到质量更高的膝关节摆动信号，提出了一种改进的集成经验模态分解和小波阈值联合的去噪方法，并对小波阈值去噪函数做了改进，得到了较高的信噪比并减小了均方误差，使VAG信号的质量得到改善。

为了达到上述目的，本发明的技术方案为：

一种改进的小波阈值膝关节摆动信号去噪算法，包括以下步骤：

第一步，对膝关节摆动信号进行集成经验模态分解(EEMD)

第二步，对分解后得到的固有模态分量做去趋势波动分析(DFA)

2.1)用DFA方法确定各阶模态分量的标度指数；

2.2)通过相关分析确定k值。

第三步，提出改进的小波阈值函数

3.1)对各阶固有模态分量做小波变换，得到相应的小波变换系数；

3.2)通过改进的小波阈值函数对小波系数做阈值计算，得到去除噪声后估计的小波系数；

第四步，对第三步得到的小波阈值处理后的前k阶模态分量，和EEMD分解的后k阶模态分量进行重构，得到最终去噪后的膝关节摆动信号。

与现有技术相比，本发明的有益效果为：本发明提出的改进的小波阈值函数和其他几种经典的小波阈值函数相比，有较好的性能，对仿真信号用不同的阈值函数对其去噪，本发明改进的小波阈值函数由最高的信噪比和最小的均方误差。本发明提出的去噪方法用在膝关节摆动信号上，可以增加膝关节摆动信号的信噪比和减小均方误差，并改善膝关节摆动信号的质量。

附图说明

图1原始信号和集经验模态分解的前3阶的固有模态分量(IMF)；

图2是集成经验模态EMD分解后的IMF4-IMF8分量；

图3是集成经验模态EMD分解后的IMF9-IMF12及剩余分量；

图4是不同小波阈值函数的波形图；

图5使用本发明提出的方法和其他方法对膝关节摆动信号去噪后的图；(a)原始的含噪声VAG信号，(b)使用本文方法对VAG信号去噪，(c)EMD对VAG信号去噪。

具体实施方式

为使本发明实施的目的、技术方案和优点更加清晰，以下结合本发明的技术方案和附图作进一步详细描述：

一种改进的小波阈值膝关节摆动信号去噪算法，该方法主要可分为四个环节，分别为对VAG信号做EEMD分解；对分解后得到的固有模态分量做去趋势波动分析(DFA)；提出改进的小波阈值函数，采用这个阈值函数对各阶模态分量做小波阈值去噪；对第三步得到的小波阈值处理后的前k阶模态分量和EEMD分解的后k阶模态分量进行重构，得到最终去噪后的膝关节摆动信号。其具体步骤如下：

步骤A.对原始的膝关节摆动信号进行集成经验模态分解(EEMD)，得到各阶固有模态分量，包括如下子步骤：

A1.在第n次试验中，信号x(t)加入一个白噪声系列u_n(t)，获得一个新的时间序列：

y_n(t)＝x(t)+u_n(t),n＝1,2,...,N

式中，N表示集成的数量。

A2.通过三次样条函数获得y_n(t)的极大值包络线u_n(t)和极小值包络线v_n(t)，其均值m_n(t)表示为：

A3.令h_n(t)＝y_n(t)-m_n(t)，如果h_n(t)不满足IMF条件(a)和(b)，令h_n(t)＝y_n(t)，重复第二个步骤，直到h_n(t)为IMF函数；

(a)信号的零点和极值点个数相同或相差一个；

(b)信号的上下包络线关于时间轴对称；

A4.将信号

作为y_n(t)重复以上两个步骤，直到r_n(t)不能继续分解，原信号分解为：

式中，i表示每次分解过程中IMFs的个数，

表示第j个IMF，r_i ⁿ表示y_n(t)在第n次分解过程的剩余量；

A5.重复上述A1～A4四个步骤n次，每次给信号加入不同的白噪声；

A6.对n次分解过程中的第j个IMF取平均值，得到的值

就是EEMD最终的固有模态分量IMF；

步骤B.对分解后得到的固有模态分量做去趋势波动分析(DFA)

采用去趋势波动分析方法确定各阶固有模态分量的标度指数，通过相关分析确定k值，使得C₁-C_k为噪声主导分量，k之后的IMF为信号主导分量。DFA的具体算法如下：

B1.加入给定长度为L分解后的

其平均值为ω_j，完整的IMF时间序列s(m)定义为：

B2.完整的IMF时间序列s(m)被分成几个等长且长度为k1的窗口段，最小二乘拟合得到的趋势定义为s_k1(m)，并适应所选择的窗口样本，局部线性趋势波动的计算是在每一个窗口段中通过从完整的时序信号s(m)减去局部线性趋势波动s_k1(m)，平均波动F(k1)是用局部线性趋势波动的均方根来表示：

采用分形标度指数α定义log₁₀F(k1)和log₁₀k1的线性关系，可以表示为F(k1)～k1^α。当0.5＜α＜1时，时间序列具有长程相关性；当0＜α＜0.5时，时间序列具有反相关性，此时的时间序列被认为是噪声。因为随机噪声不是远程相关的，膝关节摆动信号可以用远程相关的IMFs(α＞0.5)来构造重建，而具有反相关的IMFs(α＞0.5)和最终的单调余量(基线漂移)被认为是噪声。所以分解点k确定为6。

步骤C.提出改进的小波阈值函数，采用这个阈值函数对各阶模态分量做小波阈值去噪，具体步骤如下：

C1.对前k阶模态分量做小波变换，得到相应的小波变换系数；

C2.通过改进的小波阈值函数对小波系数做阈值计算，得到去除噪声后估计的小波系数；假设d_j,k2为小波变换后的系数，d'_j,k2为做阈值处理后的小波系数，j为小波分解的尺度，k2∈Z，t1为阈值，t1取为

σ为噪声的标准差，可以由每一层小波系数估计，其中，

表示第二层小波系数中间值的估计值，Y表示原始信号的长度；Z表示整数，k2是第j层的第k2个数。

硬阈值函数的d'_j,k2在d_j,k2＝t1时不连续，这样由d'_j,k2得到的重构信号会产生附加的振荡；软阈值函数的d'_j,k2在d_j,k2＝t1时是连续的。但当|d_j,k2|≥t1时，做阈值处理后的小波估计系数d'_j,k2和阈值处理前每个尺度下的小波分解系数d_j,k2存在总体的偏差，会影响信号的去噪效果。基于这些函数的缺点，本文提出新的小波阈值函数：

式中，

u的范围是0≤u≤1，可通过m进行参数调节，β为可调参数，令β的范围为[0,1]，通过调节β的值，减少阈值处理后的小波系数d'_j,k2和阈值处理前的小波系数d_j,k2存在的恒定偏差，克服了指数阈值函数和软阈值函数存在着恒定偏差的缺点，也比硬阈值函数的连续性好。表1为同一仿真信号使用四种不同的小波阈值函数去噪得到的信噪比和均方根误差值。从表中可以看出本发明提出的新的小波阈值函数相对于其他三种阈值函数，获得了最高的信噪比和最小的均方误差，说明改进的阈值函数能取得良好的去噪效果。

表1不同小波阈值函数去噪后的信噪比和均方根误差

C3.得到前k阶模态分量估计的小波系数，将这些估计的小波系数相加，得到一个小波阈值处理后的前k阶模态分量。

步骤D.对第三步得到的小波阈值处理后的前k阶模态分量，和EEMD分解的后k阶模态分量进行重构，得到最终去噪后的膝关节摆动信号。表2为为使用本发明提出的方法和其他方法对膝关节摆动信号去噪后的信噪比、均方根误差，平均绝对误差和峰值信噪比值。从表2中可以看出，本发明提出的方法获得较高的信噪比和峰值信噪比值以及较小的均方根误差，平均绝对误差。这些指标为衡量信号去噪效果的指标，因此本发明提出的方法有良好的去噪效果。

表2本文提出的方法和EMD方法的去噪性能评价指标

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，根据本发明的技术方案及其发明构思加以等同替换或改变，都应涵盖在本发明的保护范围之内。

Claims (5)

1.一种改进的小波阈值膝关节摆动信号去噪算法，其特征在于以下步骤：

步骤A.对原始的膝关节摆动信号进行集成经验模态分解EEMD，得到各阶固有模态分量IMF，如下式所示：

其中：N表示集成的数量；n表示信号分解过程的次数；

表示第j个固有模态分量；

步骤B.对分解后得到的固有模态分量做去趋势波动分析DFA

采用去趋势波动分析方法确定各阶固有模态分量的标度指数，通过相关分析确定k值，使得C₁-C_k为噪声主导分量，k之后的IMF为信号主导分量；

步骤C.提出改进的小波阈值函数，采用这个阈值函数对各阶模态分量做小波阈值去噪，具体步骤如下：

C1.对前k阶模态分量做小波变换，得到相应的小波变换系数；

C2.通过改进的小波阈值函数对小波系数做阈值计算，得到去除噪声后估计的小波系数；假设d_j,k2为小波变换后的系数，d'_j,k2为做阈值处理后的小波系数，k2∈Z，t1为阈值，小波阈值函数如下：

其中，u＝1-exp[-p(|d_j,k2|-t1)²]，u的范围是0≤u≤1；β为可调参数，令β的范围为[0,1]；

C3.得到前k阶模态分量估计的小波系数，将这些估计的小波系数相加，得到一个小波阈值处理后的前k阶模态分量；

步骤D.对第三步得到的小波阈值处理后的前k阶模态分量，和EEMD分解的后k阶模态分量进行重构，得到最终去噪后的膝关节摆动信号。

2.根据权利要求1所述的一种改进的小波阈值膝关节摆动信号去噪算法，其特征在于，所述的步骤A包括如下子步骤：

A1.在第n次试验中，信号x(t)加入一个白噪声系列u_n(t)，获得一个新的时间序列：

y_n(t)＝x(t)+u_n(t),n＝1,2,...,N

式中，N表示集成的数量；

A2.通过三次样条函数获得y_n(t)的极大值包络线u_n(t)和极小值包络线v_n(t)，其均值m_n(t)表示为：

A3.令h_n(t)＝y_n(t)-m_n(t)，如果h_n(t)不满足IMF条件(a)和(b)，令h_n(t)＝y_n(t)，重复第二个步骤，直到h_n(t)为IMF函数；

(a)信号的零点和极值点个数相同或相差一个；

(b)信号的上下包络线关于时间轴对称；

A4.将信号

作为y_n(t)重复以上两个步骤，直到

不能继续分解，原信号分解为：

式中，i表示每次分解过程中IMFs的个数，

表示第j个IMF，

表示y_n(t)在第n次分解过程的剩余量；

A5.重复上述A1～A4四个步骤n次，每次给信号加入不同的白噪声；

A6.对n次分解过程中的第j个IMF取平均值，得到的值

就是EEMD最终的固有模态分量。

3.根据权利要求1或2所述的一种改进的小波阈值膝关节摆动信号去噪算法，其特征在于，所述的步骤B中DFA的具体算法如下：

B1.加入给定长度为L分解后的

其平均值为ω_j，完整的IMF时间序列s(m)定义为：

B2.完整的IMF时间序列s(m)被分成几个等长且长度为k1的窗口段，最小二乘拟合得到的趋势定义为s_k1(m)，并适应所选择的窗口样本，局部线性趋势波动的计算是在每一个窗口段中通过从完整的时序信号s(m)减去局部线性趋势波动s_k1(m)，平均波动F(k1)是用局部线性趋势波动的均方根来表示：

采用分形标度指数α定义log₁₀F(k1)和log₁₀k1的线性关系，表示为F(k1):k1^α；当0.5＜α＜1时，时间序列具有长程相关性；当0＜α＜0.5时，时间序列具有反相关性，此时的时间序列被认为是噪声；

当α＞0.5时，膝关节摆动信号用远程相关的IMFs构造重建，而当α＜0.5时，具有反相关的IMFs和最终的单调余量被认为是噪声。

4.根据权利要求1或2所述的一种改进的小波阈值膝关节摆动信号去噪算法，其特征在于，所述的步骤C中阈值t1取为

其中，σ为噪声的标准差，由每一层小波系数估计，

表示第二层小波系数中间值的估计值；Y表示原始信号的长度；w为小波分解的尺度。

5.根据权利要求3所述的一种改进的小波阈值膝关节摆动信号去噪算法，其特征在于，所述的步骤C中阈值t1取为

其中，σ为噪声的标准差，由每一层小波系数估计，

表示第二层小波系数中间值的估计值；Y表示原始信号的长度；w为小波分解的尺度。

Images