采用 HHT (Hibert-Huang Transformation) 时间序列分析方法处理从人体采集到的脉搏波信号。方法 通过 经验模态分解(EMD)技术将一非线性、非稳态过程的原始离散数据序列分解为一组内在模态函数(IMFs),然后对每一个 IMF 进行 HT 变换,这样得到的信号幅度和瞬时频率都是时间的函数,即获得脉搏波信号幅度和频率的时间分布。再 根据已获得的 HH 谱,进而得到边际谱。这是一种更具适应性的、新型的、基于模态分解的时间序列数据处理方法。 结果 首先对一系列由标准的周期函数构造而成的时间序列信号进行了 EMD处理,验证HHT方法分解的可行性、有效 性;然后分别对一例正常人脉搏波信号和一例典型的冠心病人脉搏波信号进行分解处理,对得到结果进行了比较。结 论 HHT 方法在生物医学信号处理领域将会有广阔的应用前景。
在中医四诊中,脉象诊断占有非常重要的位置,从手腕的寸、关、尺三部按不同的轻重(沉、中、 浮)可获知人体五脏六腑的病理变化。但是长期以来中医理论对脉象的描述都是根据医师的个人判断,容 易造成脉象概念模糊、笼统。这种凭借主观经验结合其他表象诊断的检测手段,难以给出量化的、较 为确定的诊断结果,某种程度上成为了中医脉诊应用、发展和交流的制约因素。将信号处理技术运用于 脉象分析,提取其中与人体生理和病理变化密切相关的特征参量,将有助于脉象研究的客观化,为中医 诊断学的深入研究提供新的方法。在广泛使用的信号分析处理方法中,傅立叶(Fourier)变换是一种最常用的方法,它采用零到无穷大的各周期分量的叠加来刻画原始信号的时间特性,非常适用于确定性的平 稳信号,但在对非线性、非平稳过程的处理上它不仅缺乏应有的分辨率,还会因强制拟合原数据序列而产生 许多噪声信号。
本文尝试采用 HHT(Hilbert-Huang Transformation)方法[4,5]对脉搏波信号进行分析研究。HHT方法是一 种新的信号分析方法,它不受傅立叶分析的局限,能够进行非线性、非平稳信号线性化和平稳化处理,保留 信号本身特性。该方法是通过经验模态分解(Empirical Mode Decomposition,EMD)方法对一段时间内的 信号序列直接分解,产生一组内在模态函数(Intrinsic Mode Function,IMF),随后使用 HT 变换(Hilbert Transform)得出瞬时的频率和振幅,进而得到边际谱。它是一种更具适应性的时频局部化分析方法。其 特点如下:(1)精确性:更好的过滤性和更精确的数据分析能力;(2)灵活性:处理线性、非线性信号以 及平稳、非平稳信号;(3) 正确性:保留信号的内在性质。
经验模态分解(EMD)与 HH 谱分析方法简介
HHT 方法的原理
在信号分析中,频率是一个很重要的概念。实践证明,与时域表示相比,信号的频域表示往往更能体现信 号的本质特征[6]。HT 变换与其它变换不同,它把信号从时间域仍然变换到时间域。设 X(t)为一时间序列, Y(t)是它的 HT 变换,即
X(t)和 Y(t)形成一复数,可以得到 X(t)对应的解析信号 Z(t)
将信号的瞬时幅值 a(t)、瞬时相位 θ(t)表示如下
方程(1)定义的 HT 变换可以看成是 X(t)关于 1/t 的卷积。其瞬时频率定义式为
对于实信号,利用 HT 变换求出其共轭正交分量,然后对实信号进行解析表示,可以求出该信号的三个 瞬时特征参数,即瞬时幅度、瞬时相位和瞬时频率,从而实现真正意义上瞬时参数的提取。 由瞬时频率的物理意义可知,并不是任意的信号都能用瞬时频率来讨论。当信号满足只含一种振动 模态,没有复杂叠加波的情况才可行。也就是说对于含有复杂波系的信号序列,直接 HT 变换在一定程度 上失去了方法上的有效性。为此,要采用一种新型的信号处理技术——经验模态分解(EMD),其本质就 是对采集到的信号序列进行分模态处理,得到一系列平稳的内在模态函数(IMFs)分量,再对每一个分量 进行 HT 变换得到 HH 谱,进而得到边际谱。
经验模态分解 该分解过程基于一个最基本的假设,即采集的数据是由许多基本的内在模态叠加而成,每一种模态 对应于一种物理过程,它们或线性或非线性,并且具有相同数目的极值点与过零点,即要求在横坐标轴上下对称分布。 不同时间尺度的各种模态根据其特征尺度进行分离。对任意给定时间段,可能同时存在许多运动模 态,它们互相叠加得到原始获得的复杂信号。分离之后每种模态是相互独立的,在连续的过零点之间不 存在其他的极值点。 内在模态函数(IMF)所要满足的判断条件:(1)整组数据极值点和过零点的数目相同或者最多相差一 个;(2)局部极大值包络线和局部极小值包络线的平均值为 0。 在实际信号的处理过程中,完全满足第二个条件是不现实的,所以只要二者的平均值小于一个预先 确定的小量即可。根据定义,可以采用如下方法分解函数:(1)寻找到信号所有的局部极大值并用三次样 条函数插值连接获得上包络线;(2)同样的方法连接局部极小值点作为下包络线。两条包络线应包含了所 有的数据,他们的中心线记为 m1,它们的余项记为 h1,即
X (t) – m1 = h1 (5)
经过这次分离,突起的值成为新的局部极大值和局部极小值。此过程有两个目的:消除载波和让波形更加对称。 这样获得的 h1 可以恢复原来埋没在初始数据中的其他固有模态,也就是说不一定满足内在模态函数 的两个条件,随后的步骤将 h1 看作原始数据,找到包络的中心线 m11 继续分解
h1 – m11 = h11 (6)
重复以上过程,经过 k 次,直到满足以上两个限定条件,结束分解,使得 h1k 成为第一个 IMF 项,即
令 C1 = h1k (8)
C1 为从该数据分离出的第一个 IMF 分量。它应包括了信号的最小时间尺度即最短周期的模态,令原始信 号与C1的差值为剩余信号r1
X(t) -C1 = r1 由于 r1 包含了较长周期的组份,可将其视为原始数据并重复以上的过程获得 C2。其新的差值为
当 rn 为单调序列或者相对原始信号幅度极小可忽略不计时,认为完成提取信号内在模态的过程。综合以 上方程(9)和(10),最后获得:
原始信号可表示为 n 个内在模态函数 Cj和一剩余信号 rn 之和,rn 或是常数量或是一组单 调数据。所有的 IMF 分量经过逆向叠加,最后可以还原回原始数据 X(t)。
HH 谱及边际谱 对式(11)的每一个 IMF 作 HT变换后累加得
信号幅度 a(t)与瞬时频率 ω(t)都是时间的函数,因此可把幅度显示在频率 – 时间平面上,即 构成了 HH 幅度谱,HH 谱精确地描述了信号的幅值在整段上随频率和时间变化的规律。由于能量可用振幅的平方来描述,因此 H (ω, t)也在一定程度上反映了信号能量在频率(或时间)各种尺度上的分布规律。HH 谱 H (ω, t)确定以后,就可以利用下式对时间积分得到 HH 边际谱(marginal Hilbert spectrum):
HH边际谱提供了每一个频率值所对应的总幅度值,在统计意义上表征了整个时间跨度内信号在每个频率点 上能量累积的分布情况。
仿真算例
对已知信号进行EMD分解产生的一系列内在模态函数(IMFs)经过HT变换后再进行分析,是检验HHT 方法实际效果的一种最简单的方法。为此我们设计了含有 4 个分量的时间序列作为模拟信号,其解析表 达式为:
取采样频率为 500Hz 并生成一较长的时间序列。在实际的信号处理中,数据序列都是有限长的,两 端的数据无法保证恰好是极值,因而包络线在信号两端容易发散,并且随着分解的进行误差会向内传播 进而“污染”到整个数据序列,称之为“端点效应”[7]。本文采用分步骤去掉两边数据的方法以消除 “端点效应”,最终截取中间一段进行分析。( 见图 1 ) 如图 1 所示,横坐标表示时间(t),单位秒(s); 纵坐标表示波幅(α),单位:克力(g)。图 1a 代表原 始信号,图 1b 的 C1~C4 代表 4 个 IMF 分量,图 1c 代表重构图和原始图的比较。在图 1c 中,实线代表原始数据,虚线代表重构数据。可以看到经过逆向叠加后的重构数据和原始数据吻合得很好,这说明 EMD 分解技术能把叠加在信号里的 4 个信号真实地分解出来。各内在模态函数的频率和幅度特性可以通 过图 1b 清楚地体现出来。C1~C4 相应于原始信号的 4 个模态,频率依次降低,这是由分解本身的特性所 决定的,EMD 方法总是先把最高频信号分离出来,也就是按特征时间尺度从小到大的顺序。分解得到 的 IMFs 分量与原始序列相应分量符合得相当好。
根据式(13)得出边际谱图 2:
从边际谱图可以清楚地看到在整个时间区域内的频率分布以及每一种内在模态所对应的幅度,进而可 以得出系统能量在整个频域内的分布情况。 对时频分布图的分析可以看出,该方法对于时间域和频率域分辨率都很高,对于能量突变的定位与 检测能力较强,具有较好的非稳态动态变化的时频刻画能力。
讨论
通过上述分析,可看到 HHT方法的高精度分解优势与非线性动态数据时频刻画能力,但要真正有效 地用好HHT方法,还应注意以下几点:
(1) 波动现象。HHT方法多次用到三次样条函数,这可能会引起过冲或欠冲现象。如高阶 IMF 序列(低 频部分)可能会出现不应有的波动现象,对于先分解出来的 IMF,其瞬时频率也有一定程度的波动。选择一 定的 IMF 序列进行 HT 变换,可在一定程度上避免样条函数造成的波动影响。
(2) 端点效应。端点效应的解决程度直接影响 HHT 方法的应用效果,解决得好,时频分析效果就好, 否则结果可能不理想甚至是完全错误的。消除端点效应的影响是 HHT 方法一个重要研究课题。 上面介绍了 HHT 方法、仿真算例和脉搏波形的时频估计与分析等初步结果,可得如下结论:
(1) 经EMD分解变换得到的IMF序列是直接从原始时序数据中分离出来的,事先无需确定分解次数,不 受人为因素影响,因此 IMF 序列能更好反映原始数据固有的物理特性,其分解是客观的,内在的和自适应的。 每阶IMF序列都代表了某种特定意义的频带信息,而这些信息的具体含义需要与中、西医合作研究来确定。
(2) 原始信号经过EMD分解获得的IMF序列具有稳态性,不但可用HHT变换得出三维时频谱值,而且也 可用于多尺度建模、谱分析等工作,这是其他方法所不具备的优势。
(3) 将HHT方法用于对非平稳动态变化的脉搏波信号的分析,预计能够敏感地捕捉到脉搏波形信号随时 间和频率动态变化的不同分量的主要特征,可得到清晰的图像;能够反应能量突变点信息,其时频局部定位 能力较强。HHT时频谱的大部分能量主要集中在一定的时间和频率范围内,而不是整个时频空间,这可能真 实地反应了非线性序列的本质特征。 HHT方法得到的结果不受窗函数影响与时频测不准原理的限制,具有完全的时频局域化特性,可准确 描述信号的时变特征。这些已有的工作表明:HHT 方法在生物医学领域将会有广阔的应用前景。