高劍波

(北京師范大學(xué)地理科學(xué)學(xué)部地理數(shù)據(jù)與應(yīng)用分析中心,北京 100875)

引言

20 世紀(jì)科學(xué)巨匠馮·諾伊曼臨謝世時,對死亡后的人不再有思維,感到迷惑和恐懼.自2021 年8 月25 日吾師鄭哲敏先生逝世至今,我有時也有這樣的感覺,因不能與先生電話或面談而得到先生對某件事的看法.但有時又覺得好像不全如此,因為先生的思維方式至少部分傳下來了.在這篇紀(jì)念先生的文章里,我想在談學(xué)術(shù)研究前,先提一下先生對我個人的影響.這些影響基本決定了我在非線性領(lǐng)域和其他領(lǐng)域研究問題的偏好.本文主要介紹先生和我在非線性領(lǐng)域的一些工作,及后來這些工作如何被拓展并演變成一些通用的非線性問題的分析方法.這些方法已被廣泛用于自然科學(xué)、工程技術(shù)和社會科學(xué)的諸多領(lǐng)域.它們特別適用于各領(lǐng)域(包括運維)的故障診斷,生物醫(yī)學(xué)數(shù)據(jù)的分析,及不確定性的度量.大家要是覺得有些方法對研究有用,請隨時聯(lián)系我,索要相關(guān)程序和使用說明.

1 教誨

1988 年,我從浙江大學(xué)工業(yè)自動化專業(yè)畢業(yè).在浙大時,跟數(shù)學(xué)老師周先意先生,學(xué)了不少數(shù)學(xué).周老師是陳建功和蘇步青先生的學(xué)生,是數(shù)學(xué)家王元的同班同學(xué).周老師在1958 年被發(fā)配到江西農(nóng)村.空余時間都在思考教數(shù)學(xué)時如何把分析、代數(shù)和幾何融合起來.1978 年蘇步青先生將其調(diào)回浙大.由周老師教我們數(shù)學(xué),是我們的幸運.周老師對學(xué)生極好,跟我們聊天時經(jīng)常提到4 個他最得意的學(xué)生,曹銀和博士便是其中之一,其時他已在北京力學(xué)所,師從白以龍院士,跟鄭先生也熟.來京前,我沒見過曹師兄,但經(jīng)常與其通信交流.

1987 年大三時,在浙大電機系,我聽了一個外校教授關(guān)于混沌動力學(xué)的報告,覺得這個領(lǐng)域太迷人了,因此跟曹師兄提及,并表示以后很想做這方面的研究.其時力學(xué)所正籌劃于次年6 月建立非線性連續(xù)介質(zhì)力學(xué)開放實驗室(LNM),并由鄭先生掛帥,所以曹師兄問我是否愿意去力學(xué)所跟鄭先生讀研.我自然很激動.之后曹師兄又帶來一個好消息,說鄭先生會幫我辦理保送手續(xù),這樣就不用擔(dān)心因沒學(xué)過力學(xué)而考不上(不像現(xiàn)在考研,當(dāng)時只考專業(yè)課),這些溝通發(fā)生在1987 年10 月.至11 月下旬,出了些問題,曹師兄跟我說,鄭先生因太忙,忘了去給我申請免試讀研,而申請的截止期是11 月1 號.所以師兄問我,是否等一年再說.我想了想,說還是考吧.考試在2 月初.我花了2 個月,學(xué)了4 門力學(xué)課,考上了力學(xué)所.

我在力學(xué)所待了6 年.第一年在玉泉路的研究生院,與先生幾乎沒有接觸.后5 年,幾乎隨時可以找先生聊.在LNM 剛開始時,真正做非線性的年輕人就兩三位,包括何國威院士、趙力師兄和我.可惜,我因資質(zhì)有限,當(dāng)時學(xué)得很有限.下面我講五個方面.

第一點與情懷有關(guān).在追思先生時,很多人提到先生說過的一句話,“不要做錦上添花的事,要做雪中送炭的事”.這個情懷的精髓是堅持不懈地努力,做出最好的工作.這包含兩個方面.一是對自己的要求.聽說先生90 歲高齡時還自己推導(dǎo)公式.二是幫助學(xué)生和年輕人,做出最好的工作.記得2007 年,我寄給先生我的英文著作[1]后,先生很喜歡,隨即建議我寄一本給馮元禎先生.馮先生是生物醫(yī)學(xué)工程奠基人之一,生物力學(xué)之父,是鄭先生最好的朋友.馮先生當(dāng)時已是88 歲高齡.開會時經(jīng)常碰到加州大學(xué)的教授說,馮先生如此高齡,每年還會在PNAS 上發(fā)表高質(zhì)量的文章,太令人敬佩.在書寄給馮先生后,近一個月,沒收到任何回音,我忐忑不安,心里默念是不是自己書寫得不好或馮先生的專業(yè)與我相差較多.此時,馮先生的親筆書信到了,信中特別強調(diào),他在每一頁都學(xué)到了一些新東西.我特別激動,告知了鄭先生信的內(nèi)容,先生也很高興.馮先生的鼓勵對我影響很大.馮先生在信中還建議我看看他的好友,美國工程院院士黃鍔(Norden Huang)關(guān)于經(jīng)驗?zāi)B(tài)分解(empirical mode decomposition,EMD)的工作[2].馮先生太有洞見了！他看出了我們數(shù)據(jù)分析的體系還缺類似于EMD 的一塊.所幸,當(dāng)時我們已經(jīng)發(fā)展了一個新方法.這在第4 節(jié)會講到.總的來說,馮先生的提醒,很大地促進了我們發(fā)展自適應(yīng)去趨勢、去噪、多尺度分解和分形分析.

另一件事,好像是2010 年初秋,先生還在美國探訪兒子.我回國探親,回美時,路過北京.遠在美國的先生安排我去拜訪林家翹先生.因郭永懷先生和李佩先生是林先生的好友,鄭先生托李佩先生帶我去清華拜訪林先生.當(dāng)時林先生已是95 歲的高齡了.到了林先生的別墅,聊了一個多小時,只聊一個主題,如何理解蛋白質(zhì)結(jié)構(gòu).林先生最出名的工作是流體力學(xué)的流動穩(wěn)定性和星系螺旋結(jié)構(gòu)的密度波理論,但關(guān)于此卻一字未提.不囿于自己熟悉的領(lǐng)域,隨著時代發(fā)展不斷拓展新的領(lǐng)域,這就是大師們的情懷.身處百年未遇之大變局時期,我們必須發(fā)揚光大大師們的這種求真務(wù)實和開拓創(chuàng)新的精神.

第二點與心態(tài)有關(guān).2002 年8 月,我在佛羅里達大學(xué)當(dāng)助理教授,其時離開力學(xué)所和鄭先生已有8 年,自己首次帶研究生,我想,怎么教能好些呢？腦袋里突然冒出來一句先生跟我講過的話:100%的相信,同時又要120%的懷疑.先生的高明之處在于,不用另外解釋,我自然覺得這句話很有道理.在我指導(dǎo)學(xué)生時,還展開解釋了,所以落了下乘:前半句,100%的相信,表示對所選的領(lǐng)域和問題的重要性和前瞻性沒有任何懷疑,因此會全身心、充滿激情地投入到學(xué)習(xí)和研究中,只有這樣,才能將天賦真正變成才能,并發(fā)揮出來.后半句,表示在所選的領(lǐng)域中,一些最基本的東西都可能有錯.說的俗一些,就是不要被權(quán)威所束縛.其實,樸實的心態(tài)應(yīng)該是,充分欣賞領(lǐng)域里那些漂亮的工作,但不該有誰是權(quán)威、能否超越等想法,而應(yīng)該努力找到所選問題的最佳解,并成為最好的專家之一.雖然我不如先生高明,但這個教學(xué)方法還是有效的,如我在佛羅里達大學(xué)的博士生胡靜,深得此法的精髓,對后面第3，4 節(jié)描述的方法,有很大的貢獻.

第三點與思考有關(guān).先生與我聊天時,強調(diào)最多的是,要自由思考.剛到力學(xué)所時,有朋友勸我一定要多讀鄭先生出名的那些文章.但因我缺乏基礎(chǔ),同時當(dāng)時看不到那些工作與非線性科學(xué)的關(guān)系,且鄭先生也沒有要求我讀,所以在力學(xué)所的幾年,基本沒讀過鄭先生的文章.2014 年,在慶祝鄭先生90 歲生日的會議上,與我很熟的一位老教授,談慶明老師,跟我說,當(dāng)年我進力學(xué)所時,鄭先生特意囑咐他們,不要限制小高,讓他完全自由發(fā)展.談老師特別強調(diào),鄭先生這樣關(guān)照的,我是絕無僅有的一個.聽談老師這么說,我很感動,也很慚愧.我想主要是我沒力學(xué)背景,所以鄭先生特別關(guān)照我.當(dāng)然,由此也可看出鄭先生對非線性科學(xué)的重視.可惜,1994 年出國后,有近10 年時間,我只斷斷續(xù)續(xù)地作為業(yè)余愛好者繼續(xù)做些非線性問題的研究,所以成績很有限.

第四點是寬容.事情得從1989 年秋天我在LNM 做的報告說起.鄭先生請了很多人來聽,包括北京大學(xué)的朱照宣先生.那是我做的第一個報告,講得自然很不好,很多人聽后有些不滿意.不過,會上鄭先生沒有批評我,會后更沒約束我的思維.

后來我沒有轉(zhuǎn)行做量子計算,理由是,大尺度、宏觀的量子糾纏態(tài)幾十年內(nèi)是沒法實現(xiàn)的,因此不可能發(fā)展出能用于通用計算的量子計算機.這個推理不一定嚴(yán)謹(jǐn),但決定是正確的,21 世紀(jì)量子計算、量子通信需要的是實驗方面的俊才,但于此我不擅長.另外,1993 年,MIT 的Cuomo 和 Oppenheim[3-4]用Pecora 和 Carroll[5]發(fā)展的混沌同步開創(chuàng)了用混沌做安全通訊的新領(lǐng)域.在1998 年前后,美國海軍研究部的Pecora 博士組織了一個多學(xué)科大學(xué)研究計劃(Multidisciplinary University Research Initiative,MURI,類似 DARPA 的項目) 以推動混沌安全通訊的軍用和民用的研究.參加項目的有加州大學(xué)圣地亞哥分校(UCSD)、加州大學(xué)洛杉磯分校(UCLA)和斯坦福大學(xué).其時,我在UCLA 電子工程系跟隨Rubin 教授讀博.系里的兩位老師,Jia-Ming Liu 教授和Kung Yao 教授,參與了MURI 項目.姚教授做通訊,劉教授用半導(dǎo)體激光實現(xiàn)通訊.我對此很感興趣,參與了劉教授團隊的工作,因此與幾位同學(xué),Steven Huang、Howard Chen、Shuo Tang 和 C.C.Chen 成了好朋友.前三位是劉教授的學(xué)生,最后一位是姚教授的學(xué)生.除了Shuo Tang 去英屬哥倫比亞大學(xué)當(dāng)教授,其他三位都去了臺灣當(dāng)教授.與他們討論時,我告知C.C.Chen,混沌通訊的功耗應(yīng)該比Cuomo等[5]所說的大,因Cuomo 等在模擬隨機微分方程時,把隨機項處理錯了(成正比,但Cuomo 等用dt代替了 dt—當(dāng) dt?1 時,?dt,相當(dāng)于將實際噪聲縮小了幾個數(shù)量級,因而得出功耗小的結(jié)論).C.C.Chen 和姚教授發(fā)現(xiàn)我是對的,因而證明,混沌通訊也許能用于軍用,但在民用領(lǐng)域,價值很有限[6-7].之后,這方面研究的熱度在美國大降,但在香港和歐洲仍持續(xù)了很久.

第五點一定要重視科學(xué)實驗和實驗數(shù)據(jù)分析,盡量從重大的現(xiàn)實問題中,抽象出具有根本科學(xué)意義的問題,然后同時解決現(xiàn)實和科學(xué)問題.我想這是從普朗特傳至馮·卡門、馮·卡門傳至錢學(xué)森、再由錢學(xué)森傳至鄭哲敏的技術(shù)科學(xué)的精髓.

2 混沌時間序列分析

為了不讓我成為無根的浮萍,談慶明老師建議鄭先生讓我跟一個研究圓柱尾流混沌現(xiàn)象的實驗團隊一起做研究.這個項目的立項,白老師和曹師兄起了很大作用.團隊由林貞彬和鄂學(xué)全兩位老師領(lǐng)銜,成員中有位動手能力極強的年輕人,李東輝.團隊研究的目的是通過測量和分析圓柱尾流固定點速度分量的時間序列來厘清圓柱尾流的初始轉(zhuǎn)捩區(qū)域是否存在混沌現(xiàn)象,若是,再繼續(xù)研究混沌在湍流發(fā)展過程中的作用.

自1963 年Lorenz[8]發(fā)表Lorenz 混沌方程至20 世紀(jì)八九十年代,混沌現(xiàn)象的基本特性及由規(guī)則運動向混沌運動轉(zhuǎn)化的機制都已清楚.前者包括分?jǐn)?shù)維特征和系統(tǒng)對初值的敏感性,后者包括周期倍化至混沌[9]、準(zhǔn)周期至混沌[10]及間歇性至混沌[11]等路徑.這里只粗略介紹一下最基本的概念和現(xiàn)象,更詳細的請見文獻[1,12].系統(tǒng)對初值的敏感性指的是小擾動隨時間呈指數(shù)增長.令d(0)為0 時兩條任意軌跡之間的小距離,d(t)為時間t時它們之間的距離,對于真正的低維確定性混沌,我們有

其中,λ1被稱作最大的Lyapunov 指數(shù).混沌吸引子的軌跡在相空間中有界.相鄰軌跡的指數(shù)發(fā)散導(dǎo)致的不斷拉伸,加上吸引子的有界性引起的不時折疊,使混沌吸引子經(jīng)常成為分形,其特征是

其中,N(ε) 表示在相空間中完全覆蓋吸引子所需的線性長度不大于 ε 的(最小)盒子數(shù).D稱為吸引子的計盒維數(shù).通常它是一個非整數(shù).對于混沌 Lorenz吸引子來說,D為2.06.值得指出,方程(2)也適用于度量無人經(jīng)營過的山路或海岸線的長度.這種情況,D通常大于1 但小于2.容易看出,當(dāng) ε 越來越小時,山路或海岸線的長度越來越長.特別地,當(dāng) ε →0 時,山路或海岸線的長度也變成無限長.

當(dāng)時,混沌時間序列分析的基本框架也已完善.給定一個時間序列x(1),x(2),···,x(n),首先需構(gòu)造一個相空間,其中的向量為

其中,m稱作嵌入維,L是延遲時間.然后,可以用Wolf 等[13]的算法計算最大Lyapunov 指數(shù),并用Procaccia 和Grassberger[14]的算法計算關(guān)聯(lián)維(計盒維數(shù)D的一個很緊湊的下限)和關(guān)聯(lián)熵(λ1的一個很好估計).一般假設(shè),當(dāng)關(guān)聯(lián)維為非整數(shù)、Lyapunov 指數(shù)和關(guān)聯(lián)熵為正時,所研究的時間序列是混沌的.當(dāng)我們研究圓柱尾流有無混沌現(xiàn)象時,也用到了這個假設(shè)[15].但是這樣做是不對的.有諸多原因.其一是,用Wolf 等[13]的算法計算最大Lyapunov指數(shù)時,小擾動的指數(shù)增長特性,只是被假設(shè),但沒被證明.事實上可以證明,從任意一組隨機的數(shù)據(jù),都可以計算出一個正的Lyapunov 指數(shù).一方面,這與混沌來自決定性系統(tǒng)這個原理完全不一致.另一方面,這并不表示隨機系統(tǒng)里小擾動有指數(shù)發(fā)散的特征.事實上,證明一個復(fù)雜時間序列的小擾動有指數(shù)發(fā)散的特征,是一個極富挑戰(zhàn)性的問題.第二個原因是,一類很常見的隨機過程,1/fα過程,也擁有非整數(shù)的關(guān)聯(lián)維和正的關(guān)聯(lián)熵,因而會被誤判成混沌[16-17].第三個重大原因是,實驗數(shù)據(jù)總有噪聲.當(dāng)時的分析方法沒法幫助研究者撇開噪聲,把數(shù)據(jù)里最重要的特性刻畫出來.第四個原因是技術(shù)上而非原理性的:實驗數(shù)據(jù)長度有限,因此式(3) 里m和L的選擇就變得很重要.這個問題被稱作相空間的最優(yōu)重構(gòu).1993 至1994 年間,我和鄭先生發(fā)表了三篇文章,把這些問題都解決了.下面解釋一下基本原理.我們從相空間的最優(yōu)重構(gòu)講起.

從原理上來說,式(3) 里的m只要夠大就行,L可以任意.m夠大指的是不小于2 倍的計盒維數(shù)D[18].問題在于D需要估計而非已知.所以相空間的最優(yōu)重構(gòu)這個問題,本質(zhì)上是在數(shù)據(jù)有限這個約束條件下,估計最小可用的嵌入維m,然后選擇L使得動力學(xué)的演化能被最清楚地觀察到.m不能太小,可以通過考慮將南北極投影到赤道來理解-在2 維的赤道上,沒法區(qū)分南北極,所以重構(gòu)的相空間必須夠大,使得地球看起來確實像個地球.L的選擇與相空間里運動的速率有關(guān)-若速率變化忽大忽小非常極端,則動力學(xué)的演化很難被觀察清楚.當(dāng)時已有一個很好的方法,假的最近鄰(false nearest neighbor)方法,來最優(yōu)地重構(gòu)相空間[19-20].這是一個幾何的方法,其基本想法是,當(dāng)m增大至最小可用的嵌入維時,假的最近鄰的個數(shù)會急劇減少.嵌入維選定后,延遲時間可以通過讓假的最近鄰個數(shù)最小而得到.

我和鄭先生構(gòu)造的方法是動態(tài)的[21-23].基本想法是假定最近鄰?fù)ㄟ^動態(tài)迭代后的發(fā)散遠快于真的最近鄰的發(fā)散.用方程表述是

其中,Vi和Vj等是重構(gòu)的向量,k是演化時間,尖括號表示滿足下述條件的所有可能的 (Vi,Vj) 對的系綜平均

其中,εk和 Δ εk是任意選擇的小距離.這里的每一個不等式對應(yīng)一個球殼,引入球殼的根本作用是為了刻畫噪聲的影響.這一點在后面討論混沌的直接動力學(xué)判據(jù)時,會更清楚.為最優(yōu)地決定m和L,可以固定k至一個較小的值.以Rossler 吸引子的重構(gòu)為例,見圖1,我們發(fā)現(xiàn),m=3,L=8 是一個優(yōu)化解.更多的細節(jié)請見文獻[21].

圖1 Rossler 吸引子的重構(gòu)Fig.1 Reconstruction of Rossler's attractor

值得指出,這個動態(tài)的方法和假的最近鄰法,當(dāng)時是,現(xiàn)在仍是解決相空間重構(gòu)的兩個最系統(tǒng)的方法.我們的動態(tài)方法提供更多的信息.以Lorenz 吸引子Lyapunov 指數(shù)的估計為例,見圖2.當(dāng)向量對(Vi,Vj)無所限制,即下面不等式(6)里的w為1 時,Λ(k)隨時間變化的斜率比Lyapunov 指數(shù)小.原因是,對應(yīng)Lyapunov 指數(shù)為0 的軌道上的切向運動也被計算在里邊了.讓w大一些,如選取54,則 Λ (k) 隨時間的變化的斜率與最大Lyapunov 指數(shù)完全一致

圖2 Lorenz 吸引子Lyapunov 指數(shù)的估計Fig.2 Estimation of the Lyapunov exponent for the Lorenz attractor

更重要的是,這個動態(tài)的方法給出了一個混沌的非?？煽康闹苯觿恿W(xué)判據(jù).如圖3(a)所示,對真正的混沌系統(tǒng)來說,來自不同球殼的隨時間變化的指數(shù)曲線會形成一個包絡(luò)線,其斜率給出最大Lyapunov 指數(shù)的一個非常好的估計.對隨機數(shù)來說,如圖3(b)所示,來自不同球殼的隨時間變化的指數(shù)曲線,雖然都能給出一個正的Lyapunov 指數(shù)(相當(dāng)于L(k)/k),但因為是分散的,用不同的球殼去估計最大Lyapunov 指數(shù)時,得到的結(jié)果都不一樣.這是隨機系統(tǒng)的本質(zhì)特征.特別地,當(dāng)數(shù)據(jù)集的大小趨向于無限,而球殼的大小趨向于0 時,最大Lyapunov 指數(shù)將變得無限大.這也是隨機系統(tǒng)的一個根本特征——熵?zé)o限大的一個很好的體現(xiàn).

圖3 隨時間變化的指數(shù)曲線作為混沌的直接動力學(xué)判據(jù):(a) Lorenz 系統(tǒng),(b) 隨機系統(tǒng).數(shù)字從小到大對應(yīng)球殼從大到小[21]Fig.3 Time-dependent exponent curves for the chaotic:(a) Lorenz data and (b) IID random variables,where the curves,from bottom up,correspond to shells from larger to smaller[21]

有了這個混沌的直接動力學(xué)判據(jù),我們可以搞清楚圓柱近尾流低雷諾數(shù)的流動是否有混沌運動[24].圖4(a)顯示的是Re=136 時的相圖.類似的相圖在文獻里經(jīng)常出現(xiàn),如有節(jié)律的化學(xué)反應(yīng)和原子力顯微鏡[25],且被當(dāng)作混沌運動存在的證據(jù).這里,小擾動以冪律而非指數(shù)方式增長,l nεt～t1/2,見圖4(b),所以這是極限環(huán)上的布朗運動而不是混沌運動.

圖4 圓柱近尾流低雷諾數(shù)(Re=136)流動的混沌行為:(a) 相圖,(b)lnt=ln||Vi+k-Vj+k||隨時間的變化[24]Fig.4 Test of chaotic behavior in the near wake cylinder of low Reynolds number (Re=136):(a) phase diagram and (b) temporal evolution of l nt=ln||Vi+k-Vj+k|| [24]

懂了混沌的直接動力學(xué)判據(jù)后,就能很方便地理解噪聲對動力系統(tǒng)的影響.對混沌系統(tǒng)來說,噪聲的作用是破壞圖3(a)所示的包絡(luò)線,隨著噪聲越來越大,包絡(luò)線從上到下,對應(yīng)球殼從小到大,將逐漸被破壞,直至包絡(luò)線完全消失,系統(tǒng)小擾動的指數(shù)發(fā)散特征不再能被分辨出來.這個現(xiàn)象可以用來估計噪聲的大小[26].總的來說,這套方法能非常好地刻畫噪聲對動力系統(tǒng),包括半導(dǎo)體激光系統(tǒng)的動力學(xué)的影響[27],也能非常好地研究一個相反的問題——噪聲引發(fā)的混沌[28-30].細節(jié)就不再贅述了.這里只強調(diào)一下,這也是我在上一節(jié)提到的,能在沒有任何基礎(chǔ)的情況下,簡單融入劉教授的團隊,研究用半導(dǎo)體激光做混沌安全通訊的原因.

3 基于依賴指數(shù)的李雅普諾夫指數(shù)(SDLE)的多尺度分析

混沌理論對20 世紀(jì)科學(xué)思想的最大沖擊是,隨機性可以不由無限維的隨機系統(tǒng)產(chǎn)生,而由有限維的確定性系統(tǒng)產(chǎn)生.這個將隨機系統(tǒng)與確定性系統(tǒng)對立的表述,產(chǎn)生的一個很嚴(yán)重的后果是,早期的一些研究者,包括我自己,在選擇了混沌研究后,會忽略隨機系統(tǒng)的研究.1997 年在UCLA 的電子工程系研究互聯(lián)網(wǎng)上的交通流時,我意識到這個偏見很糟糕.任一時間、任一地點的互聯(lián)網(wǎng)牽涉到的用戶數(shù),都遠多于地面交通牽涉到的車輛數(shù).所以正常情況下互聯(lián)網(wǎng)上的交通流必然是隨機的,這種情況是相對于網(wǎng)絡(luò)受攻擊時的狀態(tài)而言,當(dāng)網(wǎng)絡(luò)受攻擊時,其上的交通流會有一個成分,由程序(即規(guī)則)產(chǎn)生,因此會不再是完全隨機的.互聯(lián)網(wǎng)上交通流的一個基本特征是長程相關(guān)性[31].雖然離散的混沌映射能模擬這個特征[32],但兩個著名的隨機分形模型,1/fα過程和乘法級聯(lián)多分形模型,能更全面地模擬互聯(lián)網(wǎng)上交通流的特征[33-36].值得注意,乘法級聯(lián)多分形模型是為理解湍流而被發(fā)展起來的[37-41],且至今仍是湍流最好的模型之一.早期復(fù)雜性科學(xué)研究的一大目標(biāo)是用低維混沌理論研究湍流.所以這里一個深刻的問題是如何將低維混沌與隨機分形理論有機地結(jié)合起來.SDLE 是為了實現(xiàn)這個目標(biāo)的一個很有成效的嘗試.為了充分理解SDLE 的涵義,下面先介紹一下 1/fα過程.隨機系統(tǒng)的其他主要模型,包括ARIMA 模型,Levy 過程,和乘法級聯(lián)多分形模型,就不在此贅述了(可見文獻[1]).

在表征復(fù)雜系統(tǒng)的活動類型中,一類非常普遍且令人費解的行為是 1/fα過程,其功率譜密度隨頻率(時域)或波長(空間域)以冪律方式衰減[1],因此其維數(shù)不能通過主成分分析等方法來有效降低[42].這種過程的一個子類,通常表示為 1/f2H+1過程,具有長程相關(guān)性,其中,H稱作赫斯特指數(shù)(Hurst Parameter).取決于 0＜ H ＜1/2,H=1/2,或 1/2 ＜ H ＜1,它們分別被稱作具有反持久相關(guān)性、無記憶或僅短程相關(guān)性,及持久長程相關(guān)性(長記憶)[43].此類過程的突出例子包括視覺[44]、金融[45-46]、DNA 序列[47-51]、人類認知[52]、全球恐怖主義[53]、協(xié)調(diào)[54]、姿勢[55-56]、心臟動力學(xué)[57-61]以及素數(shù)的分布[62]等.

具有長記憶過程的精確定義如下:均值為 μ、方差為 σ2,且協(xié)方差平穩(wěn)的隨機過程X={X(t),t=0,1,2,···},當(dāng)其自相關(guān)函數(shù)r(w),w≥0 以冪律方式衰減時[63]

則有長程相關(guān)性,其中 0 ＜H＜ 1.當(dāng) 1/2 ＜H＜1,這是為什么這種過程被稱作具有長記憶.過程X的功率譜密度 (PSD) 為它的積分被稱為隨機游走過程,其PSD 為作為1/f過程,它們不能被馬爾可夫過程或ARIMA 模型恰當(dāng)?shù)亟64],因為這些過程的 PSD 與1/f明顯不同.為了充分模擬1/f過程,必須使用分?jǐn)?shù)階過程.最流行的模型是分?jǐn)?shù)布朗運動模型[43].

值得強調(diào)的是,充分發(fā)展的湍流的Kolmogorov 能譜對應(yīng)H=1/3,其效應(yīng)在用雷達觀察到的海平面的海雜波上也有反映[65-66].

在上述基礎(chǔ),我們可以討論SDLE.

假設(shè)一個適當(dāng)?shù)南嗫臻g已經(jīng)重構(gòu)好了,我們考慮一個軌跡的集合.將兩個相鄰軌跡之間的初始分離表示為 ε0,它們在時刻t和t+Δt的平均距離為εt和 εt+Δt.相空間中軌跡分離的一個示意如圖5 所示.我們考察 εt和 εt+Δt的關(guān)系.當(dāng) Δt→0,我們有

圖5 相空間中軌跡分離的示意圖Fig.5 Schematic of evolution of separation of trajectories in phase space

其中,λ (εt) 就是 SDLE,由下式給出

也可以把它表示成一個微分方程

有趣的是,λ (εt) 的計算,完全可以借用上一節(jié)描述的我和鄭先生發(fā)展的方法.具體是,還是基于不等式(5)描述的球殼,我們直接得到

其中,t和 Δt是采樣時間的整數(shù)倍,尖括號表示球殼內(nèi)所有下標(biāo)為i,j的平均值.注意,SDLE 蘊含Lyapunov 指數(shù)的概念,但比后者的概念寬廣,因SDLE 是一個函數(shù),而Lyapunov 指數(shù)只是一個數(shù).

粗略地說,當(dāng)時我和鄭先生發(fā)展的方法可以稱作是一個積分的形式,SDLE 的表述的是一個微分的形式.由積分形式轉(zhuǎn)至微分形式帶來的好處是,在SDLE 的框架下,對所有已知的時間序列模型,我們都可以證明它們有其獨特的標(biāo)度律.下面列舉一些,更多的細節(jié)請見文獻[1,61-67].

性質(zhì)1:對確定性混沌,可以證明,在小尺度上

性質(zhì)2:對于有噪聲的混沌和噪聲引起的混沌,在小尺度上,有

其中 γ ＞0 決定信息衰減的速率;

性質(zhì)3:對1/f2H+1過程,可以證明

性質(zhì)4:對 α-stable Levy 過程,可以證明

性質(zhì)5:對隨機振動來說,依賴于相空間的重構(gòu),λ(ε)≈-γlnε 和 λ (ε)≈-Hε-1/H都可被觀察到;

性質(zhì)6:對于具有多尺度行為的復(fù)雜運動,上述不同的標(biāo)度律有可能在不同的 ε 范圍內(nèi)都被觀察到.

上述6 個性質(zhì),涵蓋了幾乎所有已知的時間序列模型.因而我們可以期望,SDLE 能統(tǒng)一所有已知的復(fù)雜性的度量.事實上也確實如此.比如在刻畫腦電動力學(xué)時,我們發(fā)現(xiàn),常用的復(fù)雜性的度量,包括Lyapunov 指數(shù),關(guān)聯(lián)維,關(guān)聯(lián)熵,赫斯特指數(shù),等都可以被對應(yīng)到SDLE 在某個尺度的值.細節(jié)請見文獻[1](15.7.1 節(jié))和文獻[68].因這個原因,SDLE 在臨床腦電的分析中有極大的應(yīng)用前景[69].

因性質(zhì)6 有重大意義,我們用一個例子使其具體化.考察一個動力系統(tǒng)

其中,[xn] 表示xn的整數(shù)部分,ηn是在區(qū)間 [-1,1] 上均勻分布的隨機數(shù),σ 刻畫噪聲的強度,F(y) 是個由下述方程給出的函數(shù)

函數(shù)F(y) 描述的是一個混沌動力學(xué),其Lyapunov指數(shù)為 l n(2+Δ).所以,由方程(16)描述的動力系統(tǒng)是一個整數(shù)格點上的隨機運動加上整數(shù)格點內(nèi)的混沌運動.因此,當(dāng)沒有噪聲時,我們應(yīng)該在小尺度上觀察到性質(zhì)1,而在較大尺度上觀察到性質(zhì)3.確實如此,如圖6 所示.當(dāng)有噪聲時,性質(zhì)1 被性質(zhì)2 替代.這里,圖6(a)用了5000 點,圖6(b)用了500 點.后一種情形,因點太少,整數(shù)格點上的布朗運動沒能被刻畫好(事實上,真正牽涉整數(shù)格點的個數(shù)只有幾十).總的來說,用這么少的點刻畫這么豐富的動力學(xué),很令人驚訝.

圖6 小尺度上的混沌運動與大尺度上的布朗運動Fig.6 Chaotic motion on small scales and Brownian motion on large scales

性質(zhì)6 的一個直接后果是,我們可以把一個混沌或分形時間序列拆成很多段,中間插入周期性運動,則原來的混沌或分形特征還是能被SDLE 很好地刻畫出來.事實上,前者對應(yīng)間歇性混沌.一個例子由圖7 給出.分形+振蕩是長時間監(jiān)測的心跳動力學(xué)的一個基本特征.這里,振蕩由呼吸引起.這樣的分形很難被SDLE 以外的任何方法刻畫.細節(jié)就不再此贅述了,請見文獻[61].

圖7 (a) 由Logistic 映射 xn+1=axn (1-xn) 生成的間歇性混沌時間序列,a=3.8284.(b) 由10000個點的時間序列計算出的 SDLE 曲線,m=4,L=1,最大的球殼尺寸為 (2 -13.5,2 -13).由箭頭 A 指示的混沌標(biāo)度區(qū)清晰可見.箭頭 B 和 C 表示的區(qū)域?qū)?yīng)周期性運動和混沌運動之間的過渡,箭頭 D 表示大尺度的振蕩運動Fig.7 (a) Intermittent chaos generated by Logistic map xn+1=axn (1-xn),a=3.8284.(b) SDLE curve calculated from a time series of 10000 points,m=4,L=1,and the largest shell is (2 -13.5,2 -13).The chaotic scale region indicated by arrow A is clearly visible.The regions represented by arrows B and C correspond to the transition between periodic and chaotic motions,while arrow D represents large-scale oscillatory motion

思考一下,我們就能夠明白,SDLE 的提出,拓展了傳統(tǒng)混沌動力學(xué)的研究范疇,因后者一般只研究由不變測度刻畫的含有混沌吸引子的動力系統(tǒng).在大數(shù)據(jù)時代,數(shù)據(jù)一般都是被連續(xù)觀察的.真正讓人激動的數(shù)據(jù)是那些含有狀態(tài)突變的數(shù)據(jù).這種情形,雖不易由大多數(shù)已知的方法來處理,但對SDLE 并沒造成任何挑戰(zhàn).這里,值得提一下河流徑流動力學(xué).許多研究者猜測河流動力學(xué)是混沌的,但分析結(jié)果不能讓人信服.原因是,河流動力學(xué)是強間隙性的—旱季的枯水期和雨季的豐水期,其動力學(xué)是完全不一樣的.因此,河流的流動形態(tài)中若存在混沌,該混沌必須是間隙性的.所幸,這里的混沌很容易通過SDLE 刻畫.一個例子如圖8 所示.

圖8 美國Colorado 河的間隙性混沌.這里,紅色的曲線代表由下節(jié)的非線性自適應(yīng)濾波器處理過的數(shù)據(jù)Fig.8 Intermittent chaos in the Colorado River,USA.Here,the red curves represent the data processed by the nonlinear adaptive filter to be explained in the next section

4 自適應(yīng)去趨勢、去噪、多尺度分解和分形分析

2006 年,美國因在伊拉克和阿富汗缺兵力,由陸軍研究部負責(zé),請一家公司設(shè)計了一款能用于野外很方便采集腦電的頭盔,然后請杜克大學(xué)醫(yī)學(xué)院的一個團隊采集了一組腦電,最后找到我們,希望我們能設(shè)計出一套通過腦電分析,評估受炮彈沖擊波沖擊而腦震蕩的士兵,是否已恢復(fù)可以再上前線了.因這個經(jīng)費來之不易,我們做的非常認真.做了幾個月,發(fā)現(xiàn)實驗得到的腦電信號,其實沒有任何真正腦電的信息.但我們不能就這樣告訴陸軍研究部.所以我們設(shè)計了一個新的分析方法,證明給我們的腦電,只有頭部搖晃的信號,加上電源60 Hz 的信號(中國是50 Hz),及一些噪音.陸軍研究部雖然對我們的分析很滿意,但畢竟很失望,因為設(shè)計的頭盔完全無效.我自己得到的一個深刻體會是,實驗數(shù)據(jù)采集的團隊,必須包含能做分析的人員,不然采集的數(shù)據(jù)可能都是垃圾.

我們的方法是這樣的[70-75].先把一個時間序列分成長度為奇數(shù)w=2n+1 的小段,相鄰兩段重疊n+1個點.這樣做引入了一個時間尺度=(n+1)τ,其中 τ 為采用時間間隔.對于每一段,我們擬合一個M階的最佳多項式.請注意,M=0 和 1 分別對應(yīng)于分段常數(shù)和分段線性擬合.分別用yi(l1) 和yi+1(l2)表示擬合的第i段和第(i+1)段的多項式,其中l(wèi)1,l2=1,2,···,2n+1.注意最后一段的長度可能小于2n+1.我們將重疊區(qū)域的擬合定義為

其中,w1=1-(l-1)/n,w2=(l-1)/n可以重寫為1-dj/n,j=1,2,其中dj表示那一點與y(i) 和y(i+1)的中心的距離.這意味著權(quán)重隨點和線段中心之間的距離線性減小.這種加權(quán)確保了對稱性并有效地消除了相鄰段邊界周圍的任何跳躍或不連續(xù)性.實際上,該方案確保擬合處處連續(xù),在非邊界點處光滑,在邊界處具有右導(dǎo)數(shù)和左導(dǎo)數(shù).該方法能有效找出任何一種信號的趨勢,假如存在的話.

回到杜克大學(xué)的腦電數(shù)據(jù).兩個例子如圖9 所示.紅色的曲線由我們的自適應(yīng)濾波器得到.它們對應(yīng)頭部的晃動.較細的黑色趨勢線由另一個很有名的濾波器Loess[76]得到.Loess 等價于Savitzky-Golay 濾波器[77].可以看出,紅色的趨勢線比黑色的更精確.

圖9 腦電因頭部晃動引起的趨勢Fig.9 Trend of EEG caused by head shaking

值得比較這個濾波器與前面提到EMD[2].圖10(a)顯示的是全球年海表溫度 (SST)及線性、拋物線、非線性擬合,右邊顯示的是相應(yīng)的殘差序列.圖10(b)與文獻[2]的結(jié)果完全一致.值得注意,EMD 是個二分法,時間尺度每次減半.這里的自適應(yīng)算法牽涉到的時間尺度是連續(xù)的,因此,在許多應(yīng)用場合,自適應(yīng)算法應(yīng)該比EMD 更準(zhǔn)確.

圖10 (a)全球年海表溫度(SST)及線性、拋物線、非線性擬合和(b)相應(yīng)的殘差序列Fig.10 (a) Global annual sea surface temperature (SST) and linear,parabolic and nonlinear fitting and (b) corresponding residual series

自適應(yīng)濾波器有非常好的去噪功能.一個例子如圖11 所示.可以看出,自適應(yīng)濾波器比基于混沌的方法和小波變換都好.細節(jié)請見文獻[71].

圖11 混沌洛倫茲信號的去噪Fig.11 Denoising of a chaotic Lorenz signal

自適應(yīng)濾波器也有極好的信號提取功能.一個例子如圖12 所示.其中,圖12(a)和12(b) 是以采樣頻率為250 Hz 采集的 EEG 和 ECG 信號;圖12(c)和12(d)是通過使用自適應(yīng)和小波變換從 EEG 中提取的ECG 成分.圖12(a)中的箭頭突出顯示了 EEG信號中包含的 ECG 分量.在呼吸綜合征的監(jiān)測和治療中,腦電EEG 的分析起到很重要的作用.但混雜了心電ECG 成分的腦電,使得通過分析判斷是否有呼吸綜合征變得很困難.因此,一個重要步驟是先去掉心電ECG 的成分.如圖所示,自適應(yīng)方法比小波變換更優(yōu).

圖12 (a,b)以采樣頻率為250 Hz 采集的 EEG 和 ECG 信號,單位為mV;(c,d)使用自適應(yīng)和小波變換從 EEG 中提取的ECG 成分.(a)中的箭頭突出顯示 EEG 信號中包含的 ECG 分量Fig.12 (a,b) EEG and ECG signals acquired at a sampling frequency of 250 Hz in mV;(c,d) ECG components extracted from EEG using adaptive and wavelet transforms.Arrows in (a) highlight ECG components included in the EEG signal

理解了自適應(yīng)濾波器上述功能后,就很容易理解它的多尺度分解的功能.這相當(dāng)于用一系列窗口w1,w2,...,wk,得到一連串趨勢線,任意兩個趨勢線的差,都代表原來時間序列在對應(yīng)那兩個窗口的頻段的成分.對具體例子感興趣的,請見文獻[73].

最后,我們解釋如何用自適應(yīng)濾波器作分形和多分形分析.事實上,當(dāng)數(shù)據(jù)有趨勢時,它是所有估計赫斯特指數(shù)及廣義的赫斯特指數(shù)譜的方法中最好的.這里只給出一些最基本的結(jié)果.細節(jié)請見文獻[73].

然后,我們用窗口大小為w對u(n) 擬合一個全局趨勢v(n).殘差u(i)-v(i) 表征圍繞全局趨勢的波動,其方差是赫斯特指數(shù)H的一個非常好的估計

將上式的平方與開根號稍作變化,我們就得到估計廣義的赫斯特指數(shù)譜的多分形分析

將式(20)應(yīng)用于兩個隨機游走過程u(n) 和z(n),我們可以計算它們間互相關(guān)的長程相關(guān)指數(shù)Huz

圖13(a)顯示的是一個典型的腦電信號的分形分析結(jié)果.從中我們可以定義三個參數(shù),短時間和長時間的赫斯特指數(shù)(記為HS和HL),及飽和度(saturation level),就可以把三類腦電數(shù)據(jù),正常(Group H),癲癇發(fā)作(Group S),癲癇病人未發(fā)作(Group E),幾乎完全區(qū)分開.近期,我們發(fā)現(xiàn)這個方法有很大的臨床應(yīng)用前景[78].

圖13 (a)一個典型腦電信號的分形分析結(jié)果;(b) 三類數(shù)據(jù),正常(Group H)、癲癇發(fā)作(Group S)、癲癇病人未發(fā)作(Group E)的分類,精度幾乎100%Fig.13 (a) Typical adaptive fractal analysis of a an EEG signal;(b)Classification of three types of EEG data,normal (Group H),seizures(Group S),and seizure patients not having seizure (Group E),with an accuracy of almost 100%

5 結(jié)語

5 G 時代的到來及 6 G 時代的快速臨近,科學(xué)、工程和社會中積累的大數(shù)據(jù)將很快變得更加巨大.任何人都不得不抓住這個前所未有的機會.雖然計算機科學(xué)家們正在努力開發(fā)更強大的數(shù)據(jù)庫管理和機器學(xué)習(xí)方法,但現(xiàn)在是進入下一階段的時候了.這個階段將以捕獲大數(shù)據(jù)背后的動力學(xué)及深刻理解人腦的工作機制為基本目的.就數(shù)據(jù)分析而言,我們可以洞見,機器學(xué)習(xí)和復(fù)雜性科學(xué)的主流方法在未來會相互補充并互相促進.為了幫助加速這個結(jié)合,我們提倡協(xié)同使用基于主流機器學(xué)習(xí)的方法和復(fù)雜性科學(xué)中的多尺度方法.

最后我們強調(diào),在科學(xué)研究中,最重要的是一直保有一顆好奇、不斷探索的心.只有這樣,才能發(fā)現(xiàn)真正重大的科學(xué)問題并切實解決實際問題.這種普世的精神,鄭先生等老一輩科學(xué)家們一直在言傳身教.望讀者們能和我們共勉！