基于Granger模型的黄土沟壑区典型人工林绿水流研究

符 冲, 宋孝玉, 李蓝君, 赵新凯, 李怀有

(1.西安理工大学 省部共建西北旱区生态水利国家重点实验室,陕西 西安 710048; 2.黄河水利委员会 西峰水土保持科学试验站, 甘肃 庆阳 745000)

绿水是通过蒸散过程而流向大气圈的水汽流，绿水资源中的绿水流即为实际蒸散发的总量[1]。蒸散发作为能量交换和水文循环的参与者，是土壤—植被—大气连续体(SPAC)系统中必不可少的一个环节[2]。蒸散发包括土壤蒸发和植物蒸腾，其中，土壤蒸发被认为是绿水资源中的低效绿水，植物蒸腾则被认为是绿水资源中的高效绿水[3]。黄土高原地区多为干旱或半干旱气候，干旱作为主要的气象灾害，长期严重制约该区域内农业生产和生态恢复[4]。同时，受气候特征和下垫面条件综合作用，水土流失也一直是该地区亟待解决的关键问题。针对此种情况，国家大力推广植树造林、退耕还林工程，并取得了显著的生态效益[3]，但在造林过程中，由于未切合实际的考虑到当地水资源以及树种耗水的匹配特点，亦出现了林木存活率低或树干呈现“半生半死”等不良现象，导致植树造林效率大大降低[5-7]。因此，有必要对黄土区典型人工林地绿水流变化及其影响因素进行研究，这对于充分利用绿水资源以及促进人工林地的有效管理等均具有重要作用[5-6]。目前，模拟绿水流变化的模型较多应用于大尺度上[8-10]，对小尺度方面关注较少，而且对干旱、半干旱地区高、低效绿水进行分离的相关报道较少。现有方法多基于水量平衡、能量平衡等原理发展得来，主要有Penman公式[11]、修正Penman公式[12]、Penman-Monteith公式[13]、Shuttleworth-Wallace模型[8]、双作物系数法[14]、互补相关模型等[15]。其中，Penman类模型对数据要求较为严格，除了常规气象资料外，还需要较多的植被参数，这些植被参数往往难以获取或者取值精度不高，使该类模型在资料缺乏区域难以推广应用；相对而言，基于互补相关原理的Granger模型仅需要气象资料即可计算出绿水流，应用相对简便，方便于小区域尺度的蒸散发变化和绿水研究。在绿水流动态性研究方面，目前研究对象多为小麦(Triticumaestivum)[16]、大豆(Glycinemax)、玉米(Zeamays)、水稻(Oryzasativa)[17]，对黄土沟壑区典型人工林刺槐(Robiniapseudoacacia)、侧柏(Platycladusorientalis)、油松(Pinustabuliformis)的研究较为缺乏；而且现有研究进行绿水流变化探讨时，时间尺度多以整年或者整月为主[18]，针对短时段内变化特征进行的研究也较少。因此，本文以黄土沟壑区典型小流域内人工刺槐、侧柏和油松林地为研究对象，在偏枯年份(2016年)和特枯年份(2017年)野外监测的基础上，首先对Granger模型在各样地内的适用性进行评价，然后结合实测的低效绿水，对各样地内高、低效绿水进行分离，基于分离结果，对无雨期各林地的高、低效绿水与太阳辐射、叶面积指数、土壤含水率的关系进行分析，以期为黄土沟壑区典型人工林地绿水资源的合理利用、植树造林生态效益的提升以及人工林地的系统管理等提供科学参考与指导。

1 材料与方法

1.1 研究区概况

南小河沟流域(35°41′—35°44′N，107°30′—107°37′E)位于蒲河下游、董志源西侧，流域内建有黄河水利委员会西峰水土保持科学试验站。流域面积36.3 km2，海拔1 050—1 423 m，沟壑密度为2.68 km/km2，沟道比降2.8‰。地形地貌主要由塬面、坡地、沟谷组成，具有典型的黄土沟壑区地貌特征。流域处于半干旱地带，多年平均降雨量为544.8 mm，多年平均气温为9.3 ℃，平均日最高气温为28.8 ℃，平均日最低气温为-22.6 ℃。流域内降雨量集中在5—9月，占全年的70%以上，并多以大雨、暴雨的形式发生。土壤类型主要以黑垆土、黄绵土为主，容重分布在1.25～1.45 g/cm3。

流域内无天然林分布，人工林主要是刺槐、侧柏和油松，草被群落主要是白羊草、艾蒿和苜蓿。经过实地考察，选择研究区内树龄基本一致，长势较为良好的刺槐、侧柏及油松林地作为试验样地,试验样地基本情况详见表1。

表1 试验样地基本情况

1.2 试验设计

试验于2016年4月15日至10月15日及2017年4月15日至10月15日在南小河沟流域各样地内进行，观测及测定项目主要包括土壤含水率、土壤田间持水量、土壤蒸发、植被叶面积指数(LAI)等。土壤含水率使用管式时域反射系统(Trime TDR， Germany)测量，测量间隔为3～5 d，测量深度为200 cm，深度0—100 cm时间隔深度为10 cm，100—200 cm时，间隔深度为20 cm。土壤田间持水量使用环刀法进行测定。土壤蒸发采用自制的微型蒸渗仪(深度20 cm，口径16 cm)进行，测量间隔为3～5 d；降雨及气象数据(太阳辐射、温度、湿度、2 m处风速等)利用WatchDog 2009系列气象观测系统进行监测。LAI则使用植物冠层分析系统(WinScanopy 2006， Canada)进行定点观测，观测间隔为6～7 d。

1.3 互补相关原理及Granger模型

互补相关理论是Bouchet[19]提出的一种计算区域绿水流的理论模型。他认为在充分湿润时，绿水流(ETa)是和潜在蒸散量(ETp)相等的，都为湿润环境的蒸散量(ETw)。但是随着水分减少，ETa会随之减小，而原来用于蒸散的能量也不会消失，被用来改变该区域温度、湿度等，即用来增加ETp，且其增加量应该和剩余的原来用于蒸散的能量相等。公式(1)表示充分湿润时ETa和ETp之间的关系：

ETa+ETp=2ETw

(1)

在此基础上，Granger等引入“相对蒸散发”的概念并重新定义了可能蒸散发和湿润环境蒸散发[20]，最终推导出估算绿水流的关系式，即：

(2)

式中：ETa为绿水流(mm)； Δ为温度—饱和水气压曲线斜率(mbar/℃)；Rn为太阳净辐射(W/m2) ；G为土壤热通量(W/m2)；R为相对蒸散发；γ为干湿计常数(mbar/℃)，取值为0.658(mbar/℃)；Ea为干燥力项(W/m2)，它们与干燥力D之间的关系式为：

(3)

(4)

式中：b1，b2为经验系数，Granger根据不同下垫面观测资料[20]，分析认为b1和b2取值分别为8.045,0.028。

1.4 水量平衡法

利用实测的刺槐、油松、侧柏林地的土壤含水率，根据水量平衡原理计算试验期间样地的绿水流：

ETa=W1-W2+P-I-R

(5)

式中：W1为前一次测量的土壤含水量；W2为后一次测量时的土壤含水量；P为前后两次测量期间的降雨量；I为截留量；R为径流量，该式中变量单位均为mm。

根据课题组在南小河沟流域进行的相关研究[3]，刺槐、侧柏及油松冠层截流比例大致保持在13%，20%以及18%。由于本研究试验样地和其试验样地在立地条件、植被密度、土壤性质等方面十分接近，故本研究直接采用其成果进行林地截流量的计算。同理，试验期内径流量采用李淼[21]所拟定的南小河沟流域林地产流方程进行计算。

1.5 模型评价标准

以2016年为率定期，将Granger模型计算的绿水流与水量平衡法计算结果进行对比，通过手工试错法调整模型参数，使模拟结果与实测结果接近。采用率定后的模型参数，以2017年为验证期对模型进行验证。率定与验证模型时，以纳什效率系数(NSE)[22]对模型模拟精度进行评价：

(6)

式中：Pi为Granger模型模拟的绿水流；Oi为试验测定的绿水流；Oavg为试验测定绿水流的平均值。NSE越接近于1，说明模拟值与实测值的拟合程度越好，即适用性越好。一般而言，NSE大于0.5则表示模型模拟精度良好，处于可接受的范围；NSE大于0.7表示模型模拟精度较高[23]。

1.6 绿水流分离

使用经过验证后的Granger模型计算绿水流，减去刺槐、侧柏、油松林地实测的低效绿水，即可计算得到刺槐、侧柏、油松林地的高效绿水。

2 结果与分析

2.1 模型参数率定

模型参数率定后的结果详见表2。参数率定后，模型模拟效果如图1—2所示，可以看出，率定期模型在刺槐、侧柏、油松林地内的NSE分别为0.72，0.70和0.67，在验证期则分别为0.56，0.52，0.58。整体来看，率定参数后的Granger模型对各样地内绿水流的模拟效果较好，可用于进行日尺度的绿水流计算。

表2 试验样地Granger模型绿水流模拟参数率定结果

图1 率定期(2016年)Granger模型模拟不同林地绿水流效果

图2 验证期(2017年)Granger模型模拟不同林地绿水流效果

2.2 样地尺度高低效绿水关系分析

利用多元非线性回归模型进行绿水流及其影响因子关系的拟合，即将低效绿水(E)或高效绿水(T)分别与太阳净辐射(Rn)、土壤含水率(θ)、叶面积指数(LAI)建立函数关系。建立方程时，不仅考虑各单因子对E和T影响，同时也考虑各因子之间的交互影响作用。方程具体形式为：

diθ×LAI+eiRn×LAI+

fiRn×θ+giRn×θ·LAI+hi

(7)

式中：ai,bi,ci,di,ei,fi,gi,hi,n1,n2,n3均为待定系数。考虑到E主要发生在表层土壤,因此,建立函数关系时θ分别采用0,10以及20 cm土层含水率分别进行拟合。对T而言,则采用0,10,20,30和40 cm土层的θ来分别建立函数关系。

经多次试算，最后确定n1=1，n2=2，n3=1。采用不同深度含水率建立的回归方程决定系数(表3—4)。可以看出，在其他因素确定的情况下，低效绿水对10 cm土壤含水量敏感，而高效绿水对20 cm土壤含水量敏感。采用以上两深度进行方程拟合时，低效绿水方程的决定系数在刺槐、侧柏及油松林地内为0.44，0.64，0.57；对高效绿水而言，则分别为0.62，0.57，0.63。

表3 各林地低效绿水方程决定系数R2

表4 各林地高效绿水方程决定系数R2

2.3 各因子对绿水流的影响分析

由于降雨量大小会影响各因子对绿水流的贡献度，因此，在对各影响因素进行分析前，需要确定2016和2017年的丰枯程度。根据气象科学研究院拟定的降雨量5级判别标准，本研究采用降雨量距平百分率(Pa)来进行水文年份的判别(表5)[24]。结合2016与2017年降雨量(397.6与217.8 mm)得到，研究区域2016年降雨量距平百分率Pa2016=-27%，2017年降雨量距平百分率Pa2017=-60%，因此，2016年是偏枯年份，2017年是特枯年份。

表5 水文年份判别标准

2.3.1 叶面积指数(LAI) 图3表示不同林地在不同年份无雨期内(一般为4 d)低效绿水(E)和高效绿水(T)随叶面积指数(LAI)变化。由图3可以看出，刺槐、侧柏、油松林地的E值在偏枯年(2016年)分别分布在在0.13～7.45，0.08～4.22 mm以及0.18～7.63 mm，日平均E值分别为1.86，1.91以及1.06 mm/d；在特枯年(2017年)，E值则分布在0.47～6.31 mm，0.36～7.36 mm以及0.01～9.88 mm，日均E分别为1.59，1.84，2.47 mm/d。对比分析可以看出，不同水文年刺槐林地低效绿水变化幅度较小，其次是油松，而侧柏变化最敏感，其最大E在不同水文年变化比例超过了70%。

对T而言，刺槐、侧柏、油松林地的T在偏枯年(2016年)分别分布在0.04～16.07，0.11～18.48 mm以及0.05～16.49 mm，日均T分别为4.01，4.62以及4.12 mm/d左右；在特枯年(2017年)，T则分布在0.63～14.85，0.27～9.59 mm以及0.09～10.44 mm，日均T为3.71，2.40，2.61 mm/d。可以看出，不同水文年对侧柏的高效绿水影响最大，次之油松，对刺槐的影响是最小的。

就绿水流而言，刺槐、侧柏、油松林地在偏枯年份和特枯年份的绿水流分别为181.9，193.1，178.3和159.4，152.3，179.2 mm；各样地不同水文年无雨期间绿水流变化表现为：侧柏(21.13%)>刺槐(12.37%)>油松(0.50%)，结合高低效绿水的变化情况可知侧柏是对水分最敏感的树种，油松是不同水文年绿水流变化最小的树种。

纵观三样地叶面积指数的变化可知，随着叶面积指数的增大，刺槐林地E在偏枯年份和特枯年份均呈下降的趋势，T则呈现增大的趋势；相比之下，在特枯年份，叶面积指数较小时T即会达到峰值。侧柏E随着叶面积指数的增大呈下降趋势，T则呈现波动变化，出现了几个大小不一的峰值。油松E随着叶面积指数的增大同样呈下降趋势，T则呈现稳定的上升趋势。树种间对比可以看出，叶面积指数对侧柏的影响最不稳定，这说明侧柏对水分较为敏感，因此，其LAI的增加会导致E，T间比例关系的不稳定所致。

图3 不同林地高、低效绿水流随叶面积指数变化过程

2.3.2 太阳净辐射(Rn) 图4为不同林地低效绿水(E)和高效绿水(T)随太阳净辐射(Rn)变化图。可以看出，偏枯年(2016年)Rn相对比较充足，随其增大，各林地的E，T均呈现出明显的波动上升趋势。但是在特枯年(2017年)，Rn较小时，即可使T达到峰值。对各林地进行对比，可以看出，在偏枯年份，侧柏E的最大值出现所对应的Rn最小，而刺槐和油松则保持相对一致；在特枯年份，刺槐E的最大值所对应的Rn最小，油松次之，侧柏最大。对T而言，偏枯年份和特枯年份刺槐T最大值对应的Rn都较小，侧柏和油松保持一致。对比分析可以看出，Rn对侧柏林地E的影响最大，对刺槐林地的T影响最大。

图4 不同林地高、低效绿水流随太阳净辐射变化

2.3.3 土壤含水率(θ)

(1) 对低效绿水的影响。图5是不同林地低效绿水(E)随10 cm土壤含水率(θ10 cm)变化图。由图5可以看出，E整体随着θ10 cm的增大而增大。在偏枯年份(2016年)，E总量表现为刺槐(75.6 mm)>油松(64.1 mm)>侧柏(59.2 mm)，在特枯年份(2017年)，E总量表现则为油松(97.9 mm)>侧柏(79.9 mm)>刺槐(59.6 mm)。对比分析可知，不同水文年油松低效绿水的变化最大(52.73%)，侧柏次之(34.97%)，对刺槐的影响最小(21.16%)。根据田间持水量试验中测定的0—40 cm土壤的田间持水量，设定4个土壤水分梯度，分别为：充分供水(>75%田间持水量)，轻度缺水(50%～75%田间持水量)，中度缺水(25%～50%田间持水量)，重度缺水(<25%田间持水量)。相对含水量指土壤含水量占田间持水量的百分比，3块林地对比可知，偏枯年份刺槐林地10 cm相对含水量集中保持在26.5%～35.2%之间，特枯年份保持在35.5%～54.4%之间；偏枯年份侧柏相对含水量集中保持38.8%～56.5%之间，特枯年份保持在26.4%～56.4%之间；偏枯年份油松10 cm相对含水量集中保持37.4%～52.7%之间，特枯年份保持在37.2%～59.9%之间。可以看出，在枯水年不同时期林地土壤水分大多都保持在田间持水量的60%以下，但是幼苗存活重要条件之一是土壤含水量至少要保持在田间持水量的50%以上[25]，侧柏幼苗甚至要保持在70%～80%才是最适宜的水分环境[26]。试验期中的无雨期间，两种水文年份下，刺槐林地相对含水量低于50%的占90.1%和56.8%；侧柏为64.0%和68.4%；油松是70.1%和71.2%，大部分时期人工林地都处于中度缺水状态，即存在着水分不足进而限制树木生长的现象，这是导致树苗成活率低的一个重要原因。

图5 不同林地低效绿水流随10 cm土壤含水率变化过程

(2) 对高效绿水的影响。图6是不同林地高效绿水(T)随20 cm土壤含水率(θ20 cm)变化图。由图6可以看出，随着θ20 cm的增大，各样地不同年份林木的T都呈下降趋势。在偏枯年份(2016年)无雨期内，T总量表现为：侧柏(133.9 mm)>油松(114.2 mm)>刺槐(106.3 mm)，在特枯年(2017年)，则表现为刺槐(99.8 mm)>油松(81.3 mm)>侧柏(72.4 mm)。对比分析可以看出，不同水文年侧柏林地内T的变化幅度最大(45.93%)，油松次之(28.81%)，刺槐最小(6.11%)。3块林地对比可知，偏枯年份刺槐林地20 cm相对含水量集中保持在26.9%～39.9%之间，特枯年份保持在43.3%～65.1%之间；偏枯年份侧柏林地相对含水量集中保持40.4%～59.4%之间，特枯年份保持在30.8%～52.3%之间；偏枯年份油松林地相对含水量集中保持42.2%～52.4%之间，特枯年份保持在41.7%～56.4%之间。两种水文年份下，刺槐林地20 cm相对含水量低于50%的占90.1%和61.0%；侧柏为60.1%和80.0%；油松是66.7%和50.5%。三片林地大多处于中度缺水状态，只有极少一部分时间处于轻度缺水状态。

图6 不同林地高效绿水流随20 cm土壤含水率变化过程

3 讨论与结论

水分胁迫是植物正常生长和发育过程中最容易受到的非生物环境胁迫因子，水分的供给情况往往成为限制植物生长以及成林的重要因素[27-28]。土壤水分缺乏的情况最容易出现在无雨期阶段，在黄土高原，如何保证人工林在无雨期间的存活率是保证人工林存活的关键。对主要造林树种在无雨期间绿水流变化研究的首要问题就是进行绿水流的分离。本研究结果表明，率定后的基于互补相关原理的Granger模型在南小河沟的适用性良好，在率定期纳什效率系数均在0.7左右，在验证期也都超过了0.5，表明模拟结果较好，在数据不充足或者缺失的小流域可用来模拟绿水流的变化情况。Xu等[29]指出互补相关模型在率定参数后可以在各气候区得到理想的蒸散量，这与本研究结果一致。

本研究的结果表明，E对10 cm土壤含水率敏感，T对20 cm土壤含水率敏感。宋孝玉等[24]研究表明，E与土壤表层含水率有较大关系，这与本研究结果一致。王健等[30]认为，T与40 cm土壤含水率相关性较好，与本研究不一致，这是由于在南小河沟流域，降雨多以大雨、暴雨形式发生，但持续时间短，因此，更多的降雨来自小于10 mm的小雨，这种小雨一般只能影响到表层(0—10 cm)或次表层(10—20 cm)的土壤含水率[31-32]，对40 cm土壤含水率的影响较小，所以本研究中T表现出对40 cm土壤含水率不如对20 cm土壤含水率更为敏感。

本研究表明，侧柏的T受LAI影响较小，而刺槐、油松的E，T以及侧柏的E都受LAI影响较大。相对而言，Rn对侧柏林地的E影响最大，对刺槐林地的T影响最小。不同水文年油松E的变化最大，侧柏次之，刺槐最小。不同水文年侧柏T变化幅度最大，刺槐次之，油松最小。对比E和T的变化可知，不同水文年对侧柏T的影响大于E，对刺槐和油松则正好相反。各树种中，侧柏对水分最敏感。这与以往的结果不一致[33]，这是因为从幼苗开始，侧柏就需要土壤水分保持在田间持水量的70%～80%才能正常生长[26]，而且侧柏根系分布在0—90 cm，主要分布在0—40 cm，刺槐根系分布在0—120 cm且呈均匀分布[33]，本研究所用的水文年份为偏枯年份和特枯年份，降雨量距平百分率Pa2016=-27%，Pa2017=-60%，差异较大，在特枯年份几场大暴雨的存在就决定了当年大部分的降雨量，大暴雨可以补充刺槐林地的深层土壤水分，对刺槐而言，由于根系均匀分布，其可以将深层土壤水分转化成高低效绿水；而侧柏由于根系的集中分布导致其对深层水的利用相对困难，因此，其绿水流少于刺槐。由绿水流总量可知，偏枯年份到特枯年份，油松林地绿水流总量相对稳定，这与以往的研究基本一致[5,34]。同时，刺槐相比侧柏也较为稳定，这主要是因为刺槐林地穿透雨量多于侧柏[3]，降雨对刺槐林地土壤水分的补充较为充足，使得刺槐根系所受的水分胁迫程度低于侧柏。

对土壤水分析可知，3种林地大部分阶段都处于中度缺水状态，存在着水分不足进而限制树木生长的现象，这意味着绿水资源对树种生长的支持程度有所欠缺，因此，对这3种树种而言，植树造林时，要将人工林密度控制在土壤水分对应的植被承载能力内。同时，也应通过改善林下植被层次以及枯落物覆盖等的影响，来尽量减少E，促使E朝着T转化，从而促进绿水资源更有效的利用。