应用时间序列模型分析长武塬区地下水水位的变化特征*

张盼,刘文兆,2

(1.西北农林科技大学资源环境学院,陕西杨陵712100;2.中国科学院水利部水土保持研究所,陕西杨陵712100)

地下水是水资源的重要组成部分,目前我国大部分地区的地下水形势不容乐观。长武塬区地下水资源贫乏,长期以来存在的盲目开发和超量开采造成了地下水位持续下降,进而造成水资源供需矛盾加剧。因此,研究地下水动态变化规律,预测未来变化趋势,对水资源利用管理具有重要的意义。

时间序列分析是应用比较广泛的一种科学分析方法。由于这种方法比较简单,而且可以定时定量地推测事物的未来发展,所以对预测、决策具有较高的实用价值[1]。近20年来,时序分析理论已应用于地下水资源评价、预报和管理之中。常用于地下水位预报的时序模型有自回归模型(AR),滑动平均模型(MA)、自回归滑动平均模型(AR-MA)。近年来随着计算机技术的发展和广泛应用又发展了多种模型,如求和自回归滑动平均模型(ARIMA)、乘积型季节性模型、组合模型、混合模型以及门限模型、混沌时间序列模型等[2-5]。这些模型应用灵活方便且具有较高的模拟精度,给大区域地下水动态预报分析带来了极大的便利。

1 研究区概况

1.1 自然地理概况

研究区所在的长武县位于陕西省咸阳市辖区西北角,西、北、南三面分别与甘肃省泾川县、宁县、正县、灵台县接壤,全县土地面积为583 km2。地势由西向东倾斜,海拔高度847～1 274 m,该地区系典型的黄土高塬沟壑区[6],属暖温带半湿润大陆性季风气候,降水集中在7-9月,占全年降水总量的55%以上。全县有泾河、黑河、南河3条河流。泾河沿县北向东后折向东南流去,境内流长56 km,年平均流量42.2m3。黑河、南河年平均流量为10m3/s以下。三条河流将长武县分割为巨家塬、枣元塬和长武塬三部分[7-8]。其中长武塬土地面积为297.4 km2,塬面积为101.6 km2,塬面和沟坡面积分别占34.2%和65.8%,地下水埋深30～80 m[7,9]。

1.2 水文地质概况

长武县属祁、吕、贺山子型构造的联合复合区。基岩大部分被黄土覆盖,基底是巨厚的中生界沉积砂页岩,平铺缓斜,褶皱轻微,自下而上依次为三迭系、侏罗系、白垩系[6]。研究区主要由第四系新老黄土组成,自上而下为马兰黄土、离石黄土、午城黄土。塬区潜水含水层是以风积黄土状亚砂土为主的孔隙-裂隙含水层,有良好的储水空间;且每层古土壤均含有钙质结核,含水层底部钙质结核胶结致密,作层状分布,形成较好的隔水层,构成一定储水条件,是黄土塬区主要的潜水含水层,含水深度为20～100 m。在长武北塬区中心部位,塬面平坦开阔,地下水埋深25～60 m,含水层厚度 60～80 m,为中等富水潜水层;在塬区周边,由于沟谷对地下水的疏干作用,地下水埋深则在60～85m以上,含水层厚度也相应变薄。在南部塬区,由于地形破碎,地下水赋水条件很差,地下水埋深较深,人畜饮用水困难。

塬区地下水埋深多为40～100 m,四周均被切割至基岩的沟谷所包围。黄土中地下水位高于沟谷中地表水流,更高于基岩中的承压水,因此无地下侧流补给,降水垂直入渗是塬面地下水补给几乎唯一的来源[7,10]。

2 材料与方法

2.1 数据获取

选取研究区当地水利局地下水监管站3口水位埋深监测井长期观测资料进行分析。其中,551井位于研究区西段洪家乡风口村,552井位于研究区中段地掌乡代苓村,555井位于研究区东段彭公乡丰头村。每5天(即每月1日,6日,11日,16日,21日,26日)用测绳量测监测井水位。

2.2 时间序列分析法

时间序列,也叫时间数列、历史复数或动态数列。它是将某种统计指标的数值,按时间先后顺序排到所形成的数列。时间序列预测法就是通过编制和分析时间序列,根据时间序列所反映出来的发展过程、方向和趋势,进行类推或延伸,借以预测下一段时间或以后若干年内可能达到的水平。在生产和科学研究中,对某一个或一组变量 x(t)进行观察测量,将在一系列时刻t1,t2,…,tn(t为自变量且t1＜t2＜…＜tn)所得到的离散数字组成序列集合x(t1),x(t2),…,x(tn),称之为时间序列,这种有时间意义的序列也称为动态数据。这样的动态数据在自然、经济及社会等领域都很常见。

时间序列分析是根据系统观测得到的时间序列数据,通过曲线拟合和参数估计来建立数学模型的理论和方法。它把系统看作一个暗箱,不考虑外界因素的影响,假设预测对象的变化仅与时间有关,根据客观事物发展变化过程中的内在延续性进行预测。这种内在延续性反映了外部影响因素综合作用下对象的变化过程。时间序列预测过程只依赖于历史观测数据及其数据模式,从而使预测研究更为直接和简洁。时间序列分析主要用于:系统描述、系统分析、预测未来、决策和控制。

在对地下水系统的研究中 ,把实体系统作为直接研究对象极其困难。一般通过系统模型来描述实体系统的特性及其变化规律[11]。应用时间序列(简称“时序”)分析法进行水位预报,正是将地下水系统视为“黑箱”或“灰箱”,根据地下水位动态观测资料,提取和分析历史资料本身所蕴涵的信息,找出其规律,并利用这些规律,达到预报未来的目的,无须再进行专门的试验来获取其它参数,这给大区域地下水动态预报分析带来了极大的便利,并且,该方法易于掌握,计算工作量小,易于应用推广[12]。

如今,SAS的丰富数据采集、数据管理、数据分析和信息展现的能力,使之成为决策支持的最好工具。其中,SAS/ETS提供的基本时间序列预测方法包括:逐步自回归模型、指数平滑模型和季节性模型[13]。用SAS软件所提供的这三种方法对研究区地下水位建立模型,经过比较,选取了适合研究区的最优模型——指数平滑模型。

2.2.1 时间序列模型基本组成 由于地下水水位数据具有趋势性,周期性和一个静态分量,一个非稳定的时间序列模型可以适用于地下水位数据的分析。根据地下水水位数据的特征,一个地下水非稳定的地下水时间序列模型可以用以下叠加形式表示:

式中:H(t)——地下水位埋深;R(t)——趋势项;P(t)——周期项 ;ε(t)——随机分量 。

从已知序列(观测值)H(t)(t=1,2,3,…,n)中提取各分量。提取的顺序为:趋势分量R(t)、周期分量P(t)和随机分量ε(t)。各分量的数学模型建立后,再将其线性叠加,就得到了地下水水位埋深预测模型H(t)。其中,t为地下水水位监测周期,一个单位的t表示一个地下水水位监测周期,取决地下水水位监测频率。

(1)趋势分量R(t)的确定。对于趋势分量R(t)可以用多项式逼近,即可得到估计值R'(t)。

采用多元回归方法确定待定系数b0,b1,b2,b3,…,bk和阶数k。同时需要计算趋势曲线拟合的相关系数R。

(2)周期分量P(t)的确定。趋势函数确定后,对去除趋势分量后的部分p(t)进行周期项分析,即:

采取谐波分析法进行周期分量的分析提取。对序列P(t),可以用L个谐波叠加的形式表示其估计值为

式中:P'(t)——序列 P(t)的估计值;L——谐波个数,取n/2的整数部分;k——谐波序号,k=1,2,…,L;Tk=n/k,也就是第k个谐波的周期;ak,bk——傅立叶系数,计算式为

(3)随机分量ε(t)的确定。去除趋势分量和周期分量的部分既是随机分量部分ε(t)。

式中:ε(t)——一平稳随机系列项,可直接对其用自回归模型求解。其自回归模型为

式中:p——模型阶数;φi——模型自回归系数,i=0,1,2,…,p。

对某一阶数的自回归模型,类似多元回归计算可求得自回归系数φi。本文采用了AIC准则[4]来确定模型的阶数,即:

使上式值最小时所对应的 p值为最佳阶数。

式中:n——序列数据总个数;б2——式(7)阶数为 p时残差的方差。(4)时间序列模型。将上述趋势分量,周期分量,随机分量线性叠加,即可得到地下水位的总预测模型:

2.2.2 模型的检验 应用SAS软件自带的时间序列模型对研究区地下水位进行模拟后,对其精度进行检验,进而确定其是否合理。本文采用后验差方法[3,6]对模型进行精度检验。通过后验差比值C和小误差频率P来检验(见公式10)。

式中:S1——原始数据(实测值)方差的开方;S2——残差(绝对误差)的标准差;ξ(k)——残差;¯ξ— —残差均值。

上述两个指标C,P把预测的等级划分为四等,如表1所示,通常用其作为检验预测模型的精度。

表1 预报精度评价标准

如果P,C值都在允许范围内,则模型可用于计算预报值,否则需要对模型进行检查、分析和重新调参。

图1 研究区地下水位埋深概况

2.3 所选监测井水位埋深概况

分别选取551井,552井和555井1976-1993年,1976-2008年和1978-2003年的水位观测资料进行分析。(其中552井监测情况最好,1976年至今资料保存完备;551井的监测资料自1993年后部分丢失,故只选连续较好的前18年监测资料进行分析;555井由于干枯而于2003年停测,期间1996年资料丢失。)

研究区三口井水位埋深情况见图2。由图可以看出这三口监测井中552井水位较浅水位埋深不超过35 m,555井水位较深,水位埋深大于50 m,551井水位埋深则处于35～40 m。

3 时间序列分析在研究区地下水位预报中的应用

3.1 模型的应用

利用SAS/ETS中的指数平滑模型对研究区552井(1976-2008年)、551井(1976-1992年)及555井(1978-2002年)月平均地下水位数据进行模拟,拟合结果见图2。并分别得到观测井552井2009年各月(1-8月)水位埋深预测值,551井1993年各月(1-9月)水位埋深预测值,以及555井2003年各月水位埋深预测值,与实测值进行对比,用以评价模型的预测效果。其决定系数分别为:R2=0.999(552井),R2=0.962(551井),R2=0.997(555井)。拟合结果及误差分析见表3,预测结果及误差分析见表4。

图2 模型计算值与实测值拟合曲线

由552监测井2003-2008年的计算值与实测值的拟合曲线(图2)可以更清楚地观察模型对该井的计算值与模拟值的拟合情况。表1可明显看出,该模型的拟合精度很高,其最大绝对误差(残差)为0.06 m,最小绝对误差为5.67×10-5m。相对误差大多小于0.1%,个别稍微大点的值也不超过1%。由该模型对研究区551监测井和555监测井的模拟预测拟合曲线图(图2)可以看出,两口井的水位也是逐年下降的。

3.2 模型的检验

运用公式(10)对模型的预测精度进行计算得到P值和C值(见表2),对照表1可知,SAS软件自带的时间序列模型对研究区这三口井的预测精度均达到理想的预测等级。

同时,由表3可以看出,指数平滑模型对于551井的模拟预测效果较理想,拟合精度好。其绝对误差最大值为0.14 m,最小值可达0.001m。相对误差一般都小于0.4%,最大也不超过1%。该模型对555井的井水水位拟合程度同样较高。其中,绝对误差(残差)最大为0.03 m,最小甚至不足0.001 m。相对误差大多数在0.55%之下,个别较大的值也不超过1%。

表2 精度检验结果表

在把已建立的模拟模型用于预报前,也要进行精度检验,本文用未参加建模的水位资料(551井的1993年,552井的2009年,555井的2003年)进行后验预测检验。可得到551井1993年、552井2009年以及555井2003年各月水位埋深预测值。其后验预测结果见表4。

表4中可以看出,552井预测的绝对误差(残差)最大值为0.12 m,绝对误差的最小值可达1.23×10-3m。相对误差一般都不超过1%。该模型对551井的预测精度较之552井稍好些。其绝对误差(残差)最大值为0.13 m,最小值可达9.59×10-4m。相对误差均不超过0.4%。对555井的预测其绝对误差(残差)最大值为0.15 m,最小值达到3.38×10-5m。其相对误差最小值为6.39×10-5%,最大值不超过0.3%。可以看出,该模型对研究区监测井的拟合及预测精度都较高。

表3 模型计算结果与误差分析

表4 模型模拟结果与误差分析

3.3 分析与讨论

时间序列分析是把系统看作一个暗箱,不考虑外界因素的影响,假设预测对象的变化仅与时间有关,根据客观事物发展变化过程中的内在延续性进行预测。但是利用地下水系统过去的资料进行建模是通过系统过去的规律预测未来,当系统的某些因素发生较大变化时,模型预报误差将增大,因此建模时需要考虑到系统因素的变化。影响研究区地下水变化的主要因素如下。

(1)黄土区地下水埋深大部分为40～100m,且降水几乎是唯一的补给源[7,14-15]。研究区地下水补给来源为当地降水,降水一方面转化为土壤水缓慢补给给地下水,年补给量约为10mm;另一方面通过降雨形成径流,在塬面的黄土碟、黄土胡同或者居民开挖的涝池等负地形地段汇集后补给地下水,补给滞后时间多为30～40 d,也可能有较短的为10 d左右,受降雨量、降雨强度以及补给通道湿度等因素影响[7]。研究区50 a中(利用长武县气象局1957-2006年的逐日降水数据分析),丰水、平水和枯水年分别为18 a、12 a和20 a,丰水年和枯水年数比较接近[16]。且降水年际分布不均,年平均降水量为578.7mm(1957-2006年),年最大降水量为954.3 mm(2003年),年最小降水量为 296.0 mm(1995年),相差657.3 mm,表现出年际变化大的特点。

(2)研究区塬区的机井数量由1993年101眼增加到2008年181眼,由于研究区塬面积有限,以后塬面上的机井数量不再会有较大变化。与此同时,地下水源开采量由1999年的1.84×106m3增加到2008年的3.31×106m3。1998年至2007年10 a间年均降水量为589.39 mm,除2004年年降水量(493.7 mm)不足500mm外,其余年份年降水量均超过500mm,且2003年达到954.3 mm。而研究区地下水位埋深却整体处于下降趋势,是由于人为开采量增加导致的。由于塬面上的机井数量今后不会发生较大变化结合当地部门对地下水资源保护的重视,日后将会有计划的对地下水进行开采。但是随着社会的发展,人口的增长,工农业用水依然在不断增长,故该因素为影响研究区地下水变化的又一重要因素。

上述这些因素在时序建模中没有涉及到,因此,未来某些因素发生较大变化时,模型预报误差会增大。所以,应及时把新得到的实测数据加到原序列中,实时更新模型,以提高预测的精度。

由上述分析可知,SAS/ETS中指数平滑模型可以用于研究区地下水位的研究。其拟合精度高,且需要的数据不多,获取也简单,可以考虑用于研究区地下水位动态变化规律的研究。

4 结论

(1)运用SAS软件自带的时间序列分析法中的指数平滑模型分别对长武塬区三口监测井建立地下水位埋深模型。结果表明,拟合精度和预测精度均较高。该模型较全面地反眏了地下水位动态变化规律,且计算简单,所需资料较少,易于获得,是一种较好的模拟预测模型,能很好地应用于分析地下水位动态的趋势性和周期性以及进行地下水埋深预测,为区域地下水的合理开发利用及水资源规划提供可靠依据。

(2)通过对趋势项的分析可知,552监测井32 a间水位埋深由25.3 m(水位标高1 163.12 m)上升到31.84 m(水位标高1 156.58m),即水位降幅达到6.54 m,年均下降为0.2 m。由图可知,1976-1995年间,该井水位下降缓慢,这20 a间的水位年均下降为0.1 m。1996-2008年间,呈明显下降趋势,年均下降0.4 m。尤其是近几年,地下水位下降幅度越来越大,从模拟结果可看出,如果不进行科学管理,控制开采,还有继续下降的趋势。

(3)通过对模型中的周期项分析可知,该区地下水位动态具有两个主要的周期成分。一个周期长度为一年,反映了一年内地下水位的季节性周期变化。另一个周期长度为7.4～10 a,并且随着时间的增长周期也逐渐增长,这与该地区气候(主要是降水量)多年的周期性变化情况大体一致。

(4)利用地下水系统过去的资料进行建模是通过了解过去的规律对未来进行预测,当系统的某些因素在未来时刻有较大变化时,模型预报误差将增大。因此,要及时把新得到的实测数据加到原序列中,实时更新模型,以提高预测的精度。

鉴于研究区特殊地形及气候的影响,地下水的开发利用在当地工业和生态用水中占据重要地位,有的地方甚至成为唯一的水源。因此建议当地相关部门加强对地下水位变化规律的研究,掌握预测水位未来变化趋势的方法,合理规划和开发利用地下水资源。

[1] 董殿伟,刘久荣,林沛,等.时间序列分析法在地下水水位预测中的应用[J].城市地质,2007,2(4):29-32.

[2] 丁涛,周惠成,黄健辉.混沌水文时间序列区间预测研究[J].水利学报,2004,35(12):15-20.

[3] 杨位钦,顾岚.时间序列分析与动态数据建模[M].北京:北京工业学院出版社,1986:8-14.

[4] 杨志霞.时间序列模型在深层地下水水位预测中的应用[J].河北工程技术高等专科学校学报,2000(3):34-38.

[5] 杨忠平,卢文喜,李平.时间序列模型在吉林西部地下水动态变化预测中的应用[J].水利学报,2005,36(12):1475-1479.

[6] 中国科学院水利部水土保持研究所.土地资源及生产力研究[M].北京:科学技术文献出版社,1991:4-6.

[7] 王锐.基于环境同位素的黄土塬区降水-土壤水-地下水转换关系研究[D].北京:中国科学院研究生院,2007.

[8] 长武县农业区划委员会.陕西省长武县农业资源调查和农业区划报告集[M].长武:长武县印刷厂,1985.

[9] Liu Wenzhao,Li Zhi,Wang Rui,et al.Dynamic change of groundw ater level in Changw u tableland region on the Loess Plateau of China[C]//Proceedings of 3rd Internationna lWorkshop on Yellow River Studies.Japan:Kyoto,2007.

[10] 咸阳市水政水资办,长武县水政水资办.陕西省长武县水资源评价及开发利用现状分析[M].咸阳:咸阳市地下水管理监测站印制,1993.

[11] 卢晓燕,施鑫源,眭克仁.地下水资源系统时序管理模型[J].工程勘察,1999(4):35-38.

[12] 杨金忠,蔡树英.一种区域地下水预报的时间序列分析方法[J].武汉水利电力大学学报,1996,29(4):6-10.

[13] 张小娟,蒋云钟,秦长海,等.时间序列分析在地下水位预报中的应用[J].南水北调与水利科技,2007,5(4):40-42.

[14] 杨文治,赵沛伦,张启元.不同湿度条件下土壤水分的蒸发性能和移动规律[J].土壤学报,1981,18(1):24-37.

[15] 李玉山.黄土区土壤水分循环特征及其对陆地水分循环的影响[J].生态学报,1983,3(2):24-37.

[16] 陈杰,刘文兆,王文龙,等.长武黄土高塬沟壑区降水及侵蚀性降雨特征[J].中国水土保持科学,2009,7(1):27-31.