基于系统冲击响应的时间序列建模方法

徐洪涛,王跃钢,崔洪亮,蒲 源

(1.第二炮兵工程学院,陕西西安 710025;2.第二炮兵士官学校,山东青州 262500)

0 引言

工程实践中常常要对复杂的结构进行系统辨识,然而在实际中许多时候无法获得系统的输入信息,一些大型结构更是难以采用人工激励的方式去进行动态模拟试验,也很难从它们的实际工作状态中获得完整的激励信号。这些结构在工作中的载荷虽难以测量,可是响应信号却容易测得。基于系统的响应信号,直接对系统参数识别的方法有很多:如传统系统辨识的Prony法、ITD法等[1],以及通过时序建模的ARMA(Auto-regressive Moving-average)模型法[2]。对于建立ARMA模型的时间序列,需采用平稳、正态、0均值的随机过程。但实际上许多时间序列并不满足0均值和平稳性。虽经0化处理和剔除趋势项等方法可以得到0均值、近乎平稳的随机过程。但处理后的时间序列用于建模并识别参数则会使识别结果与实际值有很大的偏差。在工程实际中很难得到系统准确的白噪声响应信号,因而给系统ARMA模型辨识带来一定的不便。

针对上述问题,本文提出一种利用系统的冲击响应信号进行ARMA时序建模的新方法,设计相应的实验对提出的算法进行仿真验证,并给出了相应的实验结果,证明本文方法的可行性。

1 ARMA时间序列

时间序列是一组按时间顺序排列的动态数据,在这些数据之间有着内在必然的联系及规律性,通过建模分析、研究这些数据便可以辨识出系统的特性参数。ARMA时间序列法是利用参数模型对这些有序随机数据进行分析、研究、处理的一种方法。其数学描述如下[3]:

若离散随机过程 x(n)服从线性差分方程:

式中,e(n)是一离散白噪声,则称随机过程x(n)为一ARMA过程。而式(1)所示的差分方程称为ARMA模型。系数a1,…,ap和b1,…,bq分别称为自回归(AR:autoregressive)参数和滑动平均(MA:moving average)参数,而p和q分别叫做AR阶数和MA阶数,以符号ARMA(p,q)表示。该方法的建模本质是将待识别的系统描述成一个以白噪声e(n)为激励,x(n)为输出的离散系统,通过输出x(n)的信息(相关函数阵)来估计模型的AR参数和MA参数。然而在实际工程应用中,获取系统的白噪声输入响应比较困难,从而无法得到输出x(n)准确的相关函数信息,因而利用常规Yule-walker方程无法解决模型参数的辨识问题。

2 方法原理

对任意ARMA模型来说,白噪声输入响应信号与冲击输入响应信号的相关函数R(τ)之间有一定的联系(成比例),而且在计算ARMA模型的参数时,仅利用了响应信号相关函数R(τ)的有关信息,因此完全可以用冲击响应信号的有关信息来替代白噪声输出来进行模型参数的计算,从而完整建立系统的ARMA模型。工程中也很容易得到冲击响应数据,中小型结构可以通过锤击法获得冲击响应信号,大型结构可以通过爆破方式激励或随机减量技术从结构的响应中提取出自由衰减信号。本文基于这一原理利用系统冲击响应进行参数辨识,避免了求取白噪声响应的不便。

3 方法推导

3.1 模型AR参数辨识

系统激励记为白噪声e(n)和冲击信号δ(n),记x(n)为白噪声输出,h(n)为系统冲击响应输出,对应输出的相关函数分别记为Rx(τ)和Rh(τ)。

由ARMA模型的Wold分解定理[3]可知,因果ARMA过程x(n)具有唯一平稳解

其相关函数可以表示为:

由于实际冲击响应信号总是有限长度,设为N,故可以得到冲击响应信号的相关函数表达式为:

故式(3)可以表示为:

从式(5)可以发现相关函数 Rx(τ)和Rh(τ)之间成正比关系,分析发现正是由于这一关系给参数的估计带来了便利,从而避免了求取系统的白噪声输出,直接使用冲击响应输出即可。

接下来给出常规模型AR参数的修正 Yule-Walker方程如下:

从方程形式可以发现将方程两边同时除以Nσ2即可得到式(7):

式(7)中用冲击响应相关函数 Rh(τ)代替Rx(τ)同样可以准确求得模型的AR参数,从而避免了求取系统白噪声输出的不便。

3.2 模型MA参数辨识

目前关于MA参数求解的方法大致可以分为三种:一是谱分解法,二是利用高阶AR模型来近似MA模型,三是最大似然估计法,其中最为有效的是第二种方法。

由Wold分解定理可知,对信号 x(n)可以建立一个无穷阶的AR模型,即:

那么可以用它表示一个q阶的MA模型,即

于是有A∞(z)B(z)=1,将此式两边取Z反变换,左边应该是a(k),k=1,…,∞和b(k),k=1,…,q的卷积 ,即

在实际工作中,只可能建立一个有限阶次的AR模型,如p,p≫q,其AR模型的系数为(1),(2),…(p),利用该组参数近似MA模型,将其代入到式(10)中必定会产生误差,即

综上所述,运用AR参数求解方法,MA参数求解步骤归纳如下:

1)由N点数据x(n)建立一个p阶的AR模型,要求p≫q,可以采用上文介绍的方法来求解p阶的AR参数(k),k=1,…,p。

2)利用^ap(k),k=1,…,p建立一个等效的q阶的AR模型,再一次利用求解AR参数的方法,即可得到MA参数b(k),k=1,…,q。

3.3 模型ARMA参数辨识

根据上文的推导,可以将常规的ARMA建模方法作一改进,得出如下计算步骤[4]:

1)采用修正 Yule-Walker方程来求解模型的AR部分参数,在求解过程中,用冲击响应函数的相关函数Rh(τ)取代 Rx(τ),求得 AR 参数a1,a2,…,ap;

3)采用上文中提到的求解MA参数的方法求出b1,b2,…,bq,从而实现整个ARMA(p,q)参数的估计。

4 仿真算例

根据上文讨论,对文中研究的辨识算法进行计算机仿真实验。建模方案采用ARMA(2n,2n-1)方案,该方案是S.M.WU和S.M.Pandit提出的一种比较适合工程应用的一种建模方案,建立的是ARMA(2n,2n-1)模型[5]。从n=1开始拟合模型,并按照n=n+1的原则进行递推,直至n=50模拟结束。为了检验每次拟合模型的准确性,考虑到用相似性系数这一标准来衡量。

相似系数:相似系数是描述估计信号与原始信号相似性的参数,其定义如下:

当 y i=csj(c为常数)时,ζij=1;当y i与sj相互独立时,ζij=0。

通过50次模型拟合,辨识出了50种阶次逐步增加的ARMA模型。对每一模型,我们在辨识出其全部参数后,利用估计参数和系统阶次来估计ARMA(2n,2n-1)模型,然后求出其估计冲击响应输出,联合我们先前得到的系统实际冲击响应输出进行相似性计算,求出相似性系数。得到50个不同的相似性系数,选择其中最大的,其相应的模型应该与实际系统更加接近,更能体现系统的特性。

算法仿真实验:

现给出一因果、稳定的ARMA(30,29)模型,其参数如下:

图1 系统理论冲击响应特性曲线Fig.1 Theoretical impulse response of the system

图2 相似系数曲线Fig.2 Curve of the similarity coefficient

采用上述算法进行参数估计,首先给出系统的冲击响应特性曲线,如图1。对冲击响应信号进行50次递推仿真,得出相似系数与n的关系曲线,如图2。从图中可以发现大约n

=13时,相似系数达到了最大值,几乎接近于1,因此可以建立ARMA(26,25)模型来描述原始系统,

通过上文中的方法来求解各部分参数:

通过估计参数来估计ARMA(26,25)模型,然后再次计算估计模型的冲击响应特性,得到估计曲线,如图3。

图3 系统估计冲击响应曲线图Fig.3 Estimated Impulse response of the system

比较系统理论和系统估计的冲击响应特性曲线可以发现,估计系统特性不仅在总体上与真实系统很接近,而且在某些局部特性(采样点50左右)也能够较完整的复现真实系统的特性。众所周知,一个系统的冲击响应特性决定了系统的本质特性,也就是说两个系统如果具有相同的冲击响应特性,那么可以说是两系统等效的。因此,仿真试验表明文中论述的建模方法能够准确的对真实系统进行建模,效果较好,达到了冲击响应特性的完整复现。

5 结论

常规ARMA时间序列建模方法需要获取系统的白噪声响应信号,而在工程实际中却很难准确得到系统的白噪声响应,本文针对该问题提出了一种基于系统冲击响应的时间序列建模新方法。该方法利用系统白噪声输入响应信号与冲击输入响应信号的相关函数R(τ)比例关系,对常规方法进行改进从而实现系统的ARMA建模。计算机仿真实验表明:本文方法利用系统冲击响应输出信号既可有效地实现系统的ARMA建模,同时又避免了常规建模方法中获取系统白噪声响应的不便。由于在工程实际中可以方便的获取系统的冲击响应信号,从而提高了该建模方法的工程实用性。

[1]李德葆,陆秋海.试验模态分析及其应用[M].北京:科学出版社,2001.

[2]郭永刚,许亮华,水小平.基于冲击响应数据的ARMA法建模以及模态参数识别[J].地震工程与工程振动,2006,26(5):167-171.GUO Yonggang,XU Lianghua,SHUI Xiaoping.ARMA model building and system indentification by using impulse time response[J].Journal of Earthquake Engineering and Engineering Vibration,2006,26(5):167-171.

[3]张贤达.现代信号处理(第二版)[M].北京:清华大学出版社,2002.

[4]胡广书.数字信号处理——理论、算法与实现[M].北京:清华大学出版社,2003.

[5]杨叔子,吴雅.时间序列分析的工程应用(上册)[M].武汉:华中理工大学出版社,1991.