CN111652404A - 一种全天气地表温度反演方法和系统 - Google Patents

Abstract

本发明公开了一种全天气地表温度反演方法和系统，该方法包括：获取与地表温度相关的数据，包括：热红外地表温度数据、被动微波数据、表面温度再分析数据、归一化植被指数和数字高程模型数据；如果被动微波数据有值，则基于随机森林模型预测地表温度；如果被动微波数据缺失，则基于温度年循环模型预测地表温度。本发明的技术方案克服了由于地表温度的被动微波数据缺失导致的地表温度反演数据不准确的问题，能够高效、准确地预测全天气地表温度。

Description

一种全天气地表温度反演方法和系统

技术领域

本发明涉及遥感数据分析和处理领域，具体地，涉及一种全天气地表温度反演方法和系统，更具体地，涉及一种基于随机森林模型和温度年循环模型(RF-ATC)的全天气地表温度反演方法和系统。

背景技术

地表温度(Land surface temperature,LST)是通过遥感技术获取的地表热力学表面温度，它与地球系统表面能收支和水平衡紧密相关。其中，全天气(All-weather)地表温度因其能完整够描述温度变化的能力，已被广泛应用于气候变化、农业监测、城市热岛、火灾探测、公共卫生和环境问题等多个研究领域。由于太阳辐射的时空变化，大气条件的不稳定性和土地表面的异质性引起地表温度的快速变化，准确的全天气地表温度数据获取有一定难度。

基于地面观测很难获得大规模持续性地表温度观测，然而遥感技术尤其是热红外(thermal infrared,TIR)卫星的发射，使得大规模高时空分辨率地表温度获取成为可能。然而，由于热红外遥感其极易受到云层污染和大气条件的影响，导致地表温度缺失超过一半，严重影响地表温度数据的实际应用。被动微波(Passive microwave,PMW)遥感因其穿透云雾的能力，被应用于全天气地表温度的反演。

尽管已有研究尝试结合热红外与被动微波遥感重建高时空分辨率全天气地表温度，目前全天气温度数据生产仍然存在以下问题：(1)辅助数据获取困难、复杂计算过程以及地表温度和地下温度的差异；(2)大量的缺失数据及其显著时空异质性；(3)微波数据由于轨道而产生的数据缺失，这三点局限严重影响了大规模全天气地表温度重建方法。

发明内容

本发明所要解决的技术问题是解决由于地表温度的被动微波数据的数据缺失而影响大规模全天气地表温度重建的技术问题。

根据本发明的一个方面，提供一种全天气地表温度反演方法，包括以下步骤：

S1:获取与地表温度相关的数据，包括：热红外地表温度数据、被动微波数据、表面温度再分析数据、归一化植被指数和数字高程模型数据；

S2:如果被动微波数据有值，则基于随机森林模型预测地表温度；

S3:如果被动微波数据缺失，则基于温度年循环模型预测地表温度。

可选地，在步骤S1之后还包括数据预处理步骤S11，包括：将来自不同坐标系的所述与地表温度相关的数据转换为相同的投影，将所述被动微波数据和热红外地表温度数据重采样到所需分辨率。

可选地，数据预处理步骤S11还包括：剔除所述热红外地表温度数据中平均温度误差>3K的数据。

可选地，在步骤S2中，采用在晴空下的热红外地表温度数据、表面温度再分析数据、归一化植被指数和数字高程模型数据对随机森林模型进行训练，随机森林模型的公式如下：

LST_i,j

＝f_RF-i(PMW_i,j(06V,06H,23V,23H,89V,89H),SKT_i,j,NDVI_i/16,j,DEM_j)

i指的是第i天(DOY)；j表示第j个像元；f_RF-i是第i天的随机森林模型；LST是晴空下热红外地表温度数据MODIS LST；PMW代表被动微波数据的六个通道，包括06V，06H，23V，23H，89V和89H；SKT，NDVI和DEM分别表示表面温度再分析数据，归一化植被指数和数字高程模型数据；

采用训练后的随机森林模型预测地表温度。

可选地，所述温度年循环模型的公式如下：

j是第j个像元；t表示一年的第t天；φ表示温度年循环模型；T_0,j,A_j和θ_j分别表示温度年循环模型的年均值，振幅和相位；d表示一年的天数；b_j是ΔT_SKT(t)的系数；ΔT_SKT(t)是观测到的表面温度和初始温度年循环模型预测的表面温度之差。

可选地，所述ΔT_SKT(t)用以下公式表示：

ΔT_SKT(t)＝T_SKT(t)-φ_SKT(T_{SKT_0,j},A_{SKT_j},θ_{SKT_j},t)

＝T_SKT(t)-(T_{SKT_0,j}+A_{SKT_j}sin(2πd^-1t+θ_{SKT_j}))

T_SKT是表面温度再分析数据；φ_SKT是基于表面温度的初始温度循环模型；T_{SKT_0,j},A_{SKT_j}和θ_{SKT_j}分别表示基于表面温度的初始温度循环模型的年均值、振幅和相位。

可选地，采用Sigma剔除法，用于剔除由于有偏训练而产生的随机森林模型预测异常值，剔除的随机森林模型预测异常值满足下列公式：

ΔT_LST＜(x_m-2S)orΔT_LST＞(x_m+2S)

ΔT_LST表示随机森林模型预测的地表温度与温度年循环模型预测的地表温度的差；x_m表还中位数；S是用中位数代替的标准偏差。

根据本发明的第二方面，提供一种全天气地表温度反演系统，包括：

数据获取模块，用于获取与地表温度相关的数据，包括：热红外地表温度数据，被动微波数据，表面温度再分析数据，归一化植被指数和数字高程模型数据；

随机森林模型预测模块，用于如果被动微波数据有值，则基于随机森林模型预测地表温度；

温度年循环模型预测模块，用于如果被动微波数据缺失，则基于温度年循环模型预测地表温度。

可选地，在所述数据获取模块之后，还包括数据预处理模块，用于将来自不同坐标系的所述与地表温度相关的数据转换为相同的投影，将所述被动微波数据和热红外地表温度数据重采样到所需分辨率。

本发明的有益技术效果包括：克服了由于地表温度的被动微波数据缺失导致的地表温度反演数据不准确的问题，能够高效、准确地预测全天气地表温度。

应当理解的是，以上的一般描述和后文的细节描述仅是示例性和解释性的，并不能限制本公开。

附图说明

为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施方式，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是根据本发明的全天气地表温度反演方法的实施例的流程图；

图2是根据本发明的全天气地表温度反演方法的实施例的一个可选实施方式的流程图。

图3表示在白天基于随机森林模型预测的结果以及本发明的方法预测的结果。

图4表示在夜晚基于随机森林模型预测的结果以及本发明的方法预测的结果。

图5是根据本发明的全天气地表温度反演系统的实施例的结构示意图。

图6表示本发明的全天气地表温度反演方法和系统的精度验证结果。

具体实施方式

这里将详细地对示例性实施例进行说明，其示例表示在附图中。下面的描述涉及附图时，除非另有表示，不同附图中的相同数字表示相同或相似的要素。以下示例性实施例中所描述的实施方式并不代表与本发明相一致的所有实施方式。相反，它们仅是与如所附权利要求书中所详述的、本发明的一些方面相一致的系统的例子。

参见图1和2，根据本发明的一个方面，提供一种全天气地表温度反演方法，包括以下步骤：

S1:获取与地表温度相关的数据，包括：热红外地表温度数据、被动微波数据、表面温度再分析数据、归一化植被指数和数字高程模型数据。

S2:如果被动微波数据有值，则基于随机森林模型预测地表温度；

S3:如果被动微波数据缺失，则基于温度年循环模型预测地表温度。

本发明的技术方案克服了由于地表温度的被动微波数据缺失导致的地表温度反演数据不准确的问题，能够高效、准确地预测全天气地表温度。

在步骤S1中，本发明中采用了以下研究数据,详情见表1和表2。

1、热红外地表温度数据

应用Aqua MODIS每日地表温度产品MYD11A1(collection 6)作为本发明中的热红外地表温度数据。MODIS LST产品基于分裂窗算法(generalization split-windows,GSW)应用(当地时间1:30和13:30)无云观测数据观测反演地表温度。已有研究证明MYD11A1的温度反演误差小于1K，均方根误差(root mean square errors,RMSE)小于0.5K。

2、被动微波数据

本发明应用GCOM-W1卫星上搭载AMSR-2传感器获取的Level 2A每日被动微波数据以补充天气信息，该卫星与MODIS Aqua同时过境(大约在当地时间1:30和13:30)。AMSR-2被动微波数据从G-Portal中下载，G-Portal是一个允许用户搜索和下载由日本航天探索局(Japan Aerospace Exploration Agency’s,JAXA’s)地球观测卫星产品的门户系统。AMSR-2上有六个水平和垂直极化通道，它们分别是6.925、10.65、18.4、23.8、36.5和89GHz，考虑到这些通道的特征和算法效率，选择了06、23和89GHz这三个通道应用于本发明中。

3、表面温度再分析数据

表面温度(Skin temperature,SKT)表征表层温度及表面通量传输，本发明应用ERA5-Land小时0.1°x 0.1°气候再分析数据集中的表面温度提供天气变化信息。表面温度从欧洲天气预报中心(European Centre for Medium-Range Weather Forecasts,ECMWF)下载。研究表明ERA-5比先前的版本(ERA-Interim)有更高的时空分辨率，适用于各种应用。

4、其他数据

归一化植被指数(Normalized Difference Vegetation Index,NDVI)和数字高程模型(Digital Elevation Model,DEM)用于补充物候，地形和空间信息。归一化植被指数选用的是EARTHDATA提供的Aqua MODIS MYD12A2，数据高程模型选用的是美国地质调查局(United States Geological Survey,USGS)提供的GTOPO30。

表1

本发明使用数据的详细信息

在步骤S1之后和步骤S2之前，还需要对本发明中所有的数据进行预处理使其空间匹配并为后续过程做好准备的数据预处理步骤S11,数据预处理包括：将来自不同坐标系的数据都转换为相同的投影然后应用最近邻方法将被动微波数据AMSR-2和表面温度SKT重采样到所需分辨率；该分辨率可以是1千米。

根据本发明实施例的一个可选实施方式，预处理步骤S11还包括剔除MODIS地表温度“平均温度误差>3K”的数据以尽可能避免云雾的影响。

在步骤S2中，由于被动微波数据中存在轨道缺失值，需要对其掩膜区中：MODISLST和AMSR-2微波数据均有值的像素标记为随机森林的训练数据；MODIS LST缺失而AMSR-2微波数据有值的像素标记为由随机森林预测的地表温度；MODIS LST和AMSR-2都缺失的像素标记为由温度年循环模型预测的地表温度。

研究表明基于统计方法重建地表温度的基本假设是：晴朗天气条件下地表温度与选定变量之间的关系类与多云条件下的关系是一致的。具体来说，通过晴朗条件下的训练样本对随机森林模型进行训练，然后应用该模型预测云遮盖情况下缺失的地表温度数据。随机森林模型主要可分为以下四个步骤：数据准备、采样和增长、袋外验证和预测。由于随机森林算法已在scikit-learn python库中实现，因此本发明不再赘述其算法具体内容。晴空条件下随机森林模型由PMW，SKT，NDVI和DEM进行训练(公式(1))，然后由训练后的模型预测地表温度缺失值。

i指的是第i天(DOY)；j表示第j个像元；f_RF-i是第i天的随机森林模型；LST是晴空下MODIS LST；PMW代表微波数据的六个通道，包括06V，06H，23V，23H，89V和89H；SKT，NDVI和DEM分别表示表明温度、归一化植被指数和数字高程模型。其中，NDVI的下标为i/16,j表示是每16天一景。

在步骤S3中，基于温度年循环模型进一步预测地表温度可以补充由AMSR-2轨道缺失导致的随机森林预测结果的缺失部分。此外，温度年循环模型因其“半物理”模型的稳健特征，可以弥补随机森林模型受训练数据时空分布影响而产生的异常预测值。研究表明，地表特征和气象条件可以用来表征初始温度年循环模式的天气波动。因此基于Liu(2019)提出的温度年循环综合框架，选用表面温度作为辅助数据由于其被证明是提高温度年循环模型精度的主要辅助数据，构建温度年循环模型如公式(2)。

j是第j个像元；t表示一年的第t天；φ表示温度年循环模型T_0,j,A_j和θ_j分别表示温度年循环模型的年均值、振幅和相位；d表示一年的天数(2018年有365天)；b_j是ΔT_SKT(t)的系数。

ΔT_SKT(t)是观测到的表面温度和初始温度年循环模型预测的表面温度之差(公式(3))。之前的研究已经证明观测的表面温度和预测的表面温度的差异与卫星反演的地表温度和模型预测的地表温度的差异是一致的。因此，应用表面温度差异补偿初始温度年循环模型的天气波动信息。

ΔT_SKT(t)＝T_SKT(t)-φ_SKT(T_{SKT_0,j},A_{SKT_j},θ_{SKT_j},t)

＝T_SKT(t)-(T_SKT__0,j+A_{SKT_j}sin(2πd^-1t+θ_{SKT_j})) (3)

T_SKT是表面温度再分析数据；φ_SKT是基于表面温度的初始温度循环模型；T_{SKT_0,j},A_{SKT_j}和θ_{SKT_j}分别表示基于表面温度的初始温度循环模型的年均值、振幅和相位。

温度年循环模型凭借其表征地表温度时间动态变化的能力，可用于剔除由于有偏训练而产生的随机森林模型预测异常值。

根据本发明实施例的一个可选实施方式，Sigma剔除法用于筛选异常值(公式(4)和(5))。基于sigma法剔除异常值后，重新构建温度年循环模型，并预测缺失地表温度。

ΔT_LST＜(x_m-2S)orΔT_LST＞(x_m+2S) (4)

ΔT_LST表示随机森林模型预测的地表温度与温度年循环模型预测的地表温度的差；x_m表还中位数；S是用中位数代替的标准偏差。

图3和图4分别展示了在白天和夜晚的MODIS LST基于随机森林模型预测的结果以及本发明提出的方法预测的结果。在白天，尽管MODIS LST数据被云覆盖严重影响，但是当AMSR-2微波数据有值时，随机森林模型生成的全天气地表温度具备刻画温度时空变化的能力。另外，通过比较随机森林的结果和本发明提出的完整算法的结果，温度年循环模型不仅可以在AMSR-2微波数据缺失时生成全天气地表温度，并且能够剔除有偏训练(例如，图3第四行)的随机森林模型生成的异常值(参见图4(b4)和(c4))。之前的研究对AMSR-2微波数据缺失的处理是直接采用时间聚合，这种解决方法会导致缺失数据处有明显的接缝。本发明提出的方法，弥补AMSR-2缺失并生成完整的全天气地表温度时，并没有引入明显的接缝。通过目视解译，本发明所提出的方法生成的全天气地表温度能够刻画温度时空变化。

参见图5，根据本发明的另一方面，本发明提供了一种全天气地表温度反演系统，包括：

数据获取模块10，用于获取与地表温度相关的数据，包括：热红外地表温度数据，被动微波数据，表面温度再分析数据，归一化植被指数和数字高程模型数据；随机森林模型预测模块12，用于如果被动微波数据有值，则基于随机森林模型预测地表温度；温度年循环模型预测模块13，用于如果被动微波数据缺失，则基于温度年循环模型预测地表温度。

根据本发明实施例的一个可选实施方式，在数据获取模块之后，还包括数据预处理模块11，用于将来自不同坐标系的所述与地表温度相关的数据转换为相同的投影，将所述被动微波数据和热红外地表温度数据重采样到所需分辨率。

精度验证

本发明应用黑河流域联合遥测实验研究(Heihe Watershed Allied TelemetryExperimental Research,HiWATER)的地表温度测量值验证本发明的方法的准确性。选用的站点包括：大沙龙站(DSL),花寨子站(HZZ)，垭口站(YK)，阿柔站(AR)。HiWATER的目的是在跨学科的基础上应用遥感技术提高水文和生态过程的可观测性。由于HiWATER已被广泛应用于MODIS LST产品和地表温度时空插补方法的精度验证中，所以本发明也采用其全天气地表温度观测开展精度验证。表2描述了选定站点的地理位置和详细信息。

表2

选用HiWATER站点的信息

图6表明，观测的地表温度与本发明预测的地表温度之前具有强相关性。在晴空下，MODIS LST与地表观测温度的平均RMSE白天为3.06K，夜晚为3.34K。在多云天气下，本发明预测的地表温度和地表观测温度的平均RMSE白天为4.12K，夜间为3.88K。夜间的RMSE相对较低，可能是由于夜间地表温度波动较小，更易模拟。MODIS LST与地表观测的相关性比本发明预测的温度与地表观测的相关性更高，这是因为预测模型的误差来自于卫星数据与地表观测温度之间的系统误差以及模型误差两方面。由于以前的研究表明MODIS LST与地表观测温度之间存在较大系统偏差(-3.44K)，所以本发明的精度需要将MODIS LST的RMSE作为系统误差后考虑模型误差。经过精度验证，本发明被证明是一种生成全天气地表温度的高效、准确的方法，在未来大规模应用中有很大潜力。

表3

精度验证结果

本发明的特征和益处通过参考实施例进行说明。相应地，本发明明确地不应局限于这些说明一些可能的非限制性特征的组合的示例性的实施例，这些特征可单独或者以特征的其它组合的形式存在。

以上所述实施例，仅为本发明的具体实施方式，用以说明本发明的技术方案，而非对其限制，本发明的保护范围并不局限于此，尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，其依然可以对前述实施例所记载的技术方案进行修改或可轻易想到变化，或者对其中部分技术特征进行等同替换；而这些修改、变化或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以权利要求的保护范围为准。

Claims (9)

1.一种全天气地表温度反演方法，其特征在于，包括以下步骤：

S1:获取与地表温度相关的数据，包括：热红外地表温度数据、被动微波数据、表面温度再分析数据、归一化植被指数和数字高程模型数据；

S2:如果所述被动微波数据有值，则基于随机森林模型预测地表温度；

S3:如果所述被动微波数据缺失，则基于温度年循环模型预测地表温度。

2.根据权利要求1所述的全天气地表温度反演方法，其特征在于，在所述步骤S1之后还包括数据预处理步骤S11，包括：将来自不同坐标系的所述与地表温度相关的数据转换为相同的投影，将所述被动微波数据和所述热红外地表温度数据重采样到所需分辨率。

3.根据权利要求2所述的全天气地表温度反演方法，其特征在于，所述数据预处理步骤S11还包括：剔除所述热红外地表温度数据中平均温度误差>3K的数据。

4.根据权利要求1或2所述的全天气地表温度反演方法，其特征在于：

在所述步骤S2中，采用在晴空下的所述热红外地表温度数据、表面温度再分析数据、归一化植被指数和数字高程模型数据对随机森林模型进行训练，随机森林模型的公式如下：

LST_i,j＝f_RF-i(PWW_i,j(06V,06H,23V,23H,89V,89H),SKT_i,j,NDVI_i/16,j,DEM_j)

i指的是第i天；j表示第j个像元；f_RF-i是第i天的随机森林模型；LST是晴空下热红外地表温度数据MODIS LST；PMW代表被动微波数据的六个通道，包括06V，06H，23V，23H，89V和89H；SKT，NDVI和DEM分别表示表面温度再分析数据、归一化植被指数和数字高程模型数据；

采用训练后的随机森林模型预测地表温度。

5.根据权利要求1或2所述的全天气地表温度反演方法，其特征在于：所述温度年循环模型的公式如下：

j是第j个像元；t表示一年的第t天；φ表示温度年循环模式；T_0,j,A_j和θ_j分别表示温度年循环模型的年均值、振幅和相位；d表示一年的天数；b_j是ΔT_SKT(t)的系数；ΔT_SKT(t)是观测到的表面温度和初始温度年循环模型预测的表面温度之差。

6.根据权利要求5所述的全天气地表温度反演方法，其特征在于：所述ΔT_SKT(t)用以下公式表示：

T_SKT是表面温度再分析数据；φ_SKT是基于表面温度的初始温度循环模型；T_{SKT_0,j},A_{SKT_j}和θ_{SKT_j}分别表示基于表面温度的初始温度循环模型的年均值，振幅和相位。

7.根据权利要求1或2所述的全天气地表温度反演方法，其特征在于：采用Sigma剔除法，用于剔除由于有偏训练而产生的随机森林模型预测异常值，剔除的随机森林模型预测异常值满足下列公式：

ΔT_LST＜(x_m-2S)或ΔT_LST＞(x_m+2S)

S＝1.4826×median{|ΔT_LSTi-x_m|}

ΔT_LST表示随机森林模型预测的地表温度与温度年循环模型预测的地表温度的差；x_m表还中位数；S是用中位数代替的标准偏差。

8.一种全天气地表温度反演系统，其特征在于，包括：

数据获取模块，用于获取与地表温度相关的数据，包括：热红外地表温度数据、被动微波数据、表面温度再分析数据、归一化植被指数和数字高程模型数据；

随机森林模型预测模块，用于如果所述被动微波数据有值，则基于随机森林模型预测地表温度；

温度年循环模型预测模块，用于如果所述被动微波数据缺失，则基于温度年循环模型预测地表温度。

9.根据权利要求8所述的全天气地表温度反演系统，其特征在于，在所述数据获取模块之后，还包括数据预处理模块，用于将来自不同坐标系的所述与地表温度相关的数据转换为相同的投影，将所述被动微波数据和热红外地表温度数据重采样到所需分辨率。

Images