环境科学  2024, Vol. 45 Issue (6): 3375-3388   PDF    
基于可解释机器学习的青藏高原草地物候变化多因素影响分析
刘慧文1, 刘欢2, 胡鹏2, 彭辉1,3,4, 王硕1     
1. 中国海洋大学环境科学与工程学院,青岛 266100;
2. 中国水利水电科学研究院流域水循环模拟与调控国家重点实验室,北京 100038;
3. 中国海洋大学海洋环境与生态教育部重点实验室,青岛 266100;
4. 中国海洋大学海洋环境地质工程山东省重点实验室,青岛 266100
摘要: 气候变化背景下, 青藏高原植被物候发生显著改变. 然而, 影响物候的水热因素众多, 目前较少有研究关注多因素对青藏高原物候的影响效应, 导致对青藏高原物候变化机制认识不足. 为此, 研究通过遥感数据解译, 在对2002~2021年青藏高原草地物候时空变化特征分析的基础上, 聚焦降水、气温、海拔和土壤等多方面, 利用可解释机器学习方法(SHAP)揭示物候变化的主导因素, 并量化分析多因素对物候的交互影响. 结果表明:①青藏高原分别有56.32%、67.65%和65.50%的草地表现出生长季开始时间(SOS)提前、生长季结束时间(EOS)延迟和生长季长度(LOS)延长趋势;②青藏高原草地SOS和LOS主要受水分条件影响, 3月0~10 cm土壤水分对SOS提前和LOS延长起促进作用的范围分别在10~25 kg·m-2和15~25 kg·m-2之间, 峰值分别在20 kg·m-2和18 kg·m-2左右;EOS则主要受温度影响, 9月和10月温度越高对EOS延迟促进作用越强, 并分别在高于8℃和-0.5℃时达到峰值;③水热等因素对物候的影响存在非线性交互效应, 3月0~10 cm土壤水分达到20 kg·m-2后, 更有利于低降水和低海拔地区SOS提前;10月温度高于0℃后较好的水分条件更有利于EOS延迟;3月0~10 cm土壤水分在12~22 kg·m-2之间时, 高降水地区LOS更长. 研究表明, 可解释机器学习方法可为物候变化的多因素影响定量分析提供一种新的方法.
关键词: 青藏高原      草地物候      时空变化      极限梯度提升(XGBoost)      SHAP      交互影响     
Multi-factor Impact Analysis of Grassland Phenology Changes on the Qinghai-Xizang Plateau Based on Interpretable Machine Learning
LIU Hui-wen1 , LIU Huan2 , HU Peng2 , PENG Hui1,3,4 , WANG Shuo1     
1. College of Environmental Science and Engineering, Ocean University of China, Qingdao 266100, China;
2. State Key Laboratory of Simulation and Regulation of Water Cycle in River Basin, China Institute of Water Resources and Hydropower Research, Beijing 100038, China;
3. Key Laboratory of Marine Environment and Ecology, Ministry of Education, Ocean University of China, Qingdao 266100, China;
4. Shandong Provincial Key Laboratory of Marine Environment and Geological Engineering, Ocean University of China, Qingdao 266100, China
Abstract: The vegetation phenology of the Qinghai-Xizang Plateau is changing significantly in the context of climate change. However, there are many hydrothermal factors affecting the phenology, and few studies have focused on the effects of multiple factors on the phenology of the Qinghai-Xizang Plateau, resulting in a lack of understanding of the mechanisms underlying phenological changes on the Qinghai-Xizang Plateau. In this study, we used remote sensing data interpretation to analyze the spatial and temporal variability of grassland phenology on the Qinghai-Xizang Plateau from 2002 to 2021, focusing on precipitation, temperature, altitude, soil, and other aspects to reveal the dominant factors of phenological variability using an interpretable machine learning method (SHAP) and to quantify the interactive effects of multiple factors on phenology. The results showed that: ① The growing season start (SOS) of grasslands on the Qinghai-Xizang Plateau mostly ranged from 110 to 150 d, with 56.32% of grasslands showing an early SOS trend; the growing season end (EOS) mostly ranged from 290-320 d, with 67.65% of grasslands showing a delayed EOS trend; and the growing season length (LOS) mostly ranged from 120 to 210 d, with 65.50% of the grasslands showing a trend towards longer growing season lengths. ② SOS in grasslands on the Qinghai-Xizang Plateau was mainly influenced by moisture conditions, in which soil moisture between 10 and 25 kg·m-2 in the 0-10 cm soil layer in March promoted the advancement of SOS and peaked at approximately 20 kg·m-2. EOS was mainly influenced by temperature, with higher temperatures in September and October having a stronger effect on EOS latency promotion and peaking at over 8℃ and -0.5℃, respectively. The main influencing factors of LOS were more consistent with SOS, in which soil moisture between 15 and 25 kg·m-2 in the 0-10 cm soil layer in March promoted the prolongation of LOS and peaked at approximately 18 kg·m-2. ③ There was an obvious interactive effect of water and heat and other factors on phenology; after soil moisture reached 20 kg·m-2 in the 0-10 cm soil layer in March, SOS was more advanced in low-precipitation and low-altitude areas. Better moisture conditions were more conducive to EOS delay at temperatures above 0℃ in October, and soil moisture in high precipitation areas promoted LOS prolongation more when soil moisture was between 12 and 22 kg·m-2 in 0-10 cm in March. The results also demonstrated that interpretable machine learning methods could provide a new approach to the analysis of the multifactorial effects of phenological change.
Key words: Qinghai-Xizang Plateau      grassland phenology      spatial and temporal variation      extreme gradient boosting (XGBoost)      shapley additive explanations (SHAP)      interactive effects     

青藏高原对气候变化高度敏感, 近几十年其变暖速度几乎是全球平均水平的两倍[1]. 青藏高原生态系统以高寒草甸和草原为主, 作为气候变化最敏感最直接的响应指标之一, 青藏高原高寒草地物候变化被广泛关注[2, 3]. 基于原位观测和卫星遥感的数据表明, 快速变暖使青藏高原草地物候发生了显著变化[4], 例如春季物候的提前[5, 6]以及秋季物候的延迟[7, 8]等. 植被物候的变化会对生态系统功能和碳循环过程产生强烈影响[9], 且可能通过改变陆地与大气之间的水和能量交换对区域气候形成反馈作用[10, 11]. 揭示植被物候变化控制机制不仅有助于改进生态水文建模过程, 而且有助于提高未来气候变化情景下对区域碳收支和生态系统变化的预测能力[4].

过去研究多关注物候对气候变化的响应, 通常认为气温和降水是青藏高原草地物候变化的主要影响因素[7, 8, 12], 青藏高原的变暖和降水增加推动了春季物候提前以及秋季物候延迟[13]. 水分可用性是限制植被生长的关键因素, Xu等[14]研究发现青藏高原草地生态系统正由热约束转向水约束, 生长季开始时间受降水影响更大, 且干旱地区相较于湿润地区对降水变化更为敏感[15, 16]. 生长季结束时间受降水和温度协同影响, 充足的降水和变暖共同促进了青藏高原生长季的延长[17]. 气候变化主导了青藏高原植被的改善[18], Wu等[19]研究认为温度是影响青藏高原植被变化最重要的因素, 水分可用性对草地产生了积极影响. 现有研究在揭示青藏高原草地物候变化机制中尚存在一定不足. 一方面, 研究多聚焦单一因素变化对物候的影响[7, 8, 20~22], 忽视了该过程中单因素与其他因素之间的交互影响. 另一方面, 研究多关注气温和降水而较少研究土壤水分变化和海拔高差等对青藏高原草地物候的影响. 土壤水分是植被最直接的水源, 而海拔会通过改变气候格局间接影响植被物候[23, 24]. 除此之外, 物候变化是多因素共同作用的结果, 其与影响因素之间通常存在非线性关系, 而目前研究使用的相关性分析、线性回归等传统的统计学方法难以描述非线性变化, 无法得到各因素对物候变化的相对贡献并且定量解释物候对影响因素变化的响应机制[11].

本文基于2002~2021年的MODIS NDVI数据提取青藏高原草地物候分析其时空变化特征, 进一步, 以海拔以及月尺度的气温、降水和土壤水分作为影响要素, 采用可解释机器学习方法量化了各要素变化对物候的影响及其交互效应.

1 材料与方法 1.1 研究区概况

青藏高原(26°~39°N, 73°~104°E)素有“世界屋脊”之称, 平均海拔4 000 m以上, 总面积约2.6 × 106 km2. 空间上, 青藏高原海拔高差较大, 低处不足2 000 m, 最高则达到8 844 m, 如图 1所示. 受地形和季风影响, 青藏高原水热条件存在较强的空间异质性. 气温呈现出中部较低, 东南部和西北部较高的分布形态[24], 且季节差异明显, 昼夜温差大[25]. 降水空间分配不均, 年平均降水量从东南部高于1 000 mm减少至西南部的低于50 mm[26], 年内有明显的干湿两季[27]. 植被类型以高寒草地为主(图 1), 面积约为1.5 × 106 km2, 占青藏高原总面积的57.7%[28].

图 1 研究区示意 Fig. 1 Schematic of the study area

1.2 数据来源与预处理

本文使用归一化植被指数(NDVI)数据、植被类型数据、高程数据、气象数据以及土壤水分数据, 数据来源见表 1. 除数字高程模型(DEM)和植被图外, 其余数据均采用2002~2021年的.

表 1 数据名称及来源 Table 1 Data names and sources

采用ArcGIS软件将NDVI数据进行格式和投影转换、拼接和裁剪, 并使用TIMESAT软件提取物候数据. 根据1∶100万植被图提取青藏高原草地区域, 并将所有数据统一分辨率为10 km × 10 km, 作为机器学习模型输入.

1.3 研究方法 1.3.1 物候提取方法

为减少云和积雪等对遥感影像的干扰, 首先通过Savitzky-Golay(SG)滤波[29]对NDVI时间序列数据进行平滑处理. 然后, 采用动态阈值法[30]提取植被物候指标, 包括生长季开始时间(SOS)、生长季结束时间(EOS)以及生长季长度(LOS). 该方法已在全球得到了广泛地使用和验证[31~36], 经反复试验将SOS和EOS的阈值设置为0.2, 采用儒略历(1月1日为1 a中的第1 d)来表示上述3个物候指标[33]. 动态阈值法的计算公式为:

(1)

式中, NDVIratio为输出比值, 其第一次超过阈值的时间定义为SOS, 第一次低于阈值的时间定义为EOS[22]. NDVIt为日期t的NDVI值, NDVImin和NDVImax分别为年度NDVI时间序列的最小值和最大值. 为消除非植被区域的影响, 本研究提取了年均NDVI大于0.1的像元[32, 34]. 物候的提取主要在TIMESAT3.3中完成.

1.3.2 趋势分析方法

采用Sen's斜率估计法计算2002~2021年青藏高原草地物候的变化趋势, 并通过Mann-Kendall检验方法判断其变化是否显著. Sen's斜率估计法优点在于可以限制异常值对斜率的影响[37], 计算公式如下:

(2)

式中, median为取中值函数, β为物候指标变化趋势斜率. 当β > 0时, 序列值呈增大趋势;当β = 0时, 序列值变化趋势不明显;当β < 0时, 序列值呈减小趋势. ij为年份, 2002 ≤ i < j ≤ 2021, xixj分别为ij时刻的物候指标值.

Mann-Kendall是一种检验数据是否具有单调趋势的非参数方法, 与Sen's斜率估计法结合可以很大程度上避免离群值的干扰[38]. 本研究采用90%置信度, 即|Z| ≥ 1.64时表示通过了显著性检验. 趋势分析利用R v4.2.1中的“trend”包实现.

1.3.3 物候变化多因素影响分析方法

以往研究表明, 水热条件包括气温、降水以及土壤水分等会直接影响植被物候[13, 23]. 考虑到青藏高原海拔高差较大, 导致水热条件的空间异质性强, 会间接影响植被物候[39]. 本研究选择聚焦气温、降水和土壤水分以及海拔等方面, 探究其对青藏高原草地物候的影响. 同时, 为进一步明确年内水热条件变化对植被物候的影响, 气温、降水以及土壤水分采用月尺度数据, 具体见表 2. SOS和EOS的影响因素为季前水热条件和海拔, 其中SOS的季前为1~5月, EOS的季前为6~11月, 二者的影响因素都作为LOS的影响因素. 对于SOS, 由于青藏高原从前一年9月开始不稳定冻结, 到当年3月开始不稳定融化[40], 期间降水多以积雪状态贮存且土壤水分基本保持稳定, 故将该阶段累计降水作为一个预测变量, 土壤水分从3月开始作为预测变量.

表 2 青藏高原草地物候影响因子1) Table 2 Factors influencing grassland phenology on the Qinghai-Xizang Plateau

采用可解释机器学习方法, 通过构建极限梯度提升树(extreme gradient boosting, XGBoost)模型并结合SHAP(Shapley additive explanations)重要性算法, 量化各影响因素对青藏高原草地物候变化的贡献大小, 并解释不同因素对物候变化的交互影响. XGBoost是一种兼具线性模型和树模型特点的优化模型, 相较于简单线性统计分析方法, XGBoost模型可以建立多因素与目标变量之间的非线性关系, 具有训练速度快, 不易过拟合等优点, 被广泛应用于解决各类回归和分类问题[41]. 但是, 它尚存在可解释性较差的缺点, 虽然可以显示预测变量对目标变量的整体影响程度, 但却不能解释预测变量的变化如何影响目标变量[42]. 此时, SHAP值[43]可以通过分离每个预测变量对目标变量的边际贡献弥补XGBoost模型的缺点[44]. 二者相结合, 不仅可以量化各个预测变量对目标变量的贡献, 还可以分析不同预测变量对目标变量的交互作用. 假设第i个目标变量的第j个预测变量为xij, 第i个目标变量的模型预测值为yi, 所有目标变量预测值的平均值为ybase, 则SHAP值服从如下公式[41]

(3)

式中, fxij)为第i个目标变量的第j个预测变量的SHAP值, 表示该预测变量对目标变量预测值的边际贡献. SHAP的绝对值反映了预测变量对模型预测的影响程度, 可用于影响因素的重要性排名, 绝对值越大, 说明该影响因素对植被物候的影响越大. 当fxij) > 0时, 说明预测变量对模型预测起正向作用, 即影响因素值增大会促进物候指标值增大;反之, 当fxij) < 0时, 说明该预测变量对模型预测有负面影响, 即影响因素值增大会导致物候指标值减小.

将3个物候指标作为目标变量, 相应的影响因素作为预测变量分别训练XGBoost模型, 然后使用SHAP值量化每个因素对物候的边际贡献. 通过SHAP重要性算法(每个预测变量的SHAP值的绝对平均值)对影响因素的重要性进行排序, 确定主要影响因素后利用SHAP值分析其对物候的影响及其交互效应. XGBoost模型利用python3.7中Scikit-learn机器学习库的“xgboost”函数来构建. SHAP分析通过python3.7中“shap”包计算.

2 结果与分析 2.1 青藏高原草地物候时空变化特征 2.1.1 草地生长期开始时间时空变化

2002~2021年青藏高原草地SOS的时空变化结果如图 2所示. 结果显示, 草地SOS多年均值在第127 d, 其中71.18%草地SOS分布在第110~150 d之间, 说明其多在4月下旬至5月末开始进入生长季. 空间上, SOS差异显著, 自东向西大体呈逐渐推迟的态势, 多年均值在东部地区多在第120 d之前, 而在西部则推迟到第140 d之后, 见图 2(a). 近20年间, 超半数(56.32%)的草地SOS呈提前趋势[图 2(b)], 提前速率均值为0.47 d·a-1, 主要发生在青藏高原的东部. 其中, 有11.14%的草地SOS出现显著提前[图 2(c), P < 0.1]. 43.68%的草地SOS出现推迟, 平均每年延迟0.48 d, 主要位于青藏高原西南地区, 其中8.71%表现出显著性(P < 0.1).

图 2 2002~2021年青藏高原草地SOS时空变化趋势 Fig. 2 Spatial and temporal trends in grassland SOS on the Qinghai-Xizang Plateau from 2002 to 2021

2.1.2 草地生长期结束时间时空变化

2002~2021年青藏高原草地EOS的时空变化计算结果见图 3. 青藏高原草地EOS多年均值在第303 d, 其中81.58%草地EOS分布在第290~320 d之间, 说明青藏高原大部分草地到10月下旬至11月中旬生长季结束. 空间上, 表现出从中部向四周逐渐推迟的趋势, 中部地区的EOS多年均值在第280 d左右, 西部和东部地区推迟到第310 d之后, 见图 3(a). 20年间有67.65%的草地表现出EOS延迟的趋势[图 3(b)], 平均延迟速率约为0.58 d·a-1, 其中18.96%的区域显著延迟[图 3(c), P < 0.1]. EOS提前的草地中仅有7.33%表现出显著性, 平均每年提前0.43 d.

图 3 2002~2021年青藏高原草地EOS时空变化趋势 Fig. 3 Spatial and temporal trends in grassland EOS on the Qinghai-Xizang Plateau from 2002 to 2021

2.1.3 草地生长季长度时空变化

2002~2021年青藏高原草地LOS时空变化计算结果见图 4. 青藏高原草地LOS多年均值为176 d, 其中92.73%的草地多年均值在120~210 d之间. 空间上, LOS呈现出自东向西逐渐缩短的趋势, 在东部地区多在190 d之上, 而在西部则缩短到不足150 d, 见图 4(a). 青藏高原有65.50%的草地LOS表现出延长趋势[图 4(b)], 平均延长速率为0.94 d·a-1, 其中18.15%的区域表现出显著延长[图 4(c), P < 0.1].

图 4 2002~2021年青藏高原草地LOS时空变化趋势 Fig. 4 Spatial and temporal trends in grassland LOS on the Qinghai-Xizang Plateau from 2002 to 2021

2.2 青藏高原草地物候影响因素分析

本研究利用SHAP值对XGBoost模型的结果进行了解释, 图 5显示了影响物候最重要的15个因素. 进一步, 研究选择对物候影响重要性排名前25%的因素作为主导因素进行重点分析. 3月的表层土壤水分对青藏高原草地SOS影响最大[图 5(a)], 其次是5月的降水量、海拔以及冻结期的降水量, 如图 6(a)所示, 3月表层土壤水分高、5月降水量大、海拔低以及冻结期降水量大都会促进青藏高原草地SOS的提前(SHAP值为负);EOS主要是由温度决定的[图 5(b)], 对其影响最大的因素是9月和10月的温度, 其次是海拔、7月表层(0~10 cm)土壤水分、9月10~40 cm土壤水分以及6月降水量, 如图 6(b)所示, SHAP值为正代表影响因素会促进EOS延迟, 温度(9月和10月)高、海拔低、生长季前期水分(7月表层土壤水分、6月降水)低、以及9月10~40 cm水分高都会促进EOS延迟. LOS的前3位影响因素与SOS较为一致[图 5(c)], 影响较大的依次是3月表层(0~10 cm)土壤水分、海拔、冻结期降水、11月温度、11月降水、8月降水、10月温度、2月温度、5月降水以及3月温度, 如图 6(c)所示, SHAP值为正表示会促进LOS延长, 海拔高和生长季前温度(2月和3月)低会导致LOS缩短, 其他影响因素值高会促进LOS延长.

y轴表示影响因素重要性排名;x轴表示各影响因素SHAP值的平均值 图 5 根据均值(|SHAP值|)对前15个因素的重要性排序 Fig. 5 Importance ranking of the top 15 variables according to the mean (|SHAP value|)

特征排名(y轴)表示影响因素的重要性;SHAP值(x轴)是表征模型中某个因素影响的统一索引;每行中的颜色条显示有关每个影响因素如何影响物候的更多详细信息, 红(蓝)点表示影响因素的值更高(更低) 图 6 输入变量对物候的影响 Fig. 6 Influence of input variables on phenology

笔者观察到影响因素对SOS的影响存在明显的阈值效应, 如图 7(a)所示, 当3月的表层土壤水分大于10 kg·m-2时, 对SOS的提前由抑制转为促进作用(SHAP值由正转负), 并在土壤水分为20 kg·m-2左右时促进作用达到最大, 而随着土壤水分进一步增加, 其促进作用开始减弱;相似地, 5月降水量在0~12 mm时, 降水越多对SOS的推迟作用越强, 而当降水超过12 mm, 则该推迟作用开始减弱并在降水达到35 mm时开始促进SOS提前, 降水量达到65 mm时促进作用达到了最大[图 7(b)]. 海拔高于4 500 m时会使得SOS推迟, 在5 000 m时这种影响到达峰值[图 7(c)]. 冻结期累计降水高于110 mm时, 开始促进SOS的提前, 但这种促进作用在累计降水超过250 mm时变得不确定[图 7(d)].

散点表示横坐标上影响因素对应的SHAP值;虚线为SHAP值为零的线 图 7 XGBoost(SOS)的SHAP依赖关系 Fig. 7 SHAP dependency plot for XGBoost (SOS)

主导因素变化对EOS的影响结果见图 8, 整体来看, 生长季晚期(9月和10月)温度高, 会促进EOS的延迟[图 8(a)8(b)]. 温度在9月高于5℃以及在10月高于-3℃都会促进EOS延迟, 但这种促进作用是有限的, 当9月温度达到8℃或10月温度达到-0.5℃时, 其对EOS延迟的促进作用达到最大且继续升温会保持不变. 海拔对EOS的影响不确定性较大, 存在明显的垂直分散[图 8(c)]. 7月表层土壤水分低(< 22 kg·m-2)以及9月10~40 cm深度土壤水分高(> 85 kg·m-2)会促进EOS延迟[图 8(d)8(e)], 9月土壤水分在85~100 kg·m-2之间时对EOS的影响变化剧烈. 生长季前期6月降水高于130 mm时会促进EOS延迟[图 8(f)].

散点表示横坐标上影响因素对应的SHAP值;虚线为SHAP值为零的线 图 8 XGBoost(EOS)的SHAP依赖关系 Fig. 8 SHAP dependency plot for XGBoost (EOS)

主导因素对LOS的影响结果见图 9, 3月表层土壤水分在15~25 kg·m-2之间时会促进LOS延迟, 这种促进作用在18 kg·m-2时达到峰值[图 9(a)]. 海拔低于4 000 m会促进LOS的延长[图 9(b)]. 通常来说, 较高的降水会促进LOS的延长[图 9(c)~9(f)]. 季前冻结期(上一年9月至当年3月)累计降水达到100 mm, 5月降水达到50 mm, 8月降水达到80 mm以及生长季结束时(11月)的降水达到5 mm都会促进LOS的延长. 生长季结束时(10~11月)温度越高越有利于LOS延长[图 9(g)9(h)], 10月温度达到-2℃或者11月温度达到-7℃由抑制LOS延长转为促进. 而生长季开始前(2~3月)温度低反而有利于LOS的延长[图 9(i)9(j)].

散点表示横坐标上影响因素对应的SHAP值;虚线为SHAP值为零的线 图 9 XGBoost(LOS)的SHAP依赖关系 Fig. 9 SHAP dependency plot for XGBoost (LOS)

2.3 青藏高原草地物候主导因素的交互影响

本研究认为降水、气温、土壤水分和地形等因素对青藏高原草地物候存在交互影响, 共同控制物候的改变. 为此, 以对物候变化贡献较高的主导因素为基础因子, 同时在其它主导因素中任选一个作为交互因子, 组成因子组, 探讨双主导因素对物候的交互影响. 本研究过程中发现并非所有影响因素之间都存在交互作用, 故以下仅展示分析存在明显交互作用的因素对物候的影响.

多因素影响SOS的交互效应图见图 10. 可以看出, 当3月表层土壤水分达到20 kg·m-2左右时, 其与其它主导因素的交互作用会对SOS变化起到显著影响. 其中, 5月降水和冻结期降水与3月表层土壤水分对SOS的交互影响相似[图 10(a)10(b)], 即当土壤水分达到20 kg·m-2左右时, 低降水量反而会促进SOS提前, 也就是说, 过高和过低的水分条件都不利于SOS提前. 另外当3月土壤水分超过20 kg·m-2时, 低海拔地区相较于高海拔地区土壤水分对SOS提前的促进作用更强[图 10(c)]. 如图 10(d)~10(f)所示, 当冻结期降水在50~110 mm时, 3月表层土壤水分高、5月降水量大以及海拔低的地区SOS推迟效应更加显著, 而当冻结期降水高于250 mm后, 较好的降水条件和较低的海拔会增强其对SOS提前的促进作用.

(a)3月0~10 cm土壤水分与5月降水对SOS的交互影响, (b)3月0~10 cm土壤水分与季前累计降水对SOS的交互影响, (c)3月0~10 cm土壤水分与海拔对SOS的交互影响, (d)季前累计降水与3月0~10 cm土壤水分对SOS的交互影响, (e)季前累计降水与5月降水对SOS的交互影响, (f)季前累计降水与海拔对SOS的交互影响, 散点表示横坐标上影响因素在产生交互作用的其他因素影响下对应的SHAP值;色柱表示与横坐标上影响因素产生交互作用的因素取值 图 10 多因素影响SOS的交互效应 Fig. 10 Interaction effects of multiple factors affecting SOS

图 11(a)~11(c)所示, 温度对EOS的影响明显受到水分条件的调控. 10月温度低于-3℃时会抑制EOS的延迟, 且水分较高时抑制作用更强;当10月温度高于0℃, 6月降水量大、7月表层土壤水分高以及9月10~40 cm土壤水分高的地区对EOS延迟的促进作用更强. 同样地, 海拔对EOS的影响表现出明显的垂直分散[图 11(d)~11(f)], 在3 000~4 500 m之间, 6月降水、7月表层土壤水分高以及9月10~40 cm土壤水分高会促进EOS延迟. 6月降水对EOS的影响受到温度调节[图 11(g)11(h)], 6月降水量低于75 mm时, 9月和10月温度越高, 对EOS延迟的促进作用越强;6月降水量超过75 mm后, 9、10月温度较低的地区EOS更延迟. 除此之外, 6月降水量较大时, 7月表层土壤水分对EOS几乎无影响[图 11(i)].

(a)10月温度与6月降水对EOS的交互影响, (b)10月温度与7月0~10 cm土壤水分对EOS的交互影响, (c)10月温度与9月10~40 cm土壤水分对EOS的交互影响, (d)海拔与6月降水对EOS的交互影响, (e)海拔与7月0~10 cm土壤水分对EOS的交互影响, (f)海拔与9月10~40 cm土壤水分对EOS的交互影响, (g)6月降水与9月温度对EOS的交互影响, (h)6月降水与10月温度对EOS的交互影响, (i)7月0~10 cm土壤水分与6月降水对EOS的交互影响, 散点表示横坐标上影响因素在产生交互作用的其他因素影响下对应的SHAP值;色柱表示与横坐标上影响因素产生交互作用的因素取值 图 11 多因素影响EOS的交互效应 Fig. 11 Interaction effects of multiple factors affecting EOS

降水、海拔和温度都会影响LOS对土壤水分的响应[图 12(a)~12(d)]. 3月表层土壤水分在12~22 kg·m-2之间时, 季前冻结期和生长季开始时的高降水更有利于LOS的延长, 但在土壤水分超过25 kg·m-2之后, 过高的水分条件导致了LOS缩短[图 12(a)12(b)]. 海拔和温度的影响较为一致, 土壤水分在17~22 kg·m-2之间时, 高温和低海拔条件下其对LOS延长的促进作用更大[图 12(c)12(d)]. 冻结期降水满足促进LOS延长的条件后(高于100 mm), 海拔和水分条件会对其促进作用大小产生影响[图 12(e)12(h)]. 海拔低的地区季前降水对LOS延长的促进作用更大[图 12(e)]. 生长季开始时降水(5月)和结束时降水(11月)较高更有利于LOS的延长[图 12(f)12(g)]. 3月表层土壤水分高也会增强季前降水对LOS延长的促进作用[图 12(h)].

(a)3月0~10 cm土壤水分与季前累计降水对LOS的交互影响, (b)3月0~10 cm土壤水分与5月降水对LOS的交互影响, (c)3月0~10 cm土壤水分与海拔对LOS的交互影响, (d)3月0~10 cm土壤水分与3月温度对LOS的交互影响, (e)季前累计降水与海拔对LOS的交互影响, (f)季前累计降水与5月降水对LOS的交互影响, (g)季前累计降水与11月降水对LOS的交互影响, (h)季前累计降水与3月0~10 cm土壤水分对LOS的交互影响, 散点表示横坐标上影响因素在产生交互作用的其他因素影响下对应的SHAP值;色柱表示与横坐标上影响因素产生交互作用的因素取值 图 12 多因素影响LOS的交互效应 Fig. 12 Interaction effects of multiple factors affecting LOS

3 讨论

基于2002~2021年青藏高原草地物候时空变化趋势结果, 整体上表现出SOS提前, EOS延迟以及LOS延长的趋势, 这与已有的研究结果基本一致[5, 13, 34, 45~49]. 本研究以海拔、气候条件和土壤水分条件作为预测变量, 物候作为目标变量构建XGBoost模型并使用SHAP法对模型结果进行了解释, 以探究影响青藏高原草地物候变化的主导因素以及其交互效应. 模型结果表明, 青藏高原草地SOS主要是由水分条件驱动, 3月表层土壤水分对SOS影响最大, 这可能是因为3月表层冻土已经开始不稳定融化[40], 底层冻土尚且处于冻结状态, 表层充足的水分保证了植被根系的生长活动. Wang等[45]的研究证明了冻土融化对物候的影响, 同时表明季前降水是影响物候的重要因素, 这与本研究的结果一致. 季前冻结期充足的降水为植被萌发提供了保障, 但过量的水分会对植被生长产生反作用. 青藏高原草地SOS主要是在5月, 该月充足的降水(≥ 35 mm)保证了植被的提前萌发. 水分的驱动作用还受到海拔的限制, 高海拔地区即使水分充足SOS也会较为延后, 这可能与高海拔地区温度较低有关[31, 50]. EOS的影响因素较多, 整体是由温度条件驱动, 受到水分条件限制, 这与先前的研究结果基本一致[4, 8, 23, 47]. 季前温度是青藏高原草地EOS的主要影响因素, 较高的温度有利于EOS的延迟, 但这种作用在温度高于阈值时受到了限制, 可能是因为植被绿化造成生长季晚期的水分胁迫[51], 水分条件已经无法满足植被生长的需求. Li等[47]的研究发现夏季植被生长旺盛往往会诱发更早的植被休眠, 所以在植被生产力较高的区域应加强管理, 适度放牧[46, 52]. LOS的主要影响因素与SOS较为一致, 整体而言水分条件对生长季的影响更大, 当水分满足条件后, 海拔和温度会对其影响程度产生调控作用.

本研究考虑了海拔、气温、降水以及土壤水分等因素, 定量解释了其对青藏高原草地物候变化的影响, 为进一步探究草地变化对气候的反馈作用以及对水碳循环的影响奠定了基础, 并为青藏高原生态恢复和放牧管理措施制定提供了理论支撑. 青藏高原范围广, 水文气候地貌条件具有空间异质性, 因此不同地区物候变化特征和主导因素可能有所差别, 且当前考虑的因素可能依然未能完整解释青藏高原草地物候的变化[48, 50, 53]. 未来应考虑纳入更多的影响因素, 根据物候变化的空间差异将青藏高原划分为多个特征相似的区域开展分区研究, 提高成果的针对性和精细度.

4 结论

(1)青藏高原草地SOS多在第110~150 d之间, 有56.32%的草地表现出SOS提前趋势, 提前速率约为0.47 d·a-1;EOS多在第290~320 d之间, 有67.65%的草地表现出EOS延迟趋势, 平均延迟速率约为0.58 d·a-1;LOS多在120~210 d之间, 有65.50%的草地表现出LOS延长趋势, 平均延长速率为0.94 d·a-1

(2)水热等因素对物候指标的影响存在明显阈值和峰值效应. 青藏高原草地SOS主要受水分条件影响, 3月表层土壤水分在10~25 kg·m-2之间时会促进SOS提前并在20 kg·m-2左右时促进作用达到峰值;EOS主要受温度影响, 生长季晚期9月和10月温度分别达到5℃和-3℃后其对EOS延迟的抑制作用转为促进作用, 分别高于8℃和-0.5℃时, 其促进作用达到峰值;LOS的主要影响因素与SOS较为一致, 3月表层土壤水分在15~25 kg·m-2之间会促进LOS延长, 在18 kg·m-2左右达到峰值, 通常降水量越大越促进LOS延长.

(3)影响青藏高原草地物候的因素之间存在明显的交互效应. 3月表层土壤水分达到20 kg·m-2后低降水和低海拔地区其对SOS提前的促进作用更强;10月温度高于0℃后较好的水分条件更有利于EOS延迟;3月表层土壤水分在12~22 kg·m-2之间时, 高降水相对于低降水更有利于LOS延长, 季前降水量达到100 mm后, 在海拔低、水分条件充足的区域其对LOS延长的促进作用更大.

参考文献
[1] Yao T D, Xue Y K, Chen D L, et al. Recent Third Pole's rapid warming accompanies cryospheric melt and water cycle intensification and interactions between monsoon and environment: multidisciplinary approach with observations, modeling, and analysis[J]. Bulletin of the American Meteorological Society, 2019, 100(3): 423-444. DOI:10.1175/BAMS-D-17-0057.1
[2] Li P, Liu Z L, Zhou X L, et al. Combined control of multiple extreme climate stressors on autumn vegetation phenology on the Tibetan Plateau under past and future climate change[J]. Agricultural and Forest Meteorology, 2021, 308-309. DOI:10.1016/j.agrformet.2021.108571
[3] Shen M G, Piao S L, Dorji T, et al. Plant phenological responses to climate change on the Tibetan Plateau: research status and challenges[J]. National Science Review, 2015, 2(4): 454-467. DOI:10.1093/nsr/nwv058
[4] An S, Chen X Q, Zhang X Y, et al. Precipitation and minimum temperature are primary climatic controls of alpine grassland autumn phenology on the Qinghai-Tibet Plateau[J]. Remote Sensing, 2020, 12(3): 431. DOI:10.3390/rs12030431
[5] Zheng Z T, Zhu W Q, Chen G S, et al. Continuous but diverse advancement of spring-summer phenology in response to climate warming across the Qinghai-Tibetan Plateau[J]. Agricultural and Forest Meteorology, 2016, 223: 194-202. DOI:10.1016/j.agrformet.2016.04.012
[6] Zhang G L, Zhang Y J, Dong J W, et al. Green-up dates in the Tibetan Plateau have continuously advanced from 1982 to 2011[J]. Proceedings of the National Academy of Sciences of the United States of America, 2013, 110(11): 4309-4314.
[7] Zhu W Q, Jiang N, Chen G S, et al. Divergent shifts and responses of plant autumn phenology to climate change on the Qinghai-Tibetan Plateau[J]. Agricultural and Forest Meteorology, 2017, 239: 166-175. DOI:10.1016/j.agrformet.2017.03.013
[8] Sun Q L, Chen G C, Li B L, et al. Multiple factors induce temporal shifts in herbaceous plant leaf senescence dates on the Qinghai-Tibetan Plateau[J]. Agricultural and Forest Meteorology, 2023, 334. DOI:10.1016/j.agrformet.2023.109441
[9] Piao S L, Liu Z, Wang T, et al. Weakening temperature control on the interannual variations of spring carbon uptake across northern lands[J]. Nature Climate Change, 2017, 7(5): 359-363. DOI:10.1038/nclimate3277
[10] Richardson A D, Keenan T F, Migliavacca M, et al. Climate change, phenology, and phenological control of vegetation feedbacks to the climate system[J]. Agricultural and Forest Meteorology, 2013, 169: 156-173. DOI:10.1016/j.agrformet.2012.09.012
[11] Piao S L, Liu Q, Chen A P, et al. Plant phenology and global climate change: current progresses and challenges[J]. Global Change Biology, 2019, 25(6): 1922-1940. DOI:10.1111/gcb.14619
[12] Cheng M, Wang Y, Zhu J X, et al. Precipitation dominates the relative contributions of climate factors to grasslands spring phenology on the Tibetan Plateau[J]. Remote Sensing, 2022, 14(3). DOI:10.3390/rs14030517
[13] Shen M G, Wang S P, Jiang N, et al. Plant phenology changes and drivers on the Qinghai-Tibetan Plateau[J]. Nature Reviews Earth & Environment, 2022, 3(10): 633-651.
[14] Xu C Y, Liu D, Wang X Y, et al. Shifting from a thermal-constrained to water-constrained ecosystem over the Tibetan Plateau[J]. Frontiers in Plant Science, 2023, 14. DOI:10.3389/fpls.2023.1125288
[15] Shen M G, Piao S L, Cong N, et al. Precipitation impacts on vegetation spring phenology on the Tibetan Plateau[J]. Global Change Biology, 2015, 21(10): 3647-3656. DOI:10.1111/gcb.12961
[16] Ma P F, Zhao J X, Zhang H Z, et al. Increased precipitation leads to earlier green-up and later senescence in Tibetan alpine grassland regardless of warming[J]. Science of the Total Environment, 2023, 871. DOI:10.1016/j.scitotenv.2023.162000
[17] Chen Z G, Shen M G, Jiang N, et al. Daytime warming strengthened delaying effect of precipitation on end of the vegetation growing season on the Tibetan Plateau[J]. Science of the Total Environment, 2023, 892. DOI:10.1016/j.scitotenv.2023.164382
[18] 徐勇, 戴强玉, 黄雯婷, 等. 2000~2020年西南地区植被NDVI时空变化及驱动机制探究[J]. 环境科学, 2023, 44(1): 323-335.
Xu Y, Dai Q Y, Huang W T, et al. Spatio-temporal variation in vegetation cover and its driving mechanism exploration in Southwest China from 2000 to 2020[J]. Environmental Science, 2023, 44(1): 323-335.
[19] Wu K, Chen J H, Yang H, et al. Spatiotemporal variations in the sensitivity of vegetation growth to typical climate factors on the Qinghai-Tibet Plateau[J]. Remote Sensing, 2023, 15(9). DOI:10.3390/rs15092355
[20] Wang C Y, Wang J N, Wang X F, et al. Phenological changes in alpine grasslands and their influencing factors in seasonally frozen ground regions across the Three Parallel Rivers Region, Qinghai-Tibet Plateau[J]. Frontiers in Earth Science, 2022, 9. DOI:10.3389/feart.2021.797928
[21] Che M, Chen B, Innes J L, et al. Spatial and temporal variations in the end date of the vegetation growing season throughout the Qinghai-Tibetan Plateau from 1982 to 2011[J]. Agricultural and Forest Meteorology, 2014, 189-190: 81-90. DOI:10.1016/j.agrformet.2014.01.004
[22] 孔冬冬, 张强, 黄文琳, 等. 1982-2013年青藏高原植被物候变化及气象因素影响[J]. 地理学报, 2017, 72(1): 39-52.
Kong D D, Zhang Q, Huang W L, et al. Vegetation phenology change in Tibetan Plateau from 1982 to 2013 and its related meteorological factors[J]. Acta Geographica Sinica, 2017, 72(1): 39-52.
[23] Wang X Y, Wu C Y, Liu Y, et al. Earlier leaf senescence dates are constrained by soil moisture[J]. Global Change Biology, 2023, 29(6): 1557-1573. DOI:10.1111/gcb.16569
[24] 徐丽娇, 胡泽勇, 赵亚楠, 等. 1961—2010年青藏高原气候变化特征分析[J]. 高原气象, 2019, 38(5): 911-919.
Xu L J, Hu Z Y, Zhao Y N, et al. Climate change characteristics in Qinghai-Tibetan Plateau during 1961-2010[J]. Plateau Meteorology, 2019, 38(5): 911-919.
[25] 管琪卉, 丁明军, 张华敏. 青藏地区高寒草地春季物候时空变化及其对气候变化的响应[J]. 山地学报, 2019, 37(5): 639-648.
Guan Q H, Ding M J, Zhang H M. Spatiotemporal variation of spring phenology in alpine grassland and response to climate changes on the Qinghai-Tibet, China[J]. Mountain Research, 2019, 37(5): 639-648.
[26] Chen X Q, An S, Inouye D W, et al. Temperature and snowfall trigger alpine vegetation green-up on the world's roof[J]. Global Change Biology, 2015, 21(10): 3635-3646. DOI:10.1111/gcb.12954
[27] 杨元合, 朴世龙. 青藏高原草地植被覆盖变化及其与气候因子的关系[J]. 植物生态学报, 2006, 30(1): 1-8.
Yang Y H, Piao S L. Variations in grassland vegetation cover in relation to climatic factors on the Tibetan Plateau[J]. Journal of Plant Ecology, 2006, 30(1): 1-8. DOI:10.3321/j.issn:1005-264X.2006.01.001
[28] 王立景, 肖燚, 孔令桥, 等. 青藏高原草地承载力空间演变特征及其预警[J]. 生态学报, 2022, 42(16): 6684-6694.
Wang L J, Xiao Y, Kong L Q, et al. Spatiotemporal patterns and early-warning of grassland carrying capacity in the Qinghai-Tibet Plateau[J]. Acta Ecologica Sinica, 2022, 42(16): 6684-6694.
[29] Chen J, Jönsson P, Tamura M, et al. A simple method for reconstructing a high-quality NDVI time-series data set based on the Savitzky-Golay filter[J]. Remote Sensing of Environment, 2004, 91(3-4): 332-344. DOI:10.1016/j.rse.2004.03.014
[30] Jönsson P, Eklundh L. TIMESAT - a program for analyzing time-series of satellite sensor data[J]. Computers & Geosciences, 2004, 30(8): 833-845.
[31] 王泰华, 杨大文. 2001—2020年三江源冻土区植被物候变化特征分析[J]. 冰川冻土, 2023, 45(2): 711-723.
Wang T H, Yang D W. Variation characteristics of vegetation phenology during 2001-2020 on frozen ground zone of the Three Rivers Source Region[J]. Journal of Glaciology and Geocryology, 2023, 45(2): 711-723.
[32] 赵心睿, 刘冀, 杨少康, 等. 北方地区典型林草地物候时空变化特征及其对气象因子的响应[J]. 生态学报, 2023, 43(9): 3744-3755.
Zhao X R, Liu J, Yang S K, et al. Spatio-temporal variations of typical woodland and grassland phenology and its response to meteorological factors in Northern China[J]. Acta Ecologica Sinica, 2023, 43(9): 3744-3755.
[33] 何宝忠, 丁建丽, 李焕, 等. 新疆植被物候时空变化特征[J]. 生态学报, 2018, 38(6): 2139-2155.
He B Z, Ding J L, Li H, et al. Spatiotemporal variation of vegetation phenology in Xinjiang from 2001 to 2016[J]. Acta Ecologica Sinica, 2018, 38(6): 2139-2155.
[34] Liu X G, Chen Y N, Li Z, et al. Driving forces of the changes in vegetation phenology in the Qinghai-Tibet Plateau[J]. Remote Sensing, 2021, 13(23). DOI:10.3390/rs13234952
[35] Wang J, Liu D S. Vegetation green-up date is more sensitive to permafrost degradation than climate change in spring across the northern permafrost region[J]. Global Change Biology, 2022, 28(4): 1569-1582. DOI:10.1111/gcb.16011
[36] 常清, 王思远, 孙云晓, 等. 青藏高原典型植被生长季遥感模型提取分析[J]. 地球信息科学学报, 2014, 16(5): 815-823.
Chang Q, Wang S Y, Sun Y X, et al. The remote sensing monitoring model of the typical vegegtation phenology in the Qinghai-Tibetan Plateau[J]. Journal of Geo-Information Science, 2014, 16(5): 815-823.
[37] 唐志光, 邓刚, 胡国杰, 等. 亚洲高山区积雪物候时空动态及其对气候变化的响应[J]. 冰川冻土, 2021, 43(5): 1400-1411.
Tang Z G, Deng G, Hu G J, et al. Spatiotemporal dynamics of snow phenology in the High Mountain Asia and its response to climate change[J]. Journal of Glaciology and Geocryology, 2021, 43(5): 1400-1411.
[38] 雷茜, 胡忠文, 王敬哲, 等. 1985-2015年中国不同生态系统NDVI时空变化及其对气候因子的响应[J]. 生态学报, 2023, 43(15): 6378-6391.
Lei Q, Hu Z W, Wang J Z, et al. Spatiotemporal dynamics of NDVI in China from 1985 to 2015: ecosystem variation, regional differences, and response to climatic factors[J]. Acta Ecologica Sinica, 2023, 43(15): 6378-6391.
[39] Shen M G, Zhang G X, Cong N, et al. Increasing altitudinal gradient of spring vegetation phenology during the last decade on the Qinghai-Tibetan Plateau[J]. Agricultural and Forest Meteorology, 2014, 189-190: 71-80. DOI:10.1016/j.agrformet.2014.01.003
[40] 彭小清, 田伟伟, 李璇佳, 等. 青藏高原和环北极冻土变化研究进展[J]. 冰川冻土, 2023, 45(2): 521-534.
Peng X Q, Tian W W, Li X J, et al. Research progress on changes in frozen ground on the Qinghai-Tibet Plateau and in the circum-Arctic region[J]. Journal of Glaciology and Geocryology, 2023, 45(2): 521-534.
[41] Wang S, Peng H, Hu Q, et al. Analysis of runoff generation driving factors based on hydrological model and interpretable machine learning method[J]. Journal of Hydrology: Regional Studies, 2022, 42. DOI:10.1016/j.ejrh.2022.10113
[42] Wang S, Peng H, Liang S K. Prediction of estuarine water quality using interpretable machine learning approach[J]. Journal of Hydrology, 2022, 605. DOI:10.1016/j.jhydrol.2021.127320
[43] Lundberg S M, Lee S L. A unified approach to interpreting model predictions[A]. In: Proceedings of the 31st International Conference on Neural Information Processing Systems[C]. Long Beach: Curran Associates Inc., 2017. 4768-4777.
[44] Zhang W Q, Luo G P, Hamdi R, et al. Can gross primary productivity products be effectively evaluated in regions with few observation data?[J]. GIScience & Remote Sensing, 2023, 60(1). DOI:10.1080/15481603.2023.2213489
[45] Wang T H, Yang D W, Zheng G H, et al. Possible negative effects of earlier thaw onset and longer thaw duration on vegetation greenness over the Tibetan Plateau[J]. Agricultural and Forest Meteorology, 2022, 326. DOI:10.1016/j.agrformet.2022.109192
[46] Liu H, Liu S L, Wang F F, et al. Management practices should be strengthened in high potential vegetation productivity areas based on vegetation phenology assessment on the Qinghai-Tibet Plateau[J]. Ecological Indicators, 2022, 140. DOI:10.1016/j.ecolind.2022.108991
[47] Li P, Peng C H, Wang M, et al. Dynamics of vegetation autumn phenology and its response to multiple environmental factors from 1982 to 2012 on Qinghai-Tibetan Plateau in China[J]. Science of the Total Environment, 2018, 637-638: 855-864. DOI:10.1016/j.scitotenv.2018.05.031
[48] 吉珍霞, 裴婷婷, 陈英, 等. 青藏高原草地物候动态及其对驱动因子的响应[J]. 草业科学, 2023, 40(1): 4-14.
Ji Z X, Pei T T, Chen Y, et al. Grassland phenological dynamics and its response to driving factors on the Qinghai-Tibet Plateau[J]. Pratacultural Science, 2023, 40(1): 4-14.
[49] 李佳, 唐志光, 邓刚, 等. 2001-2020年青藏高原草地物候变化遥感监测[J]. 水土保持研究, 2023, 30(4): 265-274.
Li J, Tang Z G, Deng G, et al. Remote sensing monitoring of grassland phenological changes in the Qinghai-Tibetan Plateau during 2001-2020[J]. Research of Soil and Water Conservation, 2023, 30(4): 265-274.
[50] Teng H F, Luo Z K, Chang J F, et al. Climate change-induced greening on the Tibetan Plateau modulated by mountainous characteristics[J]. Environmental Research Letters, 2021, 16(6). DOI:10.1088/1748-9326/abfeeb
[51] Lian X, Piao S L, Li L Z X, et al. Summer soil drying exacerbated by earlier spring greening of northern vegetation[J]. Science Advances, 2020, 6(1). DOI:10.1126/sciadv.aax0255
[52] Dong S K, Shang Z H, Gao J X, et al. Enhancing sustainability of grassland ecosystems through ecological restoration and grazing management in an era of climate change on Qinghai-Tibetan Plateau[J]. Agriculture, 2020, 287. DOI:10.1016/j.agee.2019.106684
[53] Lian X, Piao S L, Chen A P, et al. Seasonal biological carryover dominates northern vegetation growth[J]. Nature Communications, 2021, 12(1). DOI:10.1038/s41467-021-21223-2