R/S分析法在石家庄市栾城区土壤含水量动态分析中的应用

高雄飞，刘元会，郭建青，王媛英

(1.长安大学理学院，陕西西安 710064;2.长安大学环境科学与工程学院，陕西西安 710051)

土壤含水量是土壤的基本物理性质指标之一。准确掌握土壤水分的动态变化是确定灌溉制度并适时进行灌溉的重要依据［1］。对土壤含水量时间序列的时空演变过程进行研究在土壤含水量预测预报中具有重要的作用。目前，处理这类长期观测资料的方法主要有回归分析法［2］、频谱分析法［3］、随机微分方程法［4］等。这些方法的结果能够在一定程度上反映土壤动态的内在物理机制，但是，它们均未能很好地揭示土壤含水量动态过程的分形特征［5］。土壤含水量由于受到降雨、蒸发、日照时间等因素的影响，具有高度的非线性，对于这类非线性的具有统计特性的数据系列，采用重标极差分析法［6］(rescaled range analysis，简称R/S分析法)，可以较好地揭示土壤含水量动态变化过程的内在规律及其分形特征，从而有效预报或控制其变化规律。

Hurst［7］在总结尼罗河的多年水文观测资料时提出了R/S分析法，其后Mandelbrot［8］对该方法在理论上进行了补充和完善。Wallis等［9］曾利用此法对河流流量、泥浆沉积量、树木年轮、降雨量等自然现象进行研究。该方法在其他领域(如物理、生物学等)也得到了广泛的应用［10］。本文以栾城水文站2002—2008年深度为10 cm、80 cm和150 cm土壤含水量时间序列为例，采用R/S分析法对土壤含水量的动态过程进行分析，进而借助分形理论对土壤含水量的分形特征进行定量描述，期望揭示当地农田土壤水分的动态变化特征。

1 R/S分析法

1.1 R/S分析的基本原理

将时间序列{Xi=1，2，3，…，N}划分为m个时间长度为n的独立时间序列，其中m=N/n(m为整数)，子序列记为Da，其中下标a表示第a个子序列［11］。

a.每个时间长度为n的子序列Da的均值计算公式为

式中:X¯——子序列平均值;j——子序列下标。

b.每个子序列Da的累积偏差为

式中:Xj，a——子序列的第j个观测值;X¯a——第a个子序列的平均值。

c.每个子序列Da的极差为

d.每个子序列Da的标准差为

e.每个子序列的重标极差定义为

f.对每个子序列重复第b步到第e步的计算，得到一个重标极差序列，计算该序列的均值:

g.将n增加到下一个更大的值，重复第a步至第f步的计算，直到，得到一个序列。假设该序列累积偏差的极差与观测数的一个待定幂次方成正比［12］，即

其中

式中:τ——时间步长;C——待定常数，可以理解为反映自然现象固有特性的常数;H——时间序列的Hurst指数。

对式(7)两边取对数，有

式(8)利用最小二乘法进行求解，其斜率H即为Hurst指数。

H表示了时间变换的标度因子，因此时间序列具有自仿射性，经过严格的数学推导［13］，时间序列的分形维数与H的关系［14-15］如下:

式中:D——分形维数，D定量地描述了土壤含水量动态所具有的分形特征，可以表征不同深度土壤含水量变化过程的复杂程度，D越大说明变化过程越复杂［16］。

1.2 Hurst指数的意义

自相似性和长程依赖性已成为时间序列的重要特征［17-18］，定量描述长程依赖性的主要方法就是估计Hurst指数。Hurst指数是刻画时间序列的重要参数，依据分形理论，时间序列的Hurst指数介于0～1之间。以0.5为间隔，时间序列在不同的区间表现出不同的特性［19］:(a)如果H=0.5，则该时间序列为独立同分布的随机序列，时间序列在各个时间尺度上都是相互独立的;如果H≠0.5，则可判断该时间序列为非随机序列。(b)如果0.5＜H＜1，表明该时间序列是持续性的，或称其具有长期记忆性(也就是说，若时间序列在过去是上升(下降)的，现在或将来时期将继续上升(下降))。H越接近1，时间序列的持续性程度就越强，反之就越弱。如果0＜H＜0.5，表示时间序列存在反持续性，也就是说，该时间序列在过去是下降趋势就意味着将来是上升趋势，反持续行为的强度依赖于H接近0的程度。

2 实例应用

2.1 基本数据

本例数据引自文献［20］。该文献给出了2002—2008年中国科学院栾城农业生态系统试验站内试验农田15个不同深度层次的土壤含水量实测数据、观测点的地理位置及观测方法。为了探讨不同深度土壤含水量的动态变化特征，笔者选用了深度为10 cm、80 cm和150cm共3层土壤含水量的实测数据进行计算与分析。深度10 cm代表土壤表层，深度80 cm代表作物根系的主要发育区，深度150 cm代表作物根系发育区底部。

2.2 土壤含水量的动态分析

将土壤含水量看成是随机变量，即t～X(t)，利用Hurst经验公式(即式(8))对其进行计算，并将计算结果表示在双对数坐标中，再利用最小二乘法对点进行拟合，得到的直线斜率即是H的估计值。序列的D可以通过式(9)间接求得。该时间序列长度N=379，为了研究不同层次土壤含水量时间序列的动态变化特征，分别取不同的时间尺度N=100、200、379，通过Matlab软件将所得的一系列lgτ、lgF(τ)在双对数坐标上进行拟合，从而确定直线的斜率，结果见图1、图2和图3。

图1 深度为10 cm不同时间尺度土壤含水量R/S分析结果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，由此可以判断当地不同深度土壤含水量的时间序列是非随机时间序列，时间序列存在显著的长程依赖性，而且该时间序列是持续性的，即同一深度含水量时间序列总体呈现下降或上升的趋势。从图4对3个土壤深度层次变化曲线拟合的直线可以看出，深度为10 cm处含水量呈现上升趋势，深度为80 cm和150 cm处含水量呈现下降趋势。产生这一现象的原因可能是随着当地经济发展速度加快，过量开采地下水，导致地下水持续下降，引起深层土壤含水量呈现出减小的趋势。然而，由于受到降雨补充或灌溉蒸发的消耗，表层土壤含水量的变化幅度较大，呈现上升的趋势，这也与图4所示的结果基本相符。

深度为10cm、80 cm和150 cm土壤含水量时间序列的平均分形维数D分别为1.262、1.123、1.047，D随着深度的增加越来越小，由此可知，随着土壤深度的增加，土壤含水量波动幅度越来越小，噪声更少。产生这一现象的原因可能是表层土壤含水量的大小受降水和蒸发的交互影响较强，土壤含水量变化幅度较大，波动程度较高，故分形维数也较大;对深层土壤含水量而言，由于受降雨与蒸发的影响较弱，土壤含水量变化幅度较小，波动程度较低，故分形维数减小。另一方面，土壤含水量时间序列的长程依赖性随着深度的增加也越来越弱，并得出当地土壤含水量的时间序列为分形时间序列的结论。

图2 深度为80 cm不同时间尺度土壤含水量R/S分析结果Fig.2 Results of R/S analysis of soil moisture on different time scales at depth of 80 cm

图3 深度为150 cm不同时间尺度土壤含水量R/S分析结果Fig.3 Results of R/S analysis of soil moisture on different time scales at depth of 150 cm

表1 不同土壤深度在不同时间尺度的相关参数Table 1 Parameters at different depths on different time scales

3 结 语

借助分形理论中的R/S分析法，对栾城县2002—2008年深度为10 cm、80 cm和150 cm土壤含水量时间序列进行动态分析，得到如下结论:

a.在深度为10cm、80cm和150cm处土壤含水量时间序列的Hurst指数H分别为0.733、0.877、0.953，均大于0.5。不同层次的土壤含水量表现出长程的依赖性，也就是说，同一土壤层含水量在整个序列上总体呈现下降或上升的趋势，而且该序列是非随机序列。

图4 不同深度土壤含水量的变化趋势Fig.4 Trend of soil moisture dynamics at different depths

b.在10 cm、80 cm和150 cm深度，土壤含水量序列的平均分形维数D分别为1.262、1.123、1.047，其分形维数逐渐变小，说明土壤的含水量波动幅度也越来越小，噪声更少，其长程依赖性也越来越弱，另一方面也可以判断出栾城区土壤含水量的时间序列为分形时间序列。

［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逐步回归分析法构建北疆棉蚜年发生程度的预报模型［J］.安徽农业科学，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］于福荣，王友贺.频谱分析法在郑州市降水量预报中的应用［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］.地理科学，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］.地理科学进展，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］蔡江东.重标极差分析法预测径流量问题探究［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］刘志，姚运生.鄂西及邻区地震时空丛集的物理分形及R/S分析［J］.大地测量与地球动力学，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］.天津:天津大学，2005.

［12］王斌.基于分形理论的我国股票市场价格波动研究［D］.西安:西安电子科技大学，2010.

［13］屈世显，张建华.复杂系统的分形理论与应用［M］.西安:陕西人民出版社，1995:344-358.

［14］黄勇，周志芳，王锦国，等.R/S分析在地下水动态分析中的应用［J］.河海大学学报:自然科学版，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］.干旱区地理，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］王小军，张强，古璇清.基于分形理论的灌溉水有效利用系数空间尺度变异［J］.地理学报，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］孙鸿烈，胡春胜.农田生态系统卷河北栾城站(1988—2008)［M］.北京:中国农业出版社，2011:120-128.