CN111274535B - 一种地面站点空间代表性动态评价方法 - Google Patents

Abstract

本发明具体涉及一种地面站点空间代表性动态评价方法，用于评价卫星遥感地表温度验证过程中站点的空间代表性。本发明基于地表温度的时间尺度分解模型，提出一种地面站点空间代表性动态评价方法，解决基于地面站点的地表温度验证过程中的以下两个问题(1)地面站点观测地表温度能不能代表卫星像元尺度的地表温度真值(2)如果能代表，在什么时间段能代表。本发明可以对站点的空间代表性进行评价，能够回答地表温度验证中哪些站点具有好的代表性，在哪些时间段具有代表性，这对于充分利用现有观测站点数据，降低新建站点的成本具有显著的经济意义。

Description

一种地面站点空间代表性动态评价方法

技术领域

本发明具体涉及一种地面站点空间代表性动态评价方法，用于评价卫星遥感地表温度验证过程中站点的空间代表性。

背景技术

作为地表与大气之间能量交互的主要指示因子，地表温度在诸如气候变化、水文过程模型等领域扮演着重要的角色。对于从卫星遥感获取的地表温度数据，开展严格的验证评价工作有利于促进相关产品的发展和应用。目前对遥感地表温度数据的验证评价工作主要包括基于地面站点观测数据的验证、基于辐射亮度的验证、交叉比较和时序分析等四种方法，其中基于地面站点观测数据的验证方法是被广泛接受的一种验证方法。然而，尺度不匹配的影响，直接使用地面观测数据用于地表温度数据的验证会引入空间代表性误差。因此，相关学者开展了站点空间代表性的评价工作。

目前，针对地面站点空间代表性评价方法主要有两种类型：(1)比较点面特征：此类方法主要是借助高分辨率专题数据比较站点所在位置值与卫星像元尺度内的值的差异，如在足迹模型，标准差、变异系数等统计指标、两个尺度的差值等。这些方法操作上比较简单，但是长时间序列的高分辨率数据难以获得。(2)分析区域异质性：此类方法同样借助高分辨率专题数据分析区域的异质性程度，如利用分形几何和半变异函数等。除了长时序高分辨率数据难以获取以外，最根本的问题在于区域异质性不能与空间代表性完全划等号，两者有根本上的定义差异。

相对于诸如地表覆盖类型、生物量等随时间变化较慢的地表参数，地表温度无论是瞬时、日内还是年内都具有很大的变化。这主要是地表温度不仅与地表属性如地表类型有关，还与局地气象状况有关。而受地表类型的影响，地表的热属性也存在着差异，站点与像元尺度内地表在对瞬时气象状况的响应上会出现不同的响应速率。因此，简单的利用基于少量样本评价出的站点空间代表性在时间序列上不具备代表性。

发明内容

本发明的目的，就是针对上述存在的问题及不足，提出的站点空间代表性动态评价方法将基于代表性的定义结合地表属性参数、气象参数对站点空间代表性进行时间序列拓展，达到空间代表性动态评价的目的，从而对遥感地表温度进行更客观的验证。

本发明提供一种站点空间代表性评价方法，动态的评价卫星遥感地表温度验证过程中站点的空间代表性，具体包括：

代表性指在一定的时空维度上测得的数据，能够反映相同或不同时空尺度上真实情况的程度。基于此，本发明将遥感地表温度验证过程中站点的空间代表性定义为：地面传感器观测地表温度真值与卫星像元尺度地表温度真值的差异，

ΔT(t)＝T_F(t)-T_P(t) (1)

式中，ΔT为代表性评价指标，亦写作SRI，单位为K；t为时间，T_F和T_P分别表示地面站点传感器观测温度真值和卫星像元尺度地表温度真值，单位为K。因此，只要知道两个尺度的地表温度真值即可评价站点在空间上的代表性。

地表温度受到地球自转和局地天气状况的双重影响，在时间维度上快速变化。因此，代表性评价指标是随时间变化的函数。Zhan等(2014)将地表温度分解成三个组分，

式中，ATC表示地表温度的季节性变化成分；DTC表示地表温度在日内的变化情况；USC表示地表温度受局地瞬时天气因素等情况影响下的波动成分。考虑到针对同一颗卫星过境时间相近，假设天与天之间同一卫星过境时间差异可以忽略不计，即DTC成分可以忽略。因此，方程(2)可以改写为，

结合方程(3)，方程(1)可以改写为，

式中ΔATC和ΔUSC分别为两个尺度地表温度季节性变化成分的差值和天气因素影响成分差值。

在地表温度季节性变化方面，本发明采用的参数化方案如式(5)所示，

式中，

和A₀分别是季节性变化成分的均值和振幅，θ₀是该曲线的相位。因此根据该参数化方案，在少量离散的地表温度的情况下，可以得到地表温度的时间序列。因此，对于ΔATC部分，可以进一步写成，

式中，

是两个尺度ATC中

的差值；A和θ是该曲线的振幅和相位，与两个尺度ATC中的振幅和相位关系可以表示为，

其中，k为常数，取值为0，±1。

进一步地，在式(4)和(6)的基础上，可以得到离散的ΔUSC。根据前面所述，USC受到局地气象状况的影响。在区域尺度范围内，由于不同地物类型的热属性存在差异，因此对气象状况的响应存在差异，引入表征地表属性的的影响因子和局地气象因子x_i对其进行解释，i＝1,2,…,n，n为影响因子的个数。基于此，ΔUSC的时序化可以表示为，

ΔUSC(t)＝g[x₁(t),x₂(t),……,x_n(t)] (8)

式中，g是表征ΔUSC与影响因子x_i(i＝1,2,…,n)之间的参数化方程。受地理位置及参数的可获取程度及变量之间的相关性程度影响，该参数化方程形式视具体情况而定，通常可以通过多元线性回归或机器学习方法获取。

联立式(1)，(4)，(6)，(8)可以得到时序化的代表性指标ΔT，

至此，根据站点代表性的定义及温度的时间尺度分解模型，得到了时序化的代表性评价指标。基于该指标，可以动态的评价站点观测值在像元尺度上的代表性程度，实现站点空间代表性的动态评价。

本发明的有益效果为，现有的地表温度验证中通常使用地面站点观测地表温度与卫星像元地表温度进行直接比较，由于两者之间存在尺度不匹配等原因，使得验证结果中存在站点空间代表性误差。而常规的站点空间代表性评价方法对于地表温度站点的评价存在欠缺，无法解释地表温度的快速变化导致的空间代表性动态变化。结合本方法，一方面可以对站点的空间代表性进行评价，能够回答地表温度验证中哪些站点具有好的代表性，在哪些时间段具有代表性，这对于充分利用现有观测站点数据，降低新建站点的成本具有显著的经济意义；另一方面，基于该方法的站点代表性评价结果，可以量化和降低地表温度验证中的代表性误差，使得对遥感反演地表温度数据产品的评价更客观，具有明显的科学意义。

附图说明

图1站点空间代表性动态评价流程；

图2从Landsat TM/ETM+获取的MODIS像元内和地面站点视场内的TF，TP以及拟合得到的ATCF、ATCP；

图3分解得到的非稳态成分ΔUSC以及拟合得到的连续ΔUSC；

图4针对MODIS像元的代表性评价指标，图中灰色区域代表SRI绝对值小于1K的区域；

图5考虑与不考虑空间代表性的验证结果对比。

具体实施方式

下面结合附图和实施例，更详细的对本发明中的评价方法进行说明。

实施例

本实施示例中选取了分布在中国不同气候区的16个地面观测站点，如表1所示。这些站点地表温度的观测仪器均为长波辐射计，具有相对较大的视场范围。同时，下垫面类型在站点视场范围内和站点对应像元内也存在不同程度的差异。此外，每个观测站还有风速、气温、相对湿度、短波辐射、大气压等常规气象要素。用于模拟站点视场和卫星像元视场的中高分辨率数据源为Landsat TM/ETM+，空间分辨率为可见光波段30米，热红外波段120米/60米，被重采样到30米，时间分辨率为16天一次，其中热红外波段亮温用于计算地表温度，可见光波段用于计算植被指数以描述地表状况。选取的的待验证地表温度产品为MODIS地表温度产品(MOD11A1)，空间分辨率为1km，每日白天和夜间各过境一次。考虑到选取的Landsat TM/ETM+数据仅有白天，因此，仅验证白天的MODIS地表温度。

表1实施例中用到的地面站点数据基本情况

本例的处理流程如图1所示。

流程一：数据预处理：

1、估算归一化植被指数NDVI：根据近红外波段反射率(ρ_NIR)和红光波段反射率(ρ_Red)计算归一化差值植被指数，

2、反演Landsat地表温度：辐射传输方程，

式中，下标λ表示波长；L_λ为传感器波长为λ的光谱辐亮度，单位为W/(m²·sr·μm)；T_s为地表温度，单位为K；B_λ(T_s)为地表温度为T_s时对应的黑体辐亮度，单位为W/(m²·sr·μm)；τ_λ为地表到传感器的大气透过率；

分别为大气上行、下行长波辐射，单位为W/(m²·sr·μm)，后面三种大气参数通过https://atmcorr.gsfc.nasa.gov在线获取；ε_λ为地表发射率，可以通过NDVI阈值法获取。

式中，ε_vλ与ε_sλ为分别为植被、裸土的发射率；P_v为植被覆盖度；C_λ用于纠正由于地表粗糙度导致的腔体效应，对于平坦地表取值为0；NDVI_s、NDVI_v为区分纯净裸土区域、裸土与植被混合区域、纯净植被区域的阈值。

3、实测站点真实地表温度数据的获取：根据长波辐射计测量的上下行长波辐射计算得到地表的真实温度。

4、提取MODIS地表温度产品中站点所在像元，并根据质量控制标记筛选质量较好的像元。

流程二：地表温度成分分解：

步骤1：根据流程一中计算的地表温度，分别提取MODIS像元范围和站点传感器视场范围内地表温度的平均值。得到离散的T_F和T_P，计算离散的代表性评价指标ΔT。

步骤2：根据式(3)和步骤1中的T_F和T_P，运用最小二乘法分别拟合像元范围内和传感器视场范围内的ATC曲线参数。计算两个尺度的ATC曲线差值ΔATC，该ΔATC差值为逐日连续。

步骤3：根据式(4)，计算非稳定成分差值ΔUSC，此时得到的ΔUSC成分在时间上是离散的。

流程三：非稳定成分时间序列拓展

步骤1：选取NDVI作为描述地表状态变化的指标，并利用HANTS模型将其在时间序列上拓展到每天。然后计算两个尺度的NDVI差值ΔNDVI。

步骤2：假设在公里范围内气象参数的变化小且可以忽略不计，即站点观测气象参数可以代表像元尺度值。从时间序列气象参数中选取与ΔUSC时刻最接近的值，并根据日期选择对应的ΔNDVI。至此，完成ΔUSC与相关解释因子的时间匹配。

步骤3：选择多元线性回归(MLR)和随机森林(RF)两种方法分别对概念模型g(式(8))进行求解，并选取其中精度较好的方法构建最终的求解模型。因此可以将该部分分为两个小部分：I：求解模型选择：随机从匹配的ΔUSC和影响因子对中选择90％样本利用MLR和RF分别进行模型训练，剩余的10％样本对训练模型进行验证。以RMSE和R²为评价指标，选择RMSE较小，R²较大的方法求解概念模型g。II：概念模型训练及ΔUSC预测：根据I中选取的概念模型求解方法，结合匹配的ΔUSC和影响因子对进行概念模型g训练，将训练之后的g模型结合各种因子预测ΔUSC，并得到逐日连续的ΔUSC。

流程四：站点空间代表性动态评价及地表温度验证

步骤1：根据流程二中得到的时序ΔATC和流程三中得到的时序ΔUSC可以得到时序的ΔT，即空间代表性评价指标SRI。

步骤2：根据时序空间代表性评价指标SRI明确各站点的代表性变化情况；同时结合各因子分析决定和影响不同时间段站点空间代表性存在差异的主要因素。

步骤3：将SRI与LST样本对进行匹配，根据不同的应用需求选择一定阈值(以SRI绝对值为1K为例)内的样本进行统计评价，其结果为该产品的验证结果。

结果

根据流程二中的处理，从Landsat TM/ETM+中获取得到了MODIS像元内和地面站点视场内的T_F，T_P以及拟合得到了两条ATC曲线(图2)，同时还得到了对应的差值ΔT和ΔATC。ΔATC与ΔT之间的决定系数R²为0.01—0.63，用ΔATC曲线并不能很好的描述ΔT的在时间上的波动情况。进而，根据流程三中的非稳定成分时间序列拓展方法对ΔUSC进行了时间序列拓展，结果如图3所示。从图中可以看出，连续的ΔUSC与离散的ΔUSC之间具有较好的相关性，两者之间的R²为0.77—0.93，且两者之间的偏差较小，MBE值为-0.08K—0.11K。根据流程四中的处理，生成了时间连续的代表性评价指标(图4)。从图中可以看出，相对于平滑的ΔATC曲线，连续的SRI曲线更捕捉到更多的抖动信息。比较离散的和连续的SRI，两者之间的R²为0.77—0.93，MBE为-0.08K—0.11K。根据相关影响因子与SRI的关系，及时序SRI可以进一步地开展站点观测数据的空间代表性动态变化情况及影响因素分析。

以盈科站为例进行分析，该站点周围主要地表覆盖类型为农田、建筑物和防风树等。站点和像元尺度内的NDVI差异(ΔNDVI)为-0.01—0.19，是一个异质性比较强的站点，从异质性角度来看，该站点不适合作为1-km尺度的卫星产品验证。但是，从图中可以看出，该站点在DOY＝1—150以前和DOY＝340—365以及DOY＝270—300这些时间段内，温度差值明显小于1K，根据前面的选择的阈值，在该时间段内可以用于卫星产品验证，将该时间段定义为可代表性时间段，其他时间段定义为非代表性时间段。进一步地，分析这两个时间段内影响两个尺度温度差异的主要因素。在非代表性时间段内，主要为植被生长季，影响因素主要为植被覆盖差异(ΔNDVI)，其与地表温度差异(ΔT)的相关系数(R)为0.83，呈明显的相关关系，即当植被覆盖差异越大时，地表温度差异越大，具体而言，主要与植被生长过程有关。在代表性时间段，植被覆盖程度较低，裸露农田和建筑物尽管无植被覆盖，但两者之间的热惯量存在差异，同样对外界天气因素的响应存在差异，因此，在该时间段内仍然存在代表性较弱的情况。出现代表性较好的时刻主要表现为反照率(albedo)较大的时候(albedo>0.2，R＝0.70)。结合当地实地情况，该时间段内有降雪情况，即降雪使地表的热惯量趋于一致。

总体而言，这些站点表现出了不同程度的代表性差异，其中关滩站由于两个尺度内的地表覆盖差异极大而无代表性较好的时间段；部分站点(如长白山站、鼎湖山站、馆陶站、花寨子站等)因为单一且均匀下垫面在一年内几乎所有时间段均表现出了较好的代表性；部分站点尽管下垫面存在较大差异(如大兴站)，但在所选时间段内亦表现出了较好的空间代表性；剩下站点的空间代表性则在所选时间段内具有较大的波动。

将代表性评价指标用于MODIS地表温度产品的验证中，结果如图5所示。可以看出，在仅考虑代表性较好时间段的样本情况下，MODIS地表温度与站点地表温度的散点比全部样本的散点更接近1:1线，且具有较好的相关性。这主要原因是代表性较弱时间段的卫星地表温度不仅受到反演过程中的各种因素的影响(如发射率、观测角度等)，还受到验证数据由于空间代表性的影响，使得这样的验证结果中存在空间代表性误差。通过本方法的处理，可以有效的剔除验证过程中受空间代表性影响较大的样本点，极大的降低了空间代表性对验证结果的影响，从而使验证结果更加公平公正。同时，对于验证站点的选择也有意义，如关滩站在前置阈值条件(SRI绝对值小于1K)下不适用用于MODIS地表温度产品的验证。

Claims (1)

1.一种地面站点空间代表性动态评价方法，该方法用于评价卫星遥感地表温度验证过程中站点的空间代表性，所述空间代表性的定义为在一定的时空维度上测得的数据，能够反映相同或不同时空尺度上真实情况的程度；其特征在于，所述评价方法包括以下步骤：

S1、将遥感地表温度验证过程中站点的空间代表性定义为，地面站点传感器观测地表温度真值与卫星像元尺度地表温度真值的差异：

ΔT(t)＝T_F(t)-T_P(t)

其中，ΔT为代表性评价指标，t为时间，T_F和T_P分别表示地面站点传感器观测温度真值和卫星像元尺度地表温度真值，下标F和P用于区分地面站点和卫星；

因代表性评价指标是随时间变化的函数，将地表温度分解为三个组分：

其中，ATC表示地表温度的季节性变化成分，DTC表示地表温度在日内的变化情况，USC表示地表温度受局地瞬时天气因素等情况影响下的波动成分；

忽略天与天之间同一卫星的过境时间差异，得到：

代入ΔT(t)可得：

ΔT(t)＝[ATC_F(t)+DTC_F(t)+USC_F(t)]-[ATC_P(t)+DTC_P(t)+USC_P(t)]

≈[ATC_F(t)+USC_F(t)]-[ATC_P(t)+USC_P(t)]

≈[ATC_F(t)-ATC_P(t)]+[USC_F(t)-USC_P(t)]

≈ΔATC(t)+ΔUSC(t)

其中，ΔATC和ΔUSC分别为两个尺度地表温度季节性变化成分的差值和天气因素影响成分差值；

S2、参数化地表温度的季节性变化成分：

其中，

和A₀分别是季节性变化成分的均值和振幅，θ₀是该曲线的相位；

对ΔATC，可得：

其中，

是两个尺度ATC中

的差值，A和θ是该曲线的振幅和相位，与两个尺度ATC中的振幅和相位关系为：

其中，k为常数，取值为0，±1；

S3、引入表征地表属性的影响因子和局地气象因子x_i，i＝1,2,…,n，n为影响因子的个数，建立ΔUSC的时序化方程为：

ΔUSC(t)＝g[x₁(t),x₂(t),……,x_n(t)]

其中，g是表征ΔUSC与影响因子x_i之间的参数化方程，g受地理位置及参数的可获取程度及变量之间的相关性程度影响，通过多源线性回归或机器学习方法获取；

S4、根据上述步骤，得到时序化的代表性指标ΔT：

根据时序化的代表性评价指标，实现站点空间代表性动态评价。

Images