嚴捷冰, 楊會杰
(上海理工大學管理學院,上海 200093)
行走取決于非常復雜的控制機制,它是由運動系統執行的復雜的、有節奏的任務.跨步時間變異性(stride duration variability)被視為步態平衡的標志,大多數研究都集中在自發步行速度的跨步時間變異性上.目前主要的研究方法有兩種:統計學方法和復雜性理論.前者主要考察序列漲落的量級,有研究表明,較大的標準差(SD)和變異系數(CV)與步伐不穩定[1-2]和容易摔倒[3-6]有關.后者主要研究漲落的動力學特性,其中,通過計算序列的分形指數考察序列的長程相關性質是重要的復雜數學方法之一.
Hausdorff等[7]在1999年發表了一篇研究兒童步態成熟的文章,研究了50個3~14歲孩子的行走序列,用多種方法論證了他們的觀點,即對于3~14歲的孩子來說,他們的步態有隨年齡的增長趨于更加穩定的趨勢.文獻[7]應用二階去趨勢漲落分析法(DFA-2)選取觀察窗口的范圍為10≤n≤20,計算每一組數據的分形指數,并指出11~14歲孩子的行走序列的分形指數比3~4歲和6~7歲孩子的行走序列的分形指數小得多.Bollens等[8]基于不同的實驗樣本,對不同年齡段的人,用重標極差法(R/S分析法)分析了他們的行走序列(長度為512步)的分形指數.研究結果表明,行走序列的分形指數與人的年齡并沒有必然的關系.然而,Pacheco等[9]的研究結果表明,在序列長度達到215時,像去趨勢漲落分析(DFA)這樣基于方差的方法能給出偏差較小(偏差小于0.05)的赫斯特指數,對于長度小于215的序列來說,基于方差的方法的計算結果的偏差較大(偏差大于0.05);而R/S分析法在序列長度為26~216范圍內始終有較大的偏差(偏差大于0.05).鑒于行走序列屬于短序列(29),而上述兩種方法在處理這種長度的序列時都有較大的偏差,且基于不同方法得出的結果互相矛盾.因此,使用更為可靠的、適合分析短序列的方法重新評估這50個3~14歲孩子的行走序列的分形指數是十分必要的.為了能有效地計算短序列的分形指數,文獻[10]提出了擴散熵的平衡估計(BEDE),從現有的研究結果來看,該方法能在較高的精度下計算出102的短序列的分形指數[10-11].
本文應用DFA-2在更大的尺度上計算了每一個序列的分形指數,以考察在較大的尺度上DFA-2的計算結果是否與文獻[7]中的計算結果相吻合.然后,應用BEDE重新計算出了每個樣本序列的分形指數,并將其與DFA-2的結果進行了對比.在上述兩個計算結果中,均沒有發現Hausdorff等在文獻[7]提到的現象,計算結果顯示,行走序列的分形指數與實驗對象的年齡無明顯關系,這與Bollens等[8]的最新研究結果相類似.
Scafetta等[12]在提出擴散熵分析(DEA)的文章中推導出了DEA的標度參數δ與基于方差的方法的指數H,在序列具有分數布朗運動特征的時候滿足δ-H,而當序列具有Levy行走特征時滿足δ本文在上述研究的基礎上,又對這50個樣本序列是否具有分數布朗運動或Levy行走的特征進行討論.從計算結果中可以看到,在“年輕”年齡組中只有1個樣本具有Levy行走的特征.在“中等”和“年長”年齡組中一共有10個樣本具有分數布朗運動的特征,在這2個年齡組中也有10個樣本具有Levy行走的特征,其中,又有7個樣本數據同時具有分數布朗運動和Levy行走的特征.
本文沿用了Hausdorff等在文獻[7]中使用的數據,并且按照文獻[7]中敘述的辦法對原始數據進行處理,除去了一些可能的異常值.
在文獻[7]中,Hausdorff等在400m的跑道上采集了50名3~14歲的健康兒童時長為8min的自發步行速度下的行走序列.在獲取數據之后按以下3個步驟去除數據中的一些異常值:
第1步 去除前60s以及最后5s的數據,以消除啟動和結束時可能存在的干擾;
第2步 去除“暫停”(即大于2s的步幅間隔)以及暫停前后各5s的數據;
第3步 去除較大的步幅間隔以及其它異常值,其方法是排除原始數據中最大和最小的5%的數據,計算剩余數據的均值和方差,去除原始數據中與該均值相差4倍方差以上的數據點.
去趨勢波動分析法是Peng等[13]在1994年提出的,其具體步驟如下:給定一個有界時間序列xt,t∈?,首先,計算它的累積離差〈xi〉);接著,將Xt分割成窗口長度為L的序列Yj,對Yj應用最小二乘法求出使其平方誤差最小的擬合函數;然后,計算出Yj-的均方根誤差F最后,對不同的窗口大小L重復上述過程,可得到如下式所示的關系

式中,指數α是在雙對數坐標下應用最小二乘法求得的L對F(L)作線性擬合后的斜率.
為了能有效地計算短序列的分形指數,Scafetta等[12]在香農熵的基礎上提出了擴散熵,擴散熵能有效地區分出復雜系統中Levy行走序列和分數布朗運動序列.在此,對擴散熵作簡單的推導.
對一維穩態時間序列{x1,x2,…,xN}作相空間重構

式中,N為個數.
從擴散的角度出發,可以將上述N-s+1個序列視為一個從原點出發走了s步的粒子運動軌跡的N-s+1次抽樣,則第k次抽樣的位移

式中,k指第k次抽樣.
將Yk(s)的分布范圍劃分成M(s)等分并統計落在每個小區間內的位移個數n(k,s),k=1,2,…,M(s).可用下式來估計位移的概率分布函數

上述擴散過程的熵估計可表示為

式中,SDE(s)為擴散熵.
如果擴散過程具有自相似結構,則有

式中,A是一個常數.

將修正由序列長度是“有限的”所帶來的誤差視為一個最優化問題,因此,定義

則擴散熵可表示為


這里,概率p(j,s)∈[0,1],n(j,s)∈{0,1,2,…,N-s+1}.取平均誤差在整個p(j,s)∈[0,1]范圍內的最小值

式中,w[p(j,s)]是一個由具體問題所決定的權重函數,為簡便起見,可設w[p(j,s)]=1.
因此,可得平均誤差

式中,Pn(j,s)[p(j,s)]是二項分布.

要得到最小的平均誤差,就是要求下述偏導數都等于零.

經計算可得

由此可得擴散熵的平衡估計為

仿照文獻[7],按照兒童的年齡將這50個樣本數據分成“年輕”、“中等”和“年長”這3組,即3~4歲(編號:1~11),6~7歲(編號:15~34)和11~14歲(編號:39~50).此外,還有幾個過渡年齡段的樣本數據:5歲(編號:12~14),8歲(編號:35)和10歲(編號36~38).圖1~4分別顯示了“年輕”年齡組(3~4歲)、過渡年齡段、“中等”年齡組(6~7歲)以及“年長”年齡組(11~14歲)這幾組兒童行走序列樣本在BEDE下分形指數的擬合結果.

圖1 3~4歲年齡組的線性擬合Fig.1 Linear data fitting for the group of 3~4-year-olds
如圖1~4所示,k′為斜率.就這50個樣本數據而言,一共有12個樣本在較大的尺度上不能用一條直線來擬合,它們分別是3,4,6,8,12,13,15,35,38,43,44,47號樣本,不能用直線擬合的數據中,在“年輕”組中占的比例最高,在“中等”組中占的比例最小.在這12個樣本中,8,12,13,15,35,38這6個樣本與二次或三次多項式符合得非常好,且這6個數據都不屬于“年長”組.此外,7,21,22,30,36,40,49號樣本似乎具有離散分形結構,其中,以7和22這2個樣本最為明顯.其余的樣本數據都在相對較大的尺度范圍內可以用一條直線來擬合.

圖2 5,8,10歲樣本的線性擬合Fig.2 Linear data fiting for the group of 5,8and 10-year-olds
首先應用BEDE計算處理完的數據的分形指數,并與DFA-2的結果進行對比.圖5(a)顯示了BEDE的分形指數的分布情況,從圖示的結果來看,沒有看出3個年齡段之間有什么明顯的區別,這與DFA-2的結果相類似,如圖5(b)所示.然而,將窗口范圍限定在10≤n≤20之后,應用DFA-2計算出來的分形指數則顯示出了如Hausdorff等在文獻[7]中所提到的結果,即2個較小年齡組的分形指數相類似,而最大年齡組的分形指數明顯小于前2組的,如圖5(c)所示.y為年齡,α為DFA-2的分形指數,σ為BEDE的分形指數.

圖3 6~7歲年齡組的線性擬合Fig.3 Linear data fiting for the qroup of 6~7-year-olds

圖4 11~14歲年齡組的線性擬合Fig.4 Linear data fiting for the group of 11~14-year-olds
Scafetta在提出DEA的同時,還對基于不同方法計算出的指數進行了詳細的對比[12],研究表明,當序列具有分數布朗運動特征時,基于方差的方法的分形指數與DEA計算出來的分形指數是等價的,即H=δ.而當序列具有Levy行走特征時,上述兩種方法計算出來的分形指數滿足關系

圖5 行走序列的分形指數Fig.5 Fractal scaling index of walking series
表1同時給出了這兩種方法計算出的所有樣本數據的分形指數,并且對樣本數據是否具有分數布朗運動或Levy行走的特征作出了判斷.

表1 分形指數間的關系Tab.1 Relationship between fractal scaling indices
從表1中可以看到,在“年輕”年齡組中只有1號樣本能被視為Levy行走序列.在“中等”和“年長”年齡組中一共有10個樣本可以被視為分數布朗運動,它們分別是16,17,19,24,30,31,33,39,42,46.在“中等”和“年長”年齡組中也有10個樣本可以被視為Levy行走,它們分別是16,19,24,30,31,39,40,42,48,50.在這中間,16,19,24,30,31,39,42這7組數據同時滿足分數布朗運動和Levy行走.
研究行走序列對于認識人是如何控制行走步態,以及理解由于某些疾病或者衰老導致的走路不穩和摔倒,具有重要的意義.由于行走序列屬于短序列,而現有的研究方法在研究短序列時的偏差都相對較大,得不到一個較為可靠的結果,因此,尋找新的能夠用來分析短序列的可靠方法是十分迫切的需求.
本文應用DFA-2重復了前人的工作,并使用BEDE這一能夠較為精確地計算出短序列分形指數的方法研究了50個兒童行走序列,得到了一些較為有趣的現象.一方面,從本文的計算結果來看,并未發現行走序列的分形指數與實驗對象的年齡之間存在某種關系;另一方面,本文的研究結果顯示,6~7歲和11~14歲年齡組中的部分兒童的行走序列體現出了分數布朗運動或者Levy行走特征,而這一現象在3~4歲年齡組里則較為罕見.然而,正因為在判斷樣本序列是否是分數布朗運動或者Levy行走時應用了DFA-2的計算結果,對于短序列來說,該結果的計算結果偏差相對較大,因此,嚴格來說,該結論還需要在今后進一步的研究中加以證實.
[1] Beauchet O,Annweiler C,Lecordroch Y,et al.Walking speed-related changes in stride time variability:effects of decreased speed[J].Journal of Neuro Engineering and Rehabilitation,2009,5(6):32.
[2] Beauchet O,Allali G,Annweiler C,et al.Gait variability among healthy adults:low and high strideto-stride variability are both a reflection of gait stability[J].Gerontology,2009,55(6):702-706.
[3] Hausdorff J M,Mitchell S L,Firtion R,et al.Altered fractal dynamics of gait:reduced stride interval correlations with aging and Huntington’s disease[J].Journal of Applied Physiology,1997,82(1):262-269.
[4] Hausdorff J M,Edelberg H K,Mitchell S L,et al.Increased gait unsteadiness in community-dwelling elderly fallers[J].Archives of Physical Medicine and Rehabilitation,1997,78(3):278-283.
[5] Maki B E.Gait changes in older adults:predictors of falls or indicators of fear[J].Journal of the American Geriatrics Society,1997,45(3):313-320.
[6] Verghese J,Holtzer R,Lipton R B,et al.Quantitative gait markers and incident fall risk in older adults[J].Journals of Gerontology.Series A:Biological Sciences and Medical Sciences,2009,64A(8):896-901.
[7] Hausdorff J M,Zemany L,Peng C K,et al.Maturation of gait dynamics:stride-to-stride variability and its temporal organization in children[J].Journal of Applied Physiology,1999,86(3):1040-1047.
[8] Bollens B,Crevecoeur F,Detrembleur C,et al.Effects of age and walking speed on long-range autocorrelations and fluctuation magnitude of stride duration[J].Neuroscience,2012,210(17):234-242.
[9] Pacheco J C R,Roman D T.What is the required series length for correct self-similarity analysis?[EB/OL].[2012-11-22].http:∥biblioteca.coqcyt.gob.mx/bvic/Captura/upload/.../WHAT-ISTHE-REQUIRED-ARTREV.pdf.
[10] Qi J C,Yang H J.Hurst exponents for short time series[J],Physical Review E,2011,84(6):066114.
[11] Zhang W Q,Qiu L,Xiao Q,et al.Evaluation of scale invariance in physiological signals by means of balanced estimation of diffusion entropy[J].Physical Review E,2012,86(5):056107.
[12] Scafetta N,Grigolini P.Scaling detection in time series:diffusion entropy analysis[J].Physical Review E,2002,66(3):036130.
[13] Peng C K,Buldyrev S V,Havlin S,et al.Mosaic organization of DNA nucleotides[J].Physical Review E,1994,49(2):1685-1689.
[14] Kantelhardt J W,Koscielny-Bunde E,Henio H A R,et al.Detecting long-range correlations with detrended fluctuation analysis[J].Physica A,2001,295(3/4):441-454.
[15] Bonachela J A,Hinrichsen H,Munoz M A.Entropy estimates of small data sets[J].Journal of Physics A:Mathematical and Theoretical,2008,41(20):202001.