基于XGBoost的LST空间降尺度方法

颜佳楠，陈 虹，姚光林，吴 骅

(1.中国科学院地理科学与资源研究所 资源与环境信息系统国家重点实验室，北京 100101；2.中国自然资源航空物探遥感中心，北京 100083；3.山东省第八地质矿产勘查院，山东 日照 276826)

0 引言

地表温度(Land Surface Temperature，LST)是环境监测中的重要参数，可用于探测城市热岛效应[1],监测森林火灾[2]和土地覆盖变化[3]，评估地表干旱程度[4]及管理水资源[5]。通过遥感卫星可以获取到大范围的LST，但是现有的卫星、遥感器受制于设备，无法获取兼顾时间分辨率高和空间分辨率高的LST。例如中分辨率成像光谱仪(Moderate Resolution Imaging Spectroradiometer，MODIS)一天内可以获取4次1 000 m分辨率的LST，而Landsat8卫星每16天获取一次100 m分辨率的LST。LST遥感产品时空分辨率相互制约的局限性影响LST的应用。因此，许多学者尝试针对MODIS数据，发展LST空间降尺度的方法，以期获得时空分辨率较高的LST遥感产品。

目前LST空间降尺度方法主要有2种：一种是基于融合的方法；另一种是基于回归的方法。基于融合的方法通过融合高时间分辨率和高空间分辨率的影像，预测高分辨率LST[6]。典型基于回归的方法是Kustas等提出的DisTrad算法，其针对的是植被覆盖区域，以NDVI为回归核，获取回归核与低分辨率温度之间的回归关系，并将此函数关系应用于高分辨率回归核与高分辨率LST[7]。Agam等人[8]基于TsHARP算法，将二次多项式回归改为一次多项式回归，把回归核换作植被覆盖度(FC)，精度有所提高。然而，LST-FC(NDVI)的回归关系在异质性地区表现较差[9]。针对这个问题，Essa等人[10]将回归核改为SAVI，NDBI，NDBal，NDMI等15个光谱指标，将此方法推广到异质性区域。Duan等人[11]以NDVI和DEM为回归核，利用地理加权回归进行LST降尺度。

为了提高模型的精度和适用范围，有学者利用机器学习的方法获取回归核以及LST之间的关系,例如Gao[12]建立回归树获取LST与短波光谱反射率之间的关系。Hutengs[13]提出了Basic-RF算法，用随机森林的方法，获取LST与数字高程数据、土地覆盖类型、红波段和近红外波段的反射率之间的关系。Li[14]基于BasicRF提出了RRF算法，相比BasicRF，增加了回归核种类，并对模型进行鲁棒性验证，在多个不同的研究区获得了较高的精度。

考虑目前LST空间降尺度最优回归核的组合形式以及其泛化能力尚不明确，本文尝试选取波段反射率、遥感光谱指数、地形因子、地表覆盖类型、经纬度以及大气再分析数据6类回归核，构建基于极端梯度提升树(Extreme Gradrent Boost,XGBoost)算法的LST空间降尺度模型，同时开展对应的平行试验，通过对比分析6类回归核LST空间降尺度的精度，评估构建模型的泛化能力。

1 研究方法

1.1 总体技术路线

基于XGBoost算法的LST空间降尺度总体技术流程如图1所示。

图1 总体技术流程Fig.1 Flowchart of downscaling procedure

LST降尺度包括5个步骤：① 数据预处理。对图像进行质量检测、重投影、裁剪、重采样以及空间聚合，获取低分辨率(1 km)的回归核数据和LST数据以及高分辨率(100 m)的回归核数据和LST数据。② 温度校正。由于待降尺度的MODIS的1 km分辨率LST产品是通过劈窗算法获取的，而选作参考基准的Landsat-8卫星的100 m分辨率LST由单通道算法产生，搭载2个传感器的卫星过境时间存在差异，为了更好地验证LST空间降尺度的结果，本文将以Landsat-8卫星获取的LST作为基准，对MODIS获取的LST进行校正。③ XGBoost模型训练。利用XGBoost获取低分辨率回归核与LST之间的回归关系。④ 模型预测。用在低空间分辨率上训练好的模型来预测高空间分辨率LST，加上低分辨率的残差，得到温度降尺度的结果。⑤ 降尺度结果的验证。

1.2 XGBoost算法

XGBoost是由Chen[15]提出的一种提升树(Tree Boosting)算法，该算法能够应用于多种场景，在机器学习大赛和数据挖掘比赛中，多次被优胜队伍选用[15]。Chen[15]提到XGBoost在梯度提升树(Gradient Tree Boosting)的基础上，进一步对模型做了优化。例如,在稀疏数据的处理上速度更快，并行和分布式计算等。这些优化使XGBoost能够用最小的资源处理大量的数据。

在数据集D={(xi,yi)|i=1,2,…,n}(xi∈m,yi∈)中选取一个样本(xi,yi)，用K个可加性函数对样本进行预测，结果如下：

(1)

XGBoost的目标函数由损失函数和正则化项组成，这2项分别表示模型的拟合效果和模型的复杂度。目标函数如下：

(2)

对损失函数进行二阶泰勒展开，有：

1.3 基于回归的LST降尺度方法

基于回归的LST降尺度方法的基本假设是在不同尺度下，LST回归核的回归关系不变[8]，其关系式[16]为：

(3)

(4)

(5)

本文通过分析比较前人研究，共筛选出6类不同的回归核，具体情况如表1所示。地表反射率表征地表反射的太阳短波辐射能量大小，与地面状态、太阳高度角相关。光谱指数、地表类型数据提供地表状态信息，包括植被覆盖度、干旱状态及建筑信息等，提高模型在异质性区域的适用性。地形数据反映不同地形地貌对辐射的影响，提高模型在地形起伏区域的适用性。大气再分析数据反映大气生物物理量，浅层土壤含水量受浅层土壤温度的影响，且浅层土壤温度、近地表气温与地表温度有相关性。空气流动的速率与植物蒸腾作用相关，间接对地表温度产生影响。经纬度提供空间变化信息。

表1 回归核的信息

1.4 平行分组实验

为了比较不同回归核情况下LST降尺度精度的差异，对6类回归核进行不同组合，共得到12组，具体分组情况如表2所示。

表2 分组情况

2 研究区概况与数据源

2.1 研究区概况

依照地理位置分布范围广、地表类型多样以及Landsat8图像晴空无云且质量较高3个标准，以2019年为待选年份，在中国范围内选出了7景满足条件的Landsat8数据，并将其作为降尺度的研究区(83°45′E～116°17′E，33°N～44°30′N)，研究区的时空信息如表3所示。

表3 研究区时空信息

2.2 数据源及其预处理

本研究所用的Landsat8 C2L2数据来源于美国地质调查局网站(https:∥earthexplorer.usgs.gov/)。Landsat8数据集包含30 m分辨率的可见光波段、近红外波段以及100 m分辨率的热红外波段。分别对原始Landsat8可见光-近红外反射率产品以及热红外LST产品进行裁剪、重采样(双线性插值法)、空间平均聚合，得到100 m分辨率和1 km分辨率的LST以及相应的降尺度回归核。

低分辨率的LST选择分辨率为1 km的MOD11A1数据[17](https:∥lpdaac.usgs.gov/products/mod11a1v006/)。对MOD11A1数据进行重投影、裁剪，并利用质量文件删除质量较差以及有云的点。

数字高程数据来源于航天飞机雷达地形测绘使命(Shuttle Radar Topography Mission，SRTM)[18](https:∥srtm.csi.cgiar.org/srtmdata/)，分辨率为90 m，经过重投影、重采样(双线性插值法)、裁剪及空间平均聚合后，作为100 m分辨率和1 km分辨率的回归核。

地表覆盖类型数据来源于GlobeLand30数据集[19](https:∥www.webmap.cn/commres.do?method=globeIndex)，分辨率为30 m。同样对其进行重投影、重采样(双线性插值法)、裁剪、空间平均聚合的数据预处理工作，得到100 m分辨率和1 km分辨率的数据。

大气再分析数据取自欧洲中期天气预报中心(ECMWF)全球气候大气再分析资料的第五代产品——ERA5[20](https:∥cds.climate.copernicus.eu/)，其空间分辨率为0.25°。从中选取浅层土壤温度、浅层土壤体积含水量、离地2 m附近大气温度、离地10 m附近径向风速、离地10 m附近横向风速，5个数据作为代表。首先获取与卫星成像时间最近2个时刻的大气再分析数据，对大气再分析数据进行重投影、裁剪以及双线性插值处理，得到100 m分辨率和1 km分辨率的数据；再对相邻2个时刻的数据进行线性插值，得到与卫星成像时间相同的大气再分析数据。

此外，由于Landsat8获取的LST与MOD11A1获取的LST之间存在算法和观测时间之间的差异，对MOD11A1的温度产品进行了温度校正。主要操作是先对Landsat8获取的LST进行空间平均聚合操作，得到1 km分辨率的LST(Landsat8)，再以MOD11A1-LST为自变量，Landsat8-LST(1 km)为因变量做回归，对原始待降尺度的MOD11A1-LST进行了温度的校正。7个研究区温度校正过程的散点图如图2所示(图2(a)～图2(g)分别是7个研究区的校正结果)，LST校正的均方根误差(RMSE)在0.6～1 K，表明温度校正后的结果可用于LST降尺度结果的检验。

图2 MOD11A1温度校正的散点图

3 结果与讨论

3.1 研究区降尺度结果的分析

表4 12组平行试验LST降尺度算法性能统计表

为了更好地展示各研究区降尺度的性能，本文挑选了表征降尺度精度的RMSE以及表征VIF作为评价指标，统计结果如表5所示，其中研究区1～7用A～G标识，分组用G1～G12表示，RstdV和VstdV分别表示7个研究区RMSE和VIF统计指标的标准差。

表5 各研究区的均方根误差和视觉信息保真度

由表4和表5可以看出，group1(回归核为地表反射率)的RMSE结果较差，表现为降尺度的总体精度较差，且各个研究区的离散度较大，标准差为0.91；剩余group的降尺度总体精度差别不大，在2 K左右，但group3(回归核为地表反射率、光谱指数)在各个研究区的表现略微差，体现在精度存在部分差异，标准差在0.7左右。对于视觉信息保真度VIF，group1(回归核为地表反射率)、group2、group3、group4(回归核为光谱指数、地形)、group6(回归核为光谱指数、地表类型)、group7、group12(回归核为地表反射率、经纬度)能够更好地保持图像的纹理。

为了进一步目视展示降尺度的结果，本文选取了研究区1、研究区2和研究区7，绘制了不同分组回归核降尺度前后的对比图，如图3所示。

由图3可以看出，以降尺度的清晰度和LST的空间分布作为评价指标，结合表4和表5统计结果，group2、group4、group7相对更优，LST降尺度RMSE可控制在2 K，LST回归残差控制在0.5 K，而视觉信息保真度VIF可超过0.07。由于这3个分组中都包含了光谱指数回归核，因此光谱指数对于LST的降尺度更加重要。在构建回归模型时，引入光谱指数回归核，能够有效控制LST降尺度的精度。

3.2 降尺度的交叉检验与泛化能力分析

为探究利用单个研究区进行训练的模型推广至其他研究区后降尺度精度、保真度指标以及回归模型适用性的变化，即探寻降尺度模型的泛化能力，选择3组(group2，group4，group7)，针对研究区1，2和7，分别统计了选定研究区不同分组回归核条件下的降尺度的性能指标，如表6所示。

表6 不同组别的回归核在不同研究区交叉验证的统计结果

图3 选定研究区降尺度前后的对比图Fig.3 Results of comparison for studied areas before and after LST downscaling

由表6可以看出，对于group2、group4和group7分组回归核，某一研究区数据训练得到的XGBoost模型应用于本研究区或者其他研究区时，LST降尺度的RMSE变化不大，二者差值小于0.5 K。对于视觉信息保真度VIF，某一研究区数据训练得到的模型应用于其他研究区时，VIF将降低，会造成LST降尺度结果的失真，如模糊或者细节丢失。对于LST回归模型的残差，残差变化较大，本地训练模型应用于本地时，残差值较小，然而用于其他区域的时候，残差能增大到10 K以上。虽然LST降尺度的RMSE能保持在同一水平，但残差值仍然较大，说明通过局部数据训练的XGBoost模型，回归泛化能力不够，某个研究区训练得到的模型，应用于更大范围的研究区时可能会带来较大的误差。为了提高XGBoost模型的回归泛化能力，需要保证训练数据的代表性。

4 结束语

本文选取7个研究区，并对地表反射率、光谱指数、地形因子、大气再分析数据、经纬度、地表类型6种回归核进行组合实验。利用基于XGBoost算法的LST降尺度模型，在不同回归核组合的情况下，将1 km分辨率的MODISLST空间降尺度成100 m分辨率的LST，并将估计值与Landsat-8反演的100 m分辨率的LST进行了分析比较，讨论了不同回归核组合条件下LST降尺度的性能差异。本文得到的主要结论如下：

① 通过12组回归核的实验比对发现， group2，group4，group7作为回归核时，可以得到精度低且质量好的预测图像，3组的RMSE在2 K左右，LST降尺度残差控制在0.5 K左右，VIF可超过0.07。

② 由于选定的最优回归核组合(group2，group4，group7)中均包含了光谱指数，因此光谱指数是LST降尺度的关键因素。在实际应用中，可根据其他回归核(如地形、经纬度)的获取情况选择具体的回归核组合。从当前研究结果来看，引入更多的回归核会部分提高LST降尺度的精度，但在一定程度上增大了降尺度实施的复杂度。

③ 基于XGBoost算法的LST降尺度模型的泛化能力仍旧不够，虽然LST降尺度的RMSE可以得到一定的保证，但这是通过残差纠正来实现的，实际上局部数据训练的XGBoost回归模型自身无法有效从某一区域推广到另一区域。为了构建适用大范围的LST降尺度模型，确定XGBoost回归模型需要选择更具代表性的训练数据。

在本文研究的基础上，后续可尝试对降尺度模型进行改进，减少降尺度残差对结果的影响。同时将考虑组内回归核的相关性，对组内回归核数据进行主成分分析，并对回归结果进行比较分析。此外，也将重点考虑回归窗口的影响，采用全局和局部等不同策略来完善降尺度模型，以期得到一个泛化能力更强的LST降尺度模型。