多因子地表温度降尺度研究

刘英，朱蓉，钱嘉鑫，党超亚，岳辉

(西安科技大学 测绘科学与技术学院，西安 710054)

0 引言

地表温度(land surface temperature，LST)是反映地球表面与大气之间能量平衡和影响地表水循环的重要物理参数之一[1-3]，主要用于植被监测[4-5]、土壤水分监测[6-7]、城市热岛效应[8-9]等方面研究。获取地表温度的方法有地面观测和热红外遥感反演等[10]，其中地面观测获取地表温度精度高但覆盖范围小，而热红外遥感能获得较大范围内的地表温度，且能表现LST的空间分布特征，成为获取地表温度的主要手段。但热红外传感器无法同时满足高时空分辨率的数据要求，限制了地表温度在不同区域尺度上的应用。为解决热红外传感器所获取地表温度的这一问题，相关学者提出了多种地表温度降尺度方法。

目前，卫星遥感地表温度空间降尺度的方法主要基于地表参数统计回归、数据融合、调制分配及光谱混合模型等展开[11]。其中，基于地表参数统计回归的温度降尺度法主要原理为“尺度不变性”，将低分辨率下的地表温度与光谱指数的关系应用到高分辨率下，最常用的方法是Kustas等[12]于2003年提出的DisTrad(disaggregation procedure for radiometric surface temperature)法以及Agam等[13]于2007年在前者基础上提出的TsHARP (an algorithm for sharpening thermal imagery)法。DisTrad方法基于归一化植被指数(normalized difference vegetation index，NDVI)和LST间的统计关系提出，但该方法并未充分考虑土壤水分在降尺度模型中的作用，TsHARP方法用更能表征植被状态的植被覆盖度代替NDVI与LST建立线性关系，从而达到降尺度目的。数据融合和空间滤波均能提高空间分辨率，基于这一特性，Fasbender等[14]提出了贝叶斯方法；Dominguez等[15]提出了基于高分辨率的反照率及NDVI数据的地表温度降尺度方法;Yang等[16]总结出了DisEMI(disaggregation method for subpixel temperature using the remote sensing endmember index based technique)算法。Guo等[17]提出考虑地形变化的温度降尺度方法——PBIM(pixel block intensity modulation)；Norman等[18]总结出了DisALEXI (disaggregated atmosphere land exchange inverse model)算法来研究微气候测量；Merlin 等[19]利用物理模型对被动微波反演的土壤湿度进行降尺度；Merlin等[20]基于植被覆盖度、土壤蒸发率等各种因子对MODIS农业区域的地表温度进行降尺度研究。随着国内外学者对温度降尺度方法的深入研究，Gao等[21]利用基于权重函数的时空数据融合方法提出了时空自适应反射率融合模型；邬明权等[22]提出一种改进的遥感时空数据融合方法(spatial and temporal data fusion mode，SATDFM)对MODIS数据和TM数据进行时空数据融合；Yang等[23]提出基于随机森林多遥感因子的温度降尺度方法，该方法避免了降尺度因子的单一化以及回归模型线性化给温度降尺度带来的不确定性。DisTrad方法和TsHARP方法因其计算方法简便、降尺度因子容易获取而被广泛应用，但他们均是针对地表类型单一的植被覆盖区，并未考虑建筑物、裸地、水体等地物类型复杂的情况。鉴于这一缺陷，杨英宝等[10]将NDVI、归一化水体指数(modified normalized difference water index，MNDWI)及城市不透水面指数引入DisTrad模型，建立一种新的精度较高的降尺度模型，该方法未充分考虑到研究区中存在的裸土区域；李小军等[24]根据地表覆盖类型的不同，分别选择与LST相关性更好的NDVI、归一化建造指数、MNDWI以及增强型裸土指数提出了新的转换模型，并从定性和定量2个角度评价了TsHARP法和新模型的降尺度精度。该方法分别建立各因子与LST之间的降尺度模型，根据最小残差和最大决定系数，最终取其中最优的单个模型建立与LST的降尺度模型。

本文充分考虑到植被、水体、建筑物与裸地同时对地表温度的影响，且相对于文献[10]考虑到了裸土影响因子。在DisTrad和TsHARP降尺度方法的基础上，提出一种新的反映多地表覆盖类型地表温度降尺度方法——G_DisTrad方法。该方法以表征植被状况的归一化差异植被指数NDVI、表征水体的归一化水体指数MNDWI及表征建筑物与裸地的指数(normalized difference built-up and soil index，NDBSI)[25]为自变量，同时与LST建立多元线性关系，并以地表覆盖类型复杂的天津、哈尔滨、南京、郑州、重庆、西安共6个大城市区域作为实验区，并与DisTrad法、TsHARP法进行对比，进而分析G_DisTrad法的可行性与适用性。

1 研究区与数据来源

1.1 研究区

选择华北地区的天津、东北地区的哈尔滨、华东地区的南京、中南地区的郑州、西南地的重庆及西北地区的西安共6个大城市地表覆盖类型复杂区为研究对象。各研究区均包含植被、水体、建筑物和裸地等多种地物类型。

1.2 数据来源及处理

本研究所用晴空无云的Landsat-8 OLI影像来自美国地质调查局网站(https：//earthexplorer.usgs.gov/)。Landsat-8 OLI影像包括分辨率为30 m的可见光波段和分辨率为100 m的热红外波段(下载时已重采样至30 m)。用于验证的实测地表温度数据来自中国气象数据网(http：//data.cma.cn/)。

数据预处理包括辐射定标、大气校正、影像裁剪、重采样等。本文以美国地质调查局网站下载的热红外波段经辐射传输方程法[26]计算的地表温度(LST)作为间接验证数据。在ENVI中经过Pixel Aggregate像元聚合获得180 m低分辨率下的Landsat-8可见光和热红外波段数据，分别计算各指标因子并将其作为降尺度模型的自变量和因变量。根据DisTrad法、TsHARP法和本文提出的G_DisTrad法各自的原理建立各因子之间的降尺度关系模型，从而将180 m分辨率的LST数据降尺度到30 m，并且对比分析这3种方法在不同城市区域的降尺度效果。

2 降尺度方法

2.1 DisTrad模型

DisTrad模型的具体实现过程如式(1)至式(4)所示。

(1)

LST′L(NDVI)=a*NDVIL+b

(2)

det(LST)=LSTL(NDVI)-LST′L(NDVI)

(3)

LSTH(NDVI)=a*NDVIH+b+det(LST)

(4)

式中：b1、b2分别代表近红外与红光波段的反射率；LST′L(NDVI)表示模拟的低分辨率下LST；LSTL(NDVI)为低分辨率下的LST；NDVIL表示低分辨率下的NDVI；det (LST)表示回归残差；a，b为常数；NDVIH表示高分辨率下的NDVI；LSTH(NDVI)表示DisTrad方法的降尺度结果。

2.2 TsHARP模型

Weng等[27]发现地表温度与植被覆盖度相对于与NDVI具有更强的相关关系。因此，Agam等[13]在DisTrad模型基础上，利用与地表温度有更高相关性的植被覆盖度(FVC)代替NDVI，与LST建立线性关系。FVC的计算方法如式(5)所示。

(5)

式中：NDVImax表示NDVI累积频率置信度为95%时对应的NDVI值；NDVImin表示NDVI累积频率置信度为5%的NDVI值。

2.3 G_DisTrad模型

DisTrad与TsHARP模型适用于地表覆盖类型单一区域，而实际情况中地表覆盖类型不仅有植被，还有建筑物、裸地、水体等。本文基于地表类型复杂的特点，提出一种改进的DisTrad法(G_DisTrad)。该方法以NDVI、MNDWI及代表建筑物与裸地的NDBSI为自变量，与LST建立多元线性关系，从而实现降尺度的目的。具体的实现过程如式(6)～式(12)所示。

(6)

(7)

(8)

(9)

LST′L=a1*NDVIL+a2*MNDWIL+
a3*NDBSIL

(10)

det(LST)=LSTL-LST′L

(11)

LSTH=a1*NDVIH+a2*MNDWIH+
a3*NDBSIH+b+det(LST)

(12)

式中：b1、b2、b3、b4、b5为Landsat-8多光谱影像中蓝、绿、红、近红外、短波红外波谱反射率；LST′L为模拟的低分辨率下的LST；LSTL为低分辨率下的LST；NDVIL表示低分辨率下的NDVI；MNDWIL表示低分辨率下的MNDWI；NDBSIL表示低分辨率下的NDBSI；det (LST)为回归残差；a1、a2、a3、b为常数；NDVIH、MNDWIH、NDBSIH表示高分辨率下的NDVI、MNDWI和NDBSI；LSTH为G_DisTrad方法降尺度的结果。

2.4 降尺度结果的精度评价指标

降尺度结果优劣可由决定系数R2和相对均方根误差RMSE来表示。R2表示降尺度后LST的结果图像与已知LST之间的相关性，其取值范围为[0，1]，值越接近1，说明相关性越高；RMSE用来判断降尺度前后温度图像的一致性，值越小，说明降尺度前后2幅图像的一致性越高。

3 实验结果与分析

将中国气象数据网上获得的实测地表温度数据作为直接验证数据，将辐射传输模型法反演的Landsat 30 m分辨率下的LST作为间接验证数据。将180 m分辨率下的NDVI、MNDWI及NDBSI作为自变量，180 m分辨率下LST作为因变量，建立二者之间的关系模型，分别获得3种方法在不同研究区域温度降尺度结果。

3.1 地表温度降尺度结果直接验证

本研究以天津市为例进行降尺度结果的直接验证。由于研究区只有3个实测点，为了选取较多的实测验证点，本文选取2019年4、5、7、8月等不同月份各1景天津市研究区影像分别作为实测数据点来源。由图1可知，3种降尺度方法均能较好地将180 m分辨率的LST降尺度到30 m分辨率，且保留了较多的地物分布信息，达到较好的降尺度效果；本文改进的G_DisTrad方法相对于DisTrad和TsHARP方法降尺度的图像分布更接近于直接用辐射传输方程反演的30 m分辨率的地表温度分布图，尤其是在水体与建筑物地区。由图2可知，G_DisTrad与实测地表温度的决定系数最高(R2=0.697 0)，其次为TsHARP(R2=0.679 1)，DisTrad最低(R2=0.678 0)。

图1 各降尺度方法结果图

图2 实测结果验证

3.2 各地区地表温度降尺度结果间接验证

在ArcGIS中采用随机打点的方法获取500个随机点进行6个研究区域3种温度降尺度结果的间接验证。各降尺度方法在各研究区的降尺度拟合方程以及降尺度结果的相关性如表1所示。

表1 地表温度降尺度方法对比

注：y为LST；b1为NDVI；b2为MNDWI；b3为NDBSI；b4为FVC。

由表1可知，回归模型中NDVI、MNDWI、FVC与LST呈负相关，NDBSI与LST呈正相关。3种降尺度方法在这6个研究区均得到较好的降尺度效果，G_DisTrad模型平均决定系数最高(R2=0.818 0)，其次为TsHARP模型(R2=0.805 0)，DisTrad的平均R2为0.787 0，整体来看，G_DisTrad模型的拟合效果要明显优于DisTrad模型、TsHARP模型。除重庆外，G_DisTrad在其余5个地区的R2均高于DisTrad和TsHARP的R2；除南京外，TsHARP在其余5个地区的R2均高于DisTrad的R2。在重庆地区，TsHARP降尺度精度最高(R2=0.865 8)，其次分别为G_TsHARP(R2=0.807 0)、DisTrad(R2=0.797 4)。总体上，本研究提出的温度降尺度方法G_TsHARP在城市复杂区域有较好的适用性，优于TsHARP方法和DisTrad方法。

为进一步对3种方法降尺度效果进行对比分析，本文对各研究区的降尺度效果进行局部验证，分别选取研究区的植被、建筑物与裸地、水体部分，将3种降尺度结果与Landsat-8 30 m热红外波段反演的地表温度进行对比(图3～图8)。

由图3～图8可知，天津市的水体区域、建筑物与裸地区域，G_DisTrad方法的降尺度效果与反演地表温度更接近；植被区域，TsHARP方法和G_DisTrad方法降尺度效果与反演地表温度更吻合，且TsHARP方法的降尺度效果细节更加明显。哈尔滨市的水体区域，G_DisTrad方法降尺度效果与反演地表温度更接近；建筑物与裸地、植被区域，3种降尺度结果的空间分布相近。南京市的水体区域，TsHARP方法和G_DisTrad方法降尺度的效果比DisTrad方法好，DisTrad方法降尺度结果的锯齿形状相对于TsHARP方法和G_DisTrad方法比较明显；G_DisTrad方法在建筑物地区和植被区域与反演的地表温度空间分布相似性更高。郑州市的水体区域，G_DisTrad方法的降尺度效果与反演的地表温度更接近；植被、建筑物与裸地区域3种方法的降尺度结果相近。重庆市的植被区域，TsHARP方法的降尺度结果在空间上的分布与反演的地表温度更接近；水体区域，G_DisTrad方法的降尺度结果的空间分布与反演的地表温度更吻合，DisTrad方法和TsHARP方法的降尺度结果在水体区域边缘锯齿状比较明显；建筑物与裸地地区，TsHARP方法的降尺度结果的空间分布与反演的地表温度更接近。西安市的水体区域，TsHARP方法和G_DisTrad方法的降尺度结果空间分布相似，DisTrad方法的降尺度效果较差；植被区域3种方法均达到了较好的降尺度效果；建筑物地区3种方法降尺度效果与反演的地表温度分布近似，且G_DisTrad方法降尺度结果的细节更明显。综上所述，3种方法中，G_DisTrad方法在水体、建筑物与裸地地区降尺度的结果优于DisTrad方法和TsHARP方法，但在植被区域，G_DisTrad方法的结果相比于TsHARP方法优势不明显。

图3 天津市植被、水体、建筑物与裸地区域的降尺度结果局部对比

图4 哈尔滨市植被、水体、建筑物与裸地区域的降尺度结果局部对比

图5 南京市植被、水体、建筑物与裸地区域的降尺度结果局部对比

图6 郑州市植被、水体、建筑物与裸地区域的降尺度结果局部对比

图7 重庆市植被、水体、建筑物与裸地区域的降尺度结果局部对比

图8 西安市植被、水体、建筑物与裸地区域的降尺度结果局部对比

4 结束语

本文对地表覆盖类型较复杂的城市区域采用传统的DisTrad方法、TsHARP方法以及本文提出的新的G_DisTrad方法进行地表温度降尺度研究，通过对各降尺度方法的结果验证对比得出以下结论。

1)基于天津市实测地表温度的直接精度验证结果表明，G_DisTrad方法的整体降尺度效果最优，G_DisTrad方法降尺度结果的决定系数达到0.697 0，其次是TsHARP方法，其决定系数为0.679 1，DisTrad方法的决定系数最低，为0.678 0。

2)基于Landsat 30 m的LST的间接精度验证结果表明，G_DisTrad方法的温度降尺度效果最好，平均决定系数最高(R2=0.818 0)，其次为TsHARP方法(R2=0.805 0)，最低为DisTrad方法(R2=0.787 0)；除重庆外，G_DisTrad方法的决定系数在3种温度降尺度方法中均为最高。将各研究区以地物类型划分为植被、水体、建筑物与裸地，分别进行3种方法的降尺度结果对比。结果表明,G_DisTrad方法在建筑物与裸地地区的降尺度效果优于DisTrad方法和TsHARP方法，而在高植被覆盖地区(重庆)，TsHARP方法的降尺度效果相对较优。

3)本研究提出的方法原理相对简单，各指标易获取，且精度比传统的DisTrad方法和TsHARP方法在地表覆盖类型较复杂的城市区域高。后续可以继续将植被覆盖度、地形等更多影响因子考虑到G_DisTrad方法中，以期得到精度更高的温度降尺度结果。