3.1. Analysis of Spatial Pattern of LST at Different Urban Zones
In this study, the spatial pattern of LST and its variation were analyzed based on different zones between bi-temporal images.
Figure 3c and
Figure 3d show that high LST (LST
4) and sub-high LST (LST
3) areas (areas with > LST
mean) were mainly distributed in the central (zone 1) and some areas of zone 2 in 2004. Zone 1 and part of zone 2 were assembled with ISA of buildings and roads, resulting in high and sub-high LST areas on both dates. Meanwhile, the distribution of LST
4 and LST
3 areas in zone 2 and zone 3 in 2021 was notably larger than in 2004. In zone 2 and zone 3, the extent of spatial distribution of LST
4 and LST
3 areas significantly increased from 2004 to 2021, especially in zone 2 and the western part of zone 3. This indicates that urban expansion resulted in an increase in LST primarily in zone 2 and zone 3 between the two dates.
The medium (LST2) and low (LST1) LST areas were mainly distributed in zone 2 and zone 3 of peri-urban areas in both images, especially zone 3. This was limited by terrain and topography because LST2 and LST1 areas were characterized by forested and water-covered regions in mountainous areas. Urban expansion led to a transformation in land cover from natural to ISA, mainly occurring in the south and west plain areas of zone 2 and zone 3 from 2004 to 2021. The land cover types in the northern areas of zone 2 and zone 3 changed little from 2004 to 2021 due to mountainous constraints; consequently, the medium and low LST areas experienced minimal change in these regions during the period.
To further quantify the impact of high and sub-high LST areas on the thermal environment, we need to analyze the change in areal proportion of high and sub-high LST areas with urban expansion.
Table 2 shows the areal proportions of the LST
4 and LST
3 categories in each zone and the total LST
4 and LST
3 area of the study area. As shown in
Figure 3, LST
4 and LST
3 areas dominated zone 1, accounting for 94.99% in 2004. However, the areal proportions of LST
4 and LST
3 areas in zone 1 decreased slightly to 92.98% in 2021, likely due to increased attention to urban ecology, consequently, greater vegetation cover interspersed in some areas of zone 1 after 2004. For zone 2 and zone 3, the areal proportions of LST
4 and LST
3 areas increased from 56.11% and 21.08% to 62.03% and 32.49%, respectively, from 2004 to 2021. This indicates that the high and sub-high LST areas exceeded LST
mean in this study area resulting from urban expansion primarily increased in zone 2 and zone 3 during the period.
As for the areal proportions of high and sub-high LST areas in each zone relative to the total high and sub-high LST areas of the entire study area, the areal proportions of LST
4 and LST
3 areas in zone 1 decreased from 26.58% to 20.78% during the period. In zone 2, areal proportions also decreased from 34.97% to 31.04%. However, the areal proportions of LST
4 and LST
3 areas in zone 3 relative to those of the entire study area increased from 38.45% to 48.18% in the same period. This indicates that, with urban expansion, high and sub-high LST areas mainly increased in zone 3, followed by zone 2. Although the area of high and sub-high LST also increased in zone 2 as shown in
Figure 3, the rate of increase was slower compared to that of zone 3 when considering the entire study area (
Table 2).
Obviously, zone 1 had the highest LST in general, followed by zone 2, and then zone 3. The high and sub-high level LST areas were heat sources contributing to the UHI effect. The lowest LST areas are primarily vegetation areas, which act as heat sinks. Although zone 2 has a cooling thermal environment due to its vegetation area, its cooling effect is limited due to the small area of vegetation. In zone 3, there are numerous vegetation areas that serve as major heat sinks, playing the most significant role in cooling the thermal environment in the study area. Additionally, the lowest LST areas also occur in the water bodies of all three zones.
3.2. Analysis of Thermal Environment Contribution at Each Urban Zone
The quantitative analysis of the thermal environmental impacts of different urban landscape patterns becomes important in urban ecology because of urban thermal environment problems. In this study, zone 1 obviously shows a higher LST, followed by zone 2 and zone 3. However, when considering the heat effect of each zone on the thermal environment and their changes between different dates, it needs to be further analyzed. In order to further quantify the heat contribution in detail, Eq. 6 was used to calculate the CI
i and to further analyze the thermal environmental impact of different zones in this study (
Figure 5).
Urban expansion with ISA clearly resulted in an increase in LST. As shown in
Figure 5, the CI
i of zone 1 and zone 2 were 24.84% and 33.44% in 2004 respectively. With the urban expansion resulting in LST variation, the CI
i of zone 1 and zone 2 decreased to 10.61% and 25.91% in 2021 respectively. However, the CI
i of zone 3 increased from 41.72% in 2004 to 63.49% in 2021. Though ISA increased in zone 2 and zone 3 in the study area, and the area of ISA with >LST
mean also increased in zone 2 from 2004 to 2021 (
Table 2), the increase rate of ISA with >LST
mean was higher in zone 3 than in zone 2 in the period. This resulted in the decrease of CI
i in zone 2.
Figure 5 indicates that urban expansion, coupled with areas exhibiting >LST
mean, has led to thermal environment issues, predominantly escalating in zones 2 and 3 during this period, particularly in zone 3. Zone 3 and zone 2 exhibited relatively high CI
i values on the two dates compared to zone 1, with zone 3 being the primary source of the heat effect in the study area and zone 2 in second place. Although CI
i values of zone 1 obviously decreased, the area with >LST
mean decreased only slightly between the two dates (
Table 2). The increase of ISA with >LST
mean was more pronounced in zone 3 and zone 2 than in zone 1 from 2004 to 2021.
Compared with zone 1, the CIi values of zone 2 decreased relatively little between the two dates. From 2004 to 2021, the CIi of zone 3 increased significantly and contributed 63.49% of the heat effect in the entire study area in 2021, whereas the influence of zone 3 on the thermal environment was not significant in 2004. The thermal environment impact of zone 3 increased significantly in 2021 due to urban expansion. The CIi of zones 2 and 3 reached from 75.16% in 2004 to 89.40% in 2021, indicating that the heat effect contribution was determined and dominated by zones 2 and 3 in the study area.
The different urban development zones exhibited distinct landscape patterns and LST features. To further quantify the impact and pattern of the thermal environment, WLUI and RWLUI were calculated based on Eqs. (8) and (9) to assess the thermal contribution effects of different zones in the study area (
Figure 6). WLUI represents the percentage of pixels with values exceeding the LST
mean for each zone.
Figure 6 showed that WLUI of zone 1 was high, accounting for approximately 93.07% and 93.22% from 2004 to 2021, with little change between the two dates. WLUI of zone 1 indicated that areas with high and sub-high LST (>LST
mean area) dominated in zone 1, with the areal ratio reaching over 93% in the urban core zone with high development density. Meanwhile, the WLUI values for zone 2 and zone 3, respectively decreased in sequence from zone 1 to zone 3 on both dates. Comparing WLUI values on both dates, WLUI values in zone 2 and zone 3 increased from 60.91% and 26.41% in 2004 to 65.25% and 34.25% in 2021 respectively. This reflects a higher increased rate of increase in WLUI in zone 3 compared to zone 2, indicating an amplified thermal environment impact in both zones from 2004 to 2021. The WLUI values of zone 2 and zone 3 indicate that >LST
mean areas increased in areas of medium development density and low development density with expansion in the period.
The RWLUI is the percentage of the high and sub-high LST area of each zone to the whole study area (
Figure 6). The RWLUI values of zone 1 were 9.18% in 2004 and 9.21% in 2021. From the view of the whole study area, the thermal environment impact of zone 1 was nearly not strengthened from 2004 to 2021 because urban expansion mainly existed in zone 2 and zone 3 between the two dates. The RWLUI values for zone 2 and zone 3 continued to increase from 13.16% and 18.11% in 2004 to 14.09% and 23.47% in 2021 respectively.
The increase in thermal environment contribution indicates that high and sub-high LST areas increased with the rise in high percent ISA resulting from urban expansion.
Figure 5 and
Figure 6 indicated that, though the areal proportion of high and sub-high LST areas was dominant in zone 1 and the lowest in zone 3, the contribution of zone 3 to regional heat in the study area was higher than that of zone 2 and zone 1. This is because zone 3 contributes the most high and sub-high LST area among the three zones, followed by zone 2 and zone 1. Our analysis revealed the characteristics of the thermal environment effect and their variation in different development density areas at two different time points. Through analyzing CI
i, WLUI, and RWLUI of each zone and their variations, we can quantify the thermal environmental contribution of different development density areas with urban expansion. Different from other analysis of thermal environmental contribution based on different land cover types, this analysis is significant for urban planning, UHI mitigation, and urban ecology.
3.3. Driving Force Detection of Spatial Heterogeneity of LST by Geodetector
The spatial distribution and aggregation of LST were different in different urban development areas. The spatial distribution pattern of LST is related to the degree of urban expansion and impacted by different land cover types, the combination and configuration of various land cover types, and other factors. Detection of driving forces behind urban LST heterogeneity is crucial for alleviating UHI and adapting to urban climatic change. The spatial heterogeneity of LST results from the spatial distribution of different land cover types and other factors [
45]. In this study, percent ISA, fractional vegetation cover, and elevation were selected as independent variables (X) and were all categorized into five types based on the Jenks natural breaks classification method. LST was selected as the dependent variable (Y). The influences of X on the spatial differentiation of Y in different zones in 2004 and 2021 were calculated using Geodetector and further comparatively analyzed.
In this study, the factor differentiation detector was used to reveal the impact of each factor on the spatial distribution of LST.
Table 3 shows the q-values of three single factors on the spatial distribution of LST, and all were significant (p < 0.01). In the study area, the degree of influence of the three driving factors on the spatial distribution of LST varied in different development densities of the urban area. The q-values of three influencing factors in zone 1 were relatively low (⩽0.309) on both dates, and variations were also small in zone 1 between the two dates (
Table 3). This indicates that the urban landscape pattern changes and urban expansion nearly did not exist in the urban core zone 1 from 2004 to 2021. Among the influencing factors, the FVC had a greater impact on LST in zone 1 in 2004. From 2004 to 2021, the q-value of FVC decreased from 0.309 to 0.144 in zone 1. According to the principle of Geodetector, the larger the q-value, the greater the driving force of the factor on LST.
Table 3 indicates that the three factors had a lesser degree of influence on the heterogeneity of LST, and the spatial distribution patterns of LST varied minimally in zone 1 between the two dates.
The three factors had a greater influence in zone 2 and zone 3 than in zone 1 on both dates. When comparing the q values of the three influencing factors between the two dates, the q values in zone 2 and zone 3 in 2021 were larger than those in 2004, except for the q value of FVC at zone 2 in 2004, which was slightly larger than that of zone 2 in 2021. The results from the factor detection analysis showed that the largest q values were 0.581 and 0.653 in zone 2 in 2014 and 2021 respectively, and 0.494 and 0.681 in zone 3 in 2004 and 2021 respectively. This indicated that FVC had the most significant influence on LST heterogeneity in zone 2 and zone 3 in 2004, while elevation had the greatest influence on LST heterogeneity in zone 2 and zone 3 in 2021. This is because the landscape pattern changes due to urban expansion were influenced by topography (elevation) in zone 3 in 2021, which had a certain impact on the q value of zone 3. The elevation was ranked next, with values of 0.558 and 0.432 in zone 2 and zone 3 in 2004 respectively, and FVC was ranked next to elevation, with values of 0.572 and 0.602 in zone 2 and zone 3 in 2021 respectively.
When comparing q values between zone 2 and zone 3, the q values of the three influencing factors were smaller in zone 3 than in zone 2 in 2004. However, the q values of the three influencing factors were larger in zone 3 than in zone 2 in 2021. This is related to the change in urban landscape patterns in zone 2 and zone 3 resulting from urban expansion. FVC and elevation had a higher influence degree than percent ISA in zone 2 and zone 3 in 2004. The explanatory power (q value * 100%) of FVC and elevation decreased from 58.1% and 55.8% in zone 2 to 49.4% and 43.2% in zone 3 respectively, in 2004, and increased from zone 2 to zone 3 in 2021. This indicates that urban expansion also resulted in the change of landscape patterns in some higher elevation areas of zone 3. FVC and elevation had a higher influence degree on spatial heterogeneity of LST in zone 2 and zone 3 in both dates. However, the explanatory power of percent ISA significantly increased from 34.9% to 52.7% in zone 2 and from 22.4% to 54.3% in zone 3 from 2004 to 2021. Though the q values of FVC and elevation in zone 2 and zone 3 were larger than that of percent ISA in 2004, the q values of the influencing factor percent ISA significantly increased from 2004 to 2021 with urban expansion.
According to the above analysis, the three factors have weak explanatory power for the spatial heterogeneity of LST in zone 1 compared to zone 2 and zone 3. This means that the change in landscape pattern resulting from urban expansion increased the influence degree of the three driving factors on LST significantly in zone 2 and zone 3 compared to zone 1. In different urban development areas, the land cover types and their spatial distribution are different, which influences the value and spatial distribution of LST. In zone 1 of the urban core area, a high percent ISA is distributed homogeneously, resulting in the homogeneity of spatial distribution of high LST. Therefore, the influence degree of the three factors on LST heterogeneity is low. However, in zone 3 and zone 2, the three factors have a higher influence degree on the spatial distribution of LST. This indicates that they are the main driving factors affecting the LST in zone 2 and zone 3, especially the elevation factor in zone 3 in 2021. However, no one single factor dominantly affects LST distribution in zone 2 and zone 3 of the study area. Through analyzing and comparing q values of different zones, we can understand the driving force of LST distribution in different development areas, which will help in constructing a good urban ecological environment.
The spatial distribution of urban landscape and thermal patterns is usually the result of multifactor interactions. The influence of the interaction between different factors may be stronger than that of a single factor when each factor interacts with other factors on LST. Therefore, it is necessary to analyze the influence of the interaction between any two factors on the spatial difference of LST. Through interactive calculation of different driving factors in Geodetector, we analyzed the interactive detection of various driving factors on LST distribution for both dates. The q-value is defined as q (X
1∩X
2) and the results are shown in
Table 4. Comparing with the results of single factor detection in
Table 3,
Table 4 shows that the q value increases significantly when each driving factor interacts with other factors, and its influence on LST increases, indicating that LST is a geographical pattern driven by multi-factor interaction.
The interaction detector creates a new stratum by overlaying factors X1 and X2to further compare q(X
1∩X
2) with q(X
1) and q(X
2) [
46]. The results of interaction detection can show whether the combined effects of the evaluation indicators X1and X2 enhance, bi-enhance, nonlinearly enhance, weaken, bi-weaken, or nonlinearly weaken the explanatory power of LST heterogeneity, or whether the combined effects on LST are independent of each other. In this study, the explanatory power of the interactions of different factors in
Table 4 was larger than that of the single factor in
Table 3, and all the interactive factors showed nonlinear enhancement in 2004 and 2021. This indicates that the interaction between any two factors of percent ISA, FVC, and elevation amplifies their impact on LST distribution in three zones of two dates. The q values of interaction detection were significantly larger than those of single factor detection in
Table 3 on both dates. This means that interactions of different factors have higher explanatory power on the LST heterogeneity at each zone than those of a single factor. The results from interaction detection analysis showed that the q values of interactions of three factors at zone 1 were lower than those of zone 2 and zone 3 (
Table 4). This trend is the same as that in
Table 3. This is because the landscape patterns of zone 1 are significant different from those in zone 2 and zone 3. Though the q values of interactions of three factors increased at zone 1 in both dates, the q values of interactions of three factors also significantly increased at zone 2 and zone 3 in both dates when compared with those of the single factor in
Table 3.
As shown in
Table 4, the results from interaction detection analysis showed that FVC ∩ EL and percent ISA ∩ EL at zone 2 were 0.743 and 0.749 (74.3% and 74.9% explanatory power on the LST heterogeneity) in 2004 and 2021, respectively, and were the most significant influences on LST heterogeneity at zone 2. The interactions of percent ISA ∩ EL and percent ISA ∩ FVC were ranked next at zone 2 in 2004 and 2021, with values of 0.637 and 0.591 in 2004 and 2021, respectively. As for zone 3, FVC ∩ EL was the most significant influence on LST heterogeneity in both dates; the q value was 0.621 and 0.767 in 2004 and 2021, respectively. Percent ISA ∩ FVC with a value of 0.510 and percent ISA ∩ EL with a value of 0.757 were ranked next at zone 3 in 2004 and 2021, respectively. Based on the analysis of the interaction between two factors, we can determine which kind of interaction has the strongest influence on the spatial heterogeneity of LST. This can help guide the composition and configuration of the driving factors in different development areas of urban for UHI.
In the study area, it was found that the influencing factors were different in three zones. Synthetically evaluating the results of factor and interaction detection, it is worth noting that the landscape pattern of urban expansion in zone 3 was limited by geomorphologic characterization, and elevation leads to an imbalanced distribution of LST in zone 3. When urban expansion was intense in zone 3, it is clear that most of the new urban area was distributed around the mountain. This led to an imbalance in urban expansion and also had a certain impact on the q value of zone 3. In addition, as a single factor, FVC had a high influence on LST heterogeneity in zone 2 and zone 3 of the study area.
Linear correlation cannot analyze the consistency of spatial distribution between two geographic variables. Compared with linear correlation analysis, geodetector can reveal the causal relationship between independent and dependent variables [
46]. In this study, the Geodetector can show the driving force of different factors on the heterogeneity of LST resulting from urban landscape expansion through the q value. Spatial heterogeneity of LST may come from an imbalance in urban expansion. This analysis can reveal this imbalance and further help to understand the influencing factors of LST in different development areas of urban. It can provide significant information for urban planning and coping with urban climatic change.