基于符号分析的极大联合熵延迟时间求取方法*

张淑清 李新新 张立国 胡永涛 李亮

(燕山大学电气工程学院，河北省测试计量技术及仪器重点实验室，秦皇岛 066004)

(2012年11月30日收到;2013年2月4日收到修改稿)

1 引言

混沌自从被确认是长期可控制的和短期可预测以后，其研究课题引起了大家广泛的兴趣.然而相空间重构作为混沌时间序列处理中的重要课题，其最佳延迟时间的选取至关重要.对于时间序列X(n)，过小的延迟时间τ将导致信息的冗余，使得重构的坐标之间不可分辨，过大的τ将使延迟坐标之间毫不相关，不能代表真实的动力系统[1]，因此选取最佳延迟时间需要使X(n)和X(n+τ)最大限度相互独立却又不是完全无关.

确定相空间重构最佳延迟时间的方法有互信息函数法、自相关函数法和平均位移法等.其中互信息函数方法是估计重构相空间延迟时间的一种有效方法，它在相空间重构中有很广泛的应用.Shaw首先提出互信息第一次达到最小时滞时作为相空间重构的最佳延迟时间[2]，Faster给出了互信息计算的递归算法.但其提出的划分网格计算方法过于复杂，不易于编程实现及推广使用.国内学者提出了相应的改进算法并获得很好的应用，这在一定程度上激发了人们对准确快速地选取相空间重构中最佳延迟时间的兴趣[3].文献[4]提出了确定延迟时间变量的联合熵之第一极大准则，相对互信息方法减少了计算量，但在求取联合熵时采用划分网格的方法，计算相对复杂而且误差较大;文献[5，6]提出了改进的网格标记法，该方法与符号分析法相比计算仍然相对复杂;文献[7]提出的符号分析法求取互信函数，从而确定出最佳延迟时间，可避免复杂的网格划分和标记，但需同时求取X(n)和X(n+τ)的各自信息熵及其联合熵.本文提出了基于符号分析的极大联合熵延迟时间求取方法，较好地解决了上述方法中所存在的各种问题，可以快速有效地提取混沌动力系统中有用定量信息，重构系统的相空间.

2 最佳延迟时间的极大联合熵准则

考虑时间序列X(n)及其延迟时间序列Xτ(n)=X(n+τ)，n=1,2,···,N.根据互信息函数的递推公式[8]，两组序列的互信息可表示为

其中X代表时间序列X(n)，Xτ代表其延迟序列X(n+τ)，H(X)是孤立的 X 的不定性，H(X|Xτ)是已知Xτ的X的不定性，所以Xτ的已知减少了X的不定性，则I(X,Xτ)的第一极小值处的τ即为最佳延迟时间.并且为了计算I(X,Xτ)，Fraser和Swinney提出了复杂的划分网格的方法[2].

互信息I(X,Xτ)表征X和Xτ的相关程度，相空间重构最佳延迟时间的确定是在保证重构坐标之间最大限度地相互独立的基础上提出来的.I(X,Xτ)取极小值时，X和Xτ的相关度也极小，即X和Xτ的联合整体不确定度达到极大.

另一方面，为了保证重构坐标之间最大限度地相互独立，应使X和Xτ联合整体的不确定性达到最大.由信息理论可知，联合熵H(X,Xτ)是X和Xτ的联合整体的不确定性的度量，H(X,Xτ)越大，则X和Xτ联合整体的不确定性也越大.由此推断，联合熵H(X,Xτ)的第一个极大值点即为相空间重构的最佳延迟时间点[4].

同时，混沌吸引子在其不变集具有的遍历特性表现在其状态变量X(n)上，即在稳态情况下，对任意的延迟时间τ，X(n)和X(n+τ)都具有相同的概率分布和统计特性，则H(X)和H(Xτ)的值接近于常数，根据互信息函数(1)式可推得

可见，H(X,Xτ)和 I(X,Xτ)呈近似相反的变化规律.互信息的极小值点即为联合熵的极大值点在理论上得到验证.

因此通过复杂的划分网格或者进行网格标记求取互信息第一极小值点[6]从而确定相空间重构最佳延迟时间可以转换为求取联合熵H(X,Xτ)的第一极大值点.联合熵的计算公式[9]可以表示为

由此可见，求取不同延迟时间下的联合熵，找出联合熵的第一极大值点所对应的τ，即为所求的最佳延迟时间点.

3 符号分析法求取极大联合熵

符号分析法实际上是时间序列的符号化，该方法是通过在几个可能值上将其时间序列离散化，从而将千变万化的数据序列转换为几个数值或特殊符号的符号序列，该过程称为粗粒化过程，这一过程能够从动力系统中有效快速地捕获有用定量信息，并且能够降低噪声的影响.将混沌时间序列通过粗粒化方法转化为符号序列[10-14]，则其轨道演化的信息就通过符号序列达到了编码.

时间序列符号化有二进制符号化规则、角区间符号化规则和概率统计符号化规则等，二进制划分[13,15]是最简单的一种划分规则，只需给定一个阈值P0，时间序列中大于此阈值P0的取1，否则取0.阈值P0的选取有均值，零值等.则时间序列X(n)最终被转换成符号序列{S(n)}，所得符号序列表征了两种数据元模式，即

如图1所示.

图1 二进制符号化规则

二进制分割虽然简便，但是由于其划分比较粗糙，容易导致原始信息的遗失.为了解决这个问题，本文通过对相空间的分割从而使时间序列X(n)实现粗粒化过程[10-14].以符号数取d为例，划分如下：

其中Xn为时间序列X(n)在n时刻的幅值大小，n=1,2,···,N.XC0,XC1,···,XCd为阈值即临界点，Si为符号序列的符号集，并用相应的整数 i来代替，i=0,1,···,d-1. 则时间序列{X(1),X(2),X(3),···,X(N)} 可以转换成整数集序列 {S(1),S(2),S(3),···,S(N)}.

时间序列被转化成长整数集序列{S(1),S(2),···,S(N)} 后，要提取有用信息，需选择标准分割长度L.将长整数集序列分割成长度为L的短序列，L个连续符号被编码为十进制数，形成了新的整数集序列表示短序列从长整数集序列Si的第k个元素S(k)开始[10].

按照上式进行标记和辨识[10-12]，则离散的符号替代连续数据的原始时间序列，如图2所示.

图2 符号化时间序列图

符号数d的个数可以通过使符号熵最大化来寻找[11,12]，为了方便起见，将混沌序列的临界点(d+1)的个数置为11，即d=10.且分割长度L取2.

对时间序列X(n)及其延迟时间序列X(n+τ)，n=1,2,···,N，根据上面所介绍的粗粒化符号方法将其编码成能捕获用定量信息的特殊的十进制数序列Lx(n)和Lxτ(n)，则各个特殊十进制数出现的频率为时间序列分析的指标，即为联合概率P(Lx,Lxτ)，则符号分析法求取的联合熵[16](3)式可以改写为

4 数值仿真试验

为了验证该方法求取最佳延迟时间的准确性和优越性，我们分别对Lorenz，Rossler和Duf fing三种常见的混沌时间序列求取最佳延迟时间并画出其吸引子的重构.

Lorenz系统的数值仿真实验如图3所示，Lorenz系统[17]方程

图3 Lorenz系统(σ=16，b=4，r=45.92)的数值仿真实验 (a)本文方法和互信息函数法求取最佳延迟时间的对比图;(b)Lorenz吸引子在x-y平面上的投影图;(c)τ=11时的Lorenz的重构吸引子图

其中取 σ=16，b=4，r=45.92，此系统为混沌系统，用Runge-Kutta法求解方程(6)，步长h=0.01，取变量x为研究对象，去除前8000个点，得到一个7000个点的时间序列.通过符号分析法求取联合熵H(Lx,Lxτ)的极大值点和通过网格划分求取互信息I(X,Xτ)的极小值点的对比图形如图3(a)所示，从图 3(a)可以看出，H(Lx,Lxτ)和 I(X,Xτ)的变化规律近似相反，并且都是在τ=11时取得第一个所对应的极值点，故可得由联合熵第一极大值点法求得的相空间重构最佳延迟时间达到了与互信息法相同的效果;试验过程中也大大地简化了其计算量.图3(b)是Lorenz吸引子在x-y平面的投影图;图3(c)是τ=11时的重构吸引子图.从图中可以看出τ=11时重构图能够很好的反映出Lorenz系统的双圈拓扑结构.从而验证了符号分析法求取极大联合熵进而确定最佳延迟时间的有效性和优越性.

Rossler系统的数值仿真实验如图4所示，Rossler系统[18]方程

其中取 d=0.2，e=0.2，f=0.5，此时 Rossler系统为混沌系统，用Runge-Kutta法求解方程(7)，步长h=0.05，取变量x为研究对象，去除前50000个点，得到一个8000个点的时间序列.图4(a)是Rossler系统处于混沌状态时的联合熵H(Lx,Lxτ)和互信息I(X,Xτ)的变化规律图，两者的变化成近似相反的规律，并且在互信息取得极小值点处，联合熵处于极大值点，即τ=23;图4(b)是Rossler吸引子在x-y平面的投影图;图4(c)是τ=23时的重构吸引子图.比较重构前后吸引子在x-y平面投影，发现重构吸引子很好地保持Rossler吸引子的单圈结构.从而证明了符号分析法求取联合熵极大值点即为最佳延迟时间的可行性.

图4 Rossler系统(d=0.2，e=0.2，f=0.5)的数值仿真实验 (a)本文方法和互信息函数法求取最佳延迟时间的对比图;(b)Rossler吸引子在x-y平面上的投影图;(c)τ=23时的Rossler的重构吸引子图

Duf fing系统的数值仿真实验如图5所示，Duffing系统[1]方程

其中取 σ=0.06，a=0.7，F=7.5，此时 Duf fing 系统为混沌系统，用Runge-Kutta法求解方程(8)，步长h=0.05，取变量x为研究对象，舍去开始暂态点，得到一个7000个点的时间序列.从图5(a)中能够明显看出Duf fing系统此时的联合熵H(Lx,Lxτ)和互信息I(X,Xτ)的近似相反的变化规律，且由符号分析法求得的联合熵H(Lx,Lxτ)的极大值点为τ=11，与互信息函数法得到的最佳延迟时间τ相同.图5(b)是Duf fing吸引子在x-y平面的投影图;图5(c)是τ=11时的重构吸引子图.比较重构前后吸引子在x-y平面投影发现重构效果良好.

图5 Duf fing系统(σ=0.06，a=0.7，F=7.5)的数值仿真实验 (a)本文方法和互信息函数法求取最佳延迟时间的对比图;(b)Duf fing吸引子在x-y平面上的投影图;(c)τ=11时的Duf fing系统的重构吸引子图

5 结论

采用符号分析与极大联合熵结合求取相空间重构的最佳延迟时间，不仅避免了计算互信息时的复杂的网格标记，而且无需计算X(n)和X(n+τ)的各自的概率分布及信息熵.该方法能达到与互信息方法计算最佳延迟时间同样的效果，而且极大地简化了计算，是求取最佳延迟时间的改进算法.理论分析与仿真实验证明，该方法物理意义明显，实际效果较好，能够更多地保留原系统的动力学特征.

[1]JiangWL，ZhangSQ，WangYQ2005ChaosandWaveletBasedFault Information Diagnosis(Beijing：National Defence Industry Press)p46(in Chinese)[姜万录，张淑清，王益群2005基于混沌和小波的故障信息诊断(北京：国防工业出版社)第46页]

[2]Fraser A M，Swinney H L 1986 Phys.Rev.A 33 1134

[3]Zhang S Q，Zhao Y C，Jia J，Zhang L G，Shangguan H L 2010 Chin.Phys.B 19 060514

[4]Tian Y C 1997 Systems Engineering Theory and Practice 05 48(in Chinese)[田玉楚1997系统工程理论与实践05 48]

[5]Lu X Q，Cao B，Zeng M 2006 Chinese J.Comput Phys.23 184

[6]Zhang J，Fan Y Y，Li H M，Sun H Y，Jia M 2011 Chinese J.Comput Phys.28 469(in Chinese)[张菁，樊养余，李慧敏，孙恒义，贾蒙2011计算物理28 469]

[7]Xiao F H，Yan G R，Han Y H 2005 Acta Phys.Sin.54 0550(in Chinese)[肖方红，阎桂荣，韩宇航2005物理学报54 0550]

[8]Zhang S Q，Jia J，Gao M，Han X 2010 Acta Phys.Sin.59 1576(in Chinese)[张淑清，贾健，高敏，韩叙2010物理学报59 1576]

[9]Zhu X K，Jia Y H 2005 Journal of Geomatics 05 3(in Chinese)[祝晓坤，贾永红2005测绘信息与工程05 3]

[10]Xiao F H，Yan G R，Han Y H 2004 Acta Phys.Sin.53 2877(in Chinese)[肖方红，阎桂荣，韩宇航2004物理学报53 2877]

[11]Lehrman M，Rechester A B 1997 Phys.Rev.Lett.78 54

[12]Rechester A B and White R B 1991 Phys.Lett.A 158 51

[13]Zhang Y 2004 J.Xiangtan Min.Inst.19 75(in Chinese)[张雨 2004湘潭矿业学院学报19 75]

[14]Azad R K，Rao J S，Ramaswamy R 2002 Chao Solitons and Fractals 14 633

[15]Xiang K，Jiang J P 2007 PR&AI 20 154(in Chinese)[向馗，蒋静坪2007模式识别与人工智能20 154]

[16]Shi Y S，Jiang Y，Song Y X 2011 Journal of Aerospace Power 26 670(in Chinese)[史永胜，姜颖，宋云雪2011航空动力学报26 670]

[17]Yu W B 2008 Calculate Experiment and Analysis on Chaos(Beijing：Science Press)p33(in Chinese)[于万波2008混沌的计算实验与分析(北京：科学出版社)第33页]

[18]Liu Z H 2006 Fundamentals and Applications of Chaotic Dynamics(Beijing：Higher Education Press)p71(in Chinese)[刘宗华 2006 混沌动力学基础及其应用(北京：高等教育出版社)第71页]