鲍舒琪,张成福,冯霜,贺帅,苗林
(内蒙古农业大学 沙漠治理学院,呼和浩特 010000)
地表温度(land surface temperature,LST)与气温通过陆气作用相互影响,气温变化会引起地表温度的变化,直接影响土壤水热过程,从而影响植被生长和枯枝落叶的分解等一系列生态过程[1]。目前主要通过遥感技术获取大范围地表温度,能够获取实时数据,且分辨率高成本低,但不能用来预测气候变化规律及响应[2]。李召良等[3]在地表温度遥感反演方法研究进展中表明卫星遥感是获取区域或全球尺度地表温度时空分布的主要手段,但卫星传感器只能拍摄实时数据,不能预测未来地表温度的时空变化。
深度学习(deep learning,DL)算法具有强大的非线性映射能力与学习能力,广泛应用于较复杂生态环境的模拟预测。Shatnawi等[4]采用非线性自回归外生人工神经网络模型成功模拟和预测了约旦北部2000—2016年间的地表温度变化。苏扬等[5]用多时相特征连接卷积神经网络双向重建模型(MTFC-CNN)重建轨道间隙区域的地表温度值,其效果优于传统的样条空间插值和时间线性回归方法。张义峥等[6]基于降尺度残差网络构建了美国密苏里州的地表温度深度学习模型,其降尺度效果优于经典传统方法且稳定性更强。
下垫面的地形起伏状况是影响地表温度分布特征的重要因素,其复杂的自然地理环境使得地表温度在空间与时间上存在差异性[7-8]。肖尧等[9]在复杂地表地表温度反演研究进展中表明由于山地地貌地形起伏大、空间异质性强,出现长时序数据获取难、地表温度反演难等问题,导致大部分观测站架设在平坦地区,大多数地表温度反演算法建立在地面平坦、地形因素影响小的区域。山地地表温度的反演需要将地形因素融入到算法,而深度学习方法能够学习山地地形因素与地表温度的非线性关系特征,准确模拟山地地区的地表温度。
内蒙古大青山位于半干旱区,相对高差大,不同坡向植被分布差异较大。本文的研究目的是利用深度学习方法模拟预测该山地地表温度的空间分布特征,并基于深度学习模型分析特定地形和植被条件下地表温度的变化规律,以期为研究全球气候变化对山地水热过程和植被的影响提供一种定量的分析方法。
内蒙古大青山位于阴山山脉中段,海拔高度在960~2 322 m之间。大青山山体呈现东西走向,北坡较平缓,有交错分布的低山丘陵和盆地与内蒙古高原相连,南坡陡峭,为剥蚀堆积形成的构造断裂地形[10]。大青山位于温带大陆性半干旱季风气候区,区内有山地森林、灌丛和草原植被,是土默川平原及黄河上中游地区重要的水源补给区。
数据包括ASTER DEM 30 m高程数据、2019—2021年5—9月生长季的MODIS产品数据(表1)和气象数据。经辐射定标和大气校正等预处理的MOD09A1、MOD11A2、MOD13A1产品数据来源于NASA网站(https://ladsweb.modaps.eosdis.nasa.gov/)。ASTER GDEM 30 m高程数据来源于地理空间数据云(http://www.gscloud.cn/)。气象数据来源于中国气象数据网(http://data.cma.cn/)的6个气象站点的气候资料数据集。
表1 MODIS产品计算公式
用NASA提供的MRT(MODIS Reprojection Tool)软件对MODIS数据进行分层处理,提取LST、归一化植被指数(normalized difference vegetation index,NDVI)、地表反照率等数据,用ENVI、ArcGIS软件进行计算、重投影、格式转换及图像镶嵌裁剪。
用ArcGIS软件的空间分析工具基于高程数据计算海拔、坡度、坡向等地形因子数据。用ArcGIS软件的反距离插值法将6个气象站点数据插值为与研究区像元大小一致的空间栅格图像,与研究区环境因子的各像元相对应。用ArcGIS软件的投影变换功能将所有空间数据统一坐标系为WGS 1984 UTM、空间分辨率为500 m。
选择影响LST的环境因子构建模型时,为防止数据集中存在大量冗余数据增加算法的成本,避免过度拟合,使用特征选择算法挑选与目标变量高度相关的特征变量,能有效提高算法的性能[11-12]。本文在构建LST模型之前,运用决策树(decision tree,DT)构造算法中的分类回归树(classification and regression tree,CART)[13],对可能引起LST变化的环境因子(如气温、海拔、植被盖度等)进行特征重要度分析,选取特征贡献度较大的因子进行LST模型构建。
1)深度学习方法介绍。深度学习本质是利用多层的非线性信息处理来进行监督或无监督的特征提取和转换,并进行模式分析和分类,可以通过多层学习训练获得更高层的特征信息[14]。深度学习神经网络由输入层、一个或多个隐藏层、输出层和每层中的神经元组成,每个神经元有一个权值和偏置对应上一层的神经元,输入层神经元数量为输入特征变量的数量,输出层的神经元数量为训练过程中给予神经网络的输出数量,隐藏层神经元数量和隐藏层数量的正确选择取决于输入输出关系的复杂程度[15]。
2)LST模型构建。本研究在数据科学管理软件Anaconda的Jupyter notebook平台,利用以Tensorflow 2.0为后端的深度学习库Keras进行模型的构建运行与测试。
从环境因子(如气温、海拔、NDVI等)和LST的空间栅格图像上随机选取1 200个样本点,提取2019—2021年3年数据共3 600个样本作为模型数据集。以决策树模型选取的环境因子作为输入层数据集,LST作为输出层数据集建立深度学习模型。输入层数据集的70%划分为训练集,20%划分为测试集,10%划分为预测集。模型的激活函数设置为Relu,优化器选取RMSprop。经多次训练后,确定藏层、神经元个数和迭代次数,构建最优LST模型。
3)模型评估。用回归损失函数均方误差(mean squared error,MSE)和平均绝对误差(mean absolute error,MAE)评估模型预测值与真实值之间的差距,衡量模型模拟性能和预测精度。MSE用来评价数据的变化程度,MSE值越小说明模型预测具有更好的精确度,MAE可以准确反映实际预测误差情况。
4)单个因子对LST空间变化的分析。在已构建的深度学习模型基础上,将待分析的单个环境因子设置为若干等级,同时将其他环境因子固定为区域均值(表2),以此计算LST随这一因子的变化规律。
表2 各环境因子均值表
研究区NDVI空间范围在0.14~0.91,等级间距设置为0.05;日均气温空间值范围在14.18~20.57 ℃,等级间距设置为0.25 ℃;海拔空间范围在960~2 322 m,等级间距设置为50 m;地表反照率空间范围在0.06~0.21,等级间距设置为0.01;坡度空间范围在0°~66.5°,等级间距设置为5°。
5)不同植被覆盖情景下LST空间变化对于单个因子的响应。为分析各环境因子在不同植被覆盖度下对LST分布的影响,本文根据《土壤侵蚀分级标准》中的植被覆盖度分级标准,设定4种NDVI情景模式:NDVI为0.2(低植被覆盖)、0.4(中低植被覆盖)、0.5(中等植被覆盖)和0.8(高植被覆盖)。基于构建的LST深度学习模型,分析在不同植被覆盖情景下,LST随各个影响因子的变化规律,同样设置其他特征变量为均值(表2)。
选取对LST空间分布变化有影响的7个环境因子(平均气温、高程、坡度、坡向、剖面曲率、NDVI和地表反射率),经决策树模型筛选出相对贡献率较高的因子为平均气温、坡度、地表反照率、海拔和NDVI(贡献率分别为0.008、0.087、0.09、0.227、0.588),相对贡献率较低的因子为剖面曲率和坡向因子(贡献率均为0)。环境因子重要度见图1。
图1 地表温度空间分布的环境因子贡献度
图2为模型预测模拟的LST值与观测值对比图。预测结果表明模型精度最高(MAE为0.60 ℃,MSE为0.65 ℃)的隐藏层个数为3层,训练迭代次数为30 000次。图2(a)表明预测值与观测值对比较吻合,构建的深度学习模型能很好地预测模拟LST的变化趋势,预测值与观测值的变化波动规律一致,模型预测结果较准确。由图2(b)可知深度学习模型预测的LST值与LST观测值散点分布较集中,比较贴合在相关线附近,且散点图决定系数R2为0.89,说明深度学习构建的LST模型预测模拟效果较好,拟合效果较优,拟合稳定性较高。因此使用深度学习方法对LST进行模拟具有一定的可靠性。
图2 模型预测模拟的地表温度值与观测值对比图
图3是2021年LST的观测值与预测值空间分布图。对比发现观测值与预测值的空间分布特征相近,每段LST区间分布区域都较为吻合,且地表温度最低值仅相差0.32 ℃,最高值仅相差0.83 ℃,模拟值比观测值空间分布较清晰。研究区西部、北部和南部地带LST最高(30.5~37.36 ℃),因其位于大青山山脚,海拔相对较低且植被覆盖较小(NDVI低于0.5);中部区域LST最低(20.38~28.5 ℃),因其是呼和浩特大青山主要山体部分,海拔高且植被覆盖大(NDVI达0.7);由于其余地带海拔介于上述两个地带之间,LST为28.5~30.5 ℃。
图3 地表温度预测值与观测值的空间分布对比图
NDVI与LST的决定系数R2为0.95;随NDVI增加LST呈下降趋势,NDVI每增加0.1时LST降低1.41 ℃(图4(a))。平均气温与LST的决定系数R2为0.86;随平均气温的增加LST呈上升趋势,平均气温每增加1 ℃时LST上升0.33 ℃(图4(b))。海拔与LST的决定系数R2为0.85;随海拔高度的增加LST呈下降趋势,海拔高度每增加50 m时LST减小0.18 ℃(图4(c)),在海拔1 600 m以上的LST下降速率较大。坡度与LST的决定系数R2为0.61;随坡度的增大LST整体呈下降趋势,坡度每增加1°时LST下降0.23 ℃(图4(d)),在45°以下坡度的LST下降速率较小,在45°以上坡度的LST下降速率较大。地表反照率与LST的决定系数R2为0.88;随地表反照率的增加LST呈上升趋势,地表反照率每增加0.01时LST升高0.02 ℃(图4(e))。
图4 地表温度随单一控制变量变化的规律图
图5是4种NDVI情景下模拟LST随单一环境因子的变化规律,表3是图5中每个线性曲线的函数参数信息表。图5(a)是在4种NDVI的情景下模拟LST随平均气温(14.5~20.5 ℃)单一变量下的变化情况。随NDVI增加,LST随平均气温上升而增加的增速呈下降趋势。当NDVI为0.2时,LST增加速率为0.52 ℃/℃;NDVI为0.4和0.5时,LST增加速率分别为0.09 ℃/℃和0.31 ℃/℃;NDVI为0.8时,LST增加速率为0.04 ℃/℃(图5(a)、表3)。在NDVI为0.2情景下,LST受平均气温影响较其他情景模式明显,随着平均气温的增加,LST在32.53~35.38 ℃之间逐渐升高,而在NDVI为0.4和0.5两种情景模式下,LST受平均气温影响差异性较小。
图5 地表温度在不同植被覆盖情景下随平均气温、海拔、坡度、地表反照率的变化规律图
表3 地表温度在不同植被覆盖情景下随平均气温、海拔、坡度、地表反照率变化的线性函数参数表
图5(b)是在4种NDVI情景下模拟LST随海拔(960~2 322 m)单一变量下的变化情况。随NDVI增加,LST随海拔增大而下降的降速呈先下降后增大趋势。NDVI为0.2时,LST下降速率为0.18 ℃/50 m;NDVI为0.4和0.5时,LST下降速率分别为0.12 ℃/50 m和0.14 ℃/50 m;NDVI为0.8时,LST下降速率为0.27 ℃/50 m;在海拔1 200 m高度出现最大值31.87 ℃,在1 500 m以上高度LST下降至最低,这是因为在研究区海拔较低处是比较平坦的城郊山脚地区,地表植被覆盖较小,地表吸收的太阳辐射能量较多,地表温度较大,在研究区海拔较高处则是大青山山体部分,植被覆盖度较大,地表吸收的太阳辐射能量就越少,使地表温度越低。
图5(c)是4种NDVI情景下模拟LST随坡度(0°~65°)单一变量下的变化情况。随NDVI增加,LST随坡度增大而下降的降速呈增大趋势。当NDVI为0.2时,LST下降速率为0.02 ℃/(°);NDVI为0.4和0.5时,LST下降速率分别为0.07 ℃/(°)和0.10 ℃/(°);NDVI为0.8时,LST下降速率为0.15 ℃/(°)。
图5(d)是4种NDVI情景下模拟LST随地表反照率(0.06~0.21)单一变量下的变化情况。随NDVI增加,LST随地表反照率增大而增加的增速呈先下降后增大趋势。NDVI为0.2时,每增加0.1地表反照率LST增加0.31 ℃;NDVI为0.4和0.5时,每增加0.1地表反照率LST分别增加0.26 ℃和0.28 ℃;NDVI为0.8时,每增加0.1地表反照率LST增加0.51 ℃。
本研究使用深度学习方法模拟了内蒙古大青山地表温度的空间分布,得到了很好的模拟结果,表明该方法模拟复杂地貌下地表温度空间分布具有可行性。本研究使用决策树算法对环境因子进行特征重要度分析。与汪子豪等[16]基于BP神经网模拟地表温度时不对输入层数据进行筛选的研究结果(RMSE为0.98 ℃)相比,使用决策树算法减少了输入变量的数据冗余,提高了模拟精度。使用深度学习方法构建LST模型时,需要选择合理的输入因子、适当的隐藏层节点数以及训练次数才能获得较高的模型精度[17]。
基于深度学习对地表温度的模型构建多在地形平坦区域研究[18-19],而在高山地区由于其复杂的地形情况,利用深度学习是否能准确估算模拟高山地区的地表温度的研究较少。
地表温度的时空分布特征受气象温度、植被覆盖和地形地貌等多个环境因素的共同影响,且这些环境因素之间也相互作用,使得地表温度在区域尺度上的时空分布状况具有一定的复杂性和差异性[20]。研究区东南部海拔高、坡度陡、植被覆盖度高共同导致地表温度最低;北部和西部区域海拔低、坡度平缓、植被覆盖度低共同导致其地表温度较高。在以往的研究中更多的是简单分析了地表温度与环境因子的相关性特征[21-22],未考虑到其他环境因素的潜在影响,而本文在分析地表温度随某一环境因子的变化规律时,利用深度学习方法控制其他环境因子为均值,使结果更具有准确性。
本研究发现,影响大青山地表温度的主要环境因子有NDVI、海拔、地表反照率、坡度和气象站平均气温,其中影响最大的因子是NDVI,其次是海拔和坡度,而坡向影响最小,这与罗瑶等[23]的研究结果一致。NDVI是地表温度影响最大的环境因子,LST随NDVI的增大而降低,是由于植物可以吸收和转化太阳辐射从而降低地表温度[24]。本研究中LST随地表反照率的增大而升高,表明地表反照率在15%~28%之间时,地表温度随反照率的增大而升高;地表反照率大于40%时,地表温度随反照率的增大而减小[25],本文由于大青山的地表反照率范围在0.06~0.21之间,地表反照率非常小,所以仅出现了地表温度随地表反照率增大而升高的情况。在大青山复杂的地貌环境条件下,其特殊的地形因子(海拔、坡度、坡向等)也是影响地表温度的重要环境因子,在地表温度的空间分布特征中有着直接或间接影响。不同地形因子对地表温度的影响存在差异,大青山海拔与地表温度的关系不是简单的线性关系,LST随海拔的增大而降低,在海拔1 600 m以上LST下降明显。这是由于在对流层海拔越高大气吸收的地表长波辐射能量越少,大气储存的热量越低气温就越低,从而影响地表与空气之间的热交换[26]。坡度与坡向主要因其不同坡度坡向条件下坡面吸收的太阳光照的大小不同,从而影响坡面LST的大小[27]。LST随坡度的增大而降低,在0°~45°时LST下降速度较小,在45°以上时LST下降速度显著。造成这种现象的原因是不同坡度的太阳入射强度和反射率不同,导致地表吸收的能量不同,由于险坡吸收的太阳辐射较小,导致LST较小。本研究发现,坡向因子与LST之间的关系不明显,这可能与本研究使用的MODIS数据分辨率较大不能够反映当地坡向变化有关。今后在使用深度学习方法模拟山地空间温度变化时,应考虑选择空间分辨率更高的数据。
地表温度受下垫面及气候环境中多种因素的共同影响,是多种因子相互作用的结果,本文由于数据的获取难度以及研究水平的局限性仅选取了部分因子进行训练,但在选择影响因子进行模型训练过程时是否还有更多的有效因子还未被验证,在后续研究中将充分考虑多种环境因子进行模型的构建。
本文基于环境因子数据集利用深度学习算法模拟预测了大青山LST的空间分布,并分析单一环境因子和植被变化情景下各环境因子与LST的关系,得到如下结论。
1)通过决策树模型运算发现,植被覆盖NDVI对LST的贡献度最大(0.588),其次是海拔(0.227)、坡度(0.087)、地表反照率(0.09)和平均气温(0.008),最低的是剖面曲率(0)和坡向(0)因子。
2)利用深度学习神经网络构建的LST模型精度MAE达0.60 ℃,MSE达0.65 ℃;预测模拟的决定系数R2达0.89。
3)LST随NDVI的增加而降低,NDVI每增加0.1时LST下降0.41 ℃;LST随平均气温增加而升高,增加速率为0.33 ℃/℃;LST随海拔升高而降低,下降速率为0.18 ℃/50 m;LST随坡度增大而降低,下降速率为0.23 ℃/(°);LST随地表反照率增大而增加,地表反照率每增加0.01时LST升高0.31 ℃。
4)随NDVI的增大,LST随平均气温上升而升高的速率呈下降趋势,随海拔升高而下降的速率呈先减小后增大的趋势,随坡度增加而下降的速率呈增大趋势,随地表反照率增大而升高的速率呈先减小后增大的趋势。