符号动力学信息熵参数优化方法研究

丁 闯，张兵志,2，冯辅周，吴守军

（1.陆军装甲兵学院 车辆工程系，北京 100072； 2.北方特种车辆研究所，北京 100072）

符号动力学信息熵（Symbolic Dynamic Entropy，SDE）能够有效检测振动信号时间序列的动力学特性[1]，因此，SDE在特征提取及故障诊断中的应用越来越引起人们的重视。SDE具有计算速度快、对数据量要求小且能捕捉到时间序列的非平稳特性、能够快速有效提取振动信号的运行状态等特征。

然而，SDE的计算结果受自身参数的影响比较大，参数的选择主要依靠人为经验，限制了SDE的应用。本文将对参数确定的2种观点进行阐述，对比了2种观点中参数确定的典型方法，并通过对比优化参数，提高SDE作为特征的算法效率和敏感度。

1 SDE原理

SDE最早由Kurths等提出[2]，计算方法过程如下：

设振动信号序列X={x(i),i=1,2,…N}，即时域信号由N个采样点元素组成。

（1）振动信号符号化。将时间序列X转化到符号域S={s(i),i=1,2,…,N}，即s(i)代表符号域S中的第i个元素。符号域中符号使用si表示，其中，si∈{0,1,2,3}，此处si代表符号{0,1,2,3}的具体某个符号，具体转化方法为

式中：μ为时间序列X的平均值，α为权重，一般取α=0.1[3]。

（2）相空间重构。对符号域的时间序列S=进行相空间重构，得到矩阵S0

其中：m为嵌入维数，τ为延迟时间，j=1,2,3,…,K，其中，K=N-(m-1)τ。重构矩阵中的每一行S0(j)为一个重构分量，重构矩阵中共有K个重构分量。

重构矩阵中每个重构分量为一个独立的状态模式，由于S0(j)有m个元素，每个元素的符号共有4种可能，因此共有n=4m种状态模式。

（3）计算各状态模式出现的概率。统计各状态模式出现的次数，记为Z(l)(l=1,2,3,…,n)，并计算其出现的概率，记为

（4）计算符号动力学信息熵。

根据SDE原理，相空间重构是算法中的关键一步，其延迟时间和嵌入维数的选择对SDE的计算结果有较大的影响。

2 SDE参数优化

相空间重构最早由Takens在1981年提出[4]，关于延迟时间τ和嵌入维数m的选择，相关学者进行了大量的研究[5]。

目前对于参数选择有两种观点：

(1)延迟时间和嵌入维数的选择互不相关，即两者独立确定，通常使用互信息法和伪近邻法相结合，首先使用互信息法确定延迟时间，然后在延迟时间确定后，使用伪近邻法确定嵌入维数；

(2)延迟时间和嵌入维数的选择存在联系，即两个参数的选择相互依赖，通常使用关联积分法（即CC算法），通过构造统计变量和延迟时间的关系确定最佳延迟时间和嵌入窗宽，根据嵌入窗宽确定嵌入维数。针对不同的应用需要根据其具体效果确定参数选择的具体方法。

对于设备不同运行状态下的时域振动信号，采用相同参数对其进行相空间重构，将得到不同的SDE。因此，可使用系统正常状态下的振动信号确定的参数，对系统不同状态下的振动信号进行SDE计算，将熵值作为特征用于判断系统的运行状态，从而达到特征提取的目的。

2.1 延迟时间和嵌入维数的独立确定方法

（1）基于互信息法确定延迟时间

对于振动信号序列X={x(i),i=1,2,…M}，其符号域的时间序列S={s(i),1,2,…,N}，为了研究s(i)与s(i+τ)的相关性，选取延迟序列s(i+τ)构成新的序列Q={q(i),i=1,2,…M}，将S与Q的相关性大小作为确定延迟时间τ的依据。根据信息论理论，两个时间序列的互信息可表示为

式中：Ps(si)和Pq(qi)分别为对应的符号域时间序列中事件si和qj的概率，Psq(si,qj)为事件si和qj的联合分布概率，此处事件si和qj分别对应相应序列的各个符号。

对于某一特定符号域时间序列S，其与延迟后的时间序列Q的互信息值只与延迟时间有关，令

式中：I(τ)的大小代表了在已知时间序列S与时间延迟后时间序列Q的相关性。当I(τ)=0时，表示S和Q完全不相关；当I(τ)取极小值时，表示S和Q最大可能的不相关，重构过程中，使用I(τ)的第一个极小值对应的延迟时间τ作为最优延迟时间[6]，记为τ0。

（2）基于伪近邻法确定嵌入维数

伪近邻法（False Nearest Neighbor，FNN）是确定最小嵌入维数的一种有效方法，其计算过程为：若嵌入维数为m时的重构相空间的近邻点仍然为嵌入维数为m+1的近邻点，则此点为真近邻点，反之若嵌入维数为m时的重构相空间的近邻点不是嵌入维数为m+1的近邻点，则此点为伪近邻点。当嵌入维数处于某个值m时，伪近邻点百分比（False Nearest Neighbor Percent，FNNP）将突变下降至0或接近0，且随着m的增大不再变化，此突变点对应的m值即为所求嵌入维数，记为m0。

伪近邻法确定嵌入维数的具体步骤如下：

当嵌入维度为m时，针对相空间中的任一重构分量S0(j)，在重构向量中，与S0(j)欧式距离最近的重构向量称之为S0(j)的最近邻点，记为SN0N(j)

此时，有

式中：j≠i，其与重构分量S0(j)的距离为

同理，求得嵌入维数为m+1时，其重构分量S0′(j)与其最近邻点SN0,Nm+1(j)的距离为Dm+1(j)，此时可得

引入两个伪近邻点判据

式中：为符号域时间序列Q的均值，Rtol和Atol为选取的阈值。判据（1）通常用于有限长无噪声的数据，判据（2）通常用于含噪声的数据。

由此可知，重构矩阵中的K个重构分量对应K个最近邻点，若K个最近邻点有l个伪近邻点，则此时伪近邻百分比（即伪近邻率）为

在实际计算过程中，将嵌入维数m从2开始依次增加，计算伪近邻点百分比，直到伪近邻点百分比小于5%或不再下降时，此时的嵌入维数即为最小嵌入维数[8]，记为m0。

2.2 延迟时间和嵌入维数的联合确定方法

延迟时间和嵌入维数的联合确定通常使用关联积分法（C-C算法），该方法最早由Broomhead和Kim提出，是一种同时确定最佳延迟时间和最佳嵌入维数的算法。其基本思想是通过嵌入时间序列的关联积分构造能够代表非线性时间序列相关性的统计量，根据其和延迟时间的关系确定最佳延迟时间和嵌入窗宽，进而确定嵌入维数。C-C算法是根据大量的统计规律中得到的，具有一定的优势。

对于时间序列X对应的符号域时间序列S={s(i),1,2,…,N}，其相空间中的重构分量为S0(j)，得出嵌入时间序列的关联积分方程为

式中：K为相空间矩阵中重构向量的数量，即K=N-(m-1)τ；r为参考半径，即为阈值；dij为重构向量中任一两个重构分量S0(i)和S0(j)的距离，即dij=‖S0(i)-S0(j)‖ ；H(u)为Heaviside函数，其取值为

对于任意的参考半径r，C(m,τ,N,r)代表任意两个重构分量的距离小于r的个数占总数的C2K中的比例。对于符号域时间序列S，首先将时间序列进行划分，分成τ个长度相等的时间序列，τ为延迟时间。然后构造一个与延迟时间τ相关的关联方程：

当延迟时间τ=1时为单个时间序列本身，而此时，有

当延迟时间τ=2时，有2个时间序列，为和{s(2),s(4),…,s(N)} ，长度为，而此时，有

以此类推，可得出当延迟时间τ为任意自然数时，其关联方程为F(m,τ,N,r)。

当N→∞时，则

由上式可知，当时间序列中的点为独立分布且N→∞时，所有的F(m,τ,r)均为零。然而实际序列的长度为有限值，且时间序列中的点是相关的，因此F(m,τ,r)并不是总为零。若参考半径r为任意值、F(m,τ,r)均为极小值点时，时间序列最接近均匀分布。此处，将F(m,τ,r)的最大值和最小值的差定义为

在求解过程中，Brock等通过对各种分布的研究，得出结论：当时间序列的长度N≥500时，对于嵌入维数m，一般取2≤m≤5，对于参考半径r，一般取为时间序列的标准差），可得到较好的结果。而将其应用于符号动力学信息熵特征提取时，若嵌入维数取2，得出的信息熵难以区分故障模式，因此，在计算时，可取1,2,3,4。然后构造以下3个统计量

3 试验验证

为了检验并对比参数优化方法的有效性，在SDE取不同参数条件下对行星齿轮箱振动信号进行特征提取。行星齿轮箱试验台如图1所示。

图1 行星齿轮箱试验台

试验台主要由变频电动机、行星齿轮箱、涡流测功机组成，变频电动机提供动力，由涡流测功机加载，行星齿轮箱的参数如表1所示。

表1 行星齿轮箱参数

分别对行星齿轮箱齿轮正常、太阳轮轮齿局部裂纹故障、行星轮轮齿局部裂纹故障等3种状态信号进行分析[7-9]。

此试验中2种状态的振动信号均在有负载情况下完成，其中，变频电机输入转速为2400 r/min，采样频率为5120 Hz。正常信号、太阳轮裂纹故障、行星轮裂纹故障时域波形如图2所示。

图23 种状态振动信号的时域波形

（1）SDE信息熵参数确定

首先利用互信息法计算正常状态下振动信号符号域时间序列的延迟时间τ与互信息之间的关系，如图3所示。

图3 振动仿真信号互信息与延迟时间的关系

由图3可知，互信息第一次取得极小值时，最佳延迟时间τ0=8。

当输入延迟时间τ0=8时，使用伪近邻法求解时间序列的嵌入维数，使用判据（2），选择阈值Atol=2，求出此时间序列的伪近邻率与嵌入维数的关系如图4所示，由图4可知，最佳嵌入维数m0=3。

使用联合积分法即C-C方法确定延迟时间和嵌入维数，图5为由联合积分法得出的相关参数变化曲线。

由图5可知，当图中)变化曲线取第1个零点或接近零点时，τ0=20，且Fcor(τ)变化曲线取最小值时，最佳嵌入窗宽τw=25，此时根据式（21）可得，m=3；当图中(τ)变化曲线取第一个极小值时，τ0=35，同理得m=2。由前文可知，嵌入维数至少大于2，因此，得出最佳延迟时间τ0=20，最佳嵌入维数m0=3。

图4 振动仿真信号伪近邻率与嵌入维数的关系

图5 联合积分法各参数的变化曲线

（2）特征提取

取齿轮正常、太阳轮裂纹故障、行星轮裂纹故障各20个样本，使用SDE计算3种状态下各个样本的熵值。使用独立方法和联合积分法确定参数的各状态下SDE如图6、图7所示。

图6 根据独立确定方法确定参数的各状态SDE

图7 根据联合积分法确定参数的各状态SDE

为进一步突出本文2种方法确定参数的实际应用效果，在实际计算中，增加任意选定的参数计算，任意选定延迟时间τ0=4，嵌入维数m0=3，此参数下各状态下SDE如图8所示。

图8 任意给定参数的各状态SDE

由图6、图7、图8可知，使用独立方法确定的参数计算得到的3种状态下的SDE值相差最大，即分类效果最好，且同种状态下的SDE较稳定。因此使用互信息法确定延迟时间，根据伪近邻法确定嵌入维数，具有更好的分类效果。

4 结语

本文重点研究符号动力学信息熵参数优化确定方法。根据其算法原理，使用相空间重构的方法确定算法参数，克服了人为主观或凭经验确定算法参数的不确定性。根据当前确定相空间参数的两种不同观点，介绍了2种方法中典型参数确定方法的基本原理，并利用行星齿轮箱实测信号对所提参数优化方法进行了验证。结果表明，相空间重构参数的独立确定方法比联合确定方法更优，这为符号动力学信息熵算法中参数的确定奠定了基础。