Beyond Land Surface Temperature: Explainable Spatial Machine Learning Reveals Urban Morphology Effects on Human-Centric Heat Stress Yuan Wanga,b,∗, Shengao Yid , Xiaojiang Lid , Pengyuan Liue , Zhiwei Yangf , Ronita Bardhanc , Rudi Stouffsa a Department of Architecture, National University of Singapore, Singapore 117566, Singapore b Cambridge Centre for Advanced Research and Education in Singapore (CARES), Singapore 138602, Singapore c Sustainable Design Group, Department of Architecture, University of Cambridge, Cambridge, United Kingdom d Department of City and Regional Planning, University of Pennsylvania, Philadelphia, PA 19104, USA e Urban Analytics Subject Group, Urban Studies & Social Policy Division, University of Glasgow
arXiv:2604.22433v1 [cs.LG] 24 Apr 2026
f Laboratory for Earth Surface Processes, Ministry of Education, College of Urban and Environmental Sciences,
Peking University, Beijing 100871, China
Abstract Heat exposure connects the built environment and public health, directly shaping the livability and sustainability of urban areas. Fully understanding the spatial heterogeneity of heat exposure and its driving factors is therefore a vital prerequisite for climate-adaptive urban planning. However, most planning-oriented studies rely on land surface temperature (LST), and whether LST adequately represents human heat exposure and how it differs from physiologically relevant heat stress remains insufficiently examined. Here, adopting Landsat-retrieved 30-m LST and GPU-accelerated 1-m universal thermal climate index (UTCI) in Singapore—a high-density tropical coastal city characterized by strong solar radiation—this study establishes a comprehensive “Modeling-Comparing-Assessing" framework to systematically evaluate the spatial and mechanistic discrepancies between the two metrics. We further investigate pronounced non-stationary and threshold-based quantitative relationships of the two metrics with urban factors by employing a novel geographically weighted XGBoost (GW-XGBoost) and generalized additive model (GAM) workflow. Our results demonstrate notable discrepancies in spatial patterns of LST and UTCI, along with substantial spatial heterogeneity in how 2D and 3D urban factors impact these two thermal metrics, as revealed by explainable GW-XGBoost models (global out-of-bag R2 = 0.855 for LST and 0.905 for UTCI, respectively). Crucially, spatially explicit SHAP interprets that sky view factor plays a central role in explaining UTCI variability but exhibits a comparatively marginal independent contribution to LST, indicating that LST inadequately captures shading-driven and radiative processes governing actual human heat stress. Notably, SHAP-GAM analysis indicates that higher albedo is associated with increased UTCI. These novel findings provide evidence for integrating physiologically relevant thermal indices to inform targeted heat risk management and climate-adaptive urban planning. Keywords: UTCI, Urban morphology, GeoAI, Spatial Heterogeneity, Climate-adaptive planning
∗ Correspondence to: Yuan Wang ([email protected])
1. Introduction Climate-informed urban planning increasingly depends on understanding how urban morphology shapes the local thermal environment. Previous research has established strong associations between urban form and thermal conditions through parameters such as sky view factor (Oke, 2002; Eliasson, 1996; Zakšek et al., 2011; Robinson, 2006), building density (Yao et al., 2023; Shi et al., 2019), and land cover (Chen et al., 2024; Zhang et al., 2024a; Li and Zeng, 2024; Chang et al., 2023). However, this body of work overwhelmingly uses land surface temperature (LST) derived from satellite observations as the dependent variable, owing to the fine spatial resolution and long-term temporal coverage of LST products. Such a reliance on LST constitutes a fundamental methodological limitation rather than a neutral proxy choice. Although LST offers long-term, spatially continuous coverage and computational convenience (Wang et al., 2024b), it measures radiometric surface heating of roofs and treetops and therefore fails to capture the radiative exchange, convective cooling, humidity effects, and wind modulation that jointly determine human heat stress within street canyons (Zhan et al., 2025). These omissions are particularly consequential in high-density urban environments, where shading geometry, enclosure effects, and vegetation modify the radiative and convective environment in ways that cannot be inferred from surface temperatures alone (Briegel et al., 2025; Zhan et al., 2025). More fundamentally, LST captures only the heat component and excludes the humidity, wind, and multidirectional radiation that determine physiologically relevant heat stress. Existing heatexposure models that use LST as a proxy may achieve high predictive performance while remaining poorly aligned with human physiological processes (Nazarian et al., 2022; Tuholske et al., 2021). As a result, LST-based assessments routinely mischaracterize human-centered heat stress, which reflects physiological heat strain, and can lead to biased energy demand estimates and inappropriate adaptation strategies. This persistent mismatch between surface-based temperature indicators and human-centered thermal conditions exposes a critical gap in current urban morphology–climate frameworks, which lack scalable approaches capable of representing human-centered thermal stress at the city scale (Nazarian et al., 2022). Addressing this gap requires moving beyond surface-based temperature indicators towards models that explicitly resolve the spatial complexity of the heat stress and systematically identify where—and why—LST fails as a decision-relevant metric (Zhan et al., 2025; Muse et al., 2024; Nazarian et al., 2022). Human heat stress proxies such as the Universal Thermal Climate Index (UTCI) address this gap by integrating air temperature (Ta ), humidity, the mean radiant temperature (Tmrt ), and wind speed (Bröde et al., 2012; Fiala et al., 2012; Jendritzky et al., 2012; Blazejczyk et al., 2012). UTCI is particularly valuable because of its thermophysiological basis, international standardization, and 2
stable sensitivity across climates. Its accuracy depends on the estimation of Tmrt , which captures the full shortwave and longwave radiative environment surrounding the human body. To compute this at scale, physically based models such as ENVI met (Roth and Lim, 2017), RayMan (Matzarakis et al., 2007), and SOLWEIG (SOlar and LongWave Environmental Irradiance Geometry (Lindberg et al., 2008; Lindberg and Grimmond, 2011; Lindberg et al., 2018)) have been developed. Among these, SOLWEIG offers a practical balance of physical realism and computational efficiency by explicitly representing six directional radiation components. Recent advances in GPU acceleration now make it possible to apply SOLWEIG citywide at fine spatial resolution (Li and Wang, 2021; Yi et al., 2025a). Although simulation models provide physically realistic estimates, understanding how urban morphology drives spatial variations in LST and UTCI requires complementary statistical modeling. Crucially, the necessity for advanced modeling arises from the fundamental physical properties of urban microclimates. Theoretically, the relationships between urban form and thermal environment are inherently non-linear and spatially non-stationary (Wan et al., 2025; Zhou et al., 2025; Zhao et al., 2026). However, existing studies mainly rely on linear regression (Zhang and Yuan, 2023), spatial regression (Wongsai et al., 2024; Assaf and Assaad, 2024), or machine learning (ML) models (Zhang et al., 2024b; Yang et al., 2024; Chen et al., 2023) to examine these relationships. Linear models cannot represent nonlinear interactions and threshold effects, and both linear and conventional ML models typically ignore spatial heterogeneity. Spatial regression models such as geographically weighted regression (GWR) account for local variation but are limited by their linear structure and susceptibility to multicollinearity. Recent developments in geographically weighted machine learning address these limitations by combining nonlinear modeling with spatial weighting (Fouedjio and Arya, 2024; Yang et al., 2023; Wang et al., 2024a; Yi et al., 2024). The geographically weighted XGBoost (GW-XGBoost) model is particularly up-and-coming because it captures nonlinear interactions, accommodates high dimensional datasets, and reveals local variable importance while maintaining interpretability through SHapley Additive exPlanations (SHAP) analysis. Characterizing the spatial heterogeneity of these geographical patterns is essential for mapping disparities in heat exposure impacts, identifying hotspot areas, and informing the design of targeted local and national interventions to mitigate current and potential future risks. Beyond methodological considerations, a more fundamental question remains insufficiently addressed: whether commonly used surface-based temperature metrics are sufficient to inform humancentric heat mitigation and urban planning decisions, which prioritize human-related parameters over purely physical indicators (Yang et al., 2025). LST primarily reflects surface thermal properties rather than the atmospheric and radiative conditions experienced by humans. As a result, its relevance
3
for assessing heat stress, particularly in dense and vertically complex urban environments, remains uncertain. In contrast, UTCI is explicitly designed to represent human thermal stress. However, despite its growing use, it remains unclear whether UTCI is systematically more sensitive to urban spatial structure and therefore more informative for place-based, human-centered decision-making. Addressing this gap is critical for evaluating the appropriateness and limitations of temperature metrics used in urban heat assessment. This study aims to address the fundamental gaps between surface-based and human-centered thermal indicators by applying a GW-XGBoost framework to jointly analyze LST and UTCI. This approach accounts for spatial heterogeneity and correlation, enabling a detailed comparison of where and why LST diverges from UTCI in high-density tropical coastal cities under strong solar radiation, using Singapore as a case study. The study is guided by three research questions: (1) What and where are the differences in LST and UTCI? (2) How do complex 2D and 3D urban morphological features drive these discrepancies through non-linear and spatially non-stationary mechanisms? (3) Which specific urban factors are systematically misrepresented when relying solely on LST for climate-adaptive planning?
2. Methodology This study employs a high-resolution comparative framework in Singapore, integrating Landsatretrieved LST and GPU-accelerated 1-m UTCI with multi-source 2D and 3D urban covariates. By implementing GW-XGBoost and SHAP-based interpretation, we quantify the spatial non-stationarity and divergent driving mechanisms of these two thermal metrics. 2.1. Study area Singapore, a high-density city-state in Southeast Asia (Figure S1), provides an ideal natural laboratory for examining how urban morphology shapes the spatial heterogeneity of LST and UTCI. Despite its small land area, the city exhibits pronounced morphological contrasts arising from intense vertical development, diverse land-use configurations, and fine-grained variations in building height, density, and vegetation structure (Figure 1). These characteristics are embedded within a tropical equatorial climate with persistently high temperatures and humidity and minimal seasonal variability, allowing observed thermal differences to be primarily attributed to urban form. Singapore’s urban fabric comprises a heterogeneous mix of high-rise residential estates, dense commercial cores, industrial zones, and an extensive, strategically distributed network of green infrastructure, including parks, urban forests, roadside vegetation, vertical greenery, and rooftop
4
gardens. This combination of strong 3D morphological variability and climate stability makes Singapore particularly well suited for isolating and quantifying the differential influences of 2D and 3D urban form on LST and UTCI. Moreover, the availability of high-resolution urban morphology and meteorological datasets enables the application of interpretable geographically weighted machine learning models to capture spatially varying morphology–temperature relationships.
Figure 1: Study area with 1-m land use/land cover (LULC), building height, and canopy height layers, and the workflow for Tmrt and UTCI computation using SOLWEIG and the human physiological model integrated in GPU-accelerated pipelines. The area marked in red is shown as an illustrative inset to highlight high-resolution spatial details that are imperceptible at the full island scale.
2.2. Variables and Data Processing Two target temperature variables LST and UTCI, are shown in Table 1. LST characterizes the radiometric skin temperature of the land surface, and UTCI integrates Ta , Tmrt , humidity, and wind speed to reflect human heat stress. May 2019 was selected as the study period because May is typically among the hottest months of the year in Singapore, and 2019 ranks as one of the warmest years on record (Meteorological Service Singapore, 2023, 2020). The hottest climatic conditions enhance the contrast in LST and UTCI, providing an ideal context for examining urban heat patterns and spatial variability in human heat exposure across the city.
5
Table 1: Summary of variables, formulae, unit, descriptions, and references of multi-source datasets. Category LST
Variables Landsat LST
Formula Described in the main text
Unit °C
UTCI
UTCI
Described in the main text
°C
SV F
π(2αi −1) 360 Pn 2 π SVF = i=1 S 2π sin 180 sin 2n θi )
3D urban morphology
BH BHsd
m m
m m
The standard deviation of canopy heights within the spatial unit.
– patches/100 ha
Canopy density (CD) is the proportion of land area covered by canopy. Patch Density, the density of patches of each landcover class. Contagion index, measuring patch aggregation based on adjacency probabilities. Landscape Shape Index measures shape complexity compared to a maximally compact shape.
–
F AR
–
CD PD
LSI COHESION
6
PAFRAC CONTIG_AM SHDI
−
Pm
SHEI
−
Pm
DEM
i=1 pi · ln(pi ) i=1 pi ln pi ln(m)
daytime, 11:00 a.m. SGT, hourly, 1-m Annulus-weighted SVF, ratio of visible sky hemisphere, indicating proportion of hemispheric radiation accessible to a planar surface. Average building height in the spatial unit. The standard deviation of building heights within the spatial unit. Building density (BD) is the proportion of land area covered by building footprints. Floor Area Ratio (FAR) is the ratio of the total floor area of all buildings to the total land area. Average tree canopy height in the spatial unit.
total foor area total land area Pn i=1 CHi q n P N 1 2 i=1 (CHi − CH) N −1 total canopy area total land area ni × 10, 000 i hA P p ln(p ) m Pm × 100 1 + i=1 j=1 ij2 ln mij Pm ∗ eik 0.25 k=1 √ A h Pn i−1 p∗ ij √ × 1 − √1Z 1 − Pn j=1 × 100 ∗ a∗ j=1 pij ij Pn Pn Pn n i=1 ln(pi ) ln(ai ) − ( i=1 ln(pi )) ( i=1 ln(ai )) Pn Pn 2 n i=1 [ln(ai )]2 − ( i=1 ln(ai )) Pn T IGij ·aij j=1 CON Pn j=1 aij
CONTAG
% – % – – –
Reference Wang et al. (2020) Lindberg et al. (2008); Bröde et al. (2012) Lindberg and Grimmond (2010)
Cai et al. (2025)
McGarigal (2015)
Patch Cohesion Index, cohesion of each landcover class. Perimeter-Area Fractal Dimension quantifies the complexity of patch shapes by examining how patch perimeter scales with patch area. Area-Weighted Mean Contiguity Index, measures the spatial connectedness of cells within each patch of a landcover class. Shannon Diversity Index, measures the diversity of landcover types based on their proportional abundance.
–
Evenness of class distribution (0 = dominated, 1 = equal).
– 0.2453 ρblue + 0.0508 ρgreen + 0.1804 ρred + 0.3081 ρNIR + 0.1332 ρSWIR1 + 0.0521 ρSWIR2 + 0.0011
m
Digital Elevation Model (DEM) represents the terrain elevation.
ESA (2020)
–
Ratio of reflected solar radiation to incident solar radiation on a surface.
Cunha et al. (2020)
NDVI
ρNIR −ρred ρNIR +ρred
–
NDBI
ρSWIR −ρNIR ρSWIR +ρNIR
–
Albedo
WET Dsea Socio-economic factors
N −1
2 i=1 (BHi − BH)
BD
CHsd
Environmental variables
i=1 BHi
q n P N 1
–
building footprint area total land area
CH
*2D landscape indices
Pn
Description daytime ∼11:00 a.m. SGT, 16-day repeat cycle, 30-m
0.1511ρblue + 0.1972ρgreen + 0.3283ρred + 0.3407ρNIR − 0.7117ρSWIR1 − 0.4559ρSWIR2 p D(x, y) = min (x − xw )2 + (y − yw )2
– m
PopD RD IntD
(xw ,yw )∈W total population total land area total road length total land area number of intersections total land area
persons/km2 km/km2 intersections/km2
RNC
number of intersections total road length
intersections/km
Ratio of the difference to the sum of near-infrared and red reflectance, used to indicate vegetation vigor and density. The normalized differential build-up index (NDBI) serves as a proxy for the intensity of construction pressure. WET is derived from the Tasseled Cap Transformation, reflecting land surface moisture conditions, particularly soil moisture. Euclidean distance from each pixel to the nearest coastline.
Baig et al. (2014)
Population density Road density Intersection density Road network compactness (RNC) measures the density of interconnections within the road network.
*Note: aij is the area of patch j of class i; A is the total landscape area; ni is the total number of patches of class i; e∗ik is the length of edge between classes i and k; p∗ij is the perimeter of patch j of class i; a∗ij is the area of patch j of class i; Z is the total number of cells in the landscape; cij is the contiguity value of patch j of class i; pi is the proportional abundance of class i; m is the number of classes (patch types).
2.2.1. Land surface temperature retrieval Daytime LST is retrieved from Landsat 8 Thermal Infrared Sensor (TIRS) Band 10 using the Statistical Mono-Window algorithm. Prior to LST retrieval, Landsat images are processed using the C Function of Mask (CFMASK) algorithm to remove cloud-contaminated pixels (Foga et al., 2017). LST is then estimated using the Statistical Mono-Window Model following Ermida et al. (2020) and Wang et al. (2020). We retrieved LST from Google Earth Engine (GEE), where all required input datasets were assembled and processed. Landsat 8 flies in a sun-synchronous, near-polar orbit with a descendingnode equator crossing at around 10:00 a.m. (±15 min) local solar time (U.S. Geological Survey, 2025). For Singapore (UTC+8; ∼ 103.8◦ E), where civil time runs about +1 h 05 min ahead of local solar time, this corresponds to a typical overpass near 11:15 a.m. Singapore time (SGT). We then aggregated the remaining clear-sky pixels to generate a representative monthly mean LST composite for May from 2018 to 2020 at approximately 11:15 a.m. SGT. 2.2.2. UTCI calculation The Physical-modeled UTCI serves as the humidity-heat stress index and the second target variable in this study, alongside LST. The workflow of Tmrt and UTCI calculation is shown in Figure 1. The building height model was derived from a near-official urban building dataset developed by Pei and Stouffs (2025), which is not publicly available due to proprietary restrictions. Tree canopy height was obtained from the 1-m resolution global canopy height maps produced by Meta and the World Resources Institute, which were generated using self-supervised vision transformers and convolutional decoders trained on airborne lidar data (Tolan et al., 2024). The dataset was accessed on 10 August 2025 via GEE Community Catalog (https://gee-community-catalog.org/ projects/meta_trees/) and is publicly distributed under the Creative Commons Attribution 4.0 International (CC BY 4.0) license. SOLWEIG estimates mean radiant flux (Rstr , W m−2 ) from six-directional radiation (north, south, east, west, upward, and downward), weighted by the corresponding view factors, surface orientation, and shading, as expressed in Equation 1. The resulting total radiant flux is then converted to Tmrt using the Stefan–Boltzmann relation (Equation 2):
Rstr = ζk
6 X
Ki Fi + εp
i=1
Tmrt =
6 X
Li Fi
(1)
i=1
Rstr εp σ
1/4 − 273.15
(2)
where Ki and Li denote directional shortwave and longwave radiation fluxes (W m−2 ), respectively. 7
Fi are angular view factors. The absorption coefficient for shortwave radiation (ζk ) was typically 0.7 for a standing person, and the emissivity of the human body (εp ) was set to 0.97. σ is the Stefan–Boltzmann constant (5.67 × 10−8 W m−2 K−4 ). Finally, we used the operational UTCI algorithm based on the 6th-order multivariate polynomial published by Bröde et al. (2012). Hourly meteorological data at ≈2 km resolution were obtained from the National Solar Radiation Database (NSRDB) of the National Renewable Energy Laboratory (NREL) (Sengupta et al., 2018), including Ta , relative humidity, 10-m wind speed, global horizontal (GHI), direct normal (DNI), and diffuse horizontal (DHI) irradiance (Table S1 in the Supplementary Material). Subsequently, these data were downscaled to 30 m through spatial interpolation to support city-scale analyses and provide finer inputs for Tmrt and UTCI mapping. We intentionally avoided generating interpolated high-resolution surfaces from local meteorological station data using approaches such as regression-based data fusion (Han et al., 2024; Yang et al., 2024). Such methods rely on auxiliary covariates to reconstruct spatial fields and would therefore embed predictor information (e.g., land use and urban morphology) into the target variable (i.e., UTCI), violating predictor–response independence and confounding the attribution of urban drivers of heat stress. To temporally align with the LST composite, we computed the UTCI using the hourly NSRDB meteorological data. We extracted the specific 11:00 am UTCI maps for every day in May 2019 and averaged them to create a representative monthly mean UTCI map. 2.2.3. 3D urban morphology Furthermore, based on the urban building dataset, SVF, floor area ratio (FAR), building density (BD), building height (BH), and standard deviation of building height (BHsd ), are used to comprehensively characterize 3D building morphology. SVF is a key parameter in urban climatology and environmental studies, especially for high-density urban areas. It is a dimensionless parameter (ranging from 0 to 1) that quantifies the proportion of the sky visible from a given point on the ground or a surface. It reflects the openness of the surrounding environment to the sky and is influenced by urban morphology, such as building height, density, and vegetation. We adopted the grid calculation model following Lindberg and Grimmond (2010) to calculate the SVF from buildings and tree canopy. The SVF was computed using the GPU-accelerated algorithm proposed by Li and Wang (2021), which implements the grid calculation principles of Lindberg and Grimmond (2010) to enable highly efficient, city-scale processing at a 1-m resolution. The ray-casting algorithm was configured with an angular sampling of 360 directions and a search radius of 150 m to accurately capture the shading influence of surrounding urban geometry. The tree canopy was treated as a porous and height-dependent volume. To accurately represent the shading effect of dense tropical foliage, the shortwave transmissivity of light through the tree canopy was set to 3%, which is 8
consistent with the recommended parameters in the standard SOLWEIG models. The SVF values on the rooftops of the buildings and water bodies were masked before being aggregated in each analysis unit. Similarly, tree canopy density (CD), canopy height (CH), including its standard deviation (CH_sd), are calculated using the tree canopy dataset to comprehensively characterize 3D canopy morphology. 2.2.4. 2D landscape indices Drawing on the methodologies established by previous research (Chen et al., 2024; Zhang et al., 2024b; Li and Zeng, 2024), this study selected six widely used landscape indices to assess 2D urban morphology, including percentage of landscape (PLAND), patch density (PD), landscape shape index (LSI), patch cohesion index (COHESION), mean contiguity index (CONTIG), and Shannon diversity index (SHDI), which are calculated based on the land covers that can fully describe urban structural characteristics and spatial configuration of the landscape. The mathematical methods for these calculations are delineated in Table 1. All specified landscape indices were computed using the FRAGSTATS 4.2 software (McGarigal et al., 2002; McGarigal, 2015). 2.2.5. Environmental variables Environmental variables, including DEM, albedo, NDVI, normalized difference building index (NDBI), wetness (WET) derived from Landsat 8 Operational Land Imager (OLI), along with the proximity indicator distance to the coastline (Dsea ), are commonly used as indicators for monitoring the surface environment (Ramsay et al., 2025; Tanoori et al., 2024; Xu et al., 2024). Among them, NDVI is constructed based on the absorption characteristics of vegetation leaves in the red light band and their reflectance in the near-infrared band. It can reflect plant biomass and vegetation coverage. Urbanization, characterized by replacing natural ecosystems with impervious building surfaces, fundamentally reshapes surface energy and moisture fluxes. NDBI is used to delineate built-up and bare-soil features. Similarly, the WET index, derived from the Tasseled Cap Transformation, can reflect surface moisture conditions and evapotranspiration potential, particularly soil moisture (Baig et al., 2014). Similarly, we applied the same water mask method as LST and UTCI to all environmental variables, restricting covariates to terrestrial surfaces and preventing unstable or spurious values over water. 2.2.6. Socio-economic factors In addition, population density (PopD) was derived from the Resident Population by Planning Area/Subzone of Residence, Ethnic Group and Sex dataset from the Singapore Department of 9
Statistics, based on the Census of Population 2020 (Singapore Department of Statistics, 2020). While road network density (RD) and road intersection density (IntD) were derived from OpenStreetMap (OSM) data. These variables were used to represent the intensity of human activity (Liu et al., 2025, 2026). 2.3. GW-XGBoost model 2.3.1. Moran’s I statistic Moran’s I statistic (both global and local) was computed as a diagnostic tool to confirm spatial dependence in both target variables and model residuals, thereby justifying the use of spatially explicit machine-learning models. This diagnostic step provides evidence of intrinsic spatial autocorrelation, thereby justifying the use of spatially explicit methods. The value of Moran’s I ranges from −1 to 1. A positive Moran’s I indicates positive spatial correlation, and the larger the value, the more pronounced the spatial correlation, while a negative Moran’s I indicates negative spatial correlation. Global Moran’s I index can be calculated as Equation 3. n I= · S0
Pn
i=1
Pn
j=1 wij (xi − x̄)(xj − x̄) Pn 2 i=1 (xi − x̄)
S0 =
n X n X
wij
(3) (4)
i=1 j=1
where n denotes the number of spatial units (i.e., subzones), wij is the spatial weight between locations i and j (determines the degree of influence). xi is value of the variable at location i and xj is value of the variable at location j. x̄ means the global average value of the independent variable at all locations. Local spatial clusters were identified using the Local Indicators of Spatial Association (LISA) derived from Moran’s I, enabling the identification of spatial clusters (high–high and low–low) and spatial outliers (high–low and low–high) within the study area. A hybrid spatial-weights matrix was adopted, combining Queen contiguity and a nearest-neighbor correction for isolated subzones. Specifically, the weights were primarily based on Queen contiguity, where polygons sharing either a common edge or vertex were considered neighbors. And any isolated polygons (islands) were connected symmetrically to their nearest neighbor based on centroid distance. The final weights matrix was row-standardized prior to analysis. 2.3.2. Geographically weighted regression (GWR) It is essential to consider the spatial heterogeneity of relationships between located variables in the analysis of geographical processes (Wang et al., 2022; Yang et al., 2023). GWR is introduced here as
10
a conceptual baseline to motivate geographically weighted machine-learning models, rather than as a primary analytical tool in this study. Unlike traditional regression models that assume relationships between variables are the same across all locations (i.e., spatially stationary relationships), GWR creates a separate regression equation for each specific location. In the classical GWR, a simple linear regression is applied locally to capture the relationship between the dependent (response) variable and the independent (predictor) variables for each location (Huang et al., 2010). At each specific location i, GWR determines the local regression coefficients by focusing on nearby observations within a defined spatial neighborhood (typically characterized by a bandwidth). These coefficients represent the localized relationship between the dependent variable and the independent variables within the bandwidth area. The general expression for GWR models is given by Equation 5.
yi = β0 (ui , vi ) +
K X
βk (ui , vi )xik + ϵi
(5)
k=1
where yi represents the dependent variable for the ith location, (ui , vi ) denotes the coordinates of the ith location, β0 is the constant term, βk (ui , vi ) (k = 1, . . . , K) represents the local regression coefficients of the independent variables at the location i, K is the total number of independent variables (excluding the intercept), and ϵi is the random residual which usually obeys a normal distribution (Gaussian distribution) hypothesis. The observations in the neighborhood location β(uj , vj ) of location i are collected to estimate the parameters β(ui , vi ) in a geographically weighted least square approach (Equation 6): β̂(ui , vi ) = X T W (ui , vi )X
−1
X T W (ui , vi )y
1 1 X= . .. 1
(6)
x11
x12
···
x21
x22
···
.. .
.. .
..
xn1
xn2
···
11
.
x1k x2k .. . xnk
(7)
w1 (ui , vi ) 0 W (ui , vi ) = 0 .. . 0
0
0
···
0
w2 (ui , vi )
0
···
0
0
w3 (ui , vi )
···
0
.. .
.. .
..
.
.. .
0
0
···
wn (ui , vi )
(8)
where β̂(ui , vi ) denotes the estimate of the location-specific parameter, X (Equation 7) and y are matrixes of independent variables and dependent variable for neighboring study area determined by the kernel bandwidth, respectively, n is the number of observations (rows), and the matrix has dimensions n × (k + 1) because it includes the intercept column. W (ui , vi ) ∈ RN ×N in Equation 8 indicates the diagonal spatial weights matrix generated by kernel functions for total N observations. Generally, the off-diagonal elements of W (ui , vi ) are zero, and the diagonal elements of W (ui , vi ) denote the geographical weight of N samples for the specific observation i. The coefficients of the same independent variable vary across the study area because of the respective estimates at each sample location. And such coefficients can be used in further spatial analysis. Besides, selecting an appropriate kernel bandwidth is crucial in GWR, as the model’s performance is highly sensitive to this parameter (Li et al., 2010). In this study, the optimum bandwidth (or the best neighbor size) is determined by minimizing the Akaike Information Criterion (AIC) values. A bi-square kernel with an adaptive distance approach was employed to capture microclimatic variations, as it allows the spatial influence of neighboring observations to vary with local data density in compact and highly heterogeneous urban environments such as Singapore. The weight using the bi-square kernel function is expressed as Equation 9:
wj (ui , vi ) =
2 2 1 − dij ,
if dij ≤ b,
0,
if dij > b,
b
(9)
where wj (ui , vi ) is spatial weight for observation j at the target location ui , vi , dij is the distance between observation jand the target locationi, b is the bandwidth parameter, defining the maximum distance where weights are non-zero and controlling the spatial extent of the weights.
12
2.3.3. XGBoost global model As proposed by Chen and Guestrin (2016), the XGBoost algorithm is a scalable and efficient tree boosting framework designed to address the limitations of traditional gradient boosting methods. It enhances prediction accuracy by iteratively building Classification and Regression Trees (CARTs) that address the residual errors of previously constructed trees. The foundation of the XGBoost algorithm is the concept of boosting, an ensemble learning approach that aims to integrate multiple weak learners to form a strong learner. Each weak learner in the ensemble is sequentially trained to correct the errors made by its predecessors. In other words, each new CART is trained to fit the residuals of the preceding CARTs, thereby incrementally reducing the overall error. For a given set of predictions ŷ, boosting minimizes the following objective as Equation 10 and 11:
L(t) =
n X (t−1) ℓ yi , ŷi + ft (xi ) + Ω(ft ),
(10)
i=1
ŷi =
T X
ft (xi )
(11)
t=1
where ℓ is the loss function that quantifies the difference between the observed (yi ) and predicted (t)
(ŷi ) values for the ith observation, ft (xi ) is the output of the tth tree for the ith observation, T is the total number of trees in the ensemble, and Ω(ft ) is the regularization term of the tth tree that penalizes model complexity to prevent overfitting. 2.3.4. GW-XGBoost model The GW-XGBoost, which combines the concepts of GWR and XGBoost, is employed for analysis due to its potential to address the limitations of the GWR model and improve predictive performance over a non-geographically weighted XGBoost model. Similar to the regression analysis framework of GWR, GW-XGBoost consists of multiple sub-models calibrated locally using XGBoost instead of linear regression as shown in Equation 12 and 13.
L(t) (ui , vi ) =
n X
(t−1) wj (ui , vi ) ℓ yj , ŷj + ft (xj ) + Ω(ft )
(12)
j=1
ŷj =
T X
ft (xj )
(13)
t=1
where wj (ui , vi ) is the spatial weight for observation j relative to location (ui , vi ), computed as Equation 9. Although Equations 11 and 13 share the same additive prediction form, the regression trees ft in Equations 13 are estimated locally using spatially weighted loss functions, and
13
therefore differ across target locations. Moreover, in GW-XGBoost, the intercept is not explicitly modeled, as decision trees inherently account for offsets through recursive partitioning and leaf values. Consequently, an explicit intercept term is unnecessary, and the first column of Xin Equation 7 can be removed. The resulting feature matrix Xmodified is defined in Equation 14.
x11 x21 Xmodified = . .. xn1
x12
···
x22
···
.. .
..
xn2
···
.
x1k x2k . .. . xnk
(14)
For each target location, GW-XGBoost fits a local XGBoost model by emphasizing nearby observations through distance-based sample weights, while retaining the nonlinear learning capacity and interpretability of gradient-boosted trees. Specifically, spatial dependence is incorporated through instance-weighted loss minimization, where observations closer to the target location (ui , vi ) receive stronger influence during model training. Specifically, spatial weights wj (ui , vi ) were passed as sample weights to the local XGBoost training procedure, which natively supports instance weighting. This approach preserves the integrity of tree-based split criteria while allowing geographically proximate observations to contribute more strongly to local model calibration. For each target location, GW-XGBoost fits a local XGBoost model using the original feature matrix Xmodified , with spatial proximity encoded solely through sample weights rather than feature re-scaling. In this way, the framework retains the nonlinear learning capacity of gradient-boosted trees while explicitly accounting for spatial heterogeneity. To provide an overview of the pairwise relationships among LST, UTCI, and all urban factors used in this study, we included a Pearson correlation matrix in the Supplementary Material (Figure S2). This figure helps illustrate the degree of association among predictors and provides additional context for interpreting the modeling results. Model interpretability in the GW-XGBoost framework is achieved through SHAP, which quantifies each local prediction into additive feature contributions expressed in consistent units of the response variable, thereby enabling spatially explicit interpretation of model behavior. Details of the SHAP-based analysis are presented in Section 3.3.
14
2.3.5. Model evaluation To ensure robust model evaluation and avoid overfitting during hyperparameter tuning, a nested cross-validation (Nested CV) procedure was firstly employed for the global XGBoost model (Grekousis, 2025). The framework consists of two loops: an inner loop that performs hyperparameter optimization and an outer loop that evaluates model generalization on held-out subsets. In the inner loop, combinations of three primary XGBoost hyperparameters, i.e., the number of trees (n_estimators), learning rate (learning_rate), and maximum tree depth (max_depth), were systematically tested using grid search. The configuration yielding the lowest root-mean-square error (RMSE) within each training subset was selected. In the outer loop, the optimized hyperparameters obtained from the inner cross-validation were evaluated across five independent test folds to assess the model’s generalization performance. This hierarchical procedure provides an unbiased estimate of the model’s predictive capability by ensuring that hyperparameter selection and model validation are conducted on separate data partitions. The resulting performance metrics (R2 , MAE, and RMSE) were averaged across outer folds to quantify both model accuracy and generalization stability. After completing the nested CV procedure, the model was retrained using the optimal hyperparameters on the entire training dataset and subsequently evaluated on an independent hold-out test set (not involved in any CV step). These globally optimized hyperparameters were then applied uniformly across the GW-XGBoost models to ensure cross-location comparability; however, because the present implementation did not impose strict block-based spatial CV, some residual optimism in predictive accuracy due to spatial autocorrelation may remain. The hyperparameter settings for XGBoost and GW-XGBoost are shown in Table S2. The local GW-XGBoost model was validated through spatial cross-validation during the bandwidth optimization process. Tree-level hyperparameters were inherited from the global nested CV search to prevent over-tuning at the local level, while the local component was governed by the bandwidth (bw) and kernel weighting. An adaptive bisquare kernel (as shown in Equation 9) was employed to define local neighborhoods, and the optimal bw was determined by minimizing cross-validation loss using a leave-one-out spatial CV (Wong, 2015; Pang et al., 2023; Cawley and Talbot, 2008). For each candidate bw, the model predicted each observation using all other locations, with distance-based weights applied to nearby samples. The final model was trained with the optimal bw and used to generate diagnostic maps of local R2 and standardized residuals to evaluate spatial heterogeneity in model performance. To rigorously evaluate the predictive performance of the GW-XGBoost models while preventing overfitting, we utilized a pseudo-out-of-bag (OOB) validation approach inherent to Stochastic Gradient Boosting. Unlike standard Random Forests, which utilize bootstrap aggregating (bagging), standard gradient boosted trees typically process all
15
training data sequentially. To generate an internal validation metric, the subsample hyperparameter was set to 0.8. Consequently, during the construction of each localized model—which is trained using data from the target subzone and its 94 nearest neighbors defined by the adaptive kernel—the algorithm stochastically samples exactly 80% of the available local instances to grow each individual decision tree. The remaining 20% of the instances are temporarily withheld as an ’out-of-bag’ validation set for that specific iteration. As the boosting sequence progresses, the model continuously predicts the target variable for these stochastically withheld, unseen instances. The reported global OOB R2 therefore represents the aggregated predictive accuracy calculated exclusively from these pseudo-OOB instances across the entire spatial framework.
3. Results 3.1. Spatial distribution of LST and UTCI Figure 2 compares the spatial patterns of 30-m LST (Figure 2a) and 1-m UTCI (Figure 2b) at the city scale across the main island of Singapore, highlighting substantial differences between surface radiative temperatures and human heat stress. Overall, LST retrieved from Landsat 8 TIRS demonstrates both higher maximum (with temperatures exceeding 50 ◦ C) and lower minimum values than UTCI. The LST map reveals pronounced surface overheating in industrial and residential areas within industrial and residential clusters in the western, eastern, and northern regions. In contrast, UTCI represents a more moderate spatial thermal condition experienced by humans at the pedestrian level. The 1-m UTCI map shows finer variations influenced by shading, vegetation, and local morphology. The UTCI values ranging from 31.75 ◦ C to 39.76 ◦ C predominantly fall within the categories of strong heat stress (32–38 ◦ C) and very strong heat stress (38–46 ◦ C) (Bröde et al., 2012). The zoomed-in panels in Figure 2c–g further present the diurnal dynamics of UTCI in a patch area of the Clementi neighborhood, capturing the temporal variations at 08:00, 11:00, 14:00, 17:00, and 20:00. The coverage and direction of shading are changing with time. Higher UTCI values are observed around midday (14:00) in open and impervious areas, consistent with findings reported by Park et al. (2014), while cooling effects from tree canopy and building shadows are more evident in the early morning and late afternoon. Figure 3 compares LST and human-centric heat stress across land use/land cover (LULC) types. The double violins and box plots in Figure 3a show the probability density, median and quartile values of LST and UTCI for each LULC class at 30 m resolution. Built-up and impervious areas show the highest LST, indicating strong surface heating effects. However, their corresponding UTCI values are relatively lower, particularly in build-up areas. Vegetated areas with canopy cover show markedly lower LST and UTCI values than other LULC types, and the consistently lower UTCI 16
Figure 2: Spatial distribution of 30-m LST (a) and 1-m UTCI (b), and the dynamic variations of UTCI at 8:00, 11:00, 14:00, 17:00, and 20:00 in a zoomed-in patch area (c-g).
further highlights the cooling influence of shading and evapotranspiration. Moreover, barren surface and vegetation without canopy show similar LST and UTCI magnitudes. In contrast, water bodies exhibit the lowest LST but much higher and stable UTCI, reflecting the decoupling between surface temperature and near-surface thermal conditions. Similarly, Figure 3b compares two nonparametric visualizations of the LST-UTCI relationship across LULC types. The binned median with the interquartile range (IQR, 25th–75th percentile) highlights the variability of UTCI within each LST bin, while the locally weighted scatterplot smoothing (LOWESS) regression reveals the smooth overall nonlinear trends (Cleveland, 1979; Moran, 1984). Both show that UTCI generally increases with LST, but the relationship is non-linear and saturates at higher LST values. However, water shows an opposite and negative correlation, which is most prominent during hot conditions. Furthermore, buildings and impervious surfaces exhibit large scatter and wide IQRs, providing evidence of microclimatic heterogeneity due to shading and material differences. Both temperature metrics were aggregated to subzones, which are the smallest planning units in Singapore. Before aggregation, UTCI values over rooftops were subsequently excluded, as these areas do not correspond to inhabited spaces and are therefore irrelevant for human heat exposure assessment. Water pixels were also masked in both the 30-m LST and 1-m UTCI surfaces. Figure 4 presents a side-by-side comparison of the categorical thermal risk and the quantitative mismatch. Figure 4a utilizes a bivariate choropleth map to illustrate the combined spatial distribution of LST and UTCI, while Figure 4b maps the absolute standardized difference (|z(LST) − z(UTCI)|). The 17
Figure 3: Comparison between LST and human-centric UTCI across LULC types.
spatial distribution in Figure 4b shows that the most severe discrepancies (dark orange regions) are not randomly distributed but are concentrated in specific urban areas that exhibit relatively low LST yet elevated UTCI. 3.2. The performance of spatial machine learning models Figure 5 presents the spatial autocorrelation of LST (upper row) and UTCI (lower row) mean values across subzones in Singapore, using both Moran’s I scatter plots and Local Indicators of Spatial Association (LISA) cluster maps. Both global and local spatial autocorrelation analyses reveal that LST and UTCI exhibit significant positive spatial dependence across subzones (Moran’s I = 0.43 and 0.37, p < 0.05), with high LST clusters primarily concentrated in eastern dense urban residential and western industrial areas and low LST clusters in vegetated and coastal areas, while high UTCI is mainly distributed in southern city center, eastern resident area, western industrial areas, and port/airport. Table 2 compares the predictive performance of the global XGBoost model and the GW-XGBoost for estimating LST and UTCI. For both targets, the XGBoost models exhibited strong predictive
18
Figure 4: Spatial relationship and quantitative divergence between surface temperatures and human-centric heat stress. (a) Bivariate choropleth map illustrating the coupled spatial distribution of LST and UTCI. The rotated diamond legend categorizes intersecting thermal profiles into four primary extremes: (1) Bottom corner (white): Low LST and low UTCI, representing cool and comfortable environments; (2) Top corner (dark purple): High LST and high UTCI, indicating severe, compounding thermal risk; (3) Right corner (blue): High LST but low UTCI, representing hot physical surfaces where human heat stress is successfully mitigated (e.g., via high ventilation); and (4) Left corner (red/pink): Low LST but high UTCI, representing cool surfaces that nonetheless experience severe pedestrian heat stress (e.g., due to trapped radiation or lack of shade). (b) Spatial distribution of the absolute standardized mismatch per subzone (|z(LST) − z(UTCI)|). Dark orange areas highlight severe quantitative discrepancies.
Figure 5: Global and local Moran’s I statistics for LST and UTCI. A hybrid Queen–nearest-neighbor spatial-weights matrix is used, where isolated subzones are symmetrically connected to their nearest polygon to ensure full spatial connectivity.
performance, achieving test set R2 values of 0.872 for LST and 0.831 for UTCI. The corresponding MAE of 0.563 ◦ C and 0.188 ◦ C, respectively, indicate small prediction deviations. The GW-XGBoost 19
2 model for LST yielded comparable aggregated predictive performance (ROOB = 0.855, MAE =
0.600 ◦ C, RMSE = 0.808 ◦ C), showing only a modest difference from global XGBoost, whereas enhancing local interpretability by capturing spatially varying relationships with a mean local R2 of 2 0.750 ± 0.108. For UTCI, GW-XGBoost improved fit precision (ROOB = 0.905, MAE = 0.172 ◦ C,
RMSE = 0.261 ◦ C) and exhibited a mean local R2 of 0.858 ± 0.057, highlighting its superior ability to represent spatial heterogeneity. To assess spatial autocorrelation, we calculated the Global Moran’s I for the residuals of both Global XGBoost and GW-XGBoost (Figure S3). For LST, the residual autocorrelation was weak and non-significant in both models (Moran’s I = 0.012, p = 0.306 and I = 0.034, p = 0.157, respectively), indicating both models already showed minimal residual spatial dependence for LST. For UTCI, the Global XGBoost residuals showed a marginal tendency toward clustering (I = 0.037, p = 0.087). However, the implementation of the GW-XGBoost model more effectively accounted for the remaining spatial heterogeneity in UTCI (I = -0.015, p = 0.372). Overall, although the global models performed well, the GW-XGBoost framework maintained rigorous spatial randomness in the residuals while also enabling the identification of highly localized and spatially varying feature relationships. Table 2: Comparison of model performance between XGBoost and geographically weighted XGBoost (GW-XGBoost) for predicting LST and UTCI. Target
Model XGBoost
LST
GW-XGBoost XGBoost
UTCI
GW-XGBoost
Evaluation
R2
MAE (°C)
RMSE (°C)
Main Hyperparameters
Kernel
Bandwidth
Nested CV (mean ± SD) Test set Mean local performance Global OOB Nested CV (mean ± SD) Test set Mean local performance Global OOB
0.859 ± 0.032 0.872 0.750 ± 0.108 0.855 0.898 ± 0.049 0.831 0.858 ± 0.057 0.905
0.582 ± 0.073 0.563 0.022 ± 0.021 0.600 0.176 ± 0.006 0.188 0.007 ± 0.006 0.172
0.818 ± 0.051 0.757 0.817 ± 0.091 0.808 0.273 ± 0.022 0.336 0.261 ± 0.038 0.261
nestimators = 500; learning_rate = 0.05; max_depth = 2
– – Adaptive
– – 94 (best)
nestimators = 500; learning_rate = 0.05; max_depth = 2
– – Adaptive
– – 94 (best)
*Note:CV = cross-validation; OOB = out-of-bag. Mean ± SD values indicate model stability across folds (for global) or across subzones (for local). All models used identical hyperparameter settings for comparability. The adaptive kernel bandwidth of 94 was selected as optimal.
Figure 6 visualizes the spatially explicit R2 and standardized residuals of the GW-XGBoost models for LST and UTCI, respectively. Both models achieved high local explanatory power (local R2 > 0.64 for LST and local R2 > 0.80 for UTCI, respectively) across planning subzones, with slightly lower performance in the central region for LST and in the airport for UTCI. Standardized residuals also show localized over-predictions LST in the western Tuas Port and forest area (Figure 6b), while under-predictions UTCI in the airport, Pasir Panjang Terminals Port, and the western water catchment with many impervious and barren lands and no-canopy vegetation that lack shading (Figure 6d). Overall, the GW-XGBoost results indicate that UTCI exhibits a systematically stronger sensitivity to urban spatial structure than LST. While LST spatial variation is largely driven by land cover types (Wang et al., 2024b), the superior predictive accuracy of our spatial model for UTCI suggests that human heat stress is intricately governed by 3D morphological arrangements. This heightened 20
Figure 6: Spatially local R2 and standardized residuals of the GW-XGBoost model for LST (a-b) and UTCI (c-d) across subzones of Singapore.
structural sensitivity underscores the potential for targeted urban design interventions to mitigate heat stress, a theme further explored through feature importance analysis in Section 3.3. 3.3. The difference of influencing factors in LST and UTCI revealed by explainable spatial machine learning The global feature importance and local summaries derived from the SHAP framework are visualized in Fig 7. For LST, environmental variables such as WET, NDBI, and NDVI dominate the predictions, highlighting the trade-off between vegetation and imperviousness. Whereas, morphological variables (SVF and DEM), radiative factor (Albedo), and WET are more influential for UTCI, emphasizing the role of shading and vertical openness in modulating human physiological heat exposure. Figure 8 further compares the Top-5 ranked features influencing LST and UTCI for each subzone. The ranking patterns vary considerably across subzones, indicating that the spatial heterogeneity in dominant predictors of surface and thermal comfort temperatures. Moreover, Figure 9 compares local variable dominance and SHAP-based effects for GW-XGBoost models of LST and UTCI in each subzone. Figure 9a shows the conventional tree-based primary gain importance variables in each subzone, indicating the model’s performance relied on them most to explain spatial temperature variation, where five factors out of 26 predictor variables exert the strongest local influence on improving the local R2 of the LST model across all 328 subzones. However, 3D urban morphological variables are not among the primary variables of importance, indicating a limitation of LST in reflecting 3D urban structure. In contrast, Figure 9d reveals 21
Figure 7: SHAP global feature importance and local summary (beeswarm) plots for LST (a) and UTCI (b).
that 3D metrics, including SVF, canopy height (CH), and the standard deviation of canopy height (CH_sd), have the most important influence on a higher R2 of UTCI and cover most built-up areas. Furthermore, Figure 9(b,e) measures actual contribution magnitude (mean |SHAP|) of each predictor per subzone and quantifies the marginal effect of each predictor on the target variable in consistent units of the response variable (e.g., ◦ C). Figure 9b shows that the primary SHAP-based importance features are mainly WET, NDVI, NDBI, and BD. Furthermore, Figure 9e confirms SVF as the top explanatory variable across most subzones. Moreover, the signed SHAP value in Figure 9f illustrates the spatially explicit SHAP value of SVF on the predicted UTCI across different subzones. Red areas with positive SHAP values, mostly concentrated in the city center (P1), northeast residential estates such as Seletar (P2), western industrial areas and ports (P3), and the airport (P4), indicate that areas with large SVF values increase UTCI by a large margin. In contrast, blue areas with negative SHAP value in the Simpang area (N1), the central water catchment and Bukit Timah Hill (N2), and the Tengah areas (N3) show that the lower SVF tend to reduce UTCI substantially. 3.4. The nonlinear relationships between urban factors and UTCI SHAP dependence plots coupled with Generalized Additive Model (GAM) smoothers in Figure 10 further reveal the nonlinear relationships between the six most influential predictors and UTCI. The top-ranked SVF exhibits a strong monotonic increase after a tipping point at approximately 0.51, 22
Figure 8: Spatial variation of the top-five feature importance ranks across subzones for (a) LST and (b) UTCI. Each column corresponds to a predictor, and each row to a subzone. Colors represent the local ranking of each feature (1 = Top).
indicating that enclosed or shaded areas help reduce heat stress by limiting radiative load, whereas beyond this value, greater sky openness amplifies heat stress through enhanced solar exposure. WET exhibits an inverse non-linear effect: UTCI decreases with increasing WET beyond 0.612, reflecting moisture-driven cooling and vegetation influence. DEM also shows a clear cooling trend with elevation, above 10 m the GAM curve declines, suggesting enhanced airflow or reduced surface 23
Figure 9: Spatially explicit results of the GW-XGBoost model for LST and UTCI across subzones of Singapore. (a-c) Primary gain, Primary SHAP, and SHAP value of SVF contributing to LST in each subzone; (d-f) Primary gain, Primary SHAP, and SHAP value of SVF contributing to UTCI in each subzone. All these maps illustrate spatial heterogeneity in the controlling factors of LST and UTCI and the reliability of local model fits.
heat storage at higher terrain.
Figure 10: SHAP dependence plots with GAM smoothers showing the non-linear relationships between the six most influential factors and UTCI, including sky view factor (SVF), wetness index (WET), elevation (DEM), albedo, canopy density (CD), and patch density (PD). “Transition point” denotes the first zero-crossing of the smoothed SHAP curve, indicating an approximate transition from negative to positive model-attributed effect (or vice versa).
Although high albedo could reduce surface and air temperature through increased reflectivity (Schneider et al., 2023; Akbari et al., 2012; Prado and Ferreira, 2005; Yi et al., 2025b), our SHAP–GAM analysis indicates a warming effect with higher albedo values. Stratified SHAP analysis in Figure 11 further indicates that the positive association between albedo and UTCI is concentrated
24
primarily in open/high-SVF environments. Under high-SVF conditions, limited shading increases direct solar radiation receipt, while high-albedo surfaces reflect a greater proportion of incoming shortwave radiation. Although our spatial framework does not directly measure pedestrian-level shortwave flux, these SHAP interaction patterns strongly support the physical hypothesis that highalbedo pavements intensify the combined effects of multiple shortwave radiation reflections, scattering the radiant flux directly onto pedestrians. This hypothesized mechanism aligns with established urban physics, where elevated Tmrt subsequently contributes to increased UTCI (Schneider et al., 2023; Schrijvers et al., 2016; Yi et al., 2025c). Moreover, previous studies have shown that increased albedo is associated with reduced cloud cover, especially in high urban surface fractions (Hamwey, 2007; Jacobson and Ten Hoeve, 2012), which is accompanied by increased shortwave radiation.
Figure 11: SHAP dependence of Albedo with Sky View Factor (SVF) interaction. Scatter points indicate the SHAP value of Albedo for individual observations, colored by their respective SVF to highlight localized feature interactions. The solid black line represents the overall Generalized Additive Model (GAM) trend, demonstrating a non-linear relationship between albedo and UTCI.
The SHAP–GAM dependence plot for CD reveals a pronounced non-linear relationship with UTCI. At CD values below the tipping point at approximately 0.81, SHAP values remain close to zero, partly because much of the shading effect associated with sparse or moderately dense tree cover is already captured by SVF in the model. As a result, CD contributes little additional explanatory power. Beyond the tipping point, however, the dense and more continuous canopy structure leads to substantial reductions in UTCI. The cooling effect strengthens further for CD > 0.85, reflecting the combined benefits of enhanced shading, reduced Tmrt , and increased evapotranspirative cooling from cohesive tree canopies. PD was computed at the landscape level in Fragstats, meaning that all land-cover classes were included when counting patches. Low PD (< 550 patches/100 ha) represents coarse-grained
25
landscapes dominated by large continuous units, which is associated with increased UTCI. As PD increases to moderate levels (∼550–750 patches/100 ha), the landscape becomes more spatially fragmented and finely mixed, resulting in negative SHAP values. However, very high PD (> 900 patches/100 ha) reflects highly fragmented landscapes composed of many small, disconnected patches, which contributes to lower cooling efficiency (Qi et al., 2025) and a slight increase in UTCI. Overall, the SHAP–GAM analysis highlights that urban factors exhibit non-linear and thresholddependent relationships with UTCI. GAM curves combined with SHAP values reveal complex effect directions that cannot be captured by linear models, underscoring the necessity of using interpretable nonlinear modeling to understand and design effective urban-climate mitigation strategies.
4. Discussion 4.1. Discrepancies between LST and UTCI The observed discrepancies between LST and UTCI reflect conceptual and physical differences between surface radiative temperatures and human physiological heat stress. Although LST and UTCI exhibit an overall positive but non-linear correlation, they differ systematically in both magnitude and spatial configuration. Compared to UTCI, LST exhibits higher maximum and lower minimum values, indicating a wider thermal range driven by the strong radiative heating and cooling of surface materials and characteristics (Small, 2006). LST captures the instantaneous radiative temperature of the land surface, which responds sharply to solar heating. In contrast, UTCI integrates air temperature, humidity, wind speed, and radiation to represent perceived heat stress at the pedestrian level, which is buffered by atmospheric processes such as convection, radiation exchange, and physiological factors (Bröde et al., 2012; Fiala et al., 2012; Jendritzky et al., 2012; Zhan et al., 2025), leading to a more moderated spatial and temporal pattern. These spatial mismatches occur across land cover types and subzones. In built and impervious areas, UTCI values are systematically lower than LST. This discrepancy is driven by pedestrianlevel microclimatic regulation, whereby shading from surrounding buildings substantially reduces shortwave radiation exposure and mean radiant temperature, thereby alleviating human thermal stress (Briegel et al., 2025). When combined with the canopy shading and evapotranspirative cooling provided by urban vegetation, these mechanisms effectively mitigate human thermal stress —pedestrian-level benefits that are inherently missed by satellite-derived surface temperature measurements (Zhan et al., 2025). Conversely, in open areas with sparse or no canopy cover, such as grassland or barren land, UTCI may remain elevated due to direct solar exposure and limited shading, despite relatively moderate surface temperatures. Water bodies exhibit the lowest LST due to the high specific heat capacity of water, and the water’s stable, cooler temperature during hot 26
periods provides a significant local cooling effect on the air directly above it (Oke, 2002). However, the high humidity and intense solar radiation without shading above the water result in a much higher UTCI value. At the subzone level, discrepancies between LST and UTCI are most pronounced in areas characterized by complex 3D urban morphology. These divergences arise from the differential sensitivity of the two metrics to urban morphology and microclimatic regulation. LST can be better explained by built-surface characteristics (NDBI, NDVI, BD), whereas UTCI is more sensitive to morphological and radiative factors (SVF, DEM, and albedo). Reduced SVF, increased building height, and canyon enclosure enhance shading from surrounding buildings, thereby mitigating UTCI even when surface temperatures remain high (Briegel et al., 2025). As a result, LST is considerably less sensitive to independent 3D urban morphological features, making it a poor proxy of human-centric heat stress (Nazarian et al., 2022; Tuholske et al., 2021; Muse et al., 2024; Zhan et al., 2025). These systematic discrepancies highlight the spatially heterogeneous nature of urban heat exposure. The contrast between these two indices underscores that satellite-derived LST tends to overestimate heat exposure relative to pedestrian-level heat stress in shaded or ventilated environments, particularly within built-up areas (Zaerpour et al., 2025; Zhan et al., 2025). Therefore, we must recognize the potential misleading risks of using LST alone as the primary basis for climate-informed urban planning, particularly in high-density cities. Interventions guided solely by surface temperature may overlook areas of elevated human heat stress and misallocate cooling resources. Integrating human-centric metrics such as UTCI is therefore essential for identifying exposure hotspots and designing effective heat adaptation strategies. 4.2. Strengths of the explainable GW-XGBoost framework The GW-XGBoost framework proves to be a robust tool for this study, successfully balancing high-precision spatial prediction with the ability to quantify complex local interactions. The superior performance of the GW-XGBoost model, as evidenced by the higher R2 (0.905) and lower error metrics (MAE = 0.172 ◦ C) compared to global baselines, validates the necessity of accounting for spatial non-stationarity in urban thermal analysis. Unlike traditional global XGBoost algorithms that force a singular, uniform relationship across the entire study area, this geographically weighted framework incorporates spatial heterogeneity through localized weights (Fouedjio and Arya, 2024; Grekousis, 2025; Wang et al., 2024a; Yang et al., 2023). By allowing model parameters to vary across space, GW-XGBoost significantly enhances our understanding of localized variable impacts across highly heterogeneous urban subzones.
27
Beyond predictive accuracy, the explainable GW-XGBoost framework offers critical insights into the complex, non-linear mechanisms driving the discrepancies between LST and UTCI. By integrating SHAP-based feature importance, we can interpret both global and local behavior of the model. Unlike traditional tree-based importance measures (e.g., gain or split frequency), SHAP importance derives from Shapley values and quantifies each feature’s additive contribution to the prediction in consistent units of the response variable (◦ C). The mean absolute SHAP value represents the magnitude of local feature importance in each subzone, while the signed SHAP value indicates its direction of influence (warming or cooling). This approach enables spatially explicit interpretation of how different urban morphological and environmental variables drive local variations in thermal exposure across Singapore. 4.3. Implications for climate-informed urban planning The discrepancies between LST and UTCI highlight the need to rethink the indicators used in climate-informed urban planning. While LST remains valuable for assessing surface overheating and material performance, it does not reliably represent pedestrian-level heat stress, particularly in dense urban environments where shading, enclosure, and ventilation strongly modulate thermal stress (Nazarian et al., 2022; Tuholske et al., 2021; Muse et al., 2024; Zhan et al., 2025). Therefore, relying on LST alone as the primary basis for planning may misidentify heat-risk hotspots and lead to suboptimal allocation of mitigation measures. For high-density tropical cities such as Singapore, integrating human-centric thermal indicators such as UTCI is essential for identifying where heat is experienced most acutely and for designing more effective adaptation strategies. Rather than advocating for the replacement of LST, our findings support a complementary framework in which LST informs surface-level thermal characteristics, while UTCI captures human-centric heat stress. Together, they provide a more complete representation of urban thermal environments. Importantly, the spatially explicit SHAP results further reveal substantial heterogeneity in the relative importance of heat drivers across subzones, indicating that heat mitigation strategies should be spatially differentiated rather than uniformly applied. For example, in compact residential districts, UTCI is predominantly governed by enclosure-related variables such as SVF, which regulate shortwave exposure, longwave trapping, and wind attenuation. This suggests that heat mitigation in these environments should not simply maximize either openness or greening, but instead balance shading provision with the preservation of local ventilation pathways (Li and Ratti, 2018; Li et al., 2024). Beyond identifying dominant drivers, the nonlinear response revealed by SHAP-GAM analysis in Section 3.4 offers actionable insights into the intensity of intervention required to achieve meaningful cooling benefits. UTCI in high-density tropical cities is shaped by nonlinear and context-specific 28
processes, rather than by uniformly incremental responses to individual variables. Among these, SVF emerges as the dominant driver of UTCI, with a threshold at approximately 0.51. Below this value, more enclosed or shaded urban configurations help suppress UTCI by limiting radiative exposure (Oke, 2002; Nazarian et al., 2022), whereas above this threshold, increasing sky openness is associated with a monotonic rise in pedestrian heat stress. This finding suggests that open and highly exposed urban spaces should be prioritized for shading-oriented design. Tree canopy density also shows a delayed cooling response: increases below a threshold of approximately 0.81 provide limited benefit, whereas sufficiently dense and continuous canopy yields substantially stronger cooling (Li et al., 2024). Consequently, sparse greening fails to provide meaningful thermal relief, underscoring that canopy continuity is as critical as canopy presence for providing additional cooling through combined shading, reduced Tmrt , and enhanced evapotranspiration (Oke, 2002; Shashua-Bar et al., 2009). Patch density further demonstrates a nonlinear relationship, indicating that moderate landscape fragmentation (PD of 550–750 patches/100 ha) is more favorable for cooling efficiency, while highly coarse or overly fragmented configurations are associated with higher UTCI, underscoring the importance of spatial configuration alongside land-cover composition (Qi et al., 2025). Finally, contrary to conventional expectations, albedo exhibits a net warming effect beyond 0.158, indicating that reflective materials should be applied cautiously in open pedestrian environments with high SVF, where increased shortwave reflection may elevate radiant heat exposure (Schneider et al., 2023; Schrijvers et al., 2016; Prado and Ferreira, 2005). Overall, these results suggest that effective heat mitigation in Singapore should be threshold-aware and spatially differentiated: reduce excessive sky exposure in open areas, ensure that greening is sufficiently dense to become effective, avoid both overly coarse and overly fragmented patterns of landscape configuration, and evaluate reflective materials together with pedestrian-level radiation conditions rather than surface temperature alone. 4.4. Limitations and prospects While this study provides insights into the discrepancy between LST and UTCI in Singapore, the generalizability of the findings should be interpreted within specific boundary conditions. First, although the GW-XGBoost framework is generally transferable, the derived threshold values and variable importance rankings are likely to be context-specific. The observed patterns are characteristic of tropical, high-density urban environments, where high humidity and intense solar radiation are prevalent under Singapore’s equatorial climate. In such contexts, the decoupling of LST and UTCI is likely a robust phenomenon due to the significant role of Tmrt in shaping human thermal sensation. However, in temperate or arid climates, the relationship between LST and UTCI may shift, as surface radiation might dominate the thermal environment differently across seasons. Future research should validate whether the identified driving factors maintain their relative importance in cities 29
with different climatic backgrounds and urban fabrics. Additionally, while this study focused on daytime heat exposure when shortwave radiation and pedestrian vulnerability peak, the absence of nighttime analysis remains a limitation. The physical drivers of the thermal environment shift drastically after sunset. For instance, a low SVF—which provides essential cooling shade during the day—becomes a liability at night by trapping longwave radiation and impeding convective cooling. Evaluating these diurnal contrasts, specifically comparing daytime and nighttime temperatures and their shifting relationships with complex urban morphology (Wang et al., 2025), is a key focus of our ongoing research. Moreover, meteorological data were obtained from NREL with an original spatial resolution of approximately 2 km and subsequently downscaled to 30 m through spatial interpolation to support city-wide analysis. We acknowledge potential biases relative to local meteorological station observations. On one hand, it remains challenging to downscale meteorological variables to a hyperlocal level, particularly given the sparse station network and the lack of sufficient input variables (Sun et al., 2024). On the other hand, as described in the Methodology, we deliberately refrained from generating interpolated high-resolution surfaces from station data to preserve predictor–response independence. While this decision avoids embedding auxiliary covariate information into the target variable (UTCI), it may also limit the representation of fine-scale microclimatic variability captured by in situ observations. Third, it must be noted that while our UTCI calculation incorporates spatially continuous meteorological data, the native 2-km resolution of the inputs means they function primarily as meteorological forcing fields. These fields (air temperature, relative humidity, and wind speed) provide the background atmospheric state, whereas the Tmrt was resolved using fine-scale urban geometry and radiation modeling. We do not claim that these coarse meteorological fields fully capture microclimatic heterogeneity in localized fluid dynamics, such as sensible heat advection or street-canyon wind attenuation. Rather, the modeling framework represents how relatively smooth background meteorological conditions are spatially modulated by local urban form, shading, and radiative exchange, which ultimately generates the fine-scale variation in pedestrian heat stress observed in this study. Future studies could consider integrating multi-source data, including satellite products, reanalysis datasets, dense sensor networks, or mobile measurements, with physics-guided approaches, such as novel Weather Research and Forecasting (WRF) models (Gao et al., 2025; Hong et al., 2026), or uncertainty-constrained data fusion techniques to better capture micro-scale thermal variability while preserving predictor–response independence (Wen et al., 2025; Kousis et al., 2022; Liu et al., 2017). Nevertheless, the high computational cost of SOLWEIG and physically based UTCI simulation models continues to constrain their application at the full city scale, particularly when high spatial
30
resolution, long simulation periods, or multiple scenario analyses are required. These models rely on detailed radiative transfer calculations, complex urban geometry, and fine temporal discretizations, which collectively impose substantial computational and data demands. This limitation has motivated the development of data-driven and hybrid physical-informed neural network (PINN) approaches that aim to bridge the gap between physical realism and computational scalability (Shaeri et al., 2025). By embedding physical constraints such as energy balance relationships and radiative principles within machine-learning architectures, physics-informed neural networks (PINNs) can approximate key microclimatic processes while reducing computational cost (Ren et al., 2025). Such hybrid frameworks may enable efficient large-scale simulations, uncertainty quantification, and scenario exploration (Yi et al., 2025d), while retaining interpretability and physical consistency.
5. Conclusion This study provides a morphology-aware comparison between LST and the UTCI for assessing urban thermal conditions in a high-density tropical city. By integrating high-resolution UTCI modeling with explainable global and geographically weighted machine learning approaches, we demonstrate that surface-based temperature metrics and physiologically relevant heat stress indicators convey fundamentally different information about urban thermal environments and their drivers. While LST exhibits strong associations with two-dimensional surface characteristics such as imperviousness and land cover composition, it fails to adequately represent pedestrian-level thermal exposure. In contrast, UTCI is systematically more sensitive to three-dimensional urban morphology. These findings underscore that relying on LST alone may misrepresent heat-risk hotspots and therefore provide an incomplete basis for climate-informed urban planning and heat risk assessment. Collectively, this work advances urban climate research by shifting the analytical focus from surface temperatures toward human-centric heat stress metrics. For urban planners and policymakers in Singapore, these findings suggest that effective climate-informed planning in high-density tropical cities requires spatially differentiated interventions that jointly consider shading, canopy continuity, spatial configuration, and radiative conditions at the pedestrian level. In particular, open and highly exposed urban spaces should be prioritized for shading-oriented interventions, greening strategies should consider the dense and continuous canopy rather than sparse vegetation alone, and landscape patterns should avoid both overly coarse and overly fragmented configurations. Reflective materials should also be assessed together with shading design, canyon geometry, and pedestrian-level radiation exposure. Methodologically, this study further demonstrates the value of combining explainable machine learning with geographically weighted modeling to reveal spatially heterogeneous urban climate processes and support more spatially tailored heat mitigation strategies. 31
Acknowledgement This research is supported by the National Research Foundation Singapore (NRF) under its Campus for Research Excellence and Technological Enterprise (CREATE) programme.
Declaration of generative AI and AI-assisted technologies in the writing process During the preparation of this work, the authors used ChatGPT in order to improve the readability and language of the manuscript. After using this tool, the authors reviewed and edited the content as needed and take full responsibility for the content of the publication.
References Akbari, H., Matthews, H.D., Seto, D., 2012. The long-term effect of increasing the albedo of urban areas. Environmental Research Letters 7, 024004. Assaf, G., Assaad, R.H., 2024. Modeling the impact of land use/land cover (lulc) factors on diurnal and nocturnal urban heat island (uhi) intensities using spatial regression models. Urban Climate 55, 101971. Baig, M.H.A., Zhang, L., Shuai, T., Tong, Q., 2014. Derivation of a tasselled cap transformation based on landsat 8 at-satellite reflectance. Remote Sensing Letters 5, 423–431. Blazejczyk, K., Epstein, Y., Jendritzky, G., Staiger, H., Tinz, B., 2012. Comparison of utci to selected thermal indices. International journal of biometeorology 56, 515–535. Briegel, F., Pinto, J.G., Christen, A., 2025. Is satellite land surface temperature an appropriate proxy for intra-urban variability of daytime heat stress? Remote Sensing of Environment 331, 115045. Bröde, P., Fiala, D., Błażejczyk, K., Holmér, I., Jendritzky, G., Kampmann, B., Tinz, B., Havenith, G., 2012. Deriving the operational procedure for the universal thermal climate index (utci). International journal of biometeorology 56, 481–494. Cai, C., Li, B., Zhang, Q., Wang, X., Biljecki, F., Herthogs, P., 2025. Bi-directional mapping of morphology metrics and 3d city blocks for enhanced characterisation and generation of urban form. Sustainable Cities and Society , 106441. Cawley, G.C., Talbot, N.L., 2008. Efficient approximate leave-one-out cross-validation for kernel logistic regression. Machine Learning 71, 243–264.
32
Chang, Y., Xiao, J., Li, X., Weng, Q., 2023. Monitoring diurnal dynamics of surface urban heat island for urban agglomerations using ecostress land surface temperature observations. Sustainable Cities and Society 98, 104833. Chen, G., Hua, J., Shi, Y., Ren, C., 2023. Constructing air temperature and relative humidity-based hourly thermal comfort dataset for a high-density city using machine learning. Urban Climate 47, 101400. Chen, T., Guestrin, C., 2016. Xgboost: A scalable tree boosting system, in: Proceedings of the 22nd acm sigkdd international conference on knowledge discovery and data mining, pp. 785–794. Chen, Y., Ma, W., Shao, Y., Wang, N., Yu, Z., Li, H., Hu, Q., 2024. The impacts and thresholds detection of 2d/3d urban morphology on the heat island effects at the functional zone in megacity during heatwave event. Sustainable Cities and Society , 106002. Cleveland, W.S., 1979. Robust locally weighted regression and smoothing scatterplots. Journal of the American statistical association 74, 829–836. Cunha, J., Nobrega, R.L., Rufino, I., Erasmi, S., Galvão, C., Valente, F., 2020. Surface albedo as a proxy for land-cover clearing in seasonally dry forests: Evidence from the brazilian caatinga. Remote Sensing of Environment 238, 111250. Eliasson, I., 1996. Urban nocturnal temperatures, street geometry and land use. Atmospheric environment 30, 379–392. Ermida, S.L., Soares, P., Mantas, V., Göttsche, F.M., Trigo, I.F., 2020. Google earth engine open-source code for land surface temperature estimation from the landsat series. Remote Sensing 12, 1471. ESA, 2020. Copernicus dem glo-30 and glo-90. European Space Agency. URL: https://doi.org/ 10.5270/ESA-c5d3d65, doi:10.5270/ESA-c5d3d65. accessed: 2025-08-26. Fiala, D., Havenith, G., Bröde, P., Kampmann, B., Jendritzky, G., 2012. Utci-fiala multi-node model of human heat transfer and temperature regulation. International journal of biometeorology 56, 429–441. Foga, S., Scaramuzza, P.L., Guo, S., Zhu, Z., Dilley Jr, R.D., Beckmann, T., Schmidt, G.L., Dwyer, J.L., Hughes, M.J., Laue, B., 2017. Cloud detection algorithm comparison and validation for operational landsat data products. Remote sensing of environment 194, 379–390. Fouedjio, F., Arya, E., 2024. Locally varying geostatistical machine learning for spatial prediction. Artificial Intelligence in Geosciences 5, 100081. 33
Gao, M., Li, H., Chen, F., Zhou, M., Yang, G., Zhu, D., Han, D., Li, Z., 2025. A novel WRF+AutoML framework for enhanced heat estimation in urban environments. Sustainable Cities and Society 134, 106908. doi:10.1016/j.scs.2025.106908. Grekousis, G., 2025. Geographical-xgboost: a new ensemble model for spatially local regression based on gradient-boosted trees. Journal of Geographical Systems , 1–27. Hamwey, R.M., 2007. Active amplification of the terrestrial albedo to mitigate climate change: an exploratory study. Mitigation and Adaptation Strategies for Global Change 12, 419–439. Han, J., Chong, A., Lim, J., Ramasamy, S., Wong, N.H., Biljecki, F., 2024. Microclimate spatiotemporal prediction using deep learning and land use data. Building and Environment 253, 111358. Hong, C., Dong, Z., Qu, Z., Zhang, C., Qian, J., Zhang, Y., Li, X., Li, C., Wang, Z., Gu, Z., 2026. Improving WRF-BEP+BEM performance in simulating wind-temperature-humidity in high-density urban areas: A case study of a megacity. Building and Environment 288, 113992. doi:10.1016/j.enbuild.2025.113992. Huang, B., Wu, B., Barry, M., 2010. Geographically and temporally weighted regression for modeling spatio-temporal variation in house prices. International journal of geographical information science 24, 383–401. Jacobson, M.Z., Ten Hoeve, J.E., 2012. Effects of urban surfaces and white roofs on global and regional climate. Journal of climate 25, 1028–1044. Jendritzky, G., de Dear, R., Havenith, G., 2012. Utci—why another thermal index? International journal of biometeorology 56, 421–428. Kousis, I., Manni, M., Pisello, A., 2022. Environmental mobile monitoring of urban microclimates: A review. Renewable and Sustainable Energy Reviews 169, 112847. Li, H., Zhao, Y., Wang, C., Ürge-Vorsatz, D., Carmeliet, J., Bardhan, R., 2024. Cooling efficacy of trees across cities is determined by background climate, urban morphology, and tree trait. Communications Earth & Environment 5, 754. Li, K., Zeng, H., 2024. Multidisciplinary parameters for characterizing the 3d urban morphology: An overview based on the relational perspective. Sustainable Cities and Society , 105364. Li, S., Zhao, Z., Miaomiao, X., Wang, Y., 2010. Investigating spatial non-stationary and scaledependent relationships between urban surface temperature and environmental factors using geographically weighted regression. Environmental Modelling & Software 25, 1789–1800. 34
Li, X., Ratti, C., 2018. Mapping the spatial distribution of shade provision of street trees in boston using google street view panoramas. Urban Forestry & Urban Greening 31, 109–119. Li, X., Wang, G., 2021. Gpu parallel computing for mapping urban outdoor heat exposure. Theoretical and Applied Climatology 145, 1101–1111. Lindberg, F., Grimmond, C., 2010. Continuous sky view factor maps from high resolution urban digital elevation models. Climate Research 42, 177–183. Lindberg, F., Grimmond, C., 2011. The influence of vegetation and building morphology on shadow patterns and mean radiant temperatures in urban areas: model development and evaluation. Theoretical and applied climatology 105, 311–323. Lindberg, F., Grimmond, C.S.B., Gabey, A., Huang, B., Kent, C.W., Sun, T., Theeuwes, N.E., Järvi, L., Ward, H.C., Capel-Timms, I., et al., 2018. Urban multi-scale environmental predictor (umep): An integrated tool for city-based climate services. Environmental modelling & software 99, 70–87. Lindberg, F., Holmer, B., Thorsson, S., 2008. Solweig 1.0–modelling spatial variations of 3d radiant fluxes and mean radiant temperature in complex urban settings. International journal of biometeorology 52, 697–713. Liu, L., Lin, Y., Liu, J., Wang, L., Wang, D., Shui, T., Chen, X., Wu, Q., 2017. Analysis of local-scale urban heat island characteristics using an integrated method of mobile measurement and gis-based spatial interpolation. Building and Environment 117, 191–207. Liu, P., Lei, B., Huang, W., Biljecki, F., Wang, Y., Li, S., Stouffs, R., 2025. Sensing climate justice: A multi-hyper graph approach for classifying urban heat and flood vulnerability through street view imagery. Sustainable Cities and Society 118, 106016. Liu, P., Wang, Y., De Sabbata, S., Lei, B., Biljecki, F., Tang, J., Stouffs, R., 2026. Living upon networks: A heterogeneous graph neural embedding integrating waterway and street systems for urban form understanding. Environment and Planning B: Urban Analytics and City Science 53, 453–469. Matzarakis, A., Rutz, F., Mayer, H., 2007. Modelling radiation fluxes in simple and complex environments—application of the rayman model. International journal of biometeorology 51, 323–334. McGarigal, K., 2015. Fragstats help. University of Massachusetts: Amherst, MA, USA 182.
35
McGarigal, K., Cushman, S.A., Neel, M.C., Ene, E., et al., 2002. Fragstats: spatial pattern analysis program for categorical maps. Computer software program produced by the authors at the University of Massachusetts, Amherst. Available at the following web site: www. umass. edu/landeco/research/fragstats/fragstats. html 6. Meteorological Service Singapore, 2020. Annual climate assessment report 2019. https://www. weather.gov.sg/climate-annual-climate-reports/. Accessed: 2025-09-20. Meteorological Service Singapore, 2023. Climate of singapore. https://www.weather.gov.sg/ climate-climate-of-singapore/. Accessed: 2025-09-20. Moran, G.W., 1984. Locally-Weighted-Regression Scatter-Plot Smoothing (LOWESS): a graphical exploratory data analysis technique. Ph.D. thesis. Monterey, California. Naval Postgraduate School. Muse, N., Clement, A., Mach, K.J., 2024. Daytime land surface temperature and its limits as a proxy for surface air temperature in a subtropical, seasonally wet region. PLOS climate 3, e0000278. Nazarian, N., Krayenhoff, E., Bechtel, B., Hondula, D., Paolini, R., Vanos, J., Cheung, T., Chow, W., de Dear, R., Jay, O., et al., 2022. Integrated assessment of urban overheating impacts on human life. Earth’s Future 10, e2022EF002682. Oke, T.R., 2002. Boundary layer climates. Routledge. Pang, Y., Wang, Y., Lai, X., Zhang, S., Liang, P., Song, X., 2023. Enhanced kriging leave-one-out cross-validation in improving model estimation and optimization. Computer Methods in Applied Mechanics and Engineering 414, 116194. Park, S., Tuller, S.E., Jo, M., 2014. Application of universal thermal climate index (utci) for microclimatic analysis in urban thermal environments. Landscape and Urban Planning 125, 146–155. Pei, W., Stouffs, R., 2025. Parametric archetype: An incremental learning model based on a similarity measure for building material stock aggregation. Automation in Construction 172, 106064. Prado, R.T.A., Ferreira, F.L., 2005. Measurement of albedo and analysis of its influence the surface temperature of building roof materials. Energy and Buildings 37, 295–300. Qi, J., Xiong, W., Li, J., Zheng, J., Ye, Q., Hu, M., 2025. Investigating non-linear and synergistic effects of urban functional zone morphology on land surface temperature. Sustainable Cities and Society , 106546. 36
Ramsay, E.E., Wang, Y., Masoudi, M., Chai, M.W., Yin, T., Hamel, P., 2025. Assessing a decisionsupport tool to estimate the cooling potential and economic savings from urban vegetation in singapore. Sustainable Cities and Society 125, 106337. Ren, Z., Zhou, S., Liu, D., Liu, Q., 2025. Physics-informed neural networks: A review of methodological evolution, theoretical foundations, and interdisciplinary frontiers toward next-generation scientific computing. Applied Sciences 15, 8092. Robinson, D., 2006. Urban morphology and indicators of radiation availability. Solar Energy 80, 1643–1648. Roth, M., Lim, V.H., 2017. Evaluation of canopy-layer air and mean radiant temperature simulations by a microclimate model over a tropical residential neighbourhood. Building and Environment 112, 177–189. Schneider, F.A., Ortiz, J.C., Vanos, J.K., Sailor, D.J., Middel, A., 2023. Evidence-based guidance on reflective pavement for urban heat mitigation in arizona. Nature communications 14, 1467. Schrijvers, P., Jonker, H., De Roode, S., Kenjereš, S., 2016. The effect of using a high-albedo material on the universal temperature climate index within a street canyon. Urban Climate 17, 284–303. Sengupta, M., Xie, Y., Lopez, A., Habte, A., Maclaurin, G., Shelby, J., 2018. The national solar radiation data base (nsrdb). Renewable and sustainable energy reviews 89, 51–60. Shaeri, P., AlKhaled, S., Middel, A., 2025. A multimodal physics-informed neural network approach for mean radiant temperature modeling. arXiv preprint arXiv:2503.08482 . Shashua-Bar, L., Pearlmutter, D., Erell, E., 2009. The cooling efficiency of urban landscape strategies in a hot dry climate. Landscape and urban planning 92, 179–186. Shi, Y., Ren, C., Cai, M., Lau, K.K.L., Lee, T.C., Wong, W.K., 2019. Assessing spatial variability of extreme hot weather conditions in hong kong: A land use regression approach. Environmental research 171, 403–415. Singapore Department of Statistics, 2020. Resident population by planning area/subzone of residence, ethnic group and sex. URL: https://data.gov.sg/. census of Population 2020. Small, C., 2006. Comparative analysis of urban reflectance and surface temperature. Remote Sensing of Environment 104, 168–189.
37
Sun, Y., Deng, K., Ren, K., Liu, J., Deng, C., Jin, Y., 2024. Deep learning in statistical downscaling for deriving high spatial resolution gridded meteorological data: A systematic review. ISPRS Journal of Photogrammetry and Remote Sensing 208, 14–38. Tanoori, G., Soltani, A., Modiri, A., 2024. Machine learning for urban heat island (uhi) analysis: Predicting land surface temperature (lst) in urban environments. Urban Climate 55, 101962. Tolan, J., Yang, H.I., Nosarzewski, B., Couairon, G., Vo, H.V., Brandt, J., Spore, J., Majumdar, S., Haziza, D., Amaral, J.V., Moutonannini, T., Bojanowski, P., Johns, T., White, B., Tiecke, T., Couprie, C., 2024. Very high resolution canopy height maps from rgb imagery using selfsupervised vision transformer and convolutional decoder trained on aerial lidar. Remote Sensing of Environment 300, 113888. doi:10.1016/j.rse.2023.113888. Tuholske, C., Caylor, K., Funk, C., Verdin, A., Sweeney, S., Grace, K., Peterson, P., Evans, T., 2021. Global urban population exposure to extreme heat. Proceedings of the National Academy of Sciences 118, e2024792118. U.S. Geological Survey, 2025. Landsat 8 mission — u.s. geological survey. URL: https://www.usgs. gov/landsat-missions/landsat-8. Wan, Y., Du, H., Yuan, L., Xu, X., Tang, H., Zhang, J., 2025. Exploring the influence of block environmental characteristics on land surface temperature and its spatial heterogeneity for a high-density city. Sustainable Cities and Society 118, 105973. Wang, H., Huang, Z., Yin, G., Bao, Y., Zhou, X., Gao, Y., 2022. Gwrboost: A geographically weighted gradient boosting method for explainable quantification of spatially-varying relationships. arXiv preprint arXiv:2212.05814 . Wang, H., Yi, T., Lu, Y., Wang, Y., Wu, J., 2025. Patterns of nighttime surface urban heat island patch in mega urban agglomerations: a case study in the pearl river delta, china. Sustainable Cities and Society 128, 106465. Wang, S., Gao, K., Zhang, L., Yu, B., Easa, S.M., 2024a. Geographically weighted machine learning for modeling spatial heterogeneity in traffic crash frequency and determinants in us. Accident Analysis & Prevention 199, 107528. Wang, Y., Wang, H., Yao, F., Stouffs, R., Wu, J., 2024b. An integrated framework for jointly assessing spatiotemporal dynamics of surface urban heat island intensity and footprint: China, 2003–2020. Sustainable Cities and Society , 105601.
38
Wang, Y., Zhao, Y., Wu, J., 2020. Dynamic monitoring of long time series of ecological quality in urban agglomerations using google earth engine cloud computing: A case study of the guangdonghong kong-macao greater bay area, china. Acta Ecol. Sin 40, 8461–8473. Wen, Z., Zhuo, L., Gao, M., Han, D., 2025. How can we improve data integration to enhance urban air temperature estimations? International Journal of Applied Earth Observation and Geoinformation 140, 104599. Wong, T.T., 2015. Performance evaluation of classification algorithms by k-fold and leave-one-out cross validation. Pattern recognition 48, 2839–2846. Wongsai, S., Wanishsakpong, W., Suwanprasit, C., Wongsai, N., 2024. Spatial autoregressive regression analysis of surface urban heat island intensity in the tropical industrial city of rayong, thailand. Urban Climate 55, 101980. Xu, D., Wang, Y., Zhou, D., Wang, Y., Zhang, Q., Yang, Y., 2024. Influences of urban spatial factors on surface urban heat island effect and its spatial heterogeneity: A case study of xi’an. Building and Environment 248, 111072. Yang, W., Deng, M., Tang, J., Luo, L., 2023. Geographically weighted regression with the integration of machine learning for spatial prediction. Journal of Geographical Systems 25, 213–236. Yang, Z., Peng, J., Jiang, S., Yu, X., Dong, J., Corcoran, J., 2025. Human-centered urban heat island intensity in global megacities: Nonlinear response patterns and region-specific thresholds. Environmental Science & Technology 59, 19781–19791. Yang, Z., Peng, J., Jiang, S., Yu, X., Hu, T., 2024. Optimizing building spatial morphology to alleviate human thermal stress. Sustainable Cities and Society 106, 105386. Yao, X., Zeng, X., Zhu, Z., Lan, Y., Shen, Y., Liu, Q., Yang, F., 2023. Exploring the diurnal variations of the driving factors affecting block-based lst in a “furnace city” using ecostress thermal imaging. Sustainable Cities and Society 98, 104841. Yi, S., Li, X., Li, D., Dong, X., Wang, R., Xu, Q., 2025a. Hyperlocal heat stress around bus stops in philadelphia: Insights from spatio-temporal microclimate modeling and explainable ai. Computers, Environment and Urban Systems 122, 102341. Yi, S., Li, X., Liu, Y., Dong, X., Tu, W., 2025b. A sub-meter resolution urban surface albedo dataset for 34 us cities based on deep learning. Scientific Data 12, 789.
39
Yi, S., Li, X., Ma, C., Wang, R., Zhou, Y., Xu, Q., Zhao, T., 2025c. Assessing the differential impact of vegetated and built-up areas on heat exposure environment: A case study of los angeles. Building and Environment 271, 112538. Yi, S., Li, X., Tu, W., Zhao, T., 2025d. Planning for cooler cities: A multimodal ai framework for hyperlocal spatio-temporal urban heat stress prediction and mitigation. Urban Forestry & Urban Greening , 129101. Yi, S., Li, X., Wang, R., Guo, Z., Dong, X., Liu, Y., Xu, Q., 2024. Interpretable spatial machine learning insights into urban sanitation challenges: A case study of human feces distribution in san francisco. Sustainable Cities and Society 113, 105695. Zaerpour, M., Papalexiou, S.M., Pietroniro, A., 2025. Increasing tree canopy lowers urban air temperature by up to 1.5° c in heat-prone areas. npj Urban Sustainability 5, 92. Zakšek, K., Oštir, K., Kokalj, Ž., 2011. Sky-view factor as a relief visualization technique. Remote sensing 3, 398–415. Zhan, W., Bechtel, B., Du, H., Chakraborty, T., Kotthaus, S., Krayenhoff, E.S., Martilli, A., Naserikia, M., Nazarian, N., Roth, M., et al., 2025. Satellite-derived land surface temperatures strongly mischaracterise urban heat hazard. arXiv preprint arXiv:2509.16568 . Zhang, J., Yao, X., Chen, Y., Lin, M., Lin, T., Zheng, Y., Geng, H., Zheng, Y., Wu, X., Zhang, G., et al., 2024a. Spatiotemporal dynamic mapping of heat exposure risk for different populations in city based on hourly multi-source data. Sustainable Cities and Society 107, 105454. Zhang, L., Yuan, C., 2023. Multi-scale climate-sensitive planning framework to mitigate urban heat island effect: A case study in singapore. Urban Climate 49, 101451. Zhang, Q., Yang, J., Ma, X., Xin, J., Ren, J., Yu, W., Xiao, X., Xia, J., 2024b. Influence of 2d/3d urban morphology on diurnal land surface temperature from the perspective of functional zones. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing . Zhao, M., Lei, S., Li, W., 2026. Incorporating urban thermal comfort into tod planning: Non-linear heterogeneous built environment effects. Sustainable Cities and Society , 107254. Zhou, S., Geng, X., Zhao, J., Hei, J., Wu, T., Chen, Z., Wu, Z., 2025. An lcz-based machine learning framework for revealing spatial heterogeneity of thermal comfort in high-density areas: Enhancing explainability and fine-grid scale resolution. Sustainable Cities and Society , 106873.
40