R/S分析法在石家莊市欒城區(qū)土壤含水量動態(tài)分析中的應(yīng)用

高雄飛，劉元會，郭建青，王媛英

(1.長安大學(xué)理學(xué)院，陜西西安 710064;2.長安大學(xué)環(huán)境科學(xué)與工程學(xué)院，陜西西安 710051)

土壤含水量是土壤的基本物理性質(zhì)指標之一。準確掌握土壤水分的動態(tài)變化是確定灌溉制度并適時進行灌溉的重要依據(jù)［1］。對土壤含水量時間序列的時空演變過程進行研究在土壤含水量預(yù)測預(yù)報中具有重要的作用。目前，處理這類長期觀測資料的方法主要有回歸分析法［2］、頻譜分析法［3］、隨機微分方程法［4］等。這些方法的結(jié)果能夠在一定程度上反映土壤動態(tài)的內(nèi)在物理機制，但是，它們均未能很好地揭示土壤含水量動態(tài)過程的分形特征［5］。土壤含水量由于受到降雨、蒸發(fā)、日照時間等因素的影響，具有高度的非線性，對于這類非線性的具有統(tǒng)計特性的數(shù)據(jù)系列，采用重標極差分析法［6］(rescaled range analysis，簡稱R/S分析法)，可以較好地揭示土壤含水量動態(tài)變化過程的內(nèi)在規(guī)律及其分形特征，從而有效預(yù)報或控制其變化規(guī)律。

Hurst［7］在總結(jié)尼羅河的多年水文觀測資料時提出了R/S分析法，其后Mandelbrot［8］對該方法在理論上進行了補充和完善。Wallis等［9］曾利用此法對河流流量、泥漿沉積量、樹木年輪、降雨量等自然現(xiàn)象進行研究。該方法在其他領(lǐng)域(如物理、生物學(xué)等)也得到了廣泛的應(yīng)用［10］。本文以欒城水文站2002—2008年深度為10 cm、80 cm和150 cm土壤含水量時間序列為例，采用R/S分析法對土壤含水量的動態(tài)過程進行分析，進而借助分形理論對土壤含水量的分形特征進行定量描述，期望揭示當?shù)剞r(nóng)田土壤水分的動態(tài)變化特征。

1 R/S分析法

1.1 R/S分析的基本原理

將時間序列{Xi=1，2，3，…，N}劃分為m個時間長度為n的獨立時間序列，其中m=N/n(m為整數(shù))，子序列記為Da，其中下標a表示第a個子序列［11］。

a.每個時間長度為n的子序列Da的均值計算公式為

式中:Xˉ——子序列平均值;j——子序列下標。

b.每個子序列Da的累積偏差為

式中:Xj，a——子序列的第j個觀測值;Xˉa——第a個子序列的平均值。

c.每個子序列Da的極差為

d.每個子序列Da的標準差為

e.每個子序列的重標極差定義為

f.對每個子序列重復(fù)第b步到第e步的計算，得到一個重標極差序列，計算該序列的均值:

g.將n增加到下一個更大的值，重復(fù)第a步至第f步的計算，直到，得到一個序列。假設(shè)該序列累積偏差的極差與觀測數(shù)的一個待定冪次方成正比［12］，即

其中

式中:τ——時間步長;C——待定常數(shù)，可以理解為反映自然現(xiàn)象固有特性的常數(shù);H——時間序列的Hurst指數(shù)。

對式(7)兩邊取對數(shù)，有

式(8)利用最小二乘法進行求解，其斜率H即為Hurst指數(shù)。

H表示了時間變換的標度因子，因此時間序列具有自仿射性，經(jīng)過嚴格的數(shù)學(xué)推導(dǎo)［13］，時間序列的分形維數(shù)與H的關(guān)系［14-15］如下:

式中:D——分形維數(shù)，D定量地描述了土壤含水量動態(tài)所具有的分形特征，可以表征不同深度土壤含水量變化過程的復(fù)雜程度，D越大說明變化過程越復(fù)雜［16］。

1.2 Hurst指數(shù)的意義

自相似性和長程依賴性已成為時間序列的重要特征［17-18］，定量描述長程依賴性的主要方法就是估計Hurst指數(shù)。Hurst指數(shù)是刻畫時間序列的重要參數(shù)，依據(jù)分形理論，時間序列的Hurst指數(shù)介于0～1之間。以0.5為間隔，時間序列在不同的區(qū)間表現(xiàn)出不同的特性［19］:(a)如果H=0.5，則該時間序列為獨立同分布的隨機序列，時間序列在各個時間尺度上都是相互獨立的;如果H≠0.5，則可判斷該時間序列為非隨機序列。(b)如果0.5＜H＜1，表明該時間序列是持續(xù)性的，或稱其具有長期記憶性(也就是說，若時間序列在過去是上升(下降)的，現(xiàn)在或?qū)頃r期將繼續(xù)上升(下降))。H越接近1，時間序列的持續(xù)性程度就越強，反之就越弱。如果0＜H＜0.5，表示時間序列存在反持續(xù)性，也就是說，該時間序列在過去是下降趨勢就意味著將來是上升趨勢，反持續(xù)行為的強度依賴于H接近0的程度。

2 實例應(yīng)用

2.1 基本數(shù)據(jù)

本例數(shù)據(jù)引自文獻［20］。該文獻給出了2002—2008年中國科學(xué)院欒城農(nóng)業(yè)生態(tài)系統(tǒng)試驗站內(nèi)試驗農(nóng)田15個不同深度層次的土壤含水量實測數(shù)據(jù)、觀測點的地理位置及觀測方法。為了探討不同深度土壤含水量的動態(tài)變化特征，筆者選用了深度為10 cm、80 cm和150cm共3層土壤含水量的實測數(shù)據(jù)進行計算與分析。深度10 cm代表土壤表層，深度80 cm代表作物根系的主要發(fā)育區(qū)，深度150 cm代表作物根系發(fā)育區(qū)底部。

2.2 土壤含水量的動態(tài)分析

將土壤含水量看成是隨機變量，即t～X(t)，利用Hurst經(jīng)驗公式(即式(8))對其進行計算，并將計算結(jié)果表示在雙對數(shù)坐標中，再利用最小二乘法對點進行擬合，得到的直線斜率即是H的估計值。序列的D可以通過式(9)間接求得。該時間序列長度N=379，為了研究不同層次土壤含水量時間序列的動態(tài)變化特征，分別取不同的時間尺度N=100、200、379，通過Matlab軟件將所得的一系列l(wèi)gτ、lgF(τ)在雙對數(shù)坐標上進行擬合，從而確定直線的斜率，結(jié)果見圖1、圖2和圖3。

圖1 深度為10 cm不同時間尺度土壤含水量R/S分析結(jié)果Fig.1 Results of R/S analysis of soil moisture on different time scales at depth of 10 cm

表1給出了深度為10 cm、80 cm和150 cm土壤含水量時間序列在不同時間尺度上的H和D。從表1可以看出，在不同的時間尺度下，H均大于0.5，由此可以判斷當?shù)夭煌疃韧寥篮康臅r間序列是非隨機時間序列，時間序列存在顯著的長程依賴性，而且該時間序列是持續(xù)性的，即同一深度含水量時間序列總體呈現(xiàn)下降或上升的趨勢。從圖4對3個土壤深度層次變化曲線擬合的直線可以看出，深度為10 cm處含水量呈現(xiàn)上升趨勢，深度為80 cm和150 cm處含水量呈現(xiàn)下降趨勢。產(chǎn)生這一現(xiàn)象的原因可能是隨著當?shù)亟?jīng)濟發(fā)展速度加快，過量開采地下水，導(dǎo)致地下水持續(xù)下降，引起深層土壤含水量呈現(xiàn)出減小的趨勢。然而，由于受到降雨補充或灌溉蒸發(fā)的消耗，表層土壤含水量的變化幅度較大，呈現(xiàn)上升的趨勢，這也與圖4所示的結(jié)果基本相符。

深度為10cm、80 cm和150 cm土壤含水量時間序列的平均分形維數(shù)D分別為1.262、1.123、1.047，D隨著深度的增加越來越小，由此可知，隨著土壤深度的增加，土壤含水量波動幅度越來越小，噪聲更少。產(chǎn)生這一現(xiàn)象的原因可能是表層土壤含水量的大小受降水和蒸發(fā)的交互影響較強，土壤含水量變化幅度較大，波動程度較高，故分形維數(shù)也較大;對深層土壤含水量而言，由于受降雨與蒸發(fā)的影響較弱，土壤含水量變化幅度較小，波動程度較低，故分形維數(shù)減小。另一方面，土壤含水量時間序列的長程依賴性隨著深度的增加也越來越弱，并得出當?shù)赝寥篮康臅r間序列為分形時間序列的結(jié)論。

圖2 深度為80 cm不同時間尺度土壤含水量R/S分析結(jié)果Fig.2 Results of R/S analysis of soil moisture on different time scales at depth of 80 cm

圖3 深度為150 cm不同時間尺度土壤含水量R/S分析結(jié)果Fig.3 Results of R/S analysis of soil moisture on different time scales at depth of 150 cm

表1 不同土壤深度在不同時間尺度的相關(guān)參數(shù)Table 1 Parameters at different depths on different time scales

3 結(jié) 語

借助分形理論中的R/S分析法，對欒城縣2002—2008年深度為10 cm、80 cm和150 cm土壤含水量時間序列進行動態(tài)分析，得到如下結(jié)論:

a.在深度為10cm、80cm和150cm處土壤含水量時間序列的Hurst指數(shù)H分別為0.733、0.877、0.953，均大于0.5。不同層次的土壤含水量表現(xiàn)出長程的依賴性，也就是說，同一土壤層含水量在整個序列上總體呈現(xiàn)下降或上升的趨勢，而且該序列是非隨機序列。

圖4 不同深度土壤含水量的變化趨勢Fig.4 Trend of soil moisture dynamics at different depths

b.在10 cm、80 cm和150 cm深度，土壤含水量序列的平均分形維數(shù)D分別為1.262、1.123、1.047，其分形維數(shù)逐漸變小，說明土壤的含水量波動幅度也越來越小，噪聲更少，其長程依賴性也越來越弱，另一方面也可以判斷出欒城區(qū)土壤含水量的時間序列為分形時間序列。

［1］李松敏，王仰仁.考慮下界面水分通量的棉田土壤水分模擬研究［J］.人民黃河，2010，32(11):89-91.(LI Songmin，WANG Yangren.Consider the interface moisture flux simulation study of soil moisture in cotton fields［J］.Yellow River，2010，32(11):89-91.(in Chinese))

［2］劉政，孫艷，李保成，等.SPSS逐步回歸分析法構(gòu)建北疆棉蚜年發(fā)生程度的預(yù)報模型［J］.安徽農(nóng)業(yè)科學(xué)，2009，37(6):2815-2817.(LIU Zheng，SUN Yan，LI Baocheng，et al.Stablishment of prediction model on year occurrence degree ofaphidsgossipier in North Xinjiang with SPSSstepwise regression analysis［J］.Journal of Anhui Agricultural Sciences，2009，37(6):2815-2817.(in Chinese))

［3］于福榮，王友賀.頻譜分析法在鄭州市降水量預(yù)報中的應(yīng)用［J］.人民長江，2012，43(10):50-53.(YU Furong，WANG Youhe.Application of frequency spectrum analysis in annual precipitation forecast for Zhengzhou City ［J］.Yangtze River，2012，43(10):50-53.(in Chinese))

［4］劉宇峰，孫虎，原志華.基于小波分析的汾河河津站徑流與輸沙的多時間尺度特征［J］.地理科學(xué)，2012，32(6):764-770.(LIU Yufeng，SUN Hu，YUAN Zhihua.Multiple time scale features of run off and sediment discharge on wavelet analysis at Hejin Station of Fenhe River Basin［J］.Scientia Geographical Sinica，2012，32(6):764-770.(in Chinese))

［5］桑燕芳，王中根，劉昌明.水文時間序列分析方法研究進展［J］.地理科學(xué)進展，2013，32(1):20-30.(SANG Yanfang，WANG Zhonggeng，LIU Changming.Hydrological time series analysis method is reviewed［J］.Progress in Geography，2013，32(1):20-30.(in Chinese))

［6］蔡江東.重標極差分析法預(yù)測徑流量問題探究［J］.人民黃河，2012，34(7):25-27.(CAI Jiangdong.Approach to the problems of rescaled range analysis in predicting the runoff［J］.Yellow River.2012，34(7):25-27.(in Chinese))

［7］HURST H E.Long-term storage capacity of reservoirs［J］.Transactions of the American Society of Civil Engineers，1951，116:770-808.

［8］MANDELBROT B B.Fractional brownian motion，fractional noise and Application［J］.SIAM Rev，1968(10):425-437.

［9］WALLISJ R，MATALAS N C.Small sample properties ofHandK-estimations of the Hurst coefficient［J］.Water Resource Research，1970，6(6):1583-1594.

［10］劉志，姚運生.鄂西及鄰區(qū)地震時空叢集的物理分形及R/S分析［J］.大地測量與地球動力學(xué)，2012，32(1):13-16.(LIU Zhi，YAO Yunsheng.Western Hubei and its adjacent area earthquake cluster of time and space physics and fractal R/Sanalysis［J］.Journal of Geodesy and Geodynamics，2012，32(1):13-16.(in Chinese))

［11］吳克金，汪波.基于多重分形與混沌理論的金融市場研究［D］.天津:天津大學(xué)，2005.

［12］王斌.基于分形理論的我國股票市場價格波動研究［D］.西安:西安電子科技大學(xué)，2010.

［13］屈世顯，張建華.復(fù)雜系統(tǒng)的分形理論與應(yīng)用［M］.西安:陜西人民出版社，1995:344-358.

［14］黃勇，周志芳，王錦國，等.R/S分析在地下水動態(tài)分析中的應(yīng)用［J］.河海大學(xué)學(xué)報:自然科學(xué)版，2002，30(1):83-87.(HUANGYong，ZHOU Zhifang，WANGJinguo，et al.R/Sanalysis in the application of the groundwater dynamic analysis［J］.Journal of Hohai University:Natural Sciences，2002，30(1):83-87.(in Chinese))

［15］趙晶，王乃昂.近50年來蘭州城市氣候變化的 R/S分析［J］.干旱區(qū)地理，2002，25(1):90-95.(ZHAO Jing，WANG Naiang.Nearly 50 years of Lanzhou City of R/Sanalysis of climate change［J］.Arid Land Geography，2002，25(1):90-95.(in Chinese))

［16］王小軍，張強，古璇清.基于分形理論的灌溉水有效利用系數(shù)空間尺度變異［J］.地理學(xué)報，2012，67(9):1201-1212.(WANG Xiaojun，ZHANG Qiang，GU Xuanqing.Fractal-based effective utilization coefficient of irrigation water space scale variability［J］.Acta Geographica Sinica，2012，67(9):1201-1212.(in Chinese))

［17］MESA O J，POVEDA G.The Hurst effect:the scale of fluctuation approach［J］.Water Resour Res，1993，29(12):3995-4002.

［18］SHAO Quanxi，LI Ming.A new trend analysis for seasonal time series with consideration of data dependence［J］.Journal of Hydrology，2011，17:104-112.

［19］MANDELBROT B B，TOQUE M S.Robust R/S analysis of long-rain serial correlation［R］.Manila:Bulletin of the International Statistical Institute，1979:69-99.

［20］孫鴻烈，胡春勝.農(nóng)田生態(tài)系統(tǒng)卷河北欒城站(1988—2008)［M］.北京:中國農(nóng)業(yè)出版社，2011:120-128.