结合大气水汽的Landsat TIRS分裂窗口地表温度产品不确定性工作流程

《Remote Sensing of Environment》：An uncertainty workflow for landsat TIRS split window surface temperature products incorporating atmospheric water vapor

　　

{"text":"•Landsat Collection 3 will potentially use the Split Window algorithm for Level-2 ST.•Current Landsat Level-2 Single Channel uncertainty product is prone to false positives.•New uncertainty workflow for Split Window integrates total precipitable water.•Machine learning accurately estimates TPW from Landsat TIRS bands 10 and 11.•Results show improved reliability of surface temperature retrieval uncertainty.IntroductionThe launch of Landsat 9 (L9) in 2021 marked another chapter in the Landsat program, extending its already impressive archive. The sensors onboard L9 are functional copies of those onboard Landsat 8 (L8), except for baffles that were integrated into the Thermal Infrared Sensor (TIRS)-2 to mitigate any effects of stray light, observed in the TIRS-1 imagery data (Montanaro et al., 2014). As such, the algorithms developed for L8 image data can be applied to L9 with minor modifications to account for the differences in the instruments’ spectral responses.Among the available Landsat products in Collection 2 is the Level 2 Surface Temperature (ST) product, which currently employs a Single Channel (SC) algorithm to derive surface temperature (Laraby and Schott, 2018; Laraby, 2017). However, the SC approach is susceptible to false positives in the associated uncertainty maps (Laraby and Schott, 2018). On the other hand, the split window (SW) approach, widely used in ST retrievals, offers improved accuracy by leveraging the differential absorption characteristics between two thermal bands (McMillin, 1975; Galve et al., 2008). Given its robustness, starting with Collection 3, the SW approach will be adopted as the ST algorithm for the level-2 product which will require the development of a dedicated uncertainty workflow. Proper error propagation through the SW model’s various terms enables quantification of the SW algorithm’s uncertainty. To further improve the reliability of the SW ST uncertainty map, compared to its predecessor the SC uncertainty approach, an enhanced uncertainty estimation framework can be implemented—integrating atmospheric water vapor information that has been shown to improve the accuracy of the SW algorithm’s uncertainty (Galve et al., 2008; Jiménez-Mu?oz and Sobrino, 2008; Rhziel et al., 2023).In this study, our primary objective is to explore the uncertainty components of the SW algorithm in greater depth, analyze the contributing terms, and incorporate atmospheric water vapor to improve the accuracy of the SW ST uncertainty estimation workflow. We also build upon the previously introduced SW approach (Gerace et al., 2020) for L8, by reporting updated coefficients for L9.Main contributions. We summarize our research contributions as follows:•Tailored SW approach for L9. We adapt and tailor the SW methodology for use with the TIRS-2 sensor onboard L9. We derive the SW coefficients for L9 using Moderate Resolution Atmospheric Transmission (MODTRAN Berk et al., 2014) modeling of different atmospheric profiles using the Thermodynamic Initial Guess Retrieval (TIGR Laboratoire de Météorologie Dynamique (LMD), 2019) database.•Enhanced SW ST uncertainty workflow. We introduce an SW ST uncertainty workflow for L9 that accounts for various sources of error such as algorithmic error, sensor noise, and surface emissivity—through error propagation. Using the TIGR atmospheric profile database and MODTRAN simulations across varying water vapor conditions, we derive a relationship between total precipitable water (TPW) and SW algorithmic uncertainty. This relationship allows TPW estimates to be used as inputs for predicting per-pixel SW algorithmic uncertainty and to contribute to the total SW ST uncertainty, which in turn motivates our next contribution on real-time TPW estimation.•Machine learning–based TPW estimation. To obtain the per-pixel TPW estimates required by the proposed uncertainty framework, we develop a boosting-based machine learning (ML) model using XGBoost (Chen and Guestrin, 2016). The model predicts TPW from L9 TIRS-2 thermal infrared bands. We construct the training data using coincident L8/9 observations as inputs and two independent TPW reference sources: (1) the Moderate Resolution Imaging Spectroradiometer (MODIS) and (2) the Aerosol Robotic Network (AERONET). We further validate model performance using independent ground-based TPW measurements from Atmospheric Radiation Measurement (ARM) microwave radiometers (MWR) to ensure generalization.•Robust TPW data processing and feature optimization pipeline. We construct a robust modeling pipeline for TPW estimation, starting with the collection of coincident Landsat TIRS and MODIS/AERONET TPW data, followed by comprehensive masking and filtering to ensure high-quality input/output pairs. We propose a pipeline that integrates mathematical feature engineering and feature selection using Jostar (Hassanzadeh et al., 2021). The feature engineering step expands the predictor space by generating physically and statistically informative linear and nonlinear transformations of the original variables, while feature selection identifies a subset that retains the most predictive information and removes redundancy. This improves model generalization by addressing bias-variance trade off.L8 and L9 data are available in two collection tiers—Collection 1 and Collection 2—each offering two primary product levels: Level 1 and Level 2. Collection 2 features improved radiometric and geometric calibration and updated processing algorithms compared to Collection 1. Level 1 products are radiometrically and geometrically corrected, while Level 2 products undergo atmospheric correction. Among the Level 2 products are surface reflectance (SR), ST, and vegetation indices (VIs) (Wulder et al., 2022; U.S. Geological Survey, 2024). The ST product, specifically, plays an important role in water resource management and crop mapping through evapotranspiration, disaster management, and ecological observations (Anderson et al., 2012; Deliry et al., 2020; Guha et al., 2018).The current Level 2, Collection 2 ST product relies on a SC algorithm that uses a radiative transfer model to solve the governing equation for a broad-band infrared channel (Band 10) for land surface temperature (LST) estimation. The corresponding SC uncertainty workflow primarily uses standard error propagation (SEP) to estimate per-pixel uncertainty. However, due to a perceived underestimation of the uncertainties based on empirical observations (Laraby and Schott, 2018), additional per-pixel uncertainty is introduced using an unknown error matrix: error accounting for uncertainties not explained by SEP. Essentially, additional uncertainty is introduced to the per-pixel uncertainty maps based on each pixel’s distance-to-cloud (DTC). Notionally, it can be appropriate to use a model to capture uncertainties not characterized by SEP. However, the DTC map is derived from the cloud mask generated by C Function of Mask (CFMask; Foga et al., 2017), which is approximately 95% accurate. CFMask is a decision-tree algorithm for cloud, cloud-shadow, snow, and water detection. CFMask generates final classification layers that are stored in the Quality Assurance (QA) band. CFMask is particularly susceptible to false positives, where highly reflective surfaces with cold brightness temperatures, such as sand or aluminum, may be misclassified as clouds. These false-positives are most evident on cloud-free days when conditions are most favorable for accurate surface temperature retrieval.Fig. 1 shows an example of this phenomenon for an arbitrary cloud-free scene. From left-to-right, the RGB and ST reference images are shown to illustrate the cloud-free nature of the scene. The DTC product is shown in (c) to show that several pixels are falsely identified as clouds. Finally, (d) shows the corresponding uncertainty map to illustrate how a significant fraction of pixels becomes contaminated using the current workflow.With the transition to Landsat Collection 3, the SW algorithm discussed in this paper will potentially be used operationally for ST retrieval from L8 and L9. The SW algorithm is capable of accurately estimating ST using two long-wave infrared bands and has been widely applied across various sensors and domains, including sea and land surface temperature retrieval (Becker and Li, 1990; Wan and Dozier, 1996; Coll et al., 1994; Ulivieri et al., 1994). SW was implemented during early satellite missions such as Advanced Very High Resolution Radiometer (AVHRR) (Llewellyn-Jones et al., 1984) and MODIS (Brown et al., 1999). SW has been validated in multiple studies using data from the Surface Radiation Budget Network (SURFRAD) sites and buoy-based ST measurements across the United States (Duan et al., 2021; Malakar et al., 2018a). For L8 TIRS-class sensors, Gerace et al. (2020) proposed an architecture that builds on prior work to enable an operational ST product (Gerace et al., 2020). More recently, Eon et al. (2023) conducted an under-flight experiment as part of a cross-validation study between L8 SW-derived ST and unmanned aerial system (UAS)-based thermal imagery (Eon et al., 2023). Their findings reported maximum temperature differences of 1 K over water and 2 K over sand. This highly-utilized algorithm has been adapted and modified for a variety of sensors and can be trained for specific applications (Gerace et al., 2020; Mao et al., 2005; Kerr et al., 1992; Jiang and Li, 2008).In this study, we apply the generalized SW algorithm originally proposed by Wan and Dozier (1996) to derive ST from imagery acquired by TIRS sensors, as shown in Eq. (1), where,–ST is the desired surface temperature [K]–bk (for k = 0,1,…,7) are sensor and algorithmic coefficients derived through a training process using MODTRAN–i and j correspond to the two thermal bands (Bands 10 and 11 for TIRS)–Δ?=?i??j is the difference in band-effective emissivities–?=(?i+?j)/2 is the average of the band-effective emissivities–Ti,Tj are the brightness temperatures in the two thermal bandsST=b0+(b1+b21???+b3Δ??2)Ti+Tj2+(b4+b51???+b6Δ??2)Ti?Tj2+b7(Ti?Tj)2,Similar to its SC counterpart, the proposed SW approach warrants an uncertainty workflow tailored to its formulation. A substantial body of prior work has addressed uncertainty budget construction for surface temperature retrieval algorithms. Mittaz et al. (2019) applied metrological uncertainty analysis to satellite Earth observations, illustrating how principles of measurement uncertainty can be applied traceably from sensor to geophysical product (Mittaz et al., 2019). Other initiatives, including the Climate Change Initiative program, have emphasized rigorous uncertainty quantification for essential climate variables; for example, independent uncertainty estimation methodologies have been developed for sea ST retrievals from Along Track Scanning Radiometer (ATSR) instruments (Bulgin et al., 2016). In the context of LST, Ghent et al. (2019) introduced a comprehensive approach to defining an uncertainty budget for MODIS LST products that distinguishes different uncertainty components and quantifies them consistently (Ghent et al., 2019). In this study, we follow the principles outlined in the Guide to the Expression of Uncertainty in Measurement (GUM; ISO and OIML, 1993) to construct an uncertainty budget tailored to SW LST retrieval.To clarify the role of the proposed GUM-based uncertainty workflow, it is useful to distinguish it from two commonly used uncertainty-assessment strategies for LST retrievals. The first is validation-based uncertainty quantification, in which satellite-derived temperatures are compared against high-accuracy ground networks such as SURFRAD, AERONET or buoy estimates. While this approach provides useful estimates of bias and error, it only yields a scalar uncertainty estimate rather than per-pixel uncertainty, as uncertainty is assessed only at validation sites (Wan, 2008). The second approach relies on radiative-transfer-based uncertainty estimation, where models such as MODTRAN are used to perturb atmospheric and surface parameters. In this framework, LST uncertainty is commonly parameterized as a function of TPW, emissivity contrast, and viewing geometry, with forward-simulated radiances used to evaluate their combined influence on LST error. This approach is also effective, but it requires operational MODTRAN simulations, which can be computationally limiting.Literature supports the adoption of different variations of the SW algorithm for ST retrieval, along with uncertainty characterization across different satellite missions. Galve et al. (2008) proposed a modified SW algorithm for MODIS, incorporating surface emissivity, temperature, and atmospheric water vapor as independent variables, with emissivity and temperature terms being decoupled (Galve et al., 2008). Their results showed that the error between the estimated ST and the ground reference increased with higher atmospheric water vapor content (see Fig. 5 in Galve et al., 2008). They conducted a detailed error propagation analysis accounting for uncertainties in brightness temperature, emissivity, water vapor, and algorithmic fitting coefficients, with total error estimates around 0.4 K. For real-time TPW, the study utilized MOD07 data to obtain total precipitable water estimates to support operational ST retrievals.Motivated by these findings, studies have examined the role of atmospheric water vapor in the SW algorithm, demonstrating that explicitly accounting for atmospheric water vapor improves the accuracy of the SW algorithm and reduces the associated uncertainty (Wan and Dozier, 1996; Ottle et al., 1997; Harris and Mason, 1992; Ulivieri et al., 1994; Yu and Barton, 1994; Fran?ois and Ottlé, 2002; Galve et al., 2008; Jiménez-Mu?oz and Sobrino, 2008; Rhziel et al., 2023). Accounting for water vapor often requires estimating TPW from available spectral bands or relying on external products that provide independent TPW measurements."}

在总降水水汽（TPW）估计中，已有多种方法利用亮温的协方差与方差之比来估计整层大气水汽，从而提升分裂窗（SW）算法的性能（Ottle et al., 1997; Sobrino et al., 1994）。Ottle等人（1997）报告称，使用协方差-方差法得到的TPW估计值的均方根误差（RMSE）小于0.5 g/cm2，该结果已通过AATSR和AVHRR数据进行验证（Ottle et al., 1997）。类似地，Liu等人（2017）提出了一种广义SW算法的改进方案，通过引入一个中心波长位于6.7 μm的第三波段，专门用于MODIS的水汽估计（Liu et al., 2017）。

鉴于大气水汽知识对支持SW算法的重要性已得到充分证明，对于TIRS类双红外传感器而言，使用热红外观测进行TPW反演似乎初始上是次优的，因为最强的水汽吸收特征出现在近红外区域（940–960 nm）。然而，基于热红外的TPW反演也已在文献中被研究和实施。可见光红外成像辐射仪套件（VIIRS）的TPW反演算法采用晴空回归框架，将红外亮温（8.6、10.8和12.0 μm）与从探空数据（如美国国家海洋和大气管理局唯一组合大气处理系统NUCAPS）获得的水汽廓线建立关系，并使用具有全球代表性的训练数据集（GAO and Li, 2020）。类似地，MODIS的MOD07产品利用经过探空仪获取的大气廓线训练的回归算法，从长波红外辐射中估计TPW和大气水汽廓线（Seemann et al., 2003）。这些方法共同确立了在反演框架中使用热红外信息估计TPW的先例，尤其是在近红外观测不可用、受污染或不同位的业务化场景中。考虑到SW算法可能会在L8和L9上的TIRS传感器的Collection 3中被采用，而L8和L9缺乏直接反演TPW的专用波段，可以利用TIRS波段中的信息来估计TPW，并将其纳入SW地表温度（ST）不确定性工作流中。

在卫星TPW反演验证的各种方法中，AERONET观测、全球导航卫星系统（GNSS）和全球定位系统（GPS）测量以及地基微波辐射仪（MWR）是最常用的方法（Van Baelen et al., 2005）。由于具有连续的全天候观测能力，MWR是估计TPW最准确的方法之一（Liljegren et al., 2001; Turner et al., 2007）。MWR已被广泛用于验证MODIS等卫星TPW产品（Grody et al., 2001）。ARM网络提供了长期的、全球分布的MWR观测，为卫星TPW反演的评估和验证提供了可靠的基准（Revercomb et al., 2003）。

在本研究中，我们使用标准误差传播开发了一种用于SW算法的不确定性工作流。SW算法在估计地表温度方面尤其具有优势，因为它不需要大气知识，即在训练过程中进行的回归消除了这一需求。然而，在开发不确定性工作流时，这带来了一个独特的挑战，因为大气是温度估计不确定性的主要贡献者。因此，我们使用模拟数据开发了一个模型，该模型明确地重新引入了大气水汽的影响，以支持逐像元的适当不确定性估计，即此处开发的不确定性工作流由两个部分组成：标准误差传播（SEP）和一个由模拟数据驱动、以TPW为函数估计不确定性的模型（TPW模型），见公式(2)。该工作流的详细信息将在第2.2节中讨论。

σST=(σSEP)2+(σTPW)2

将TPW纳入SW地表温度不确定性工作流需要两个步骤：(1) 使用MODTRAN建立TPW与SW算法不确定性之间的关系；(2) 在运行时运行性地估计TPW并将其映射到此预建立的关系中以量化相关的LST不确定性贡献。第一步通过将SW地表温度估计中的残余误差与TPW相关联来考虑残余误差，TPW作为主导贡献者，其影响可以被显式参数化。在运行时，我们使用有监督机器学习（ML）模型估计TPW，并利用建立的关系将估计的TPW值映射到SW算法不确定性以推断TPW引起的SW算法不确定性。我们认为这种方法为将TPW效应纳入SW地表温度不确定性工作流提供了一种更具物理意义且更稳健的方法。

关于运行时的运行性TPW估计，文献中的大多数方法依赖于外部来源，如数值天气预报、类似MOD07的卫星反演或辐射传输模型（RTM）推导的水汽廓线。因此，与TPW相关的不确定性来源于地表温度反演之外。虽然这些方法是可行的，但它们依赖于摄入来自其他传感器或模型的辅助产品，而这些产品可能并非始终可用或时间一致。为了避免这种依赖并保持自包含框架，我们提出使用ML模型在两个不同数据集上训练，在推理时直接从Landsat观测中估计TPW，从而实现自洽的逐像元TPW估计。

我们提出两种不同的参考数据集用于训练ML模型：用于逐像元训练的MODIS MOD05L2 TPW产品和用于站点级训练的AERONET日照计地面参考站点。如前所述，L8/9不配备直接测量TPW所需的光学波段。因此，使用基于比值的方法进行TPW反演（类似于MODIS的双通道和三通道比值方法（Gao and Kaufman, 2003））是不可行的。对于L8和L9准确反演TPW可能需要依赖专门设计用于探测水汽的仪器所提供的外部产品。MODIS以1 km地面采样距离（GSD）提供MOD05 TPW产品。在2022年10月星座离轨机动之前（MODIS Characterization Support Team, 2022），MODIS和L9的过境时间几乎重合，最大时间偏移为30分钟。这种紧密的时间对齐为应用有监督ML方法提供了机会，这些方法能够检测数据中的复杂非线性模式并从相应的逐像元Landsat数据中估计TPW。与此同时，AERONET网络提供了TPW估计的长期地面参考，全球约有600个站点（Holben et al., 1998）。AERONET辐射仪采用日照计测量近红外水汽吸收波段内的直射太阳辐照度，并使用经验校准的吸收模型反演TPW。AERONET的TPW估计已被广泛用于卫星水汽数据集的定标和验证，包括MODIS TPW产品。

基于这些互补的数据源，我们使用Landsat TIRS观测（第10和第11波段）以像元级别反演TPW，其中MODIS TPW或AERONET TPW作为训练参考数据。然后我们从训练好的模型生成空间连续的TPW地图，并将这些估计值映射到SW算法不确定性框架中。这种方法使我们能够在SW衍生LST产品旁边生成逐像元的SW地表温度不确定性栅格。本研究中使用的缩写和缩略语完整列表见附录表A.1。

各章节摘要

方法

在本节中，我们概述了用于推动不确定性讨论的特定分裂窗（SW）算法。此外，我们报告了SW系数，这将使用户能够将该算法应用于Landsat 9（L9）上的热红外传感器（TIRS）-2（Gerace et al., 2020）。然后我们详细考察了两组分SW不确定性工作流（见公式(2)），重点介绍使用标准误差传播（SEP）推导的项以及那些项

结果与讨论

在本节中，我们首先展示基于中分辨率成像光谱仪（MODIS）和气溶胶机器人网络（AERONET）数据的总降水水汽（TPW）建模结果，并使用微波辐射仪（MWR）的TPW测量值验证所提出的模型。然后我们展示SW不确定性地图的结果，并将其与单通道（SC）不确定性地图在各种场景中进行了定性比较。

结论

用于地表温度（LST）估计的分裂窗（SW）算法及其各种衍生方法在文献中已被研究了三十多年（Becker and Li, 1990; Wan and Dozier, 1996; Coll et al., 1994; Ulivieri et al., 1994）。随着配备类似双长波红外热波段的Landsat 8（L8）及其后继者Landsat 9（L9）的发射，为实施用于地表温度（ST）产品的SW算法提供了机会

CRediT作者贡献声明

Amirhossein Hassanzadeh：原创手稿撰写、可视化、软件、方法论、调查、形式分析。Robert Mancini：软件、方法论、调查。Aaron Gerace：原创手稿撰写、可视化、验证、监督、项目管理、资金支持获取、形式分析、概念化。Rehman Eon：原创手稿撰写、可视化、验证、软件、方法论、形式分析、概念化。Matthew Montanaro：审阅与

利益冲突声明

作者声明以下可能被视为潜在利益冲突的财务利益/个人关系：

Amirhossein Hassanzadeh报告称获得了美国地质调查局提供的财务支持。如果有其他作者，他们声明没有已知的可能影响本文所述工作的竞争财务利益或个人关系。

致谢

本文基于美国地质调查局合作协议编号G19AC00176支持下进行的工作。本文包含的观点和结论是作者的，不应被解释为代表美国地质调查局的观点或政策。提及商品名称或商业产品不构成美国地质调查局的认可。本稿件提交发表的前提是，美国Amirhossein Hassanzadeh|Robert Mancini|Aaron Gerace|Rehman Eon|Matthew Montanaro|Elphas Khata