CN109885959A - 一种地表温度鲁棒降尺度方法 - Google Patents

Abstract

本发明公开了一种地表温度鲁棒降尺度方法，包括以下八个步骤：1、地表参数的获取；2、大气参数的获取；3、配套数据集的制备；4、机器学习算法的训练；5、机器学习误差的估计；6、地表温度偏差的估算；7、地表温度降尺度的计算；8、降尺度结果的加权融合。本发明在前端的输入部分，通过额外引入大气降尺度因子来实现对地表温度空间分布和变化进行准确的表征；在后端的输出部分，通过小尺度机器学习算法训练以及算法加权融合，来解决算法的尺度问题和过拟合问题，最终实现不同区域不同下垫面都适用且精度有保障的地表温度降尺度。

Description

一种地表温度鲁棒降尺度方法

技术领域

本发明涉及一种方法，尤其涉及一种地表温度鲁棒降尺度方法，属于热红外遥感数据处理技术领域。

背景技术

地表温度作为描述区域或者全球范围内地表与大气的相互作用、能量交换以及水循环过程的最重要的环境参数之一，一直在陆地生态系统研究、陆表能量状态平衡以及陆表物理参数求解等多个学科领域发挥着不可替代的作用。精确的地表温度不仅能够辅助评估陆表水文平衡条件，测算土壤湿度，还能作为输入参数，辅助地表土壤蒸散发参数的计算。虽然目前热红外遥感数据已经成为全球和区域地表温度产品获取的最主要的来源之一，但热红外传感器的空间分辨率普遍较低，严重制约了地表温度的广泛应用。因此，发展对大尺度低空间分辨率的热红外传感器的地表温度产品的降尺度算法，以获取兼有高时间分辨率和高空间分辨率的小尺度地表温度产品是目前热红外遥感领域的研究热点之一。

现有的地表温度降尺度方法主要是基于地表温度与地表参数之间关系的尺度不变性假设，通过在低分辨率尺度上建立地表温度与地表降尺度因子之间的函数关系，进而将其应用至高分辨率尺度上，最终得到地表温度降尺度结果。然而受不同区域地表下垫面条件、不同时相以及不同空间分辨率环境的影响，现有的大多数地表温度降尺度模型主要存在以下三方面的技术缺陷：第一，地表温度的空间分布自身具有较大的复杂性，是多种地表和大气环境因子的相互作用和影响的结果，现有技术仅仅依靠地表降尺度因子无法精确的描述其分布特征和变化，导致了地表温度降尺度的精度受限；第二，地表温度与地表参数之间关系的尺度不变性假设不合理，导致了降尺度的精度无法有效提高；第三，现有技术仅依赖于单种算法，没有融合考虑多种算法的优劣，容易产生训练数据效果佳但验证数据牵强人意的过拟合现象，适用区域和范围有限。

发明内容

为了解决上述技术所存在的不足之处，本发明提供了一种地表温度鲁棒降尺度方法。

为了解决以上技术问题，本发明采用的技术方案是：一种地表温度鲁棒降尺度方法，包括以下步骤：

步骤1、地表参数的获取；

将降尺度前的原始地表温度产品表示为T_O，其空间分辨率表示为R_O；将预期得到的降尺度后的地表温度产品表示为T_J，其空间分辨率表示为R_J；

根据原始地表温度产品T_O的获取时间，下载空间分辨率为R_J且与T_O的获取时间接近且与地表温度空间分布和变化相关的多个地表参数集S_O；所述地表参数集S_O包括地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度；

根据R_J的数值大小，下载空间分辨率近似相等，且与T_O的获取时间接近的晴空无云热红外地表温度产品T_R，其空间分辨率为R_R；同时，根据地表温度产品T_R的获取时间，下载相同时间且空间分辨率相同的与地表温度空间分布和变化相关的多个地表参数集S_R，其空间分辨率为R_R；所述地表参数集S_R包括地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度；T_R后续将作为机器学习的因变量，而S_R后续将作为机器学习的自变量之一；

步骤2、大气参数的获取；

根据原始地表温度产品T_O获取时间及其空间分辨率数值大小，下载空间分辨率接近且同一天的多个时刻的大气再分析数据，并从中挑选出与地表温度空间分布和变化相关的多个大气参数集A_O，其空间分辨率为R_A；同时，根据热红外地表温度产品T_R获取时间及其空间分辨率数值大小，下载空间分辨率接近且同一天的多个时刻的大气再分析数据，并从中挑选出与地表温度空间分布和变化相关的多个大气参数集A_R，其空间分辨率也为R_A；大气参数集A_R后续将和地表参数集S_R一块作为机器学习的自变量；大气参数集A_O以及A_R都包括近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量；

步骤3、配套数据集的制备；

制备空间分辨率为R_J的地表参数集S_O和大气参数集A_O，作为降尺度算法的输入；制备空间分辨率为R_R的地表参数集S_R和大气参数集A_R，用于机器学习算法的训练；

步骤4、机器学习算法的训练；

利用制备好的地表温度产品T_R和地表参数集S_R和大气参数集A_R完成机器学习算法的训练；为了便于机器算法的训练和后续机器学习算法的误差估计，将对应好的地表温度产品和地表和大气参数集按60％和40％进行划分，其中60％的数据用于机器学习地表温度降尺度模型的训练，另外的40％用于机器学习算法误差的估计；

地表温度产品T_R将作为机器学习的因变量，制备好的配套的地表参数集S_R和大气参数集A_R将作为机器学习的自变量；配套的因变量和自变量都作为机器学习算法的输出和输入，机器学习将自动模拟人脑的学习行为以自动获取自变量与因变量之间存在的非线性复杂关系f；

机器学习地表温度降尺度模型的映射训练过程由下式表示：

LST_R,k＝f_n(S_R,k ¹,…,S_R,k ^p,…,S_R,k ¹⁹,A_R,k ¹,…,A_R,k ^q,…,A_R,k ⁸) 公式12

其中，LST_R,k为地表温度产品T_R中第k个地表温度；S_R,k ^p为地表参数集S_R中第k个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_R,k ^q为大气参数集A_R中第k个大气参数所对应的第q个大气降尺度因子，q取值1～8；f_n是第n种机器学习算法训练好的非线性复杂关系；

步骤5、机器学习误差的估计；

由于实际机器学习过程中会不可避免的存在偏差，因此需要进行每个机器学习算法的误差估计，并以此计算融合权重系数；除掉用于机器学习算法训练的60％的数据以外，剩下的40％将用于算法误差的估计；算法误差由下式计算得：

其中，σ_n是第n种机器学习算法的算法误差；k是序列数；LST_R,k是用于算法误差估计的地表温度产品T_R中第k个地表温度产品；是用第n种机器学习算法估算的第k个地表温度产品；∑是数学上累加求和符号；N是全部用于算法误差估计的地表温度产品的个数；f_n是第n种机器学习算法训练好的非线性复杂关系；S_R,k ^p为地表参数集S_R中第k个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_R,k ^q为大气参数集A_R中第k个大气参数所对应的第q个大气降尺度因子，q取值1～8；

步骤6、地表温度偏差的估算；

由于机器学习算法所采用的地表和大气参数个数仍旧有限，并不能够解释地表温度所有的变异情况，因此需要同一空间分辨率R_O条件下，比较降尺度处理前后的地表温度产品间的差异，并以此作为未解释的地表温度的偏差；地表温度偏差的计算公式为：

其中，ΔLST_O,l ⁿ是地表温度产品T_O中第l个地表温度偏差；LST_O,l是地表温度产品T_O中第l个地表温度实际值；是用第n种机器学习算法估算的第l个地表温度产品值；∑是数学上累加求和符号；U是第l个地表温度粗像元嵌套包含细像元的总个数；u是序列数；W_O,u是细像元在粗像元内的面积百分比；f_n是第n种机器学习算法训练好的非线性复杂关系；S_O,u ^p为粗像元所包含的第u个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_O,u ^q为粗像元所包含的第u个大气参数所对应的第q个大气降尺度因子，q取值1～8；

步骤7、地表温度降尺度的计算；

利用制备好的空间分辨率为R_J的地表参数集S_O和大气参数集A_O、训练好的机器学习算法以及估算的地表温度偏差，进行地表温度的降尺度，计算公式为：

其中，LST_O,l ⁿ是采用第n种机器学习算法降尺度后的第l个地表温度产品的结果；是用第n种机器学习算法估算的第l个地表温度产品初始值；ΔLST_O,l ⁿ是对应的第l个地表温度偏差；f_n是第n种机器学习算法训练好的非线性复杂关系；S_O,l ^p为地表参数集S_O中第l个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_O,l ^q为大气参数集A_O中第l个大气参数所对应的第q个大气降尺度因子，q取值1～8；

步骤8、降尺度结果的加权融合；

由于降尺度过程可采用多种机器学习算法，因此为了达到算法间的取长补短优势互补，采用下式对各种机器学习算法降尺度的结果进行了加权融合，具体公式如下：

其中，T_J是降尺度后的第l个地表温度最终结果；LST_O,l ⁿ是采用第n种机器学习算法降尺度后的第l个地表温度产品的结果；∑是数学上累加求和符号；Z是所选择的机器学习算法的总个数；z是序列数；是第n种机器学习算法的权重系数，该权重系数可以表示为：

其中，σ_n是第n种机器学习算法的算法误差；∑是数学上累加求和符号；N是全部用于地表温度降尺度的机器学习算法的个数；n是序列数。

进一步地，步骤1中地表参数集S_O及S_R中的地表反射率产品，均包括红光反射率产品、绿光反射率产品、蓝光反射率产品、近红外反射率产品、短波红外反射率产品和短波红外反射率产品。

进一步地，步骤3的具体过程为：

第一步：降尺度输入数据集的制备，数据集包含地表参数集S_O和大气参数集A_O，空间分辨率为R_J；

首先，根据降尺度后的地表温度产品T_J的空间分辨率和投影坐标，将下载的地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度分别做几何裁剪、几何纠正和重投影，以使这些地表参数数据都能够完全配准到一块；

接着，利用重投影后的地表反射率产品，根据下式分别计算降尺度所需的7类光谱指数；这7类光谱指数将添补至地表参数集；具体的计算公式如下：

其中，ρ_R为红光反射率产品，ρ_NIR为近红外反射率产品，NDVI为归一化植被；Fv为植被覆盖度，SAVI为土壤调节植被指数；NDVI_max和NDVI_min分别是NDVI值的最大值和最小值；NMDI为归一化多波段干旱指数，ρ_SWIR1和ρ_SWIR2分别为两个短波红外反射率产品；NDDI为归一化干旱指数，ρ_B为蓝光反射率产品；MNDWI为改进的归一化水指数，ρ_G为绿光反射率产品；NDBI为归一化建筑指数；

然后，利用地表数字高程模型，分别在3*3的窗口内，依次计算坡度SLOPE和坡向ASPECT参数，具体公式如下：

其中，f_x和f_y由下式计算得到：

f_x和f_y是中间参数，表示3*3窗口中心像素点(i，j)位置处的水平和垂直坡度；H_i,j代表像素点(i，j)位置处的高程，可直接从地表数字高程模型中提取；R_J是地表数字高程模型的空间分辨率；同样，坡度SLOPE和坡向ASPECT参数也将添补至地表参数集；

根据降尺度后的地表温度产品T_J的获取时间，将下载的同一天多个时刻的近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量先进行时间插值，线性内插至地表温度获取时间；同时，依旧参考地表温度产品T_J的空间分辨率和投影坐标，将内插好的近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量分别做几何裁剪、几何纠正和重投影，以使这些大气参数数据都能够和地表温度产品完全配准到一块；

至此，所有的地表参数集S_O和大气参数集A_O都具有相同的空间分辨率和投影坐标系，能够完全配准在一起，完成了降尺度算法的输入配套数据集的制备；

第二步：机器学习训练数据集制备，数据集包含地表参数集S_R和大气参数集A_R，空间分辨率为R_R；

使用跟第一步相同的方法，但选择的空间分辨率为R_R，制备的对象为地表参数集S_R和大气参数集A_R，利用公式1～11添补地表参数集S_R，并将地表温度产品T_R和地表参数集S_R和大气参数集A_R配准在一起，完成机器学习算法训练的配套数据集的制备。

进一步地，步骤4中，所述地表参数共有19个降尺度因子，包括地表红光反射率产品、地表绿光反射率产品、地表蓝光反射率产品、地表近红外反射率产品、地表短波红外反射率产品—两个波段、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度、归一化植被、土壤调节植被指数、植被覆盖度、归一化多波段干旱指数、归一化干旱指数、改进的归一化水指数、归一化建筑指数、坡度和坡向；

所述大气参数共有8个降尺度因子，包括近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量。

进一步地，步骤4中，由于机器学习算法之间精度存在优劣，使用多种不同的机器学习来模拟非线性复杂关系f，以提高机器学习的泛化能力；多种不同的机器学习包括基于支持向量机法、神经网络法和随机森林法。

本发明提出了一种基于机器学习的地表温度鲁棒降尺度方法，除了将地表参数作为降尺度因子以外，还额外将大气参数同时也作为降尺度因子，克服了仅用地表降尺度因子无法有效表征地表温度空间分布和变化的弊端。通过直接在高分辨率数据上进行降尺度机器学习，本发明还解决了低分辨率训练的降尺度算法直接用于高分辨率数据的算法尺度不变性不合理假设导致的适用性问题。此外，考虑到单一降尺度算法各有利弊，通过对不同机器学习降尺度算法进行了加权融合，弥补了单一算法精度的不足，解决了过拟合问题，扩大了适用范围的同时还提高了降尺度的精度。

附图说明

图1为本发明的总体步骤流程图。

具体实施方式

下面结合附图和具体实施方式对本发明作进一步详细的说明。

图1所示的一种地表温度鲁棒降尺度方法，包括以下步骤：

步骤1、地表参数的获取

所谓降尺度，就是把大尺度、低分辨率的全球气候模式输出的信息转化为小尺度、高分辨率的区域地面气候变化信息的一种方法。将降尺度前的原始地表温度产品表示为T_O，其空间分辨率表示为R_O；将预期得到的降尺度后的地表温度产品表示为T_J，其空间分辨率表示为R_J(该空间分辨率是根据需要可直接设定的，为已知参数)。

根据原始地表温度产品T_O的获取时间，下载空间分辨率为R_J且与T_O的获取时间最为接近且与地表温度空间分布和变化密切相关的多个地表参数集S_O，该地表参数集S_O包括：地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度。

根据R_J的数值大小，下载空间分辨率近似相等，且与T_O的获取时间比较接近的晴空无云热红外地表温度产品T_R，其空间分辨率为R_R(该操作是因为：不同的卫星获取的地表温度产品的空间分辨率不同，因此在实际操作时很难找到一模一样的地表温度产品，只能先找近似相等的，然后后期进行配准修正)。同时，根据地表温度产品T_R的获取时间，下载相同时间且空间分辨率相同的与地表温度空间分布和变化密切相关的多个地表参数集S_R，其空间分辨率为R_R，包括：地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度。T_R后续将作为机器学习的因变量，而S_R后续将作为机器学习的自变量之一。

地表参数集S_O及S_R中的地表反射率产品，均包括红光反射率产品(中心波长介于0.63um-0.69um)、绿光反射率产品(中心波长介于0.52um-0.60um)、蓝光反射率产品(中心波长介于0.45um-0.52um)、近红外反射率产品(中心波长介于0.77um-0.90um)、短波红外反射率产品(中心波长介于1.55um-1.75um)和短波红外反射率产品(中心波长介于2.09um-2.35um)。

步骤2、大气参数的获取

根据原始地表温度产品T_O获取时间及其空间分辨率数值大小，下载空间分辨率最接近且同一天的多个时刻的大气再分析数据，并从中挑选出与地表温度空间分布和变化密切相关的多个大气参数集A_O，其空间分辨率为R_A。同时，根据热红外地表温度产品T_R获取时间及其空间分辨率数值大小，下载空间分辨率最接近且同一天的多个时刻的大气再分析数据，并从中挑选出与地表温度空间分布和变化密切相关的多个大气参数集A_R，其空间分辨率也为R_A。大气参数集A_R后续将和地表参数集S_R一块作为机器学习的自变量。无论是大气参数集A_O还会A_R都包括近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量。

步骤3、配套数据集的制备

分为两步：制备空间分辨率为R_J的地表参数集S_O和大气参数集A_O，作为降尺度算法的输入；制备空间分辨率为R_R的地表参数集S_R和大气参数集A_R，用于机器学习算法的训练。

第一步：降尺度输入数据集的制备，该数据集包含地表参数集S_O和大气参数集A_O，空间分辨率为R_J。

首先，根据降尺度后的地表温度产品T_J的空间分辨率和投影坐标，将下载的地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度分别做几何裁剪、几何纠正和重投影，以使得这些地表参数数据都能够完全配准到一块。

接着，利用重投影后的地表反射率产品，根据下式分别计算降尺度所需的7类光谱指数。这7类光谱指数将添补至地表参数集，这些光谱指数包括：归一化植被NDVI、土壤调节植被指数SAVI、植被覆盖度Fv、归一化多波段干旱指数NMDI、归一化干旱指数NDDI、改进的归一化水指数MNDWI和归一化建筑指数NDBI，其中NDVI、SAVI和Fv反映了植被变化对地表温度空间分布和变化的影响；NMDI和NDDI反映了土壤湿度变化对地表温度空间分布和变化的影响；MNDWI和NDBI分别反映了水体和人造建筑物的比例变化对地表温度空间分布和变化的影响。这些光谱指数的使用有利于进一步表示不同下垫面情况局部环境差异引起的地表温度变化，其具体的计算公式如下：

其中，ρ_R为红光反射率产品，ρ_NIR为近红外反射率产品，NDVI为归一化植被；Fv为植被覆盖度，SAVI为土壤调节植被指数；NDVI_max和NDVI_min分别是NDVI值的最大值和最小值；NMDI为归一化多波段干旱指数，ρ_SWIR1和ρ_SWIR2分别为两个短波红外反射率产品；NDDI为归一化干旱指数，ρ_B为蓝光反射率产品；MNDWI为改进的归一化水指数，ρ_G为绿光反射率产品；NDBI为归一化建筑指数。

然后，利用地表数字高程模型(Digital Elevation Model，简称DEM)，分别在3*3的窗口内，依次计算坡度SLOPE和坡向ASPECT参数，具体公式如下：

其中，SLOPE表示坡度，ASPECT表示坡向；f_x和f_y是中间参数，可以由下式计算得到：

其中，f_x和f_y是中间参数，表示3*3窗口中心像素点(i，j)位置处的水平和垂直坡度；H_i,j代表像素点(i，j)位置处的高程，可直接从地表数字高程模型中提取；R_J是地表数字高程模型的空间分辨率。同样，坡度SLOPE和坡向ASPECT参数也将添补至地表参数集。

最后，考虑到大气参数随时空变化相对较为平稳，因此根据降尺度后的地表温度产品T_J的获取时间，将下载的同一天多个时刻的近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量先进行时间插值，线性内插至地表温度获取时间。同时，依旧参考地表温度产品T_J的空间分辨率和投影坐标，将内插好的近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量分别做几何裁剪、几何纠正和重投影，以使得这些大气参数数据都能够和地表温度产品完全配准到一块。

至此，所有的地表参数集S_O和大气参数集A_O都具有相同的空间分辨率和投影坐标系，能够完全配准在一起，完成了降尺度算法的输入配套数据集的制备。

第二步：机器学习训练数据集制备，该数据集包含地表参数集S_R和大气参数集A_R，空间分辨率为R_R。

使用跟第一步相同的方法，但选择的空间分辨率为R_R，制备的对象为地表参数集S_R和大气参数集A_R，利用公式1～11添补地表参数集S_R，并将地表温度产品T_R和地表参数集S_R和大气参数集A_R配准在一起，完成机器学习算法训练的配套数据集的制备。

步骤4、机器学习算法的训练

利用制备好的地表温度产品T_R和地表参数集S_R和大气参数集A_R完成机器学习算法的训练。为了便于机器算法的训练和后续机器学习算法的误差估计，这里将对应好的地表温度产品和地表和大气参数集按60％和40％进行划分，其中60％的数据用于机器学习地表温度降尺度模型的训练，另外的40％用于机器学习算法误差的估计。

地表温度产品T_R将作为机器学习的因变量(输出变量)，制备好的配套的地表参数集S_R和大气参数集A_R将作为机器学习的自变量(输入变量)。其中，地表参数共有19个降尺度因子，包括地表红光反射率产品、地表绿光反射率产品、地表蓝光反射率产品、地表近红外反射率产品、地表短波红外反射率产品(两个波段)、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度、归一化植被、土壤调节植被指数、植被覆盖度、归一化多波段干旱指数、归一化干旱指数、改进的归一化水指数、归一化建筑指数、坡度和坡向；大气参数共有8个降尺度因子，包括近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量。配套的因变量和自变量都作为机器学习算法的输出和输入，机器学习将自动模拟人脑的学习行为以自动获取自变量与因变量之间存在的非线性复杂关系f。这个过程中，由于机器学习算法之间精度存在优劣，因此可使用多种不同的机器学习来模拟非线性复杂关系f，例如包括基于支持向量机法、神经网络法和随机森林法等，但不限于这3种方法，以提高机器学习的泛化能力。

机器学习地表温度降尺度模型的映射训练过程可以由下式表示，即：

LST_R,k＝f_n(S_R,k ¹,…,S_R,k ^p,…,S_R,k ¹⁹,A_R,k ¹,…,A_R,k ^q,…,A_R,k ⁸) 公式12

其中，LST_R,k为地表温度产品T_R中第k个地表温度；S_R,k ^p为地表参数集S_R中第k个地表参数所对应的第p个地表降尺度因子(p在1到19中取值，代表制备好的19个地表降尺度因子其中之一)；A_R,k ^q为大气参数集A_R中第k个大气参数所对应的第q个大气降尺度因子(q在1到8中取值，代表制备好的8个大气降尺度因子其中之一)；f_n是第n种机器学习算法训练好的非线性复杂关系。

步骤5、机器学习误差的估计

由于实际机器学习过程中会不可避免的存在偏差，因此需要进行每个机器学习算法的误差估计，并以此计算融合权重系数。除掉用于机器学习算法训练的60％的数据以外，剩下的40％将用于算法误差的估计。算法误差可以由下式计算得：

其中，σ_n是第n种机器学习算法的算法误差；k是序列数；LST_R,k是用于算法误差估计的地表温度产品T_R中第k个地表温度产品；是用第n种机器学习算法估算的第k个地表温度产品；∑是数学上累加求和符号；N是全部用于算法误差估计的地表温度产品的个数；f_n是第n种机器学习算法训练好的非线性复杂关系；S_R,k ^p为地表参数集S_R中第k个地表参数所对应的第p个地表降尺度因子(p在1到19中取值，代表制备好的19个地表降尺度因子其中之一)；A_R,k ^q为大气参数集A_R中第k个大气参数所对应的第q个大气降尺度因子(q在1到8中取值，代表制备好的8个大气降尺度因子其中之一)。

步骤6、地表温度偏差的估算

由于机器学习算法所采用的地表和大气参数个数仍旧有限，并不能够解释地表温度所有的变异情况，因此需要同一空间分辨率R_O条件下，比较降尺度处理前后的地表温度产品间的差异，并以此作为未解释的地表温度的偏差。地表温度偏差的计算公式为：

其中，ΔLST_O,l ⁿ是地表温度产品T_O中第l个地表温度偏差；LST_O,l是地表温度产品T_O中第l个地表温度实际值；是用第n种机器学习算法估算的第l个地表温度产品值；∑是数学上累加求和符号；U是第l个地表温度粗像元(空间分辨率为R_O)嵌套包含细像元(空间分辨率为R_J)的总个数；u是序列数；W_O,u是细像元在粗像元内的面积百分比；f_n是第n种机器学习算法训练好的非线性复杂关系；S_O,u ^p为粗像元所包含的第u个地表参数所对应的第p个地表降尺度因子(p在1到19中取值，代表制备好的19个地表降尺度因子其中之一)；A_O,u ^q为粗像元所包含的第u个大气参数所对应的第q个大气降尺度因子(q在1到8中取值，代表制备好的8个大气降尺度因子其中之一)。

步骤7、地表温度降尺度的计算

利用制备好的空间分辨率为R_J的地表参数集S_O和大气参数集A_O、训练好的机器学习算法以及估算的地表温度偏差，进行地表温度的降尺度，计算公式为：

其中，LST_O,l ⁿ是采用第n种机器学习算法降尺度后的第l个地表温度产品的结果；是用第n种机器学习算法估算的第l个地表温度产品初始值；ΔLST_O,l ⁿ是对应的地表温度产品T_O中第l个地表温度偏差；f_n是第n种机器学习算法训练好的非线性复杂关系；S_O,l ^p为地表参数集S_O中第l个地表参数所对应的第p个地表降尺度因子(p在1到19中取值，代表制备好的19个地表降尺度因子其中之一)；A_O,l ^q为大气参数集A_O中第l个大气参数所对应的第q个大气降尺度因子(q在1到8中取值，代表制备好的8个大气降尺度因子其中之一)。

步骤8、降尺度结果的加权融合

由于降尺度过程采用了多种机器学习算法，因此为了达到算法间的取长补短优势互补，这里采用下式对各种机器学习算法降尺度的结果进行了加权融合，具体公式如下：

其中，T_J是降尺度后的第l个地表温度最终结果；LST_O,l ⁿ是采用第n种机器学习算法降尺度后的第l个地表温度产品的结果；∑是数学上累加求和符号；Z是所选择的机器学习算法的总个数；z是序列数；是第n种机器学习算法的权重系数，该权重系数可以表示为：

其中，σ_n是第n种机器学习算法的算法误差；∑是数学上累加求和符号；N是全部用于地表温度降尺度的机器学习算法的个数；n是序列数。

本发明与现有技术相比，具有以下优势：

1、本发明利用相关性度量准则，除了地表参数以外，还额外选择与地表温度最为密切相关的大气参数作为降尺度因子，克服了仅用地表参数作为降尺度因子所导致的无法有效表征地表温度空间变化的弊端；既本发明同时考虑到了地表和大气参数对地表温度空间分布和变化的影响，将与地表温度最为密切相关的地表和大气参数同时作为降尺度因子，实现地表温度空间分布和变化的精确描述和表征；

2、本发明通过直接在高分辨率上训练机器学习算法，减少了对配套数据集进行升尺度的处理环节，摆脱了降尺度算法尺度不变性的假设束缚，解决了低分辨率训练的降尺度算法直接用于高分辨率数据的算法尺度不变性不合理假设导致的适用性问题；

3、本发明通过将算法误差作为权重系数，将多种机器学习算法进行了有效融合，扬长避短，实现了多种机器学习算法的优势互补，在合理描述了地表温度与降尺度因子之间的非线性复杂关系的同时，还解决了单一算法过拟合问题，提高了地表温度降尺度的鲁棒性，确保了地表温度降尺度的精度。

上述实施方式并非是对本发明的限制，本发明也并不仅限于上述举例，本技术领域的技术人员在本发明的技术方案范围内所做出的变化、改型、添加或替换，也均属于本发明的保护范围。

Claims (5)

1.一种地表温度鲁棒降尺度方法，其特征在于：所述方法包括以下步骤：

步骤1、地表参数的获取；

将降尺度前的原始地表温度产品表示为T_O，其空间分辨率表示为R_O；将预期得到的降尺度后的地表温度产品表示为T_J，其空间分辨率表示为R_J；

根据原始地表温度产品T_O的获取时间，下载空间分辨率为R_J且与T_O的获取时间接近且与地表温度空间分布和变化相关的多个地表参数集S_O；所述地表参数集S_O包括地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度；

根据R_J的数值大小，下载空间分辨率近似相等，且与T_O的获取时间接近的晴空无云热红外地表温度产品T_R，其空间分辨率为R_R；同时，根据地表温度产品T_R的获取时间，下载相同时间且空间分辨率相同的与地表温度空间分布和变化相关的多个地表参数集S_R，其空间分辨率为R_R；所述地表参数集S_R包括地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度；T_R后续将作为机器学习的因变量，而S_R后续将作为机器学习的自变量之一；

步骤2、大气参数的获取；

根据原始地表温度产品T_O获取时间及其空间分辨率数值大小，下载空间分辨率接近且同一天的多个时刻的大气再分析数据，并从中挑选出与地表温度空间分布和变化相关的多个大气参数集A_O，其空间分辨率为R_A；同时，根据热红外地表温度产品T_R获取时间及其空间分辨率数值大小，下载空间分辨率接近且同一天的多个时刻的大气再分析数据，并从中挑选出与地表温度空间分布和变化相关的多个大气参数集A_R，其空间分辨率也为R_A；大气参数集A_R后续将和地表参数集S_R一块作为机器学习的自变量；大气参数集A_O以及A_R都包括近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量；

步骤3、配套数据集的制备；

制备空间分辨率为R_J的地表参数集S_O和大气参数集A_O，作为降尺度算法的输入；制备空间分辨率为R_R的地表参数集S_R和大气参数集A_R，用于机器学习算法的训练；

步骤4、机器学习算法的训练；

利用制备好的地表温度产品T_R和地表参数集S_R和大气参数集A_R完成机器学习算法的训练；为了便于机器算法的训练和后续机器学习算法的误差估计，将对应好的地表温度产品和地表和大气参数集按60％和40％进行划分，其中60％的数据用于机器学习地表温度降尺度模型的训练，另外的40％用于机器学习算法误差的估计；

地表温度产品T_R将作为机器学习的因变量，制备好的配套的地表参数集S_R和大气参数集A_R将作为机器学习的自变量；配套的因变量和自变量都作为机器学习算法的输出和输入，机器学习将自动模拟人脑的学习行为以自动获取自变量与因变量之间存在的非线性复杂关系f；

机器学习地表温度降尺度模型的映射训练过程由下式表示：

LST_R,k＝f_n(S_R,k ¹,…,S_R,k ^p,…,S_R,k ¹⁹,A_R,k ¹,…,A_R,k ^q,…,A_R,k ⁸) 公式12

其中，LST_R,k为地表温度产品T_R中第k个地表温度；S_R,k ^p为地表参数集S_R中第k个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_R,k ^q为大气参数集A_R中第k个大气参数所对应的第q个大气降尺度因子，q取值1～8；f_n是第n种机器学习算法训练好的非线性复杂关系；

步骤5、机器学习误差的估计；

由于实际机器学习过程中会不可避免的存在偏差，因此需要进行每个机器学习算法的误差估计，并以此计算融合权重系数；除掉用于机器学习算法训练的60％的数据以外，剩下的40％将用于算法误差的估计；算法误差由下式计算得：

其中，σ_n是第n种机器学习算法的算法误差；k是序列数；LST_R,k是用于算法误差估计的地表温度产品T_R中第k个地表温度产品；是用第n种机器学习算法估算的第k个地表温度产品；∑是数学上累加求和符号；N是全部用于算法误差估计的地表温度产品的个数；f_n是第n种机器学习算法训练好的非线性复杂关系；S_R,k ^p为地表参数集S_R中第k个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_R,k ^q为大气参数集A_R中第k个大气参数所对应的第q个大气降尺度因子，q取值1～8；

步骤6、地表温度偏差的估算；

由于机器学习算法所采用的地表和大气参数个数仍旧有限，并不能够解释地表温度所有的变异情况，因此需要同一空间分辨率R_O条件下，比较降尺度处理前后的地表温度产品间的差异，并以此作为未解释的地表温度的偏差；地表温度偏差的计算公式为：

其中，ΔLST_O,l ⁿ是地表温度产品T_O中第l个地表温度偏差；LST_O,l是地表温度产品T_O中第l个地表温度实际值；是用第n种机器学习算法估算的第l个地表温度产品值；∑是数学上累加求和符号；U是第l个地表温度粗像元嵌套包含细像元的总个数；u是序列数；W_O,u是细像元在粗像元内的面积百分比；f_n是第n种机器学习算法训练好的非线性复杂关系；S_O,u ^p为粗像元所包含的第u个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_O,u ^q为粗像元所包含的第u个大气参数所对应的第q个大气降尺度因子，q取值1～8；

步骤7、地表温度降尺度的计算；

利用制备好的空间分辨率为R_J的地表参数集S_O和大气参数集A_O、训练好的机器学习算法以及估算的地表温度偏差，进行地表温度的降尺度，计算公式为：

其中，LST_O,l ⁿ是采用第n种机器学习算法降尺度后的第l个地表温度产品的结果；是用第n种机器学习算法估算的第l个地表温度产品初始值；ΔLST_O,l ⁿ是对应的第l个地表温度偏差；f_n是第n种机器学习算法训练好的非线性复杂关系；S_O,l ^p为地表参数集S_O中第l个地表参数所对应的第p个地表降尺度因子，p取值1～19；A_O,l ^q为大气参数集A_O中第l个大气参数所对应的第q个大气降尺度因子，q取值1～8；

步骤8、降尺度结果的加权融合；

由于降尺度过程可采用多种机器学习算法，因此为了达到算法间的取长补短优势互补，采用下式对各种机器学习算法降尺度的结果进行了加权融合，具体公式如下：

其中，T_J是降尺度后的第l个地表温度最终结果；LST_O,l ⁿ是采用第n种机器学习算法降尺度后的第l个地表温度产品的结果；∑是数学上累加求和符号；Z是所选择的机器学习算法的总个数；z是序列数；是第n种机器学习算法的权重系数，该权重系数可以表示为：

其中，σ_n是第n种机器学习算法的算法误差；∑是数学上累加求和符号；N是全部用于地表温度降尺度的机器学习算法的个数；n是序列数。

2.根据权利要求1所述的地表温度鲁棒降尺度方法，其特征在于：所述步骤1中地表参数集S_O及S_R中的地表反射率产品，均包括红光反射率产品、绿光反射率产品、蓝光反射率产品、近红外反射率产品、短波红外反射率产品和短波红外反射率产品。

3.根据权利要求2所述的地表温度鲁棒降尺度方法，其特征在于：所述步骤3的具体过程为：

第一步：降尺度输入数据集的制备，数据集包含地表参数集S_O和大气参数集A_O，空间分辨率为R_J；

首先，根据降尺度后的地表温度产品T_J的空间分辨率和投影坐标，将下载的地表反射率产品、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度分别做几何裁剪、几何纠正和重投影，以使这些地表参数数据都能够完全配准到一块；

接着，利用重投影后的地表反射率产品，根据下式分别计算降尺度所需的7类光谱指数；这7类光谱指数将添补至地表参数集；具体的计算公式如下：

其中，ρ_R为红光反射率产品，ρ_NIR为近红外反射率产品，NDVI为归一化植被；Fv为植被覆盖度，SAVI为土壤调节植被指数；NDVI_max和NDVI_min分别是NDVI值的最大值和最小值；NMDI为归一化多波段干旱指数，ρ_SWIR1和ρ_SWIR2分别为两个短波红外反射率产品；NDDI为归一化干旱指数，ρ_B为蓝光反射率产品；MNDWI为改进的归一化水指数，ρ_G为绿光反射率产品；NDBI为归一化建筑指数；

然后，利用地表数字高程模型，分别在3*3的窗口内，依次计算坡度SLOPE和坡向ASPECT参数，具体公式如下：

其中，f_x和f_y由下式计算得到：

f_x和f_y是中间参数，表示3*3窗口中心像素点(i，j)位置处的水平和垂直坡度；H_i,j代表像素点(i，j)位置处的高程，可直接从地表数字高程模型中提取；R_J是地表数字高程模型的空间分辨率；同样，坡度SLOPE和坡向ASPECT参数也将添补至地表参数集；

根据降尺度后的地表温度产品T_J的获取时间，将下载的同一天多个时刻的近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量先进行时间插值，线性内插至地表温度获取时间；同时，依旧参考地表温度产品T_J的空间分辨率和投影坐标，将内插好的近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量分别做几何裁剪、几何纠正和重投影，以使这些大气参数数据都能够和地表温度产品完全配准到一块；

至此，所有的地表参数集S_O和大气参数集A_O都具有相同的空间分辨率和投影坐标系，能够完全配准在一起，完成了降尺度算法的输入配套数据集的制备；

第二步：机器学习训练数据集制备，数据集包含地表参数集S_R和大气参数集A_R，空间分辨率为R_R；

使用跟第一步相同的方法，但选择的空间分辨率为R_R，制备的对象为地表参数集S_R和大气参数集A_R，利用公式1～11添补地表参数集S_R，并将地表温度产品T_R和地表参数集S_R和大气参数集A_R配准在一起，完成机器学习算法训练的配套数据集的制备。

4.根据权利要求3所述的地表温度鲁棒降尺度方法，其特征在于：所述步骤4中，所述地表参数共有19个降尺度因子，包括地表红光反射率产品、地表绿光反射率产品、地表蓝光反射率产品、地表近红外反射率产品、地表短波红外反射率产品—两个波段、地表土地覆盖和土地利用分类图、地表数字高程模型、地表经度、地表纬度、归一化植被、土壤调节植被指数、植被覆盖度、归一化多波段干旱指数、归一化干旱指数、改进的归一化水指数、归一化建筑指数、坡度和坡向；

所述大气参数共有8个降尺度因子，包括近地表气温、大气可降水量、水平风速、垂直风速、短波下行通量、短波上行通量、长波下行通量和长波上行通量。

5.根据权利要求4所述的地表温度鲁棒降尺度方法，其特征在于：所述步骤4中，由于机器学习算法之间精度存在优劣，使用多种不同的机器学习来模拟非线性复杂关系f，以提高机器学习的泛化能力；多种不同的机器学习包括基于支持向量机法、神经网络法和随机森林法。