黄河流域农业干旱时空演变规律及其自然恢复期特征

2022-03-26薛万来

干旱地区农业研究 2022年2期

刘晔，薛万来

(1.河北工程大学水利水电学院，河北邯郸 056000；2.北京市水科学技术研究院，北京 100084)

关键字：农业干旱；干旱恢复时间；空间异质性；黄河流域

黄河流域水资源量仅占全国的2%，却承担着全国15%耕地面积和12%人口的供水任务，人均水资源量约为408 m3，仅为全国水资源平均水平的20%左右[1]。同时，灌溉农业是黄河水资源高消耗用户，2018年黄河总取水量为516.22×108m3，其中农田灌溉取水量为327.7×108m3，占总取水量的63.48%[2]。此外，受全球变暖和高强度下垫面变化的影响，地处气候湿润-干旱过渡带的黄河流域河川径流量逐年减少，农业干旱呈频发和广发的态势，导致流域粮食生产安全面临着严峻挑战[3-4]。因此，精准评估黄河流域不同气候区农业干旱发生、发展的时空演变格局及其自然恢复期，对农业抗旱、防汛工作的顺利开展以及旱涝保收水平的提高具有实际意义。

近年来，国内外学者针对黄河流域干旱的时空演变规律开展了大量研究，并取得了较好的评估结果[5-7]。在气象干旱评估方面，王璐等[8]采用Copula函数和综合干旱指数(MSDI)分析了干旱多属性概率特征及其动态变化，并运用交叉小波变换方法探究了太阳黑子、大气环流异常因子和植被覆盖对黄河干旱的影响。王飞等[9]以标准化降水蒸散指数SPEI作为干旱指标，应用极点对称模态分解方法对SPEI序列进行时频分解，并从干旱年际变化、季节变化、干旱频率和干旱强度等方面揭示了黄河流域8个水资源分区干旱的时空分布格局。在农业干旱评估方面，钱云平等[10]以土壤含水率为基础建立了作物干旱评估指标，研究了黄河流域农业干旱的地域性、季节性、持续性等特征。王飞以黄河流域为研究对象，基于MODIS遥感卫星数据的归一化植被指数NDVI与地表温度LST数据，采用S-G滤波方法重构NDVI和LST时间序列，构建了5种不同的农业干旱评估指标(VCI、TCI、VHI、MTVDI、NVSWI)，揭示了2000—2015年黄河流域农业干旱的时间演变、空间分布特征[5]。梳理过去的研究发现，多数研究侧重于气象干旱评估，对农业干旱的研究较少，且对其自然恢复期方面的研究存在明显不足，这将使得流域粮食生产安全面临着严峻挑战。

因此，本文以典型气候过渡区黄河流域为例，考虑到气象、下垫面要素的空间异质性，以水资源分区为标准将黄河流域划分为8个子区域，并利用1970—2014年气象数据集，构建了修正的帕默尔干旱指数(sc_PDSI)，揭示了不同气候区农业干旱事件的时空演变特征；其次，采用游程理论方法探究了干旱特征变量空间分布格局，并以分辨率0.25°×0.25°栅格为响应单元，耦合农业干旱检验结果，利用陆地水储量数据(TWS)计算了农业干旱的恢复期；最后，采用地理统计方法(全局和局部Moran’sI指数)揭示了不同子区域农业干旱恢复期的空间异质性分布特征。

1 研究区概况及数据来源

黄河是我国第二大河，从源头至入海途经青海、四川、甘肃、宁夏、内蒙古、山西、陕西、河南、山东共9个省份，干流长度为5 464 km，流域面积约为79.5万km2，属于典型的复杂气候过渡带(半湿润、干旱、半干旱)。降水量年内分配不均，主要集中于6—10月份，占1970—2014年多年平均降水量的69.50%左右，降水量自东南至西北逐渐递减且多年平均气温介于-4℃～14℃之间，气温西北低、东南高[11-12]。考虑到气候要素和下垫面因子的空间异质性影响，本文以二级水资源分区为例，分别为龙羊峡以上、龙羊峡至兰州、兰州至河口镇、内流区、河口镇至龙门、龙门至三门峡、三门峡至花园口以及花园口以下(分别简称为区域 I、II、III、IV、V、VI、VII和VIII)。流域地理位置及陆地水储量站点空间分布如图1所示(见 YE页)。

图1 陆地水储量站点空间分布及流域地理位置图(I、II、III、IV、V、VI、VII和VIII分别为龙羊峡以上、龙羊峡至兰州、兰州至河口镇、内流区、河口镇至龙门、龙门至三门峡、三门峡至花园口以及花园口以下区间边界；TWS代表陆地水储量)

研究数据包括：1970—2014年月尺度降水量、平均气温、潜在蒸发量、土壤体积含水量以及陆地水储量等数据，其中潜在蒸发量采用彭曼公式计算得到。气象数据集来源于英国东英吉利大学的气候研究中心(Climate Research Unit，CRU)，版本选用CRU TS 4.03，数据由气象站实测数据空间插值获得，空间分辨率为0.5°×0.5°。土壤体积含水量和陆地水储量数据来源于GLDAS(https://disc.gsfc.nasa.gov/)数据集中的CLSM模式，其为由美国航空航天局(NASA)、戈达德太空飞行中心(GSFC)和美国国家海洋和大气管理局(NOAA)国家环境预报中心(NCEP)联合开发的全球高分辨率的陆面数据，空间分辨率为0.25°×0.25°。为保证空间分辨率一致，将所有数据通过反距离权重方法(IDW)重采样为0.25°×0.25°的空间分辨率。

土壤湿度单位(kg·m-2)转化为国际单位制的体积含水量(m3·m-3)的计算公式为[13]：

(1)

式中，θ为土壤体积含水量，m3·m-3；D为土壤层厚度，mm；kg为重量单位；mm为厘米单位。值得指出的是，由于GLDAS土壤湿度产品的单位为kg·m-2,通过公式(1)将其转化为单位体积含水量。

2 研究方法

2.1 修正的帕默尔干旱指数(sc_PDSI)

(1)计算农业干旱水分异常指数Z

(2)

将水分距平d与相应月份的气候权重系数K相乘，获得Palmer Z指数[14]：

Z=dK

(3)

(2)气候权重系数K的计算公式为：

(4)

(5)

(3)PDSI指数的计算公式为：

(6)

2004年WELLS等[15]提出了修正的PDSI指数，计算出每个区域所对应的权重系数和持续因子，具有较好的空间可比性，sc_PDSI干旱等级标准见表1。

表1 sc_PDSI干旱等级划分标准

2.2 游程理论

游程理论方法是YEVJEVICH于1967年提出的广泛用于识别干旱事件的重要手段[16]。为了揭示农业干旱事件的多变量属性特征，本文利用游程理论各个响应单元的干旱特征量，其详细计算过程为：当sc_PDSI≤-1时，记作一次干旱事件，该次干旱事件从发生至结束记为干旱历时(月)，该次干旱历时相应的sc_PDSI绝对值的累加之和记为干旱烈度，其与干旱历时的比值记为干旱强度；所有干旱事件烈度的最大值记为干旱峰值，即特定干旱事件中最严重月份的干旱程度。

2.3 干旱恢复时间

基于游程理论农业干旱的识别过程和陆地水储量数据，农业干旱的自然恢复时间计算过程为：

(1)当sc_PDSI≤-1时，记作一次干旱事件，并识别出该次干旱事件的起始和终止时间；

(2)利用陆地水储量时间序列，识别出每次干旱事件对应的平均土壤含水量赤字水量，识别出农业干旱起始时刻对应的土壤含水量，并计算其差值，即农业干旱恢复其初始状态所缺的水量；

(3)以一次农业干旱终止时刻为干旱自然恢复期的起始时间，当(2)中所缺水量大于零时，相应的恢复历时记为农业干旱自然恢复期；

(4)将所有干旱事件干旱恢复时间进行平均，即为该响应单元对应的干旱恢复期。

2.4 地理统计方法(全局和局部Moran’s I指数)

Moran’sI指数是广泛用于度量空间相关性的一个重要指标，其包含全局Moran’sI(Global Moran’s index，GMI)和局部Moran’sI(Local Moran’s index，LMI)。前者是由澳大利亚统计学家PATRICK[17]于1950年开发的空间自相关性度量指标，后者则由美国ANSELIN教授[18]于1995年提出。其计算公式分别为：

z(ui),z(uj)∈δ(u,h)

(7)

(8)

式中,n为样本个数；z为空间变量的属性值，即本文的农业干旱恢复时间；h为相邻栅格之间的距离；u为空间位置；δ(u，h)表示满足||ui-uj||≈h的点对集合；ωij为空间权重矩阵；m为样本均值。

全局Moran’sI可以分析整个不同子区域的干旱恢复时间在空间上的潜在依赖性，用单一的值来反映其空间自相关程度；GMI>0(<0)代表干旱恢复时间在空间上呈正相关(负相关)，其值越大(小)空间相关性越强(弱)，GMI∈[-1,1]；局部Moran’sI指数揭示了栅格响应单元与相邻单元在干旱恢复时间的相关性，即揭示空间聚集或离散状态。

局域性关联性确定空间权重的计算公式为：

(9)

(10)

3 结果与分析

3.1 不同子区域农业干旱事件的时空演变特征

基于8个子区域1970—2014年面平均气象和土壤水数据，采用修正的帕默尔干旱指数计算了月尺度农业干旱事件的时程变化特征。分析图2发现：黄河流域不同子区域农业干旱的时程演变特征大致分为两个阶段：1970—1990年和1990—2014年，其中，1990年以前，不同子区域对应的农业干旱时程变化规律交替变化，即蓝色和红色渲染相当，而1990年以后，不同子区域间农业干旱程度差异显著，如子区域 I 多以蓝色渲染为主，而子区域 III～VII 多以红色渲染为主，即黄河源区农业干旱呈暖湿化的趋势，而其他区域呈暖干化的趋势；从sc_PDSI值角度出发，由黄河上游至黄河下游，农业干旱程度逐渐加剧(红色渲染对应的干旱面积较大)，尤其是1990年以后，这可能归因于该时期人类活动影响加剧，如黄河中下游存在四大农业灌区(宁蒙灌区、河套灌区、汾河灌区和关中灌区)，农业粮食生产大量取水导致河川径流量显著减少[19]；由此看来，黄河流域农业干旱在空间上具有异质性，且随着气候变化，这一差距愈加突显。

图2 黄河流域8个子区域农业干旱事件的时程变化特征

为进一步揭示黄河流域8个子区域农业干旱的空间变化规律，本文采用Mann-Kendall趋势检验方法对不同响应单元下的sc_PDSI序列进行趋势性检验，并采用反距离权重方法(IDW)对趋势变化检验统计量Z值进行空间插值，结果如图3所示(见218YE页)。分析图3发现，龙羊峡以上(I)、龙羊峡至兰州(II)西北部区域统计量Z值均大于0，即该区域农业干旱整体呈湿润化特征，且区域 I 大部分区域在99%置信水平下是显著的(统计量Z>2.57)，这一结果侧面印证了黄河上游农业干旱时间尺度农业干旱化程度降低的科学性和合理性；内流域(IV)农业干旱趋势检验的统计量Z值均小于-2.57(99%置信水平)，即未来内流区农业干旱化程度加剧；通过分析不同子区域统计量Z值可知，黄河流域中下游区域农业干旱愈加严重，而黄河上游未来农业干旱程度将得以缓解，该研究结果与叶培龙等[20]得到的近40 a来黄河上游全区域暖湿化趋势的气候特征基本一致。

图3 黄河流域8个子区域农业干旱事件的空间变化趋势

整体而言，黄河流域农业干旱化演变趋势在空间上干湿分明，黄河上游地处高海拔寒区，气候变暖导致冰川融水，局地降水增加，农业干旱程度减弱，而中下游存在大型农业灌区、高密集人类活动区，天然河道取用水增大，导致农业干旱程度加剧。

3.2 基于游程理论的干旱特征变量识别结果

由于农业干旱具有多变量属性(干旱次数、干旱历时、干旱强度和干旱峰值)，采用游程理论方法从各响应单元sc_PDIS序列中分离出农业干旱的多变量属性，统计结果如图4所示。

分析图4发现：由历史时间序列(1970—2014年)干旱特征统计结果箱型图中值可知，干旱特征变量属性在不同子区域之间差异明显，且干旱历时和干旱强度具有很好的长程相依性，即干旱历时长，干旱强度大；子区域I、II、III、IV、V、VI、VII和VIII对应的干旱次数分别为18、16、23、24、25、17、27、35次，相应的干旱历时分别为13、16、10、14、15、17、10、9个月(箱体中位数可知)，需要指出的是子区域VI的干旱强度和干旱峰值最大，子区域VIII最小，即前者发生极端干旱的风险较大；由干旱峰值箱体宽度可知，子区域 I、子区域 VII和子区域VIII的箱体宽度较小，而子区域II～IV的箱体宽度较大，这归因于前者多属于半湿润气候区，而后者属于干旱、半干旱过渡带，导致计算的农业干旱特征变量值在年际直接差异程度较大。

图4 黄河流域8个子区域1970—2014年农业干旱特征变量的统计结果

通过深入分析发现，黄河流域农业干旱化程度与地理位置息息相关，流域中游区域干旱化程度最高，而下游次之，上游最小，这可能归因于上游和下游均属于半湿润气候区，而中游地处干旱、半干旱过渡区。

通过上述分析发现，受气候、下垫面类型的空间异质性影响，干旱特征变量属性值在区域之间具有显著的差异性。为了更好地展示其空间差异性，基于各响应单元对应的干旱特征值，采用IDW方法进行空间插值，结果如图5所示。

分析图5发现：(1)由西至东，黄河流域农业干旱次数整体呈逐渐增加的变化趋势，最大值出现在子区域 VIII；(2)通过对比干旱历时、干旱强度和干旱峰值的空间变化趋势可以发现，三者的空间分布规律基本相同，但局部区域略有差异；(3)子区域II(龙羊峡至兰州)南部和子区域VI(龙门至三门峡)农业干旱历时长、强度和峰值均较大，这可能归因于该区域为农业高耗水区(汾河灌区和关中灌区)，导致农业干旱化程度较高；(4)纵观上述研究发现，各子区域之间，农业干旱既有差异性、又有相似性，有必要深入厘清农业干旱属性值在空间上的空间异质性特征，重点布控干旱高发区。

3.3 不同子区域农业干旱自然恢复时间的空间分布规律

由前文分析可知，不同子区域间农业干旱时空变化特征具有显著的差异性，本节通过耦合游程理论检验和陆地水储量变化过程，计算了不同子区域的农业干旱恢复时间，结果如图6和图7所示。

分析图6发现：黄河流域农业干旱恢复时间介于4～8个月之间(箱型图中位数)，各区域之间存在差异性，整体呈现为：子区域 IV的干旱恢复时间长，其中位数高达8个月左右；子区域 III次之，约为6个月；子区域 II和 VII的干旱恢复时间最小，约为4个月。值得注意的是子区域 IV的干旱修复时间最长，这主要归因于该区域属于典型的内流域，与外界不存在水资源之间的交换，仅为降水补给，从而导致其农业恢复时间较长。

图6 不同子区域农业干旱自然恢复期的差异性特征

图7为不同子区域各响应单元下农业干旱恢复时间的空间分布规律。由图可知，黄河流域北部和南部区域干旱恢复时间较短，而流域中部区域干旱恢复时间较长，局部区域干旱恢复时间高达36个月；不同子区域，农业干旱恢复时间存在明显的聚集性空间分布特征，如子区域 I的中部、子区域 II的南部等；黄河流域农业干旱恢复时间介于1～36个月，4个月以内所占的响应单元面积比重为42.68%，相应的大于4个月的响应单元面积比重为57.32%。

图7 不同子区域农业干旱自然恢复期的空间差异性

总的来说，我们发现干旱地区(北部)的生态系统从干旱中迅速恢复，这归因于该生态系统中的植物通过不同的生理、解剖和功能机制来适应长期缺水，帮助它们减少水分损失、呼吸成本、光合活性和最大限度地吸收水分。

3.4 农业干旱恢复时间的全局和局部空间自相关性

受气候、下垫面要素等影响，干旱因子特征变量属性值存在一定的空间相关性，即农业干旱恢复时间亦具有“相近相似”的空间分布格局[21]。因此，本文以不同子区域为例，采用地理统计分析方法(全局和局部Moran’sI指数)揭示农业干旱恢复时间的空间自相关性。

图8为不同子区域农业干旱自然恢复时间的全局空间自相关性检验结果。分析图8可知：除子区域 II，其他子区域的全局空间自相关性统计量Moran’sI值均大于0.20，且在99%置信水平下是显著的(T-检验)；子区域I、II、III、IV、V、VI、VII和VIII对应的Moran’sI值分别为0.44、0.02、0.35、0.41、0.26、0.26、0.66和0.26；在99%置信水平下，统计量T值分别为12.45、0.60、11.43、6.95、7.10、9.14、11.19和4.50，这意味着除了子区域 II，其他区域农业干旱恢复时间在空间上均存在显著的空间自相关性。

图8 不同子区域农业干旱自然恢复时间的全局空间自相关性

图9为不同子区域农业干旱自然恢复时间的局部空间自相关性检验结果。由图可知：黄河流域大部分区域局部空间自相关性检验结果为非显著性，即相邻区域之间农业干旱恢复时间的差异性较小，而流域中部区域存在显著的高高聚集现象；以子区域为研究对象，子区域 I中部、子区域 II北部、子区域 III西南部、子区域 IV 南部、子区域 V西南部和子区域 VI的北部均存在显著性的“高高聚集”现象，即该区域的农业干旱修复时间明显高于相邻区域，且在99%置信水平下是显著的；流域中部的子区域高高聚集位置周边存在“低高异常”现象，这意味着随着农业干旱程度的加剧，其干旱恢复时间将增加(被同化)；通过本文的研究可以精准识别出农业干旱发生时其所需要的自然干旱恢复时间，并检验出其区域性空间集聚特征，可为抗旱减灾提供依据，能够适时提高农业粮食安全生产的保障水平。

图9 不同子区域农业干旱自然恢复时间的局部空间自相关性

4 讨论

随着气候变化的加剧，农业干旱已成为农业生产和粮食安全的重要影响因素[22]。本文以黄河流域为例，通过构建修正的帕默尔干旱指数(sc_PDSI)，探究了不同气候区农业干旱的时空演变规律，结果表明：流域上游农业干旱程度最小，中游和下游农业干旱化程度加剧，这一结论与王飞采用VCI、TCI等农业干旱指标得到的干旱趋势检验结果一致[5]。同时，庞爱萍等[23]研究发现黄河上游降水量呈增加趋势，而其他区域显著减少，该结论侧面证实了本研究结论。

本研究发现，干旱区生态系统发生农业干旱可以迅速恢复，而干旱、半干旱过渡带的农业干旱时间较长，同时半湿润区(黄河源区)的干旱恢复时间也较长。该结论与ZHANG等[24]得到的不同气候区干旱恢复时间不同、且干旱气候区对干旱恢复时间较短的结论相一致。本研究对于水资源短缺的黄河流域抗旱减灾工作的顺利开展可提供一定的理论支撑。

在本研究基础上，尚需开展以下研究，如土壤需水量数据源的选择仍具有不确定性，未来将增加集合土壤水预报数据，并利用GPP或NPP等植被生态系统要素，深入探究农业干旱发生后生态系统对干旱的响应机理。此外，本研究以相邻区域进行全局和局部空间自相关性检验，为了增加检验结果的鲁棒性，未来将以0.25°×0.25°为响应单元，按照距离的方式计算农业干旱恢复时间空间自相关性。