基于Matlab的信号平稳性检验系统
摘 要:信号的平稳性检验在信号处理中起着十分重要的作用。介绍Matlab环境下设计和实现信号平稳性检验系统。该系统主要利用替代数据的平稳性特点,通过在时频域中分别计算原数据和对应替代数据的平稳度并相互比较,以实现对信号平稳性的检验。它可以求出输入数据的替代数据,并分析原始数据和替代数据的频域和时频域性质,同时还可以通过计算原始数据和替代数据各自时频域的变化程度来判断原始数据的平稳性。该系统提供了友好的用户界面。实验表明,该系统可方便地完成信号的平稳性检验,对于测试数据给出了较好的检验结果。
关键词:平稳性检验;替代数据;时频分布;随机信号
中图分类号:TP274 文献标识码:A
文章编号:1004-373X(2010)03-083-04
System of Stationarity Test Based on Matlab
SUI Ye,LI Ming
(School of Information Science & Technology,East China Normal University,Shanghai,200241,China)
Abstract:Stationarity test remains a challenge problem in the field of signal processing.Due to the importance of stationa-rity test,to find a solution to that problem is greatly desired.A stationarity test system designed and implemented based on Matlab is introduced.The system realizes stationarity test via comparing characteristics of original data and surrogate in time-frequency distribution which represent their stationarity.And the comparison on the basis of stationarization property of surrogate which is first explored in time-frequency perspective.The system can make surrogates from original data and display their characteristics in both frequency domain and time-frequency domain.The system can also compare original data with its surrogates in time-frequency domain to determine whether it is stationary or non-stationary.The system gives a friendly interface and is convenient to use.The experiment indicates our system can test stationarity of signal well.
0 引 言
信号的平稳性检验在随机信号处理中起着十分基础的作用。由于平稳信号和非平稳信号的性质差别显著,因此在处理信号之前先行判断它的平稳性就显得尤为重要。虽然信号平稳性的定义十分明确,但是实际判断过程却是复杂的,例如观察尺度对信号平稳性判断就有很大的影响。
这一领域的研究已经取得了一定的成果。一些人提出了受限和带参数的非平稳性判定方法,而另一些人则将他们的平稳性判定建立在对原始数据的一些假设上[1]。而对于更一般信号的平稳性检验的研究还没有取得太多成果。文献[2,3]中又提到了这一问题,并且提出了一种新的检验平稳性框架。这一框架混合了时频透视法和有名的替代数据法[4,5]。它的基本思想是引入“可控噪声”,即替代数据。并且由于替代数据的一些特性,它可以作为平稳性的评判标准。本文参考了文献[6]中的平稳性检验方法,设计了一个信号平稳性检验系统,并在Matlab的GUI 开发环境下实现了图形用户界面的设计。实践表明,本系统不但提供了友好的用户界面,并且可以方便地完成信号的平稳性检验。
1 平稳性检验原理
1.1 平稳性定义及其检验的重要性
假设有一个高斯过程{xl(t)}(-∞
为任意确定t时刻的全体平均。同时:
r(t1,t2)=E[xl(t1)xl(t2)](2)
被称为自相关函数(ACF)。
对于一个弱平稳过程,它的μx(t)和r(t1,t2)都是时不变的或者说与时间无关的。因此有:
μx(t)=const(3)
r(t1,t2)=E[xl(t+τ)xl(t)]=r(τ)(4)
式中:τ=t1-t2被称为时延。因此,对于平稳高斯过程{xl(t)},它的自相关函数或者它的功率谱密度函数(PSD)为:
Sxx(ω)=∫∞-∞rxx(τ)e-jωτdτ(5)
足以确定它的性质[7]。
另一方面,如果{xl(t)}是非平稳的,它的μx(t)和r(t1,t2)就是时变的或者说和时间相关的。这样它的PSD就应该放在时频域分析[8]。
由此可见,平稳性检验是任何信号处理前必不可少的一步,它决定了后续处理可以使用何种方法。
1.2 替代数据
替代数据的概念最初是由Theiler和其合作作者提出的[4],这种技术是用来产生一种所谓的“替代数据”,这种替代数据是平稳的,同时保持了原数据的一些相关的统计特性。
Theiler在文献[4]中提出了一种具体的产生替代数据的方法。由这种方法产生的替代数据是平稳的,同时保持了原数据的二阶统计特性。具体地说,替代数据保持了原数据功率谱的幅度值不变。
根据Wiener-Khintchin理论,信号的功率谱等于其傅里叶变换的幅值平方。因此保持信号的功率谱幅度值不变,就是保持其傅里叶变换的幅度值不变。因此,假设原数据为x(t),它的傅里叶变换为X(f)=∫e-i2πtfx(t)dt。则替代数据s(t)由:
s(t)=∫ei2πtf|X(f)|eiφfdf
(6)
产生。其中,φf是在[-π,π]上均匀分布的随机相位。这样就保证了s(t)和x(t)有相同的傅里叶变换幅值。在下面的例子中也可以看到,这样产生的s(t)也是平稳的。
1.3 时频分布
时频分布主要用于分析非平稳随机信号的功率谱。由于非平稳随机信号的功率谱是时变的,因此在原来功率谱的基础上再引入时间轴,成为时频分布(TFD)。TFD可以显示出信号的功率谱随时间的变化情况。
具体来说,根据文献[9]中的定义,信号x(t)的时频分布Sx,K(t,f)可以表示为:
Sx,K(t,f)=1K∑Kk=1∫+∞-∞x(s)hk(s-t)e-i2πfsds2(7)
式中:hk(t)是k阶Hermite函数,定义式为:
hk(t)=1k!2kπe-t2/2Hk(t)(8)
式中:Hk(t)是k阶Hermite多项式。
1.4 平稳性检验
平稳性可以体现在频谱随时间的波动上。具体来说,对于平稳信号,其频谱不随时间变化;而对于非平稳信号,其频谱会随时间改变。因此,可以通过比较不同时间点上频谱的相似程度来判断信号的平稳性。
按照文献[4]中的检验方法,定义不同时间点上的频谱与频谱平均值的距离c(x)n为:
c(x)n=κ(Sx,K(tn,•),
(9)
式中:符号“<>”表示求取平均。这里采用文献[10]中的距离定义:
κ(G,H)=1+∫logG(f)H(f)df•
∫(G~(f)-H~(f))logG~(f)H~(f)df(10)
式中:符号“~”表示对应函数的归一化函数,例如f~(x)=f(x)/max(f(x))。
距离c(x)n随时间的波动情况Θ1可以被定义为c(x)n的方差,即:
Θ1=1N∑Nn=1(c(x)n-
根据替代数据的定义,它是平稳的。将替代数据的c(x)n随时间的波动情况记作Θ0。通过比较Θ1和Θ0,可以确定原数据的平稳性。具体来说,将Θ0的概率密度函数记作f(Θ0),选定适当的门限γ,若f(Θ1)<γ,则判为非平稳随机信号,反之则判为平稳随机信号。
2 用户界面生成
2.1 Matlab中用户界面的生成
Matlab为用户设计图形界面提供了一个高效、方便的集成环境。在Matlab中,基本的图形对象主要包括坐标轴、控件、下拉菜单和内容菜单。用户可以通过这些对象设计出界面友好,功能强大,操作简单的图形用户界面。图形用户界面的生成主要分为以下几个步骤:
(1) 规划所设计的图形用户界面,主要包括:确定需要哪些窗口,每个窗口怎样布局,窗口中的各个对象各有什么功能,对象之间如何配合工作,以及相应的异常处理;
(2) 在Matlab提示行下输入GUIDE,载入用户界面开发环境;
(3) 利用Layout Editor,完成用户面板以及界面的制作,并对相应的按钮及控件属性进行适当的设置;
(4) 在Programme Editor中编辑各个对象的回调函数,实现各个对象的具体功能;
(5) 利用Mfile 编译器生成客户端,完成随机数据仿真系统的设计。
2.2 用户界面介绍
本文所实现的用户界面主要包括两个窗口,分别是主窗口和数据生成窗口。由于Matlab对保存绘图区域有限制,因此设计时没有在窗口中设置固定的绘图区域。窗口只相当于一个命令菜单,所有的绘图将会以独立窗口的形式根据用户需求动态产生。这样便于用户对比和保存图片。下面对主要窗口分别加以介绍。
2.2.1 主窗口介绍
主窗口如图1所示。主窗口主要用于绘制原数据和替代数据的各种波形以及显示平稳性检验结果。
图1 主窗口示意图
其中,“Create / Open Original Data”按钮用来打开数据生成窗口。
“View / Change Parameters”按钮用来查看或改变当前仿真参数,它在原始数据存在的情况下才有效。主要的仿真参数有:
“Time Scale of TFD”和”Frequency Scale of TFD”用于确定绘制TFD图片时的时间/频率轴采样周期,由于计算和显示时频分布图比较耗费时间,将采样周期设大,可以提高速度,但是相应的时频分布图的分辨率会下降。
“Max Level of Hermite Function”用于确定求TFD时所使用Hermite函数的最高阶数。最高阶数越高,则分辨率越高,但是相应的计算时间会加长。
“The Number of Surrogates”用于确定平稳性检验时所用的参考替代数据个数。个数越多,则检验结果越精确,但是会极大地延长计算时间。
“Create Surrogate”按钮用于产生替代数据,其在原始数据存在的情况下才有效。由于替代数据具有随机性,因此用户可以多次产生不同的替代数据,观察它们的性质。
右上方的下拉菜单用于选择需要绘图或者保存数据的对象,主要包括原始数据的时域、频域和时频域图,替代数据的时域、频域和时频域图,以及替代数据的平稳度分布。它在原始数据存在的情况下才有效。
“Show Selected Plot”按钮用于在新窗口中绘制下拉菜单所选图线,它在原始数据存在的情况下才有效。
“Save Selected Data”按钮用于保存下拉菜单所选图线对应的数据,它在原始数据存在的情况下才有效。
2.2.2 数据生成窗口
数据生成窗口如图2所示。数据生成窗口主要用来产生实验用数据或者打开已经存在的实验数据。
图2 数据生成窗口示意图
“Creat Original Data”按钮用来产生测试用数据。按下此按钮后会提示输入产生数据用的参数。由于数据是通过公式:
x(t)=sin[P1sin(P2t)t](12)
产生的调频信号,因此需要确定参数P1和P2,另外还要确定t的区间和采样周期。数据成功产生后会在新建窗口中显示该数据时域波形。如果当前存在数据波形,将会覆盖它。
“Open Original Data”按钮用来打开已经存在的数据文件。选择好文件后会提示输入参数。主要包括读入数据的时间起点、时间采样周期和数据长度。数据成功读入后会在新建窗口中显示该数据时域波形。如果当前存在数据波形,将会覆盖它。
“Confirm”按钮用于确认新建窗口显示的数据就是用户想要的数据,并返回主窗口。它在创建或打开的数据存在的情况下才有效。
3 数据仿真和分析
将实验数据取为调频信号x(t)=sin(sin(t/8)πt)。t起始为0,采样周期为0.1 s,数据长度为400个点,其时域波形和频域波形如图3所示。
图3 原始数据时域和频域示意图
由式(6)产生的替代数据s(t)的时域波形和频域波形如图4所示。
图4 替代数据时域和频域示意图
由图3和图4不难看出,替代数据与原数据的傅里叶变换幅值相同,但替代数据傅里叶变换的相位是随机的。
图5显示了由式(7)计算得到的原数据和替代数据的时频分布图。由图5中可见,原数据的时频分布图有明显的结构性。它表明是非平稳的,而替代数据的时频分布图的结构性较原数据有明显减弱,表明替代数据的平稳性增加。
图5 原数据与替代数据时频分布示意图
由式(11)计算得到的Θ0的概率密度函数f(Θ0)如图6所示,其中一共计算了1 000次替代数据。
由图6可见,替代数据的平稳度主要分布在0.02附近。数据的平稳度落在0~0.04之间可以认为是平稳的,而在此之外可以认为是非平稳的。
由式(11)计算得到的Θ1=0.046。位于上述区间之外,因此判为非平稳。这一结果也与图5所示的结果相吻合。
图6 替代数据平稳度分布示意图
4 结 语
利用替代数据法和时频透视法,并采用Matlab的GUI开发环境,设计了一个信号平稳性检验系统。该系统的用户界面友好。利用该系统可以观察信号及其替代数据的频域和时频域波形,检验信号的平稳性。
参考文献
[1]Hobijin B,Franses P H,Ooms M.Generalization of the KPSS-test for Stationarity[J].Statitica Neerlandica,2004,58(4):483-502.
[2]Xiao J,Borgnat P,Flandrin P.Testing Stationarity with Time-frequency Surrogates[A].Proc.EUSIPCO-07[C].2007.
[3]Xiao J,Borgnat P,Flandrin P,et al.Testing Stationarity with Surrogates a One-class SVM Approach[A].Proc.IEEE Stat.Sig.Proc.Workshop SSP-07[C].Madison,2007:720-724.
[4]Theiler J,Eubank S,Longtin A,et al.Testing for Nonlinea-rity in Time Series:The Method of Surrogate Data[J].Physica D,1992:58-77.
[5]Schreiber T,Schmitz A.Surrogate Time Series[J].Physica D,2000,142:346-382.
[6]Pierre Borgnat,Patrick Flandrin.Stationarization via Surrogates[J].Journal of Statistical Mechanics,2009(1):1-14.
[7]Papoulis A.Probability Random Variables and Stochastic Process[M].New York:McGraw-Hill,1984.
[8]Priestley M B.Evolutionary Spectra and Non-Stationary Process[J].Roy.Stat.Soc.Series B,1965,27(2):204-237.
[9]Bayram M,Baraniuk R G.Multiple Window Time-varying Spectrum Estimation[A].Nonlinear and Nonstationary Signal Processing ed W J Fitzgerald et al[C].Cambridge:Cambridge University Press,2000:292-316.
[10]Proestley M B,Rao T S.A Test of Non-stationarity of Time-series[J].Royal Statictic Society,1969:140-149.
