Explore the Mitigation Mechanism of Urban Thermal Environment by Integrating Geographic Detector and Standard Deviation Ellipse (SDE)

Abstract

The urban surface temperature is a complex integrated natural-human geographic phenomena; with the development of geostatistical methods and the application of multisource data, its research has gradually shifted from a single perspective to a study that integrates multiple factors such as nature and humanity. However, based on the context of the integration of natural and human factors and mutual constraints of each factor, the research on the mechanism of influence on urban habitat thermal environment needs to be further deepened. Therefore, this paper explores the spatial and temporal heterogeneity of urban surface temperature in Zhengzhou City during the summer of 2013–2020 from the perspective of multi-source data fusion, and uses the Geodetector model to quantitatively reveal the main influencing factors of urban surface temperature and the impact of superimposed factors on the compound effect of surface temperature. The results show that: (1) the urban thermal environment in the central of Zhengzhou city (region within the first ring) is obvious, and it is mainly concentrated in commercial and densely populated areas. (2) According to trend analysis, the northwest-southeast direction of the city continues to increase in temperature from 2013–2020, coupled with the direction of urban development. (3) Among the factors affecting urban surface temperature, normalized difference built-up index (NDBI), normalized difference vegetation index (NDVI), tasseled cap wetness (TCW), and human elements are particularly typical. NDVI and TCW are strongly negatively correlated with the urban thermal environment, while NDBI and human elements are strongly positively correlated. (4) Mitigation of the urban thermal environment can start with the interaction mechanism of positive and negative factors. This study provides new ideas for the mechanism analysis of spatial and temporal evolution patterns of the urban thermal environment under multifactorial constraints, and provides suggestions and decisions for promoting green and sustainable urban development.

1. Introduction

The urban social activities, industrial production, and other human activities consume a large amount of fuel. Moreover, due to inefficient utilization, disorderly management, and urban sprawl, urban habitat problems are becoming more prominent [1,2,3,4], which accelerates the change in surface temperature within cities. The change in urban surface temperature will cause some hazards, such as prolonged high temperature affecting human comfort, increasing the frequency of high temperature heat waves, and intensifying urban air pollution due to heat island circulation. Especially in recent years, the phenomenon of “high temperature” (with summer temperature extremes that are maintained at 35 °C) has been appearing increasingly frequently in urban areas, which has an important impact on building a healthy, comfortable, livable, and sustainable urban habitat.

Currently, the urban habitat thermal environment effect has become one of the hot spots for research in geography, environmental science, urban planning, architecture, and other disciplines. At the beginning of the study, some scholars mainly used meteorological actual measurement data for the study of the urban heat island effect. In 1972, Rao used the thermal infrared remote sensing monitoring technique for the first time to study the urban thermal environment effect [5,6]. Since then, due to the development of remote sensing satellite technology, the satellite thermal infrared band has been more widely adopted and applied in urban surface temperature studies, and has occupied an important position in urban thermal environment investigation studies. In terms of research content, the research on urban thermal environment mainly includes remote sensing monitoring, formation mechanism, simulation prediction, and ecological safety evaluation of urban thermal environment. For example, there are some studies on the relationship between land use change and urban thermal environment [7,8,9]; urban thermal environment studies based on land surface temperature and socio-ecological variables, and other studies that only list the effects of a single natural or a single human-economic perspective on urban thermal environment [10,11,12]; and at the architectural level, studies in which the main focus is on the planning of urban ventilation corridors and how to mitigate heat islands [13,14].

Recent studies have focused on surface temperature as a quantitative indicator of the impact of urban thermal environment (UTE), and its influence mechanisms [15,16]. Previous studies have shown that surface temperature is closely related to natural variables [17,18,19], but most studies have used only a set of remote sensing indices in revealing the drivers of UTE. However, these indices do not seem to be able to fully explain such a complex geographical phenomenon as the urban thermal environment, as some human activities that do not alter the physical properties of the surface may also contribute to the formation of UTE. Indeed, several empirical studies have shown that population density and social economics have a strong impact on the formation of UTE [20,21]. For example, with the same conditions, a unit area of land with different population and road densities (PD, RD) will release different amounts of heat. Scholars at home and abroad have conducted a lot of research on the spatial and temporal divergence patterns and influence mechanism of the urban habitat thermal environment, and have achieved more fruitful results. However, the study of using Geodetector models to explore the mechanism of the influence of natural and human factors on surface temperature and how to overcome the limitations of traditional mathematical statistical models of homogeneity and normality to detect the influence of the two factors after their interaction on surface temperature needs to be further deepened.

In summary, we use remote sensing images, spatial big data, and socio-economic data as the basic data sources with the support of RS and GIS technologies to analyze the spatial and temporal patterns of urban habitat thermal environment using interannual variation trends. Meanwhile, standard deviation ellipses and geographic detectors are combined to explore the influencing factors of the thermal environment effects of urban habitats in Zhengzhou. This will lead to a better understanding of the evolution patterns, intrinsic driving mechanisms, and development patterns of the urban thermal environment, and provide a theoretical basis and decision support for urban planning, urban construction, and the improvement of the quality of the habitat environment.

2. Materials and Methods

2.1. Study Area

Zhengzhou city is in central China (Figure 1), and the urban area lies between 112°42′ E–114°14′ E and 34°16′ N–34°58′ N. The area is adjacent to the middle and lower reaches of the Yellow River and the northeast side of the Funiu Mountain Range, in the transition zone to the Yellow and Huai Plain; the altitude is 40–295 m. It is a north temperate continental monsoon, with four distinct seasons. The annual precipitation is 623.3 mm. July is the hottest month with an average temperature of 27.3 °C.

It consists of five administrative districts: Jinshui District, Erqi District, Guanchenghuizi District, Zhongyuan District, Huiji District. The built-up area is 1181.51 square kilometers (Figure 1) [22]. According to the seventh national census of the People’s Republic of China in 2020, the resident population of Zhengzhou is 12.6 million.

We collected monthly temperature changes from January 2013 to July 2021 in the central city (Figure 2). Monthly temperature comparisons for 2013–2020 showed that summer temperature extremes were maintained at 35 °C. With the emergence of increasingly extreme temperatures, mitigating abnormal temperature phenomena, improving the quality of urban living environment, and achieving sustainable development of urban green spaces have become urgent issues in the development of large cities.

2.2. Datasets

The main data sources used in this study include remote sensing images, temperature data from near-surface observation points, land use data, socioeconomic data and demographic data, road network data, and city boundaries, as shown in

Table 1. Data sources.

Data Type	Data Description	Data Source	Data Display
Remote Sensing Data	NPP-VIIRS (500 m resolution NTL data, 2013–2020) [23]	http://ngdc.noaa.gov/eog, accessed on 10 November 2021
	Landsat8 TIRS (100 m resolution thermal, 2013–2020, cloud cover < 5%)	http://earthexplorer.usgs.gov, accessed on 10 November 2021
	MODIS 11 L2 (1 km resolution LST product, 2013–2020)	http://lpdaac.usgs.gov, accessed on 10 November 2021
	Digital elevation (ASTER GDEM)	https://earthexplorer.usgs.gov/, accessed on 10 November 2021
Humanities and Economic Data	Road network	GaoDe Map, accessed on 10 November 2021
	Vector boundary
	Population data	National Bureau of Statistics of China (http://www.stats.gov.cn/, accessed on 10 November 2021)
	Land use data	http://www.globallandcover.com/, accessed on 10 November 2021
	Daily temperature data of the study area from 2013 to 2021	http://data.cma.cn/, accessed on 10 November 2021

.

In this paper, we mainly used Landsat satellite for the surface temperature inversion; in order to check the accuracy of the results, since there is only one for the ground observed temperature station in Zhengzhou City, which cannot provide full coverage of the land surface temperature, the existing land temperature product (MODIS) was chosen for comparative analysis (Figure 3) [23].

Synthetic inversions were performed for surface temperatures in summer (June–August), winter (December–February), and year-round (January–December) from 2013 to 2020. After resampling and comparing with the mainstream MODIS temperature products, according to the comparison, the fit was better in summer, so the summer time (June–August) period was chosen for the follow-up operation in this paper.

In order to overcome the situation of incomplete single-day single-impact imaging caused by cloud occlusion in the study area, this paper used select thermal infrared remote sensing images (Landsat8 TIRS) with less than 5% clouds. The selected time range was summer, a total of 90 days, and the time used for analysis was about 5–6 days per year. We synthesized the complete remote sensing image map of the study area using the median synthesis method, which could meet the accuracy requirements of synthetic surface temperature, and performed pre-processing such as study area cropping; all the above processing was done in the Google Earth Engine (https://code.earthengine.google.com/, accessed on 10 November 2021) platform. In addition, human and economic data were cleaned using python 3.8, and all data in ARCGIS10.2 were unified into a personal file geodatabase with a unified coordinate system of WGS84 World Mercator.

2.3. Methods

The overall workflow is shown in Figure 4. In this study, the atmospheric correction method was mainly used to invert the surface temperature of Zhengzhou City from 2013 to 2020, and the slope trend analysis method was used to analyze the temporal trend of surface temperature from 2013 to 2020. In addition, the urban–social–economic–natural composite system was analyzed from multiple perspectives by selecting indicators related to nature and humans to explore the interaction mechanism of multiple factors and provide certain suggestions for the heat island mitigation mechanism.

2.3.1. Land Surface Temperature Inversion

For Landsat remote sensing images, the atmospheric correction method ([24,25,26]) was used to calculate the black body radiance at the same temperature and the thermal infrared radiance received by the satellite sensor. The method involves three parts: atmospheric upward radiation brightness, the satellite sensor received through the atmosphere of the Earth’s true radiation brightness value, and atmospheric radiation down to the ground after the reflection of energy. Therefore, the radiant luminance value received by the sensor was given as follows:

$${\mathrm{L}}_{\mathsf{\lambda }}=\left[{\mathsf{\epsilon }\mathrm{L}}_{\mathrm{black}}+\left(1-\mathsf{\epsilon }\right){\mathrm{L}}_{\mathrm{down}}\right]\mathsf{\tau }+{\mathrm{L}}_{\mathrm{up}}$$

In the formula, $\mathsf{\epsilon }$ is the surface-to-surface radiation rate, ${\mathrm{L}}_{\mathrm{black}}$ is the radiation brightness (K) of the black body with temperature (${\mathrm{T}}_{\mathrm{s}})$, $\mathsf{\tau }$ represents atmospheric transport rate, ${\mathrm{L}}_{\mathrm{up}}$ is the atmospheric upstream radiation brightness, and ${\mathrm{L}}_{\mathrm{down}}$ is the atmospheric downstream radiation brightness. Since the radiant luminance of the thermal radiation has been obtained, the radiant luminance of a blackbody with a temperature of ${\mathrm{T}}_{\mathrm{s}}$ can be obtained ${\mathrm{L}}_{\mathrm{black}}$ according to the following equation:

$${\mathrm{L}}_{\mathrm{black}}=\frac{{\mathrm{L}}_{\mathsf{\lambda }}-{\mathrm{L}}_{\mathrm{up}}-\mathsf{\tau }\left(1-\mathsf{\epsilon }\right){\mathrm{L}}_{\mathrm{down}}}{\mathsf{\tau }\mathsf{\epsilon }}$$

Finally, from the inverse function of Planck’s law, the true surface temperature $\mathrm{LST}$ (K) can be calculated:

$$\mathrm{LST}=\frac{{\mathrm{K}}_2}{\ln \left(\frac{{\mathrm{K}}_1}{{\mathrm{L}}_{\mathrm{black}}}+1\right)}$$

${\mathrm{K}}_1$ and ${\mathrm{K}}_2$ are constants that can be obtained from a header file. For the Landsat 8 thermal infrared band 10, ${\mathrm{L}}_{\mathrm{up}}$ is 1.25 $\mathrm{W}/\left({\mathrm{m}}^2·\mathrm{sr}·\mathsf{\mu }\mathrm{m}\right)$, ${ \mathrm{L}}_{\mathrm{down}}$ is 2.37 $\mathrm{W}/\left({\mathrm{m}}^2·\mathrm{sr}·\mathsf{\mu }\mathrm{m}\right)$, ${\mathrm{K}}_1$ is 774.89 K, and ${\mathrm{K}}_2$ is 121.08 K.

2.3.2. Research Index Features

We discussed its multi-factor influence and interaction mechanism, combined with the existing research. The research factors are the following indices (

Table 2. Factors and their calculation formulas.

Layer		Factors	Formula
Natural Factors Indicators	Greenness	NDVI	$\mathrm{NDVI}=\frac{\mathrm{NIR}-\mathrm{RED}}{\mathrm{NIR}+\mathrm{RED}}$
		SAVI	$\mathrm{SAVI}=\frac{\left(\mathrm{NIR}-\mathrm{Red}\right)\times \left(1+\mathrm{L}\right)}{\mathrm{NIR}+\mathrm{Red}+\mathrm{L}}$
	Imperviousness	NDBI	$\mathrm{NDBI}=\frac{\mathrm{SWIR}-\mathrm{NIR}}{\mathrm{SWIR}+\mathrm{NIR}}$
		IBI	$\mathrm{IBI}=\frac{\left[\mathrm{NDBI}-\left(\mathrm{SAVI}+\mathrm{MNDWI}\right)/2\right]}{\left[\mathrm{NDBI}+\left(\mathrm{SAVI}+\mathrm{MNDWI}\right)/2\right]}$
	Wetness	MNDWI	$\mathrm{MNDWI}=\frac{\mathrm{Green}-\mathrm{SWIR}1}{\mathrm{Green}+\mathrm{SWIR}1}$
		TCW	$\mathrm{TCW}=0.1511\mathrm{Blue}+0.1973\mathrm{Green}+0.3283\mathrm{Red}+0.3407\mathrm{NIR}+\left(-0.7117\mathrm{SWIR}1\right)+\left(-0.4559\mathrm{SWIR}2\right)$
Humanities and Economic Factor Indicators		NLI	$\mathrm{NLI}=\ln \left({\mathrm{VIIRS}}_i\right)$
		PD	$\mathrm{PD}=\frac{{\mathrm{N}}_{\mathrm{pop}}}{{\mathrm{Area}}_{\mathrm{circle}}}$
		RD	$\mathrm{RD}=\frac{{\mathrm{R}}_{\mathrm{length}}}{{\mathrm{Area}}_{\mathrm{circle}}}$

), and they fall into the following categories: (1) natural indices: green layer (NDVI, SAVI), wetness layer (MNDWI, TCW) and impervious layer (NDBI, IBI); and (2) Humanities and Economic indices: night light index (NLI), population density (PD), and road density (RD).

The normalized difference vegetation index (NDVI) is an important index for monitoring the land surface vegetation [27]. The wetness component (TCW) is used to detect the moisture in soil and vegetation and the structure of vegetation [28]. The soil adjusted vegetation index (SAVI) can effectively reduce the image of topsoil and improve the linear relationship between the vegetation index and leaf area index [29]. The modified normalized difference water index (MNDWI) refers to the normalized difference of processing image bands containing information about water bodies [30]. Therefore, we used this wetness component as an indicator of soil moisture [31]. NDBI represents the normalized building index [32]. The index-based construction land index (IBI) uses an index band derived from a multi-spectral band and is a ratio-based normalized difference index. It can effectively enhance the construction land information in the image and then extract a more accurate construction land index [33], where Green, Red, NIR, SWIR1, SWIR2, and TIRS1 correspond to bands 3, 4, 5, 6, 7, and 10 of Landsat 8.

This paper selected three elements—economic development, urban construction, and population distribution—for comprehensive analysis. The economic level is a comprehensive reflection of the urban population, industrial structure, geographical space, etc.; the average intensity of nighttime lights is a comprehensive reflection of the results of the interaction of these factors. Based on the above considerations, this paper selected the regional nighttime light (NIL) distribution to build a socio-economic index reflecting its economic development level [34]. Roads are a good carrier in measuring urban land expansion and urban expansion; using the length of roads in the region to calculate the road density in this region can effectively reflect the degree of urban expansion [35]. The degree of population aggregation can also identify the hotspot areas of the city and respond to the more developed areas within the city, so the population density was also selected as a measure [36].

2.3.3. Standard Deviation Ellipse (SDE)

The directional distribution (the standard deviation ellipse), proposed by Lafever [37], which can reveal the spatial distribution of various economic phenomena using the spatial statistical method and the spatial distribution characteristics of geographical elements, is widely used in geography. The standard deviation ellipse describes the distribution of the important cases of time and space, quantitatively based on the ellipse center, long axis, short axis, and azimuth angle to confirm the development direction of the leading trend of economic activities. In this study, we extracted the heat-island regions in 2013–2020 to map the heat-island trends, in order to analyze the evolution of the standard deviation ellipse of the thermal environment in the study area over four periods, and to identify the direction of urban heat island development and migration trajectories.

2.3.4. Trend Analysis Based on Image Element Scale

Based on the time series inversion of surface temperature, using a one-dimensional linear regression equation to simulate the surface temperature trend in the central of Zhengzhou city from 2013 to 2020, the calculation formula is as follows.

$$Slope=\frac{n{{\displaystyle \sum }}_{i=1}^n\left(i\times c_i\right)-{{\displaystyle \sum }}_{i=1}^{n}i\times {{\displaystyle \sum }}_{i=1}^n{c}_i}{n{{\displaystyle \sum }}_{i=1}^n{i}^2-{\left({{\displaystyle \sum }}_{i=1}^{n}i\right)}^2}$$

where $Slope$ is the slope of the image element regression equation, and a positive or negative value represents the interannual variation linear trend of the image element (when Slope > 0, it means that the image element has an increasing trend under the time series; when Slope = 0, it means that the image element is basically constant under the time series; when Slope < 0, it means that the image element has a decreasing trend under the time series). $c_i$ is the mean value of the image element in $i$ year, and n is the length of the studied time series.

The significance of the trend was tested using the F-test, and the significance represents the degree of confidence in the trend of change. The formula is as follows:

$$F=U\times \frac{n-2}{\mathrm{Q}}$$

where $U={\displaystyle \sum }_{i=1}^n{\left({\widehat{y}}_i-\overline{y}\right)}^2$ is called the error sum of squares; $Q={\displaystyle \sum }_{i=1}^n{\left(y_i-{\widehat{y}}_i\right)}^2$ is called the regression sum of squares; $y_i$ is the image element value in i year; ${\widehat{y}}_i$ is its regression value; y is the average value of image elements in the monitoring period; and n is the number of testing years.

The trends are classified into four levels according to the test results, significant decrease ($Slope$ < 0, $F$ > 4.6), non-significant decrease ($Slope$ < 0, $F$ < 4.6), non-significant increase ($Slope$ > 0, $F$ < 4.6), and a significant increase ($Slope$ > 0, $F$ > 4.6).

2.3.5. Exploration of Spatial Divergence Pattern Factors—Geographic Detector

The geographic detector (Geodetector) is a new statistical method for detecting spatial heterogeneity and revealing the driving factors behind it. The most striking advantage of the Geodetector model is that it does not require any linearity assumptions to detect the relationship between drivers and geographical phenomena, and its calculation process and results are not affected by multivariate covariance, with elegant form and clear physical implications [38]. The basic idea is that the study area is assumed to be divided into subregions, and if the sum of the variances of the sub-regions is smaller than the total variance of the area, then there is spatial heterogeneity; conversely, if the spatial distribution of the two variables converge, then there is statistical correlation between them. The results of interactive detection, ecological detection, factor detection, and risk detection in geographic detectors were analyzed to identify the relationship between the independence and interaction of factors, the relationship between interaction enhancement and attenuation, the linear and non-linear interaction, and to better detect the deep-seated causes of the heat-island phenomenon and eliminate the overlapping factors for a better understanding of the LST drive mechanism.

3. Analysis and Results

3.1. Urban Surface Temperature Classification Method

In order to visually describe the variation of surface temperature in the study area, this paper uses the standard deviation method in GIS software to divide the surface temperature into: highest temperature region (Hest-T), higher temperature region (Her-T), high temperature region (H-T), medium temperature region (M-T), low temperature region (L-T), low temperature region (Ler-T), and lowest temperature region (Lest-T). The specific grading is shown in

Table 3. Surface temperature class classification.

Temperature Level	Lowest Temperature Region	Lower Temperature Region	Low Temperature Region	Medium Temperature Region	High Temperature Region	Higher Temperature Region	Highest Temperature Region
Temperature interval	T < u − 2.5 std	u − 2.5 std ≤ T < u − 1.5 std	u − 1.5 std ≤ T < u − 0.5 std	u − 0.5 std ≤ T < u + 0.5 std	u + 0.5 std ≤ T < u + 1.5 std	u + 1.5 std ≤ T < u + 2.5 std	T ≥ u + 2.5 std

. According to the heat island classification rules, M-T, H-T, Her-T and H-est T regions are considered as urban heat island regions.

3.2. Spatial Distribution Characteristics of Urban Thermal Environment

As can be seen in Figure 5, in 2013–2020, the UTE of the research area formed a core-periphery structure in space; the surface temperature shows a gradual change from the center to the periphery. The major core areas have not changed significantly. The high temperature area of Zhengzhou was mainly concentrated in the old city within the Second Ring Road. The area has gradually evolved from the cluster to the surrounding area, and the range of the heat island has been expanding. Analysis of urban surface temperature inversion maps for these four periods from 2000 to 2020 reveals that the M-T, H-T, Her-T and Hest-T areas first expand along the main urban growth area, and then multiple heat island cores are added at the periphery following the growth of the city. The main areas in the city, such as arable land, rivers, urban water environment, mountainous areas, and the surrounding areas of the city boundary, show a “cold” phenomenon, especially in the areas along the Yellow River in the north and the riverbank in the city, etc. Their spatial features of the “cold” phenomenon are more obvious, and the changes of space position are small.

Overall, the main high temperature area in the research area presents a strong spatial diversity. While the “cold” region in 2013 is mainly concentrated in the suburbs of the city and along the Yellow River, with a minimum temperature of 23.73 °C and a difference of 19.39 °C between the high and low temperatures, the urban thermal environment is more prominent. Throughout 2013–2020, the spatial location of the cold region has changed little; it remains in the urban water body and the Yellow River area in the north of the city, and the mountainous area in the southwest. Moreover, the temperature difference increased from 19.39 °C in 2013 to 26.17 °C in 2020.

With the opening of the Zhengzhou subway system in December 2013, Zhengzhou entered the stage of rapid expansion. In 2015 (Figure 5(b-1–b-3)), 2018 (Figure 5(c-1–c-3)), and 2020 (Figure 5(d-1–d-3)), the thermal region expanded from the second ring (Old City) to the third ring road, after the large-scale urban construction, transformation of the old city, and other man-made activities. From the perspective of surface temperature performance, the location of the old city center is still a hot region, and the changes are small. From 2013 to 2020, the highest temperature increased from 43.1 °C to 44.23 °C, while the lowest temperature dropped from 23.7 °C to 18.17 °C. The temperature difference was further expanded, and the thermal environment effect was continuously strengthened. The temperature difference between the central city and the suburbs was further expanded, reaching 26.17 °C. In terms of the spatial distribution of the heat island, the standard deviation of the heat island area between regions decreased from 39.94 in 2013 to 34.15 in 2020, and Zhengzhou started to move from a single core to multi-core spatial differentiation pattern.

Regarding area change (Figure 6), the overall thermal environment area did not change much from 2013 to 2020. For the whole picture, the thermal environment area between the third and fourth rings was at a high level, while the area outside the fourth ring, due to the expanding of urban construction, was gradually increasing in recent years. The thermal environment area inside the first and second rings was more stable. This also confirms the “core-periphery” structure of the thermal environment. which means that the thermal environment in the core city of Zhengzhou was relatively stable, while the area of the peripheral thermal environment was gradually expanding. In addition, according to the changes of the area of H-T, Her-T, and Hest-T within each loop type in the city, it is known that the area of the Hest-T type had less changes between 2000 and 2020, and was in a more stable state. As for the Her-T type, the distribution was more diffuse in 2013, but with the construction of the city, the area range of the Her-T type tended to be in a stable state from 2015 to 2020. Moreover, with the construction and development of the city in recent years, the area of the H-T type was gradually expanding in areas near the fourth ring, the economic and technological development zones, and commercial areas.

3.3. Temporal Evolution Mechanism of the Urban Thermal Environment

Interannual variation trend line analysis was used to linearly fit the surface temperature from 2013 to 2020 on an image-by-image basis; as shown in Figure 7, the temperature of the central region did not change significantly on annual monitoring. With the development and expansion of the city, new thermal zones are gradually added around the city, mainly concentrated in the northwestern and northern regions, while some hot cores also appear sporadically scattered in the southeast. After 2018, the surface temperature shows a gradual growth trend in the southeast direction. Before that, urban construction had not been fully developed; however, with the construction of Zhengzhou (national) Airport Economic Comprehensive Pilot Zone in 2018, the urban hot region pattern gradually extended to the southeast. The spatial location of the overall hot region in Zhengzhou city from 2013 to 2020 basically remained unchanged; with the development of urban infrastructure and the transformation of the land surface by human activities, the spatial differentiation pattern of the large core is further expanded based on the original one, and the urban surface temperature is also increasing year by year, with the highest temperature rising to 44.3 °C. The high-temperature area is mainly around the built-up urban area, commercial area, and industrial area, while the low-temperature area is still mainly concentrated along the Yellow River and the southwest mountainous area.

According to the heat island development direction and mass center migration map of Zhengzhou City, from 2013 to 2020, the main axis of Zhengzhou UTE still developed in the northwest-south direction. In 2013, the standard deviation eccentricity of the urban heat island was 0.7304, and the ellipse eccentricity was maintained at the same level from 2013 to 2021, indicating that its distribution range and area are roughly the same, and the development intensity has been maintained in the northwest-southeast direction along the urban spine, in line with the main line of the emerging urban hotspots. It can be seen in Figure 8 that the center gravity of UTE in Zhengzhou deflected 7.3° to the west-northwest in 2013–2015, with a distance shift of 2.114 km. The center of the thermal environment in Zhengzhou deflected 43.3° to the east-northwest in 2015–2018, with a distance of 2.561 km. Regarding the thermal environment in 2018–2020, the center deflects 36.9° to the northeast again, moving a distance of 2.291 km. Additionally, the angle and distance of the second deflection was the largest in the study years 2013–2020, indicating that in 2015–2018, the urban space of Zhengzhou City turned to the northeast, moving 2.561 km compared to 2015. The movement was within the framework of the two-axis strategic development, along the fourth ring road, forming a north–south axis of development from the main urban area to the aviation city, becoming the main skeleton of the spatial expansion of the central urban area. The two axes in the central district form the east–west urban development axis.

3.4. Comprehensive Detection of Influencing Factors of Urban Thermal Environment

3.4.1. Single Factor Detection

In the process of urban development, the surface temperature in some areas increases, which may be caused by urban expansion. However, there are also some areas in the middle of the city with high population and road construction density, but the urban heat island effect is not significant. Therefore, to better study the temperature change mechanism in the urban warming and cooling areas, it is necessary to use the corresponding geostatistical methods. Exploring the factors and mechanisms affecting the urban residential thermal environment is key to analyzing the evolution of LST. Therefore, in this paper, for the 2013, 2015, 2018, and 2020 data of the study area, we carried out the sampling in a 30 m × 30 m grid, and according to the area ratio of each temperature grade. In addition, for the obtained natural data and human economic data, we extract and filter the null values.

The land surface temperature extracted from sampling points was selected as the dependent variable. The factor detection results reflect the influence importance of each factor selected on the dependent variable, as shown in Figure 9.

The largest contributions to the correlation with urban temperature are the vegetation index (NDVI) and soil conditioning factor (SAVI) (Figure 9). They always have a high negative correlation with urban temperature from 2013 to 2020, indicating that in urban thermal environments, greenness can mitigate urban temperature. For the impervious layer, the building index (NDBI) has a high positive influence on the urban thermal environment, indicating that the densely built areas promote the formation of high urban temperatures and are consistently at a high positive level on the time study scale. Conversely, in the wetness layer, TCW also has a strong negative correlation index, and similarly, urban water bodies “cool” the urban thermal environment. For the human-economic index, the road density (RD), economic index (NIL), and population density (PD), which reflect the urban construction, have a strong positive correlation with the formation of the urban thermal environment, indicating that the urban temperature is higher in areas of rapid construction and economic development.

3.4.2. Multi-factor Interaction Detection

Although correlation analysis can measure the degree of correlation between influencing factors and LST, it cannot reveal the spatially stratified heterogeneity of LST or determine the interaction effects of factors. To overcome these shortcomings, in this paper, using the Geodetector model, NDVI for greenness, TCW for wetness, NDBI for impervious surface, human indicators were selected as independent variables $X$, and surface temperature was the dependent variable $Y$. By identifying the interactions between the independent variables $X$, the influence of these drivers on surface temperature $Y$ could be determined, as well as the explanatory power of $Y$. In addition, the analysis determines whether the explanatory power of the drivers (independent variable $X$) on the urban surface temperature (dependent variable $Y$) are independent of each other.

The interaction detection results show that the interaction of the six drivers are significant for LST (all q-values are greater than 0), suggesting that not a single factor influences LST. The analysis of the interaction detection results of the influencing factors reveals that the more significant effects on LST are mainly the two-factor enhancement and the non-linear enhancement.

As can be seen from Figure 10, in 2013 and 2015, when greenness was superimposed on wetness, imperviousness, NIL, POP, and ROAD, wetness was superimposed on imperviousness, NIL, POP, and ROAD, imperviousness was superimposed on NIL, POP, and ROAD, NIL was superimposed on POP and ROAD, and POP was superimposed on ROAD, the interaction on urban LST was bivariately enhanced (P(X1 ∩ X2) > max(P(X1),P(X2))); the interaction on urban LST was nonlinearly enhanced (P(X1 ∩ X2) > P(X1) + P(X2)) with wetness superimposed on POP and ROAD in mid-2018, and with wetness superimposed on NIL, POP, and ROAD in 2020.

According to the results of the multiple linear regression and calculations of the superimposed factor correlation (

Table 4. Factor superposition correlation.

Factor	2013	2015	2018	2020
Imperviousness	0.608	0.673	0.488	0.537
Imperviousness + Greenness	−0.24186	−0.38674	−0.27084	−0.18157
Imperviousness + Wetness	−0.45666	−0.45527	−0.4473	−0.42594
Imperviousness + PD	0.261	0.27964	0.23982	0.25941
Imperviousness + NIL	0.4094	0.44785	0.34508	0.38212
Imperviousness + RD	0.40363	0.43086	0.31862	0.38718
Greenness	−0.516	−0.607	−0.351	−0.371
Greenness + Wetness	−0.55786	−0.56527	−0.44834	−0.52687
Greenness + PD	0.261	0.27963	0.23982	0.25941
Greenness + NIL	0.40852	0.448	0.34482	0.38197
Greenness + RD	0.40357	0.43099	0.31858	0.38715
Wetness	−0.439	−0.500	−0.36	−0.404
Wetness + PD	0.2608	0.26987	0.23952	0.25913
Wetness + NIL	0.12008	0.40665	0.25468	0.30644
Wetness + RD	0.38512	0.40225	0.30683	0.37447
PD	0.225	0.201	0.16	0.166
PD+NIL	0.2612	0.46001	0.24055	0.26032
PD+RD	0.26376	0.3788	0.24426	0.26454
NIL	0.355	0.043	0.248	0.234
NIL+RD	0.40855	0.40855	0.40767	0.39634
RD	0.328	0.315	0.226	0.240

), the trend and intensity of the impact of the strong impact factor superimposed on the weak impact factor are correlated with the trend and intensity of the strong impact factor. Imperviousness and greenness showed the same results as above when superimposed with wetness, POP, NIL, and RD, respectively. All are consistent with the above-mentioned superposition of two factors that are independent of each other, with the strong influence factor dominating the trend for the superimposed influence.

The superposition of two isotropic strong influence factors has an increasing trend for surface temperature in the same direction. In 2013, 2015, 2018, and 2020, POP is superimposed with NIL and RD, respectively, and NIL is superimposed with RD; the trend of these strong influencing factors on surface temperature after two super impositions is the same enhancement trend.

For the influence factors with opposite trends, the superposition of positive and negative influence factors will neutralize the influence trends, which also provides a new idea for urban governance. It is not only possible to mitigate the stronger factors and strengthen the weaker ones, but also to mitigate the urban thermal environment by finding the opposite of the factors, which are related to each other but have opposite influence trends. For example, imperviousness and greenness showed a positive driving force for urban thermal environment formation, while greenness showed a negative driving force for urban thermal environment formation; however, after superposition, there was a weakening effect.

4. Discussion

In this paper, Zhengzhou, a typical rapidly urbanizing city, is used as the research object to investigate the spatial and temporal evolution characteristics of the thermal environment in the study area by inversion of surface temperature using the atmospheric correction method. we used a geographic detector to explore the effects of multiple indicator overlays on the formation of urban thermal environments. This research provides decision support for urban planning, improvement of living environment quality, summer heat disaster management, and urban emergency response.

With the acceleration of urban construction, the hot area of the main city will not change much, but the hardening of the surrounding land and the construction of public facilities will cause the overall temperature of the city to rise and the area of the thermal environment to continue to expand. The main cold areas are in the southwestern mountainous region and the main Yellow River channel and urban inland waterways in the north, and the areas with decreasing temperatures from 2013 to 2020 are mainly located around the newly built urban water bodies. According to the result of the slope trend analysis, the main temperature growth areas are concentrated in the northwest and north areas of the city. Through the analysis, it is found that the change trend of the hot spot area obtained by using standard deviation ellipse is consistent with the change trend of the city development; as such, the future development direction of the city through reasonable planning will have a certain positive effect on alleviating the heat island effect.

The effects of single factors on the thermal environment effect are mainly focused on natural factors, among which greenness and wetness have a strong negative correlation. This also indicates that vegetation and water bodies have a mitigating effect on urban thermal environmental effects. Additionally, the strong positive correlation of impervious surface on the effect of urban thermal environment indicates that urban land hardening and urban construction have exacerbated the increase in surface temperature to some extent. In addition, densely populated economic centers or road-networked areas also have a certain influence on the urban thermal environment effect.

The interaction test results show that the interaction of the six selected drivers has significant explanatory power for LST (all q-values are greater than 0), and the superimposed effects on the explanatory power of LST are mainly two-factor enhanced and non-linear enhanced. The interactive effect of two influencing factors on the thermal environmental effect is greater than the independent effect of a single influencing factor on the thermal environmental effect; furthermore, the combination of natural and human factors all significantly increase the intensity of thermal environmental effects.

In addition, the construction of sponge cities and green China in recent years will alleviate the extreme urban heat to some extent, but the disorderly and unreasonable planning will further aggravate it, and future studies should also focus on the trend and mechanism of the urban heat phenomenon behind these unreasonable plans.

5. Conclusions

For other typical cities similar to Zhengzhou, which is a monocentric city with a large population, the study of urban habitat problems has become an important topic. In this study, we analyzed the interannual variation trend of the thermal zone of Zhengzhou city from 2013 to 2020, combined with the standard deviation ellipse, to find that the main high temperature zone of the city is located in the northwest–southeast axis and has some relationship with the construction of new cities. We also analyzed the single and superimposed factors affecting the surface temperature with the help of a geographic probe model to explore the urban surface temperature model deeply. The following is proposed:

• In the future, we should start to strengthen the protection of urban green areas, water bodies, and landscapes to avoid the over-concentration of green areas and the destruction of ecological environments such as water body sewage, and for the improvement of the cooling effect of green areas through homogenization, decentralization, irregularity, and centralization of boundaries.
• Through the analysis, it was found that urban development is related to the urban heat zone. It is important to optimize urban spatial structure and form, control urban building high-density development, in order to appropriately control urban construction intensity, reasonably control the growth scale of construction land, increase the area of non-construction land, and promote intensive, green and low-carbon land use, which can help mitigate the urban heat environment effect. This would actively guide the construction of urban ventilation corridors and sponge cities.
• According to the results of the study on the causal mechanism of urban thermal environment, to improve the urban thermal environment, the protection and maintenance of natural ecosystems such as parks, green areas, and wetlands should be strengthened, and the area of ecological land should be increased in an appropriate amount.