（运筹学与控制论专业论文）时间序列非线性分析的若干研究.pdf

中文摘要 摘要：降水在很大程度上影响着人类的工农业等生产生活。由于很多非确定性因 素的影响，降水时间序列表现出了很大的非平稳性和复杂性。本文旨在利用除趋 势波动分析( d f a ) 方法研究降水时间序列的相关性性质，并根据混沌特征量和混 沌的判别方法揭示降水时间序列的混沌特性，从而为降水时间序列基于长相关和 基于混沌的预测提供一定的基础。 文中为降水时间序列加上了与原序列相关的常见趋势，利用d f a 方法得到了 一个独立于各种非周期趋势的交叉点以及一个可以揭示时间序列混沌特性的标度 性质变化过程。由d f a 方法对周期趋势的无效性我们引入了奇异值分解( s v d ) 方法，研究中我们利用自相关法求出时间延迟，得到了s v d 方法在去除周期、准周 期趋势方面的良好效果。在降水时间序列的混沌特性研究中，我们把d f a 与混沌 方法结合起来，发现把d f a 去除趋势利用在日降水序列后，所得到的时间序列具 有更好的混沌性质。这就为我们进行后继的预测分析提供了更多的信息。上述三 部分相互依托，基于以上内容本文的结构如下： 第一章，d f a 在降水时间序列长相关性上的应用。 第二章，基于s v d 方法的除趋势分析。 第三章，降水时间序列的混沌特性分析。 文巾图2 l 幅，表8 个，参考文献2 3 篇。 关键词：除趋势波动分析( d f a ) ；降水时间序列；趋势；混沌；奇异值分解( s v d ) 分类号： a bs t r a c t a b s t r a c t ：p r e c i p i t a t i o nh a sad r a m a t i ce f f e c to na g r i c u l t u r ea n di n d u s t r yt h a ti t t a k e sas i g n i f i c a n tr o l ei n h u m a n sa c t i v i t i e s b e c a u s eo ft h ee f f e c t so fm a n y n o n d e t e r m i n i s t i cf a c t o r st h e p r e c i p i t a t i o n s e r i e ss h o wg r e a t n o n s t a t i o n a r i t ya n d c o m p l e x i t y t h i sp a p e ra i m sa tu s i n gd f at oi n v e s t i g a t et h ec o r r e l a t i o np r o p e r t i e so f t h e p r e c i p i t a t i o ns e r i e sa n ds t u d y i n gt h ec h a o sp r o p e r t i e so ft h es e r i e sb a s e do nt h e c h a r a c t e ra n dt h ee s t i m a t em e t h o d so fc h a o ss oa st op r o v i d es o m eb a s i ct h e o r yf o r f u r t h e rs t u d yo f p r o d i c t i o no nl o n g r a n g ec o r r e l a t i o na n dc h a o s i nt h i s p a p e r , b a s e do nd f aw ea d ds o m ec o r r e l a t e dt r e n d st ot h eo g i n a l p r e c i p i t a t i o ns e r i e s t h a tac r o s s o v e rw h i c hi si n d e p e n d e n to fd i v e r s i f i e dn o n p e r i o d i c t r e n d si so b t a i n e da n dw ef i n da s c a l i n gp r o p e r t yv a r y i n gp r o c e s sw h i c hm a ys h o wt h e c h a o sp r o p e r t yo ft h es e r i e sa sw e l l s v dc o m e si n t oo u rp a p e rb e c a u s ed f ai s e f f e c t l e s st op e r i o d i ca n dq u a s i - p e r i o d i ct r e n d a n dt h e nw e s t u d yt h ec h a o sp r o p e r t yo f t h ep r e c i p i t a t i o ns e i r e sw i t hd f a t h e s et h r e e p a a sc o n t a c te a c ho t h e rw i t ht h e a r i t h m e t i co fd f a b a s e do nt h et o p i c s a b o v e ，t h i sp a p e rc o n s i s t so ft h r e ep a r t sa s b e l o w ： c h a p t e r i ：a p p l i c a t i o no fd f ao nt h es t u d yo fl o n g r a n g ec o r r e l a t i o no f p e c i p i t a t i o ns e i r e s c h a p t e r 2 ：e f f e c to fs v do nr e m o v i n gp e r i o d i ca n dq u a s i p e r i o d i ct r e n d sb a s e d o nr a n d o ms i g n a l s c h a p t e r 3 ：a n a l y s i so fc h a o so f t h ep r e c i p i t a t i o ns e r i e s t h e r ea r e21p i c t u r e s ，8t a b l e sa n d2 3r e f e r e n c e s k e y w o r d s ：d e t r e n d e df l u c t u a t i o n a n a l y s i s ( d f a ) ；p r e c i p i t a t i o ns e r i e s ；t r e n d ； c h a o s ；s i n g u l a r v a l u ed e c o m p o s i t i o n ( s v d ) c i 。a s s n o ： 学位论文版权使用授权书 本学位论文作者完全了解北京交通大学有关保留、使用学位论文的规定。特 授权北京交通大学可以将学位论文的全部或部分内容编入有关数据库进行检索， 并采用影印、缩印或扫描等复制手段保存、汇编以供查阅和借阅。同意学校向国 家有关部门或机构送交论文的复印件和磁盘。 ( 保密的学位论文在解密后适用本授权说明) 靴做储獬：李锻 签字日期：污年月乞日 ， 新虢触 签字日期：口i 净易月2 _ e t 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工作和取得的研 究成果，除了文中特别加以标注和致谢之处外，论文中不包含其他人已经发表或 撰写过的研究成果，也不包含为获得北京交通大学或其他教育机构的学位或证书 而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均已在论文中作 了明确的说明并表示了谢意。 学位论文作者虢荔钍签字嗍。孑年6 月乙日 致谢 本论文的工作是在我的导师商朋见教授的悉心指导下完成的，商朋见教授严 谨的治学态度和科学的工作方法给了我极大的帮助和影响。在此衷心感谢两年来 商老师对我的关心和指导。 商朋见教授悉心指导我们完成了实验室的科研工作，在学习上和生活上都给 予了我很大的关心和帮助，在此向商老师表示衷心的谢意。 商朋见教授对于我的科研工作和论文都提出了许多的宝贵意见，在此表示衷 心的感谢。 在实验室工作及撰写论文期间，董科强、林艾静等同学对我论文的研究工作 给予了热情帮助，在此向他们表达我的感激之情。 另外也感谢我的家人，他们的理解和支持使我能够在学校专心完成我的学 匕。 1d f a 在降水时间序列长相关性上的应用 1 o 引言 很多物理系统和生物系统都表现出长程幂率相关性( 长相关性) 的复杂性质。 对于非平稳的信号，幂谱以及相关性分析等传统方法都无法量化信号的长相关性。 1 9 9 4 年f l ：l p e n g 等在研究d n a 序列时提出来的除趋势波动分析方法( d f a ) 则是一 种计算非平稳时间序列长相关性的有效方法。d f a 方法是一种标度分析方法，利 用它可以计算出一个系数标度指数a ，由这个指数来表征信号的相关性特征。 d f a 优于其它方法之处不仅在于它能检验非平稳时间序列的长相关性，它还能避 免将时间序列的非平稳性误判断为长相关性。因此，d f a 方法一经提出就迅速被 应用于诸多领域，如d n a 序列、生理学、股票市场、云结构、经济时间序列、物理 学等。 降水对于人类活动和农业生产都起着至关重要的作用，降水的研究可以应用 到诸多方面，比如水文结构设计以及洪水预防等等。降水量的多少不但直接影响 水资源的开发和利用，而且与干旱、洪涝等自然灾害的发生、发展及其强度密切 相关，因此准确地预测降水量，既可为制定水资源开发利用策略提供科学依据， 又可为防灾、抗灾、救灾提供有效的指导。然而，由于降水是一种复杂的自然现 象，其形成、发生、发展的全过程均呈现出不确定性与随机性，难以定量描述因 此，建立科学、合理的预测预报模型的同时，迫切需要迸一步探讨并揭示复杂降 水现象的不确定性的行为机制，确定其有效的时序相空间，从而为确保有效、准 确地构造其动力学预测预报模型提供前提和基础。 降水的研究方法很多，像标准统计方法中的均值，方差，变差系数等等。人 们发展了很多研究降水的方法也建立了很多降水模型，对于降水序列也有很多关 于标度性质和多重分形的研究。但是降水序列的这些传统分析方法经常由于降水 信号的非线性性质而产生诸多谬误的结果。d f a 方法的出现正可以克服此类的问 题，本章我们就重点来看d f a 方法在降水时间序列上的应用阳1 。 本章主要检验日降水时间序列在长相关性方面的性质，我们利用d f a 方法研 究降水时间序列的相关性质在加上相关趋势之后的变化情况。我们知道外部趋势 是影响时间序列相关性质的主要凶素，人们都希望消除他们来对数据有更深层次 的理解。可是在很多问题中，趋势和数据内在的长相关性质并不容易准确分辨， 所以在这里我们给序列加上与序列相关连的趋势，也就是原序列的函数构成的趋 势，来检验加趋势后序列的相关性质发生了怎样的变化，以求更深层次得探知原 序列的内在标度规律和相关性质口q 们。 本章的结构如下：1 1 节简要介绍d f a 方法的计算步骤以及它与长相关性的关 系。1 2 节介绍降水数据的来源及对其的预处理。1 3 节重点介绍原始时间序列在加 上相关趋势前后经过d f a 处理所表现出的不同标度性质，我们把所得到的结果总 结在1 4 节。 1 1 基本方法 1 1 1 长相关性 对一平稳时间序列 石( f ) ) ， f = 1 ，2 ，) ，时间尺度为s 的自相关函数为 c ( j ) = 而缶n - a ( 一酏，一i ) ( 1 ) 式中平均值x 为序列的均值。如果对于s o 而言，c ( s ) = o ，则 x ( 讲是不相关的； 如果c ( d e x p ( s s ，) ，为衰减时间，则 x ( o 是短程相关的；如果c ( s ) 一s ， ( 0 y 1 ) ，即c ( s ) 呈幂律下降，则 x ( o 是长相关的。计算自相关指数y 经常受 到 x ( f ) ) 含有噪声或某些趋势成分的影响。 实验数据经常受到非平稳性的影响。为了找到时问序列波动函数的正确标度 性质，我们需要把非平稳的趋势与系统内在的波动性严格分开，但是做到这一点 也是不容易的。如果数据中存在趋势，贝j j h u r s t 分析和其他除趋势方法可能会由于 数据的非平稳性而给出错误的结果。而且通常实际数据中的内在趋势是未知的， 甚至连内在趋势的尺度都是未知的。但是d f a 则能很好的克服这些困难，它不必知 道趋势的形态和起因就能很稳定的确定噪声数据的标度性质。 1 1 2d f a 计算步骤 2 对序列 x ( m ，我们将d f a 简要概述为下面四个步骤： ( i ) 计算时间序列 工( m 的累积离差y ( o ， 】，( d = ( - 7 0 ( 2 ) x 为序列的均值。把】，( 力等分成m 个不重叠的等时间长度s 的区间，其中m = n s 】 ( 即取整数) 。由于序列长度并不总是时间长度5 的倍数，因此有小部分序列后面的 数据信息未能被利用。因此，对该序列的逆序进行同样的操作，共有2 f 个等长 度的区间。 ( i i ) 对于每个区间 ，用最d , - 乘法拟合数据。得到局部趋势。滤去该趋势后的时 间序列记为j ：( f ) ，表示原序列与拟合值之差，即： r ( 力= y ( o - p , ( o( 3 ) 式中a ( ，) 为第v 区间的拟合多项式。如果拟合的多项式采用的是线性的、二次的、 三次的，甚至是更高阶的多项式，则分别记为d f a ( d f a i ) ，d f a 2 ，d f a 3 ，等。显然 ，l 阶的d f a 滤去了累积离差中的力阶趋势成分或原始序列中的咖j 阶趋势成分。 ( 谳) 计算每个区间滤去趋势后的方差，得： 盹s ) = y 主。y , 2 ( v - 1 ) s + i l ( 吼2 。m ) ( 4 ) 以及f ( v j ) = 三喜j ：2 【一( v m ) n 妇( 萨m “，m + 2 ， 2 m ) ( 5 ) ( i v ) 对所有等长度区间求均值并开方，计算标准d f a 波动函数( f l u c t u a t i o n f u n c t i o n ) ，得： f ( s ) = 如果波动函数f ( j ) 和时间尺度s 之间存在幂率关系即： ( 6 ) f ( s ) s 4( 7 ) 则表示存在标度性质。a 为标度指数或称为波动指数，它代表了数据的相关性 质。对于相关指数y ，有一个关于波动函数的近似估计： f ( s ) 5 1 一班( 8 ) 由上面两个式子可以得到标度指数和相关指数的关系： 口= l 一7 2( 9 ) 其中o y 1 。因此我们可以通过计算波动指数口来确定相关指数y 。所以当口 = 0 5 时，不存在相关性。数据是不相关的信号( 白噪声) 。当a o 5 时，数据是长相关的。 1 2 原始数据介绍 本文研究的降水数据取自北京的2 1 个气象站从1 9 7 3 年1 月1 日至2 0 0 4 年1 2 月3 1 日，1 1 6 8 8 天的气象纪录。我国地处大陆性季风气候区，四季降水差异很大， 夏季有些天降水十分丰富，冬季有些天降水又十分稀少。在降水稀少的天气里， 由于降水量很小，鉴于技术上或是其他的原因气象站没有详细记载准确数据，只 用一个统一的标准符号来表示，当然还有很多天是没有降水的。所以为了研究上 的方便我们将降水很少的日子的降水量统一为0 5 ( 0 1 m m ) ，没有降水的日子降水 量自然为0 。我们将2 1 个气象站的降水纪录按测量误差情况等因素取加权平均值 来作为我们进行研究的原始数据。 图1 1 原始时间序列图。数据单位为0 1 m m 从图中的数据变化情况可以清楚地看到数据的周期性趋势，而这一点从下面 的d f a 图中也可以看到。 1 3 数据分析 1 3 1 原始数据的d f a 4 首先，我们来看一下原始时间序列的相关性性质，为了得到关于时间序列更 多的信息，我们使用从1 到5 阶的除趋势波动分析。基于序列长度和研究的需要， 我们主要讨论的范围是时间尺度s o 5 对 应于长相关性。而第二个交叉点两边的区间段内，标度指数从耐帕 l 转变到 旗町 0 5 ，这意味着在大的时间尺度s 上序列是反相关的。第二个交叉点前后的标 度性质的明显变化让我们想到了正弦序列的情形。图1 2 ( b ) 中即是正弦序列的d f a 图。比较图1 2 ( a ) 和( b ) 我们可以看到他们都存在一个非常明显的交叉点。交叉点前 后标度指数发生了巨大的改变，而且标度指数都变得很小，而正弦序列的交叉点 即是其周期。由此我们可以推测原始序列的趋势可能是由周期趋势和其他趋势叠 加而成。这些发现和推测有助于我们开发新的方法去重新审视和研究降水纪录。 1 3 2 加趋势a x ( i ) y 序列的d f a 拥有两个交叉点s 譬，蟛说明原始序列有着复杂的标度性质，这就给我们对 序列的理解和预测带来了困难。所以我们换一个角度，研究序列加上相关趋势后 相关性质的改变情况，以希望通过d f a 方法发现原序列更多的相关性性质。这里所 加的都是一些简单而常见的趋势：比如最简单的是幂率函数铆z ( 明，a 是系数，p 代表次幂。也就是说我们对新序列y ( d = 工( f ) + 彳【z ( 明p 进行除趋势波动分析。我们 把4 和p 作为变量来展开研究。 先蜘为常量l ，p 为变量。研究次幂p 对r ( i ) 的d f a 的影响。 1 3 2 1p 为正整数 因为力阶d f a 能除去 1 阶的趋势，我们令p 从l 到6 来大体上看一下d f a 对序列 j ，( f ) 的影响。 6 图1 3 加幂相关函数趋势铆x ( 明p 序列：j r ( f ) = x ( f ) + 彳【x ( 明p 的d f a ，其中爿；l ，( a 驴= 1 ( b ) p = 2 ( c ) p 。3 ( d ) p 叫( e ) p 一5 ( 0p = 6 。 图1 3 中，在尺度s 1 0 3 内( a ) ，( b ) 的每阶d f a 都有两个交叉点，而( c ) e e 只有 d f a 3 ，d f a 4 和d f a 5 有两个交叉点，而且位置也变得非常靠近。( d ) ，( e ) 和( f ) 中，虽 然仍然存在交叉点，但我们几乎察觉不到它们，仔细观察会发现箭头所指处确实 有交叉点存在。虽然( d ) ，( e ) 和( f ) 中每阶d f a 都只有一个交叉点，但为了讨论的方便， 我们仍然用s 婴表示它们。图中我们用箭头指出了一些有代表性的交叉点。通过观 察箭头我们发现s 卵的位置在图1 3 的6 个子图中几乎是一致的。另外从总体上看随 着次幂p 的增加，双对数图的整体偏差情况也越来越明显。我们利用上一节中的方 法计算出这些交叉点和交叉点两边标度指数的值列于表格1 2 0 7 。 表格1 2 加趋势4 x ( 明p 序列的标度指数口枷) 和交叉点书，a = i ，p 从1 取到6 。 7 d f a ld f a 2d f a 3d f a 4d f a 5 町0 6 7 1 8 0 6 2 6 30 6 1 2 20 6 3 1 80 6 4 1 5 8 4 鸲1 4 62 3 42 7 4 p = l a f 。) 1 1 4 0 7 1 2 7 7 21 3 6 5 01 4 8 5 71 4 9 2 7 艘 4 0 65 5 67 0 48 9 29 6 5 钟o 1 5 1 2 o 1 4 1 9o 1 2 3 60 1 2 6 3o 1 2 9 6 靠0 6 1 5 3 0 5 9 2 70 5 8 5 60 5 9 6 60 6 0 7 3 1 0 6 1 1 6 1 4 62 3 42 7 4 p = 2 西o 0 8 7 3 20 9 4 5 20 9 8 9 01 0 6 4 31 0 5 9 3 螂 4 0 6 5 5 6 7 0 4 8 9 29 6 5 芷- 0 2 6 8 0 0 2 2 5 6o 1 6 8 70 1 5 8 3o 1 5 4 4 d 4 0 6 0 3 1 0 6 0 2 10 6 0 1 50 6 0 2 60 6 0 2 4 p = 3 艘 4 0 65 5 67 0 4 8 9 2 9 6 5 芷一) 0 3 3 6 1 0 2 7 7 20 2 2 2 40 2 1 5 10 2 1 9 7 掣o 5 4 9 9 0 5 5 7 00 5 5 6 50 5 6 0 00 5 6 3 3 p = 4 掣4 0 6 5 5 67 0 48 9 29 6 5 口9 0 3 7 0 90 3 1 9 8o 3 0 1 00 3 0 0 80 3 0 7 3 一。 0 5 2 7 1 0 5 3 6 70 5 3 6 50 5 4 0 70 5 4 5 4 p = 5 璎4 0 6 5 5 67 0 48 9 29 6 5 斑一 0 3 9 5 1 0 3 5 3 80 3 6 3 80 3 6 7 20 3 7 2 1 耐”0 5 1 7 70 5 2 7 60 , 5 2 7 90 5 3 2 20 5 3 7 7 p = 6 彰 4 0 6 5 5 67 0 48 9 29 6 5 硝0 4 1 3 8 o 3 7 6 00 4 0 1 20 4 0 6 00 4 1 0 6 图1 3 ( c ) 中d f a 3 ，d f a 4 ，d f a 5 _ 虽有两个交叉点，但是交叉点铅两边的标度性 质是相同的，所以我们在s 磐前拟合曲线，计算得到一个标度指数研哪即表1 2 中旷3 栏。从表1 2 我们很清晰的看到每个不同次劫的交叉点s 婴都是相同的，而且与原 始序列的交叉点也是一样的。这就说明了交叉点蠼是独立于相关趋势4 工( 明户的 幂p 。再者我们看到随着p 的增加交叉点尤前的标度指数磁m 和耐帕有着向o 5 衰减 的规律。而? 后面的标度指数脚则有向o 5 增加的趋势。为了更清楚地看清这个 现象。对：萨1 种一2 ，在交叉点嫂前我们计算一个复合的标度指数耐，把重新计算 出的科帕和趟帕列于表格1 3 中。 8 表格1 3 重新计算矿l 和旷2 的标度指数西”后，耐“和吐”的比较。 d f a ld f a 2d f a 3d f a 4d f a 5 p = j0 8 7 7 50 8 3 0 00 8 2 3 80 8 1 0 60 7 8 8 5 p = 2o 7 l1 30 6 9 2 40 6 9 0 70 6 8 6 20 6 7 7 9 硝帕 p = 3o 6 0 3 l0 6 0 2 l0 6 0 1 50 6 0 2 6 0 6 0 2 4 p = 40 5 4 9 90 5 5 7 00 5 5 6 50 5 6 0 00 5 6 3 3 p = 50 5 2 7 10 5 3 6 70 5 3 6 50 5 4 0 70 5 4 5 4 p = 6o 5l7 70 5 2 7 60 5 2 7 90 5 3 2 20 5 3 7 7 p = jo 1 5 1 20 1 4 1 90 1 2 3 60 1 2 6 30 1 2 9 6 p = 2 0 2 6 8 00 2 2 5 6 0 1 6 8 70 1 5 8 30 1 5 4 4 趟砷 p = 30 3 3 6 10 2 7 7 20 2 2 2 40 2 1 5 10 2 1 9 7 p = 4o 3 7 0 90 3 1 9 80 3 0 1 0 0 3 0 0 8 0 3 0 7 3 p = 50 3 9 5 10 3 5 3 80 3 6 3 80 3 6 7 20 3 7 2 1 p = 60 4 1 3 80 3 7 6 00 4 0 1 20 4 0 6 00 4 1 0 6 表格1 3 非常清楚地显现出了随着p 的增加，标度指数耐嚣逐渐减小到o 5 的趋 势。而”则显然是有逐渐增大到o 5 的趋势。其中畦”增加到o 5 的速度要比耐帕衰 减的速度慢一些。我们再从另一个更加简单的模式表现出这种趋势：我们计算在 全时间尺度j 上的总体标度指数口( ”并列于表格1 4 中。 表格1 4 每阶d f a 在全尺度区间j 上的标度指数。 ”)d f a ld r 气2d f a 3d f a 4 d f a 5 p = j0 6 9 5 8 0 7 4 5 80 7 7 1 60 7 8 9 50 8 0 8 6 p = 20 5 9 0 90 6 2 4 30 6 4 3 20 6 6 0 30 6 8 1 3 p = 30 5 2 3 30 5 4 7 30 5 6 2 3 0 5 7 9 80 6 0 11 p = 40 4 9 5 l 0 5 1 5 00 5 2 8 70 5 4 5 20 5 6 4 4 p = 5o 4 8 6 90 5 0 4 3 0 ，5 1 7 40 5 3 2 70 5 4 9 8 p - - 60 4 8 5 6o 5 0 1 2 0 5 1 3 90 5 2 8 6o 5 4 4 4 从一般意义上讲2 阶d f a i 鬟j 效果最好，而表1 4 中随都的增j 3 1 ：i d f a 2 也是最逼近 0 5 的。我们知道o y 1 。 图1 3 中交叉点? 独立于次幂p 但并不这样。比较图1 3 ( a ) 和( b ) ，当旷l 时，( a ) 中础 1 0 0 。图1 3 中( b ) ，( c ) 和( d ) 显示了交叉点逐 渐“消失”的过程，虽然在这个交叉点前后的标度性质相同，但随都的增加它的消 失也是一个值得注意，值得研究的现象，因为这恰恰可能就是d f a 体现出其作用的 地方。我们取p 为小数来研究这个过程。 1 d 1 矿 图1 4 加幂相关函数趋势铆x ( 例p 的序列：】，( ，) = x ( ，) + 川z ( 例9i 拘d f a ，其中 爿= 1 ，( a 庐1 8 ( b ) 矿1 9 ( c ) 旷2 8 ( d ) p - 2 9 ( e ) 旷3 8 ( d 矿3 9 。c a ) - - - ) ( b ) 是交y a 础由s 4 ( 1 0 0 ) 的过程。( c ) - - 9 ( d ) 是交叉点s 0 和的“消失过程。( e ) _ ( f ) 则是交叉点 s ( ，n ( n = 3 ，4 ，5 ) n 消失”过程。一些典型的交叉点由箭头予以标出。 图l ，4 中，可以看到这三组图都是在p = 朋+ d 8 到删+ n 9 之间( m = l ，2 ，3 ) ，非常接近于 p = m + l 。这个现象可以归结于在1 1 节中介绍的n 阶d f a 可以有效的消除r 1 阶的趋 势的结论。当我们使用5 阶d f a 时，整个韶的消失过程刚好在p 为3 9 即接近于5 1 = 4 时。因此我们可以做出如下的假设：如果我们使用d f a 到时1 阶，则在每个旷i + 0 9 附近，j = 1 ，2 珈l ，即会有一个重要的性质发生改变，而最后p 达n n - 0 1 时，即接近 于y 时改变态达到稳定态，此时p 再增加也不会再有相关性质的改变，除非继续提 高d f a 去趋势的阶数刀。 1 3 2 3p _ 2 d + l ( d - 黾混沌吸引子的维数) ，重 构相空间与实际系统的动力学特征一致。所以只要我们得到了蹴可以取到适当的 嵌入维数历了。计算关联维数硪们采用g r a s s b e r g e r 和p r o c a c c i a ( 1 9 8 3 ) 首次基于相空 间重构所发展出来的g p 算法： s t e p l 取上面所计算出来的时间延迟f ，并从嵌入维m = 2 开始对时问序列 工( f ) ) o = 1 ，2 ，丁) 重构相空间： y ( i ) = 【z ( f ) ，x ( f + 力，工( f + 2 f ) 9o o ，x ( f + ( ，”一1 ) f ) 】，中f = 1 ，2 n ，n = t - ( m - 1 ) z s t e p 2 计算重构相空间y ( j ) 的关联积分( c o r r e l a t i o ni n t e g r a l ) i 垂l 数 c ( 咖丽2 面。；丕何( ，一旷帅 ( 1 3 ) 其中：h 是h e a v i s i d e 单位函数。慷一0 表示相空中任两点r ，艺之间的距离。 b r o c k ( 1 9 8 6 ) 证明了关联维数不依赖于距离定义的形式不妨以矢量的最大分量差 作为距离 q 2 忙一- - m a x ，。l 珞一i ( 1 4 ) s t e p 3 恰当地选取，的值，关联维数鸦关联积分函数c 纠满足关系c ( r ) 一， 取对数得到：d ( m ) = i n c ( r ) i n ( r ) 。从而拟合求出对应m 的关联维数估计值。 s t e p 4 不断增加嵌入维数m ，重复以上步骤2 ，3 ，直到相应关联维数的估计值在 一定范围内不变为止，即线性部份的斜率不再随m 的增大而改变，此时维数达饱和 嵌入维数，相应的斜率即为吸引子关联维数d 。 由此我们取一适当的肫三2 舟j 进行后续的计算。 2 1 3 算法 给出t s v d 的基本定义以及时间延迟和嵌入维数的确定方法后我们给出本章 1 9 中的研究算法： s t e p i 给原幂率噪声序列 以) 加周期和准周期趋势 瓦) ，得到序列 l ) 。 s t e p 2 对 艺 进行相空间重构，嵌入维数m 和时间延迟彳按照上面2 1 2 中方 法确定a 我们有苁= ( k ，乓，+ ( m - 。) ，) ，l 七n - ( m - 0 r 。从而得到矩阵 a ：a = 兀 ： ，m = 刀一( 所一l 弦 s t e p 3 对彳进行奇异值分解，得到a = u x z x v r 。我们把中数值较大的特征 值呈i ；o o ，这些较大的特征值正是表征匕中给以所加趋势的主频。把后面的特征 值依次前移，得到矩阵z 。 s t e p 4 滤波矩阵彳= u z v 7 ，得到戎= ( ，砭，砭( + 。，) ， 1 k ，加趋势序列( 争) 以及s v d 去除趋势后得到 的回溯序列( 一) 的波动函数与时间尺度的双对数图。振幅a = l 。a ) ，b ) 和c ) 对应 标度指数喁= 0 5 5 的噪声序列，所加周期趋势的周期分别为t ：= 2 3 ，2 5 , 2 7 。d ) ，e ) 和 d 对应标度指数- - 0 7 5 的噪声序列，所加周期趋势的周期分别为t 2 3 ，2 5 , 2 7 。 表格2 1 相关参数列表。q 吒玛表示由s v d 得到的特征矩阵的前三个特征值。彳为时 间延迟，硝关联维数用以选定嵌入维数m 。 = 0 5 5 c l = 0 7 5 吒c r 2吒 f d q吒吒 fd s i n ( 2 刀i 2 3 ) 1 0 7 71 0 7 50 o o l 22 一一 1 0 7 7 i 0 7 50 0 0 1 4 2 一 1 0 0 09 7 3 l0 0 0 l l51 1 71 0 0 0 9 7 3 10 o o l l51 1 6 s i n ( 2 万i 2 ，) 3 7 8 32 9 9 90 0 0 0 51 91 1 53 7 8 32 9 9 9 0 o 0 0 6 1 9 1 0 8 s i n ( 2 万i 2 ，) 从图2 1 可以看出加趋势序列波动函数与时间尺度的双对数图的交叉点位置即 大概为周期，说明周期趋势给原噪声序列带来了谬误的交叉点。而从a ) ，b ) 和c ) 的 原始噪声序列与s v d 去除趋势回溯序列的双对数图来看，随着周期确增加，s v d 2 1 方法去趋势的效果变差。但是从a ) 与d ) ，b ) 与e ) ，c ) 与d 的比较来看，随着标度指数 的增加，s v d 的去趋势效果明显变强。所以s v d 对于具有较高频率的周期趋势和 较高标度指数的噪声序列有较好的除趋势效果。 接下来我们固定周期兀改变系麴的值，研究s v d 对趋势振幅的敏感性。取 产2 5 。锄分别为0 0 1 ，o 1 ，1 0 ，1 0 0 。对两组噪声数据做图如下： 图2 2 原始噪声序列( - 4 - ) ，加趋势序列( - - o - ) 以及s v d 去除趋势后得到的回溯序列( 一) 的波动函数与时间尺度的双对数图，周期产2 5 。a ) ，b ) ，c ) 和d ) 对应标度指数强= 0 5 5 的噪声 序列。所加周期趋势的振幅分别为爿= 0 0 1 ，0 1 ，1 0 ，1 0 0 。e ) d ，g ) 和h ) 对应标度指数a 2 = 0 7 5 的噪声序列，所加周期趋势的振幅分别为一= 0 0 1 ，o 1 ，1 0 ，1 0 0 。 由图2 2 比较a ) ，b ) ，c ) 和d ) 以及比较e ) ，f ) ，g ) 和h ) 我们得知，虽然加趋势序列的波动 函数随着振幅4 的变化而变化，但是s v d 去除趋势的效果保持不变。所以我们得到 s v d 方法对于周期趋势的振幅不敏感。或者说s v d 方法去除周期趋势的效果独立 于趋势的振幅。 2 2 ，1 2 加上个周期趋势 上面讨论了加一个周期趋势的情况，但是我们知道实际数据中可能存在着两 个以至更多的趋势成分，s v d 方法是否同样可以有效的去除他们呢? 给原始噪声 序列加埘周期趋势得到新序列：z = ， 五+ 4s i n ( 2 研互) + 4s i n ( 2 石i t 2 ) + + 以s i n ( 2 历瓦) = 五+ a js i n ( 2 7 r i t j ) j = l 为便于与2 2 1 1 中结果进行比较和讨论我们取三个周期分别为：正- 2 3 ，t = 2 5 ， ， 己= 2 7 。将a js i n ( 2 万i t j ) 分别取为s i n ( 2 n i t 1 ) + s i n ( 2 n t 2 ) ， = l s i n ( 2 x i t 2 ) + s i n ( 2 n i t 3 ) 以及s i n ( 2 m t i ) + s i n ( 2 x i t 2 ) + s i n ( 2 x i t 3 ) ，对两组 噪声序列分别做出原始噪声序列，加趋势序列以及s v d 方法去除趋势后得到的回 溯序列的波动函数与时间尺度的双对数图如下图2 3 ，并与图2 1 比较。 图2 3 原始噪声序列( + ) ，加趋势序列( 毛卜) 以及s v d 去除趋势后得到的回溯序列 ( - 一) 的波动函数与时间尺度的双对数图。a ) ，b ) 和c ) 对应标度指数珥= o 5 5 的噪声序列， 所加周期趋势的周期分别为产s i n ( 2 删五) + s i n ( 2 ，r i t 2 ) ，s i n ( 2 1 r i t 2 ) + s i n ( 2 n i t 3 ) 以及 s i n ( 2 n i t 1 ) + s i n ( 2 z r i t 2 ) + s i n ( 2 ，r i t 3 ) 。d ) ，e ) 和d 对应标度指数鸭= 0 7 5 的噪声序列，所 加周期趋势的周期同上。 表格2 2 相关参数列表。仉( f _ 1 ，2 ，7 ) 表示i 扫s v d 得到的特征矩阵的特征 值。f 为时间延迟，硝关联维数用以选定嵌入维数m 。 q = 0 5 5 o 、o l6 j6 o ，o b0 1 td t r e n d l1 11 71 0 7 61 0 5 51 0 5 40 0 0 1 20 0 0 1 20 0 0 l l21 7 5 t r e n d 29 7 6 19 1 6 69 1 5 48 8 ，1 80 。0 0 0 90 0 0 0 90 0 0 0 972 2 l t r e n d 3l1 4 31 0 8 21 0 5 11 0 3 49 9 5 89 6 5 40 0 0 1 l33 3 3 呸= 0 7 5 吒 吒q吼c r 5 吗 fd 1 r e n d i1 1 1 71 0 7 61 0 5 51 0 5 4