当前位置: 首页 > 期刊 > 《中国生物医学工程学报》 > 2000年第3期
编号:10503492
超声多普勒信号准平稳分段的研究
http://www.100md.com 《中国生物医学工程学报》 2000年第3期
     张乃强 王威琪 汪源源 余建国 吴晓峰 刘斌

    超声多普勒音频信号是一种非平稳时间序列,其功率谱是时变谱。为了能够准确有效地计算它的功率谱,要兼顾时间分辨力和时间窗截断效应。以前采用10ms固定长度时间窗进行功率谱估计虽有其合理性但并不完善。本文利用线性调频信号研究超声多普勒信号的非平稳特性,并在此基础上得出了超声多普勒信号准平稳时间段随心动时相不同而不同的结论,为超声多普勒信号谱分析的时间窗选取提供依据。

    关键词:超声多普勒;非平稳;功率谱;准平稳分段

    0 引 言

    包括超声多普勒信号在内的许多生物医学信号大多是非平稳的信号。这些信号的非平稳特性体现在它们的统计特征、频谱特征、相关函数等是随时间变化的。

    非平稳信号的谱估计方法与平稳信号有所不同。由于非平稳信号的自相关函数Rx(t,t+τ)不只与时间间隔τ有关,而且与时间起点t有关,即关于τ不对称,根据维纳-辛钦定理,其功率谱为非实,物理解释发生困难。因此严格地讲,只有计算信号的时频分布(如维格纳-维利分布等)才能有效地解决这一问题[1]。当然,在实际应用中往往可以采取近似的方法,用时窗或频窗对信号截取后进行分析[2]。但存在着由于信号时变特性所带来的时间分辨率和频率准确性的矛盾。如何解决此矛盾是非平稳信号谱估计中的关键问题。
, 百拇医药
    目前,在超声多普勒信号的临床应用研究中大都采用声谱图,而且发现一些血管方面的疾病与超声血流声谱的谱宽度有关系[3]。因此,谱估计的质量将直接影响到诊断的准确性[4]

    采用经典的傅氏变换进行谱估计,其时间分辨率与频率分辨率是一对矛盾。这是由于时间窗截断的效应。时间窗长度越短,由于时间分辨力越高,但是截断效应也越大,使得频谱展宽越严重,谱的精确度就越低。反之,时间窗长度越长,虽然截断效应越小,所得频谱结果似乎会越接近这一时间段内信号的实际情况。但此时所得的频谱仍不能准确反映信号的频率成分随时间的变化,也就是说时间分辨力越低。由于傅氏变换没有时域的局域表征能力,我们无法确定频谱中各成分在时间窗中出现的时刻。如图1所示两种不同的时域波形却具有相同的频谱。图1(a)的时域信号为:34201.gif (1100 bytes)
, 百拇医药
    经过频谱分析,显示出有f1、f2两个频率成分;而图1(b)的时域信号为:

    f(t)=cos(f1t)+cos(f2t) t∈[0,2T]

    经过频谱分析,却也显示出有f1、f2两个频率成分。t34301.gif (8235 bytes)

    图1 不同的时域信号a、b及其频谱

    在超声多普勒信号的谱估计中,目前在整个心动周期中采用固定长的时间窗(一般取10ms)进行分析,即认为信号的准平稳段长度为10ms[5]。这样做虽有一定的合理性,但在信号特性变化快时,所得频谱仍无法及时反映谱结构的变化;另一方面在信号特性变化慢时,为了减小截断效应还可以取更长的时间窗。
, 百拇医药
    准平稳分段就是寻找信号的一个时间段,在此时间段内既满足时间分辨力又满足系统准确性的要求。准平稳分段在生物医学信号的研究中有广泛运用,已有数种方法报道[6]。但这些方法的对象是脑电且多集中于时域波形的特征。超声多普勒信号的特殊性在于其准平稳分段方法可结合谱估计来考虑。本文介绍一种利用模拟的非平稳信号来研究谱估计准确性问题的方法,并进而将其用于超声多普勒信号的准平稳分段上。

    1 方法和原理

    超声多普勒频移信号可以表示为[7]

    xd(t)=x(t)ejθ(t)(1)

    其中,x(t)为平稳随机信号,θ(t)为关于t的非线性多项式(t的零次项和一次项可归入平稳随机信号x(t)中)。在一定的近似情况下可忽略二次以上项,即可表示为:
, 百拇医药
    θ(t)=πβt2(2)

    其中斜率β表示t时刻信号瞬时平均频率的变化率,也表示了非平稳的程度,单位为kHz/ms,即瞬时平均频率可表示为:

    finst(t)=2πβt(3)

    当β取一定值时xd(t)为线性调频信号。加时间窗以后的多普勒信号xwd(t)为:

    xwd(t)=w(t)xd(t)=w(t)x(t)ejπβt2=w1(t)x(t)(4)

    其中w(t)为窗函数,w1(t)为:

, http://www.100md.com     w1(t)=w(t)ejπβt2(5)

    这样,信号的时间窗、非平稳性均在w1(t)中。多普勒信号xwd(t)的功率谱密度Sxd(f)为:

    Sxd(f)=W1(f)2Sx(f)(6)

    其中,Sx(f)是x(t)的功率谱,是平稳的,W1(f)2是w1(t)的功率谱,从(5)式可以看出,w1(t)的谱估计结果W1(f)中兼有时间窗w(t)的卷积效应和信号本身的时变特性ejπβt2带入的影响。
, 百拇医药
    为了考察非平稳信号时间分辨力和频率准确性(频率分辨力)的情况,对每一给定的β值,可以产生一个线性调频信号,然后分析不同时间窗下的谱估计结果。当分析的窗长度为Δt=t2-t1时,线性调频信号的平均频率变化范围是[βt1,βt2]。设B为w1(t)经傅氏变换所得功率谱密度的-3dB宽度,则B的大小和窗长度Δt有关,如图2所示。单对w(t)而言,Δt越小,窗引起的展宽越厉害,B越大,在图2中曲线开始部分随Δt单调下降;另一方面对ejπβt2而言,Δt越大,信号的非平稳引起B的展宽越大,在图2中的直线部分随Δt单调上升。由于w1(t)是w(t)和ejπβt2的乘积,综合两方面的因素,在某一特定点Δt,谱宽度B达到极小值,此时最优的Δt可进行谱估计的时间分段,也就可以作为准平稳分段的标准。t34401.gif (1538 bytes)
, 百拇医药
    图2 功率谱宽度B与谱分析时间窗宽度

    Δt的关系,直线表示信号实际覆盖宽度βΔt

    理论上可以计算最优Δt的位置,但已有人证明除了在特殊情况下,没有一般解析式可以表示该曲线[8]。但对于高斯型时间窗,可以用以下几式可以说明这两者对谱宽的影响[8]34401.gif (526 bytes) (7)34402.gif (851 bytes) (8)

    其中34403.gif (1172 bytes)是高斯窗的RMS宽度,此时w1(t)功率谱的RMS宽度δm的变化情况是:当取定β值,时窗长度从短逐渐增长,谱宽度逐渐减小并趋向极小值;随着窗宽度继续增加,谱宽度开始增加。
, 百拇医药
    2 计算方法和结果

    我们采用不同β的线性调频信号和最常见的矩形窗对上述方法进行了模拟和验证,取得了一些结果。

    图3(a)给出了颈动脉多普勒信号的平均频率曲线。由于谱宽B只和β的绝对值有关,图3(b)给出了相应的平均频率变化率的绝对值。可以看出,颈动脉血流多普勒信号的平均频率变化率随着心动周期的不同时相而变化的。在心脏收缩期,由于是射血过程,血流速度急剧变化,相应的多普勒频移变化较大,通常可达到0.04kHz/ms~0.06kHz/ms,最大可超过0.10kHz/ms;而在舒张期,由于血流速度趋向平缓,相应频率变化也较小,基本上在0.01kHz/ms~0.02kHz/ms甚至更低。频率变化率在整个心动周期中的变化范围是比较大的。t34501.gif (5509 bytes)
, http://www.100md.com
    图3 (a)平均频率和(b)平均频率变化率的绝对值随时间的变化

    根据人体颈动脉多普勒信号的频率变化情况(图3(a)(b)),我们选取β值从0.01kHz/ms到0.20kHz/ms的线性调频信号作为模拟研究的信号。对于每一个不同的β,我们依次用不同长度的时间窗去取样数据,用傅氏变换计算其功率谱并计算谱的宽度。计算结果如图4中所示。对每一个β值,我们可以得出相应的极值点,如图5所示。并把这一时间段作为该信号的准平稳段。由于采用数值算法带来一定误差,曲线并不光滑。t34502.gif (2727 bytes)

    图4 采用不同β值的线性调频信号计算其谱宽与时间窗长度的关系所得的曲线,图中自下而上分别是β为0.01,0.02,0.05,0.10和0.20kHz/ms时的曲线t34503.gif (2016 bytes)
, 百拇医药
    图5 斜率β与时间窗长度极值点的关系

    根据上面所述,尽管超声多普勒信号的频谱特性是时变的,但是其瞬时的特性可以用频率变化的斜率β来近似。采用人体颈动脉血流的多普勒信号计算在心动周期各个时相的频率变化率。即以传统方法检测出平均频率并以此计算出多普勒信号的频率变化率β的情况。再用上述计算得到的β与Δt的关系,给出整个心动周期中各个时刻对应的准平稳段长度。结果如图6所示。它表明,在收缩期,准平稳段短至6~8ms,而在舒张期可达14ms或更长。由此可见传统的假设准平稳段长度在整个心动周期内是10ms是不尽合理的,尤其在收缩期和舒张期会有一些偏离。t34601.gif (4585 bytes)

    图6 (a)根据颈动脉平均频率数据计算出的(b)随心动时相变化的准平稳段时间长度
, 百拇医药
    3 结论和讨论

    对于谱分析的时间精度和频谱精度而言,有着如物理学中的测不准原理那样的约束:在一个方面提高了精度,在另一方面就降低了精度。因此必须在两者之间取一折衷。

    目前,虽有一些方法对生物医学信号进行准平稳分段,如在对诱发脑电的准平稳分段上各国研究者提出了大量基于自适应滤波、参数模型建模的方法,比较典型的有逆滤波法[9],自适应AR模型法[10]等。这些方法都面临一个问题,就是如何确定分段的判据。无论采用哪一种判据,都是从模型参数的角度出发,没有联系功率谱。而对超声多普勒信号的准平稳分段研究还基本没有。

    本文提出的超声多普勒信号准平稳分段的方法,与其它对非平稳信号采用自适应滤波、建模等进行准平稳分段的方法不同。它不仅考虑了信号的非平稳特性,而且还结合了时间窗的截断效应进行分析。它对信号的谱分析过程有直接的意义。分析表明,心动周期内的准平稳时间随心动时相而不同。在收缩期,准平稳段较短(6~8ms),而在舒张期则较长(可达14ms或更长)。
, 百拇医药
    将本方法应用于实际的信号处理中时,还需注意一些问题:

    (1)对于实时信号,在没有进行分析之前是无法预知其当前的平均频率变化率即它的非平稳程度的,因而不能确定这时的时间窗长度。这时可以采用对实时信号先用定长时间窗进行“预分析”,待计算得到频率变化率后再采用相应的最优时间窗进行分析的方法,这样所得的谱估计要比实时信号有延迟,或者采用自适应方法根据前几个时刻的时间窗长度确定当前的时间窗长度,待进行谱分析后再检验结果是否符合原来的假定并进行修改。

    (2)其次,即使是对离线的信号作非实时的分析,也较难得到瞬时的频率变化率β。因为多普勒信号并不是单频信号而是有一定带宽的信号。现在用平均频率来计算频率变化率β只是一种近似的方法。

    (3)对于不同的超声多普勒系统,应该有不同的Δt结果。因为不同的系统有不同的超声频率,产生的多普勒频移大小不同;血流速度和夹角变化也会导致多普勒信号的频率范围和频率变化率β不一样。
, 百拇医药
    国家自然科学基金资助项目

    教育部博士点基金资助项目

    张乃强(复旦大学电子工程系,复旦大学生物医学工程研究所,上海 200433)

    王威琪(复旦大学电子工程系,复旦大学生物医学工程研究所,上海 200433)

    汪源源(复旦大学电子工程系,复旦大学生物医学工程研究所,上海 200433)

    余建国(复旦大学电子工程系,复旦大学生物医学工程研究所,上海 200433)

    吴晓峰(复旦大学电子工程系,复旦大学生物医学工程研究所,上海 200433)

    刘斌(复旦大学电子工程系,复旦大学生物医学工程研究所,上海 200433)
, 百拇医药
    参考文献

    1,杨福生.随机信号分析.北京:清华大学出版社,1990

    2,张作生主编.生物医学工程前沿.合肥:中国科技大学出版社,1993,10~18

    3,Kassam M, Johnson KW, Cobbold RSC. Quantitative estimation of spectral broadening for the diagnosis of carotid arterial disease: Method and in vitro results. Ultrasound Med Biol,1985,11:425~433

    4,余建国,汪源源,邵谦明等.伪彩色实时声谱系统.复旦学报(自然科学版),1993,32(3):274~279
, http://www.100md.com
    5,Mo LYL, Cobbold RSC. Speckle in CW Doppler ultrasound spectra: a simulation study. IEEE Trans UFFC, 1986,33(6):747~752

    6,杨福生,高上凯.生物医学信号处理.北京:高等教育出版社,1989,543

    7,Fish PJ. Nonstationarity broadening in pulsed Doppler spectrum measurements. Ultrasound Med Biol, 1991,17(2):147~155

    8,Skolnik M. Radar handbook. New York: McGraw-Hill, 1970

    9,Bodenstein G, et al. Feature Extraction from EEG by Adaptive Segmentation. Proc IEEE, 1977,642~652

    10,杨先庆等.基于自适应AR模型的准平稳随机信号分段方法.第二届全国信号处理学术会议论文集,1986, http://www.100md.com