Lyapunov指数在混沌时间序列预测中的应用

Lyapunov指数的早期应用可以追溯到1959年Cesari的工作。1964年，Henon和Heiles首次在研究混沌系统的相空间邻近轨道的发散时，引入了Lyapunov指数的数值计算。1968年，Lyapunov指数被Oseiedec正式应用在动力系统和各态变历理论的研究中，在以后的文献中有大量关于Lyapunov指数的使用、探讨和研究报道，在当今非线性问题的研究中Lyapunov指数得到越来越多的应用。

混沌系统初值敏感性是指相空间中初始距离很近的两条轨迹会以指数速率发散，Lyapunov指数即是根据相轨迹有无扩散运动特征来判别系统的混沌特性。在相空间中，轨迹间的距离分别表现为线度、面积和体积等。对一维映射x（t+1）=f[x（t）]，假设初始位置x（t0）附近有一点x（t0）+δx（t0），则经过n次迭代后，有

所以

式中：t0与tn分别为预选择的初始时间与当前时间。

定义2.1[12]设相轨迹上两点之间的初始距离为|δx（t0）|，用|δx（tn）|表示经过n次迭代后该两点之间的距离，由式（2.13），有

则称

为系统Lyapunov指数。

当λ＜0时，系统有稳定的不动点，λ=0时，对应着分岔点或系统的周期解，即系统出现周期现象；λ＞0时，系统具有混沌特征。

对二维映射，则用面积改变量代替线度改变量，按以上思路可求得二维Lyapunov指数λ1+λ2，即两个Lyapunov指数之和，表征平均“面积”的辐散率。依此类推，三维相空间可得平均“体积”辐散率。m维相空间中有∑λi，i=1，…，m，表征m 维体积的辐散率。如果将λi按从大到小顺序排列，可得到Lyapunov指数谱，表征系统在相空间中（或沿相空间某一方向）运动的指数增长率。

当相空间是一维（单变量）时，吸引子只可能是稳定定点（不动点），其周围的点（状态）都要趋于这一点，此时必有λ＜0。反之，如果不存在定点或定点不稳定，即不存在吸引子，则λ＞0。由此可见，单变量自治系统不可能有振荡解。

对于二维情形，吸引子或者是稳定定点或者是极限环。对于稳定定点，任意方向的Δxi都要收缩，故此时两个Lyapunov指数都是负的：λ1＜0，λ2＜0。至于极限环，取Δxi沿环线的法线方向，它一定要收缩，故此法线方向的Lyapunov指数总是小于零；当取Δxi沿环的切线方向时，可以想象，它既不会增大也不会缩小，因此沿此切线方向应有λ=0。事实上，哈肯（H.Harken）证明[13，14]所有不终止于定点而又有界的轨道（或吸引子）都至少有一个Lyapunov指数等于零，它表示沿轨线的切线方向既无扩张又无收缩的趋势。所以二维极限环的Lyapunov指数是（λ1，λ2）=（0，—）。

对于三维的情况，类似的分析结果是[15]：

（1）（λ1，λ2，λ3）=（—，—，—）：稳定定点。

（2）（λ1，λ2，λ3）=（0，—，—）：极限环（周期运动）。

（3）（λ1，λ2，λ3）=（0，0，—）：二维环面（准周期振荡）。

（4）（λ1，λ2，λ3）=（+，+，0）：不稳极限环。

（5）（λ1，λ2，λ3）=（+，0，0）：不稳定二维环面。

（6）（λ1，λ2，λ3）=（+，0，—）：奇异吸引子。

上面第一种情况是明显的。对于三维相空间中的极限环，由于其他轨线在垂直于环线的两个方向上都要趋于极限环，故有两个λi之值是负的。对于二维环面，环面法线方向的Lyapunov指数自然是负的，另外两个在环面切面上互相垂直的方向上的λi则都应等于零，即为第三种情形。对于不稳定极限环和不稳定二维环面，是把第二和第三种中的λi的负号变为正号，分别是第四和第五两种情形。对于混沌运动的奇异吸引子，它是全局稳定（收缩作用才能形成吸引子）和局域不稳定（蝴蝶效应）两种因素共同作用的结果，前者要使Lyapunov指数小于零，后者却要使它大于零，于是产生了第六种情况。从上面的分析可以看出，实际的混沌运动空间必须大于二维，n（n＞2）维空间中的某一时刻，两条相邻的轨线之间的距离可以分解在n个不同的方向（各个方向是按轨迹定义的）上，并且不同方向的距离增长率是不同的，每个增长率代表一个Lyapunov指数，既可以大于零，也可以小于或等于零。对于混沌系统，至少有一个Lyapunov指数必须大于零，这是区分奇异吸引子与其他吸引子的重要特征。

可以看出Lyapunov指数的确可以表征系统运动特征：所有Lyapunov指数取值的集合（Lyapunov指数谱）决定系统在相空间的轨线（吸引子）的性质；沿某一方向Lyapunov指数λi取值的正负和大小表示长时间系统在相空间中相邻轨线沿该方向平均发散（λi＞0）或平均收敛（λi＜0）的快慢程度；最大Lyapunov指数λ1决定相邻轨线是（λ1≤0）否（λ1＞0）能靠拢形成稳定轨道或稳定定点；最小Lyapunov指数λn则表示相空间中所有轨线是（λn＜0）否（λn＞0）能收缩成稳定的吸引子。

在实际动力系统混沌识别中，通常只估计最大Lyapunov指数λ1，其主要方法有Wolf法[16，17]、Jacobian法[18]、小数据量法[8，19]等。Wolf法或其演化方法从Lyapunov指数的定义出发，直接基于相轨迹、相平面、相体积等演化来估计Lyapunov指数，在混沌研究中应用十分广泛，图2.3为Logistic映射的Lyapunov指数随参数μ变化而变化的曲线图。下面将分别介绍Wolf法与小数据量法。

图2.3　Logistic映射的Lyapunov指数

1.Wolf法

Wolf等（1985年）提出直接基于相轨线、相平面、相体积等演化来估计Lyapunov指数。这类方法统称为Wolf法，在混沌的研究和基于Lyapunov指数的混沌时间序列预测中，应用十分广泛。

设混沌时间序列x（1），x（2），…，x（t），…，嵌入维数m，时间延迟τ，则重构相空间

取初始时刻为t0，当前时刻为ti，终点时刻为tM，M=N—（m+1）τ，N 为时间序列的终点。设初始点为X（t0），设其与最近邻点X0（t0）的距离为L0，追踪这两点的时间演化，直到ti时刻，其间距超过某预先给的阈值ε（ε＞0），L′0=|X（t1）—X（t0）|＞ε，保留X（t1），并在X（t1）邻域另找一个点X1（t1），使得L1=|X（t1）—X1（t1）|＜ε，并且与之夹角尽可能地小，继续上述过程，直至x（t）到达时间序列的终点N，这时追踪演化过程总的迭代次数为tM—t0，则最大Lyapunov指数λ1为

其中，L′i=|X（ti）—X（ti—1）|，Li=|X（ti）—Xi（ti）|，Xi（ti）为ti时刻在状态X（ti）以ε为半径的邻域内的一点。

如果要计算次大的Lyapunov指数λ2，则要跟踪一个点以及邻近两个点构成的三角形，若这个三角形变得倾斜或其面积较大，需重新找到一个两边与原三角形两条边夹角最小的三角形，继续跟踪，直到终点，则次大的Lyapunov指数λ2为(www.xing528.com)

式中：Ai为三角形的面积。

由式（2.18）可以得到λ2，同理可求得λ3、λ4等。原则上可以求到最后λn，但由于实际资料的长度N 有限及噪声的影响等，只能较为可靠地估计最大Lyapunov指数λ1。

2.小数据量法

Rosenstein等[20]针对上述方法中的一些缺点提出了小数据量法，以下详细介绍小数据量法的原理与具体计算过程。设混沌时间序列x（1），x（2），…，x（t），…，嵌入维数m，时间延迟τ，则重构相空间

X（t）=﹛x（t），x（t+τ），…，x[t+（m—1）τ]﹜，t=1，2，…，M

其中M=N—（m—1）τ。在重构相空间后，寻找给定轨道上每个点的最近邻点，即

其中=1，2，…，M，且t≠，p为时间序列的平均周期，它可以通过能量光谱平均频率的倒数估计出来，X（t）为相空间中的状态点，dt（0）代表在初始时刻一对最近邻点之间的距离。

最大Lyapunov指数可以通过基本轨道上每个点最近邻近点的平均发散率估计出来。Sato等估计最大Lyapunov指数为

式中：Δt为样本周期；dt（i）为基本轨道上第t对最近邻近点对经过i个离散时间步长的距离。

后来，Sato等改进估计表达式为

式中：k为常数。最大Lyapunov指数的几何意义是量化初始闭轨道的指数发散和估计系统的总体混沌水平的量，所以，结合Sato等的估计式有：

将式（2.22）两边取对数得到

显然，最大Lyapunov指数可以近似看成式（2.23）这组直线的斜率。它可以通过最小二乘法逼近这组直线而得到，即

其中，〈·〉表示所有关于t的平均值。

小数据量法的计算步骤如下：

（1）对时间序列﹛x（i），i=1，2，…，N﹜进行快速傅立叶变换，计算出时间延迟τ和平均周期p。

（2）计算出关联维数D，再由m≥2D+1确定嵌入维数m。

（3）根据时间延迟τ和嵌入维数m 重构相空间﹛X（t），t=1，2，…，M﹜。

（4）找相空间中每个点X（t）的最近邻点，并限制短暂分离，即

其中，=1，2，…，M 且t≠。

（5）对相空间中每个点X（t），计算出该邻域点对的i个离散时间步后的距离dt（i）

（6）对每个i求出所有t的ln dt（i）平均x（i），即：

式中：q为非零dt（i）的数目，并用最小二乘法做出回归直线，该直线的斜率为最大Lyapunov指数。