蒙古高原土壤水分时空格局演变特征分析

董金义,罗敏,孟凡浩,萨楚拉,包玉海

(1.内蒙古师范大学地理科学学院,呼和浩特 010022;2.内蒙古自治区遥感与地理信息系统重点实验室,呼和浩特 010022)

土壤水分作为大气-植被-土壤-地下水系统的纽带和核心,同时也是地表植被赖以生存的水分条件基础,降水、灌溉水以及地下水只有转化成土壤水分才能被植被吸收利用[1-2]。聚焦土壤水分的时空分异特征及变化规律并实现关键驱动因素的量化对区域生态恢复及可持续发展都具有重要的现实意义[3]。蒙古高原位于典型的干旱半干旱气候过渡带,水资源稀缺是制约当地经济发展和生态稳定的关键因素[4-5]。特别是近几十年受气候变化以及经济快速发展的影响,蒙古高原土壤干化现象明显、植被生产力退化加剧,严重威胁区域生态系统稳定和已经取得的生态建设成效[6-7]。蒙古高原特殊的地理位置,导致长时间、连续准确的土壤水分监测数据缺乏,对蒙古高原土壤水分的时空分异特征及关键驱动要素很难进行全面分析[8-9]。

遥感监测、反演以及陆面过程模型模拟的土壤水分产品具有快速、宏观和动态等特点,使得土壤水分相关研究实现了由点向面的跨越,对大尺度土壤水分时空变化规律研究具有重要意义[10-11]。遥感土壤水分产品AMSR-E,SMOS,ECV 以及再分析土壤水分产品NCEP-CFSR,GLDAS,MERRA,ERA5等广泛用于全球及区域尺度土壤水分相关研究[12]。如范科科等[13]基于NCEP-CFSR 等再分析土壤水分产品,揭示了青藏高原地区土壤水分的时空变化特征;Qiu等[14]联合ECV,ERA-Interim 和站点监测数据,量化了1996—2010年中国土壤水分变化特征。然而,作为气候变化的敏感区域,目前针对蒙古高原土壤水分时空变化规律的研究相对较少,主要集中在数据精度评价方面[15]。因此,亟待基于遥感或再分析土壤水分产品开展气候变化背景下蒙古高原土壤水分时空分布特征及其变化规律的研究。

针对以上事实,为全面揭示蒙古高原土壤水分时空分异特征及变化规律,本文基于长时序、高精度的ERA5再分析土壤水分产品,研究蒙古高原不同深度土壤水分的时空分异特征及动态变化趋势;在此基础上进一步揭示蒙古高原气象、NDVI和TWSA 等关键驱动要素的影响,以期为气候变化背景下该区域水资源管理、生态系统保护及干旱等自然灾害的预防提供科学依据。

1 研究区概况与数据来源

1.1 研究区概况

蒙古高原主体包括中国内蒙古自治区(以下简称“内蒙古”)和蒙古国,地理位置介于88°43′—126°04′E,37°24′—53°23′N。蒙古高原北起萨彦岭、雅布洛诺夫山脉,南至阴山山脉,东起大兴安岭、西至阿尔泰山;地势整体西高东低,平均海拔1 580 m(图1)。研究区属温带大陆性气候,多年平均降水量约200 mm;冬季寒冷干燥且时间跨度长,夏季炎热昼夜温差大。受高原水热条件影响,植被覆盖由东北向西南依次跨越森林、森林草原、典型草原、荒漠草原、戈壁荒漠,生态环境多样且脆弱[16]。

图1 蒙古高原地形和土地覆盖类型Fig.1 Mongolian plateau terrain and land cover types

1.2 数据来源

本研究采用的数据包括土壤水分、土地利用覆盖变化(LUCC)、降水、温度、数字高程模型(DEM)、土壤类型、NDVI和TWSA 等数据。土壤水分、降水及温度数据采用2000—2020 年的ERA5 再分析数据集,时空间分辨率分别为1个月和0.1°。其中,ERA5产品包括0—7,7—28,28—100,100—289 cm 共4层不同深度土壤水分。DEM 数据源于中国地理空间数据云平台,用于对研究区坡度和坡向的计算。土壤类型数据来源于FAO 的HWSD1.2全球土壤同化数据库,该数据包括土壤中砂粒、粉粒和黏粒含量等参数。NDVI数据来源于美国地质勘察局(USGS)网站发布的MODIS每8 d的MOD13A2产品,空间分辨率为500 mm。TWSA数据采用Li等[17]通过机器学习技术重建的1979年以来全球尺度的数据集,其时空分辨率为0.5°/月。LUCC数据来自1∶100万的内蒙古和蒙古国植被类型数据。为了保证数据分辨率的统一,本研究通过最邻近插值法将所有数据均重采样到0.1°。

2 研究方法

2.1 Theil-Sen Median趋势分析和Mann-Kendall显著性检验

Theil-Sen Median趋势分析是一种非参数估计方法,该方法受数据异常值的影响较小,能够很好地揭示时间序列的变化趋势[18]。在像元尺度上,采用该方法计算2000—2020年土壤水分和其驱动因子的变化趋势。具体计算公式如下:

式中:xj与xi分别为第i和第j年的像元观测值;β为时间序列的斜率;当β＞0时,表示研究时段内时间序列呈增加趋势,反之呈下降趋势。

Mann-Kendall检验方法可用来判断Theil-Sen Median趋势的显著性。相关计算公式为:

式中:xj与xi分别为第i和第j年的像元观测值;n为时间序列长度;sgn为符号函数。在给定显著性水平α下,当统计量|Z|＞u1-α/2时,表示时间序列变化趋势显著。本研究取α=0.05,相应Z的阈值为±1.96。

2.2 重标极差(R/S)分析

R/S分析用于定量描述时间序列的长程依赖性。本研究基于R/S分析方法的Hurst指数分析蒙古高原土壤水分变化趋势的持续性特征[19]。对于蒙古高原土壤水分的时间序列SMt,t=1,2,3,…,n,对于任意正整数τ,定义该时间序列的均值序列:

累计离差:

极差:

标准差:

Hurst指数:

Hurst(H值)包含3种形式:当0.5＜Hurst≤1,表明时间序列是一个持续性序列,即未来变化与过去变化趋势一致,且H值越接近于1,则持续性越强;如果H=0.5,则说明时间序列为随机序列,不存在长期相关性;如果0≤H＜0.5,表明此时间序列具有反持续性,即未来变化趋势与过去保持相反,H值越接近于0,则相逆性越强。

2.3 地理探测器

地理探测器是探测地理要素的空间分异特性并阐明其背后驱动力的一种统计方法[20]。本研究分别利用地理探测器的因子探测、交互探测和生态探测来探究不同深度土壤水分空间分异的驱动力、多因子交互作用以及因子间作用的显著性差异。

2.3.1 因子探测 选择因子探测器计算各因子的q值,来量化不同因子对蒙古高原土壤水分空间分异的贡献,计算公式如下:

式中:因子总数记为K;N与Ng分别为全区与层g的单元数;σ2与σ2g分别为全区与层g的Y值的方差;与Nσ2分别为层内方差之和以及全区总方差;值取值范围介于0～1。选择温度(TEMP)、陆地水储量异常(TWSA)、降水(PREP)、归一化植被指数(NDVI)、土地利用/土地覆盖变化(LUCC)、砂粒(Sand)、粉粒(Silt)、黏粒(Clay)、高程(DEM)、坡度(Slope)与坡向(Aspect)共11个驱动因子进行蒙古高原土壤水分空间分异驱动力分析。

2.3.2 交互作用探测 交互作用探测用于探究两个影响因子共同作用时,对因变量Y的解释能力较单因子为增强还是减弱。首先通过因子探测计算两个自变量X1与X2别对Y的q值,得到q(X1)与q(X2),然后再计算两个因子交互作用下的q值,即q(X1∩X2),最后确定交互类型。

2.3.3 生态探测 生态探测用于分析两个影响因子对因变量空间分布的影响是否呈显著差异,用G统计量来衡量:

式中:NX1与NX2为因子样本量,其余符号含义均与因子探测相同。

2.4 回归分析

采用回归分析定量区分主导环境因子的变化对蒙古高原土壤水分的影响[21]。首先对参与回归的时间序列进行标准化处理:

式中:γ为原始时间序列;和s(γ)分别为时间序列的均值和标准差。然后对标准化之后的时间序列建立多元线性回归方程:

式中:a,b,c,d分别为降水、温度、NDVI和陆地水储量异常的回归系数,可用于量化不同因子对土壤水分影响的大小;θ为回归常数项。

3 结果与分析

3.1 蒙古高原土壤水分时空分布及变化格局

2000—2020年蒙古高原土壤水分空间分布呈现由西南向东北阶梯式增加的分布特征,且北部山区土壤水分明显高于其他地区(图2)。高原0—7,7—28,28—100以及100—289 深度土壤水分均值分别为0.21,0.22,0.21,0.22 m3/m3。不同深度土壤水分变化速率存在强烈的时空异质性,整体上随着土层深度的加深,土壤水分呈显著减少趋势的面积进一步扩张(图3)。蒙古高原0—100 cm 土壤水分呈减少趋势的面积占比分别为44.4%,42.7%,52.2%,主要分布在蒙古高原西北及东南部地区,且大部分区域减少趋势并不显著;而100—289 cm深度土壤水分呈减少趋势的面积占比高达74.7%,且大部分区域(65.2%)土壤水分减少趋势均通过显著性检验,呈增加趋势的区域集中分布在高原东北部。

图2 2000-2020年蒙古高原不同深度土壤水分的空间分布与变化趋势空间分布格局Fig.2 Spatial distribution and change trend spatial distribution pattern of soil moisture at different depths in the Mongolian Plateau from 2000 to 2020

图3 2000-2020年不同深度土壤水分时间序列变化趋势Fig.3 Change trends of soil moisture time series at different depths from 2000 to 2020

从时间尺度来看,2000—2020年蒙古高原0—7,7—28 cm土壤水分整体呈不显著增加趋势,增长率分别为0.000 15,0.000 24 m3/(m3·a);随着土层深度的增加,28—100,100—289 cm 土壤水分开始呈现减少趋势,减少速率分别为-0.000 03,-0.001 3,且100—289 cm土壤水分减少趋势通过显著性检验(p＜0.05)。该结果说明随着土层深度的增加,土壤水分呈现由增到减的变化趋势特征。

3.2 蒙古高原土壤水分变化持续性分析

为了确定土壤水分变化趋势的持续性,通过R/S分析得到蒙古高原土壤水分变化的Hurst指数空间分布图(图4)。结果发现,蒙古高原土壤水分的Hurst指数介于0.26～0.75,从表层到第4层土壤水分呈正向持续性变化(Hurst＞0.5)的区域面积占比分别为58.5%,76.7%,91.3%,98.8%,其中Hurst指数大于0.5 的区域主要分布在肯特山地区、阴山地区、杭爱山北部及大兴安岭南部,说明这些区域土壤水分处于持续稳定变化状态。随着土层深度的增加,Hurst指数大于0.5的区域面积呈现逐渐增加的趋势,即保持当前变化趋势的区域面积将会逐渐增大,说明深层土壤水分变化趋势的持续稳定性更高。

图4 蒙古高原植土壤水分Hurst指数空间分布及变化趋势Fig.4 Spatial distribution and changing trend of Hurst index of plant soil moisture on the Mongolian Plateau

通过Hurst分析结果与变化趋势结果叠加,进一步得出蒙古高原土壤水分变化趋势持续性分类图(图5)。不难看出,从表层到第3层土壤水分呈现持续增加趋势的面积占比分别为38.5%,48.4%,44.2%,主要集中在大兴安岭东部丘陵平原以及色楞格河上游地区。0—100 cm土壤水分持续减少区域主要分布在内蒙古西部和蒙古国西北部区域,占比分别为20.2%,28.4%,47.2%。除东北部和西南部部分区域外,100—289 cm土壤水分均呈现持续减少趋势,占比高达73.8%,其中呈持续显著减少趋势的面积占65.2%。随着深度的加深,呈不稳定变化趋势的地区占比呈现减少趋势(41.3%,23.2%,8.6%1.3%)且主要集中在0—28 cm,表现为前后变化相反或波动变化。

图5 蒙古高原土壤水分因子探测q 值Fig.5 Detection q value of soil moisture factors on the Mongolian Plateau

3.3 土壤水分空间分异驱动力分析

图5显示了不同环境因子对土壤水分的空间分异的驱动力(q值)。TEMP,TWSA,PREP和NDVI对不同深度土壤水分空间分异的解释力均在75%以上(q＞0.75),是蒙古高原土壤水分空间分布的主导环境因子。其中,TEMP 对不同深度土壤水分空间分异的解释度最高,q值介于0.90～0.93。TWSA 次之,q值也均在0.89以上。PREP对不同深度土壤水分的驱动力均排名第3。NDVI的驱动力位于第4,q值介于0.77～0.84。LUCC对不同深度土壤水分驱动力差异不大,q值介于0.66～0.70。在土壤质地因子中,Sand的驱动力相对最大,q值在0.5左右,其次是Silt和Clay,可能与不同土壤质地的持水能力有关。而地形因子对不同深度土壤水分驱动力最弱,q值均不足0.3。相比较而言,不同地形因子中DEM 的驱动力相对最大(q＞0.25);Aspect驱动力最小,q值在0.20 以下,与Slope二者q值均未通过显著性检验(p＞0.01)。

蒙古高原土壤水分空间分异驱动因素的交互探测与生态探测结果如图6所示。主导环境因子(除NDVI外)与土壤质地因子、LUCC以及土壤质地因子与地形因子,交互后的q值均大于单个因子的q值,呈现出双因子增强作用,说明蒙古高原土壤水分空间分异是多因子共同作用的结果。双因子交互后平均解释率最高的4组是:TEMP∩Sand,TEMP∩Silt,PREP∩Sand和PREP∩Silt,交互作用均达到0.93,说明特定土壤质地显著增加了温度、降水对土壤水分的解释力。TWSA 与土壤因子次之,解释率均在0.9以上。Aspect∩Slope交互作用最小,解释率介于0.19～0.25。NDVI与其他所有因子交互均表现为单因子减弱,可能与植物生长需要消耗土壤水分有关。

图6 蒙古高原土壤水分交互探测与生态探测Fig.6 Interactive detection and ecological detection of soil moisture in the Mongolian Plateau

在0—7,7—28,100—289 cm 土层中,TWSA 与NDVI的生态探测结果显示对蒙古高原土壤水分的空间分异特征无显著性差异(检测结果为N),其余因子间均呈现显著性差异,说明在特定深度下TWSA与NDVI对蒙古高原土壤水分的空间分异特征具有相似的机理。

3.4 环境因子时空变化格局

为了深入探讨环境变化对土壤水分动态的影响,利用Theil-Sen Median趋势分析法探究了NDVI、降水、温度与TWSA的变化趋势。如图7所示,NDVI变化趋势呈现较强的空间异质性,变化速率介于-0.01～0.01/a。NDVI呈增加趋势的面积占94.7%,呈减少趋势的面积仅占5.3%,零星分布于蒙古国西北部和内蒙古中部地区。除高原东北部部分区域(10.8%)TWSA表现为上升趋势,其他地区TWSA 均呈现不同程度的减少趋势。近21 a,研究区降水变化速率介于-4.66～12.31 mm/a,其呈现增长趋势的面积占研究区植被覆盖总面积的68.5%。降水增加较为明显的区域主要分布在蒙古高原东部,如大兴安岭、肯特山东部山地地区;而阴山山脉、杭爱山附近地区降水呈显著下降趋势。蒙古高原温度以增加趋势为主,平均增速0.013℃/a。特别是肯特山东部、大兴安岭西部以及东部平原地区的温度增长率明显高于其他地区。温度呈减少趋势的面积占比为15.1%,集中分布在库苏尔湖东南部、大兴安岭东部平原以及西部高原。

图7 2000-2020年蒙古高原影响因子变化趋势空间分布Fig.7 Spatial distribution of changing trends of influencing factors in the Mongolian Plateau from 2000 to 2020

3.5 土壤水分对环境因子变化的响应

基于回归分析法进一步量化各环境因子变化对蒙古高原土壤水分的影响,结果如图8所示。0—7,7—28 cm 土壤水分与降水呈现显著正相关性,分别有99.7%,97.5%的区域回归系数介于0.33～1.00,说明0—28 cm 土壤水分受降水正向影响显著。

通过对比不同环境因子对土壤水分影响的大小,发现降水为0—7,7—28 cm 土壤水分主控因子的区域占比高达98.7%,94.8%(图9)。随着土层深度的加深,降水的主控作用明显降低。81.2%的区域降水对28—100 cm土壤水分呈正相关关系,但降水占主控因子的区域占比下降为38.8%。当研究区土层深度为100—289 cm时,降水正向影响的区域呈现出断层式下降,表现为负相关的区域面积高达65.9%,而降水占主控因子的区域面积仅占10.5%,且主要分布于正向影响区域。TWSA与0—7 cm及7—28 cm土壤水分以正相关(51.0%,72.8%)为主,相关性较低且空间分析相似。与降水相反,随着土层深度的加深,TWSA 与土壤水分的正相关性显著增强,TWSA 与28—100 cm 及100—289 cm土壤水分呈正相关区域的面积占比分别为83.7%,91.1%。此外,TWSA 作为0—7,7—28,28—100,100—289 cm土壤水分主控影响因子的占比分别为0.6%,1.8%,38.8%,58.6%。温度与浅层土壤水分(0—28 cm)以负相关关系为主,温度与0—7 cm 和7—28 cm 土壤水分呈负相关的区域占比分别为89.2%,58.3%;呈正相关的区域集中分布在西北部山区。随着土层深度加深,温度与28—100 cm和100—289 cm土壤水分呈正相关关系的区域进一步扩大,面积占比分别为56.0%,64.6%。需要指出的是,温度作为不同深度土壤水分变化的主控因子的面积占比为0.5%～14.8%。NDVI与不同深度土壤水分的相关性差异较大。NDVI与0—7,7—28,28—100 cm土壤水分呈正相关关系的区域面积分别为51.0%,72.8%,72.0%,且主要分布在高原西北部和东南部地区。NDVI与100—289 cm 土壤水分以负相关为主,呈负相关区域面积占比为61.4%,主要分布在蒙古国东部和内蒙古的中部地区。随着土层深度的增加,NDVI占主导作用的面积比例由0.2%升高至16.2%。

图9 蒙古高原不同深度土壤水分变化主控环境因子空间分布格局Fig.9 Spatial distribution pattern of main environmental factors controlling soil moisture changes at different depths in the Mongolian Plateau

4 讨论

土壤水分作为影响陆面过程的关键因素之一,可以用于对气候和环境变化的实时监测。卫星遥感反演土壤水分产品仅能获取表层(0—5 cm)土壤水分且抗干扰能力较差,而ERA5再分析土壤水分产品可获得不同深度土壤水分数据,具有抗干扰性强和实效性可靠等优点[22-23]。该数据在干旱半干旱区得到广泛应用且精度较高。有研究指出,尽管ERA5 土壤水分数据与实测站点数据相关性随土层深度加深而减弱,但整体能够表现深层土壤水分的分布特征及变化规律[24]。因此,本研究利用ERA5不同深度土壤水分产品对2000—2020年蒙古高原土壤水分时空分异特征及其变化规律进行分析。研究结果表明,蒙古高原土壤水分空间异质性明显,不同深度土壤水分均呈现由西南向东北递增的分布格局。其空间分异主要由温度、TWSA、降水和NDVI等环境因子共同决定。随着土层深度的增加,降水和植被条件对土壤水分的解释力呈现降低趋势,可能是由于深层土壤水分的稳定性更高,而浅层土壤水分更易受到外界环境因素随机变化的影响。陈斌等[25]也指出深层土壤水分较浅层土壤水分稳定性更强。温度主要通过影响蒸散发间接影响土壤水分含量。温度对不同深度土壤水分空间分异的解释力均在90%以上(q＞0.9)。可能是由于在干旱地区土壤水分含量较低,降水量少且转化为土壤水分的能力有限,蒸散发反过来限制降水与土壤水分的转化,导致温度影响略高于降水,但二者均属于主要气候影响因子。此外,在分析温度与土壤水分的关系时,发现二者关系不均为负相关关系,这与韩新生[26]、吴芹芹[27]等的结论不同,与程文举等[28]的结果较为一致。究其原因,蒙古高原属于干旱半干旱区,在本就水分极度匮乏的地区温度升高,并不一定引起土壤水分的减少,同时温度对浅层土壤水分较深层土壤水分的影响力更强,进一步佐证了本研究结论,也反映出土壤水分空间分布的控制因素受基本气候条件影响较大。相较于气候和植被因素,地形要素和土壤质地对蒙古高原土壤水分空间分异的影响较小,但与主导环境因子进行交互作用时,主要呈现出双因子增强的作用。这种现象可能与蒙古高原特殊的土壤质地、地形地貌、水热条件以及特有的植被类型等因素有关。

从时间尺度来看,2000年以来研究区域整体呈现暖湿化趋势。降水作为浅层土壤水分变化的主控因子,降水的增加导致表层和次表层土壤水分呈不显著增加趋势,抵消了升温导致的土壤水分流失。随着土层深度的加深,降水的主控地位显著降低。Yuan等[31]指出蒙古高原低强度降水(＜5 mm)比例远高于中等强度和高强度降水比例,在部分区域该比例甚至高于80%。而低强度降水优先补给浅层土壤水分而无法渗透到深层。这可能是降水主控地位随土层深度加深而显著降低的原因。而对于28—100 cm 土壤水分,降水对其变化呈主导作用的区域主要集中在蒙古国北部和内蒙古东北部山区。可能是由于高山区降水丰沛,中、强降水事件相对较多。与降水相反,随着土层深度的加深,TWSA 对土壤水分的控制作用逐渐凸显。特别是对于100—289 cm土壤水分,TWSA 主控区域远高于其他影响因子。在干旱半干旱区,浅层地下水是深层土壤水分的重要补给来源[32]。考虑到二者高度一致的空间变化分布特征及正相关关系,TWSA 的降低可能是导致研究区深层(100—289 cm)土壤水分显著降低的主要原因。此外,由于土层深度加深导致的降水主控作用的降低,温度和NDVI的主控作用也逐渐凸显。特别是蒙古国东部植被覆盖度增加导致的耗水增加以及高原西南部和东北部部分地区强烈的升温,进一步加剧了深层土壤干化。综上所述,尽管2000年以来蒙古高原降水增加,但深层土壤水分变化与气候干湿变化并不对应,依然存在较大的潜在干旱风险。朱飙等[33]发现潜在蒸散发的增加可能是导致二者不对称的主要原因。后续研究应结合潜在蒸散发等因子进一步深入探索研究区深层土壤水分持续下降的原因。Hurst分析结果也表明,土层深度越深,未来土壤水分变化趋势的可持续性越强。特别是蒙古国西北部和内蒙古中部的0—100 cm 土壤水分以及蒙古国北部和内蒙古中部的100—289 cm 的土壤水分(图5)未来可能会持续干化,应当予以重视。

5 结论

(1)2000—2020年蒙古高原不同深度土壤水分整体呈西南低、东北高的分布格局,存在明显的空间异质性。蒙古高原土壤水分由浅至深呈现先增加后减少的变化趋势特征,仅100—289 cm 土壤水分变化趋势通过显著性检验(p＜0.05),说明深层土壤水分干化情况更加显著。0—7,7—28,28—100,100—289 cm 土壤水分呈现减少趋势的比例分别为44.38%,42.7%,52.2%,74.7%。

(2)蒙古高原土壤水分的Hurst指数介于0.26～0.75,从表层到第4层分别有58.5%,76.7%,91.3%,98.8%的区域Hurst指数大于0.5,说明未来大部分区域土壤水分变化趋势与过去相同;结合趋势分析结果,蒙古高原西北部及内蒙古中部地区土壤水分干化情况可能会进一步加重。

(3)温度、TWSA、降水和NDVI是影响土壤水分空间分布的主导环境因子,地形因子的解释力相对最弱。土壤水分空间分异是多因子共同作用的结果,主导环境因子(除NDVI外)、地形因子与土壤质地因子、LUCC交互后,均呈现出双因子增强作用。

(4)表层和次表层土壤水分的变化主要受降水的正向影响,主控区域分别占98.7%,94.8%;降水和TWSA 对28—100 cm 土壤水分的主控区域分别集中在森林和草原覆盖区,主导区域面积占比分别为38.7%,38.8%;TWSA 降低是导致100—289 cm 土壤水分干化的主要驱动力其主控区域面积占比为58.6%;植被耗水的增加以及温度的升高分别主控16.2%,14.8%区域100—289 cm 土壤水分变化。