Strengthening of surface urban heat island effect driven primarily by urban size under rapid urbanization: national evidence from China

ABSTRACT Surface urban heat islands (SUHIs) have received worldwide attention owing to their significant adverse impacts on socio-economic and environmental conditions. However, how the SUHIs vary across cities of different sizes and whether there is a threshold for city size that affect SUHIs were not explicit due to the limited selection of large cities in previous studies. Our study attempts to support comprehensive decision-making for national sustainable spatial management of SUHI in urbanized areas of different sizes. We used impervious surface area data to define the urbanized boundary for obtaining the SUHI in China. Then we analyzed the spatio-temporal patterns of SUHIs under China’s rapid urbanization. We furtherly constructed spatial autoregressive and structural equation models to verify the drivers’ direct and indirect contributions to SUHI. The SUHI effects existed extensively in most urbanized areas, not only in large or mega urbanized areas (average value: 2.08 ± 0.68°C) but also in petty urbanized areas (average value: 0.64 ± 0.19°C). We also found that the size of urbanized areas contributed the most to the increase in SUHI intensity (SUHII). Natural drivers are mediators indirectly contributing to variations in SUHII, influenced by anthropogenic drivers. In addition, we also found that the size of urbanized areas contributed the most to the increase in SUHI intensity (SUHII). The size of the urbanized areas had a positive non-linear relationship with SUHIIs. When the size was larger than 400 km2, the growth of SUHIIs maintained an equilibrium state. This study highlights the importance of impervious surface expansion for increasing SUHII. To confront SUHIs, it is necessary to perform proper urban planning.


Introduction
Urbanization is occurring at a remarkable rate worldwide. By 2050, the global urban population is expected to increase by approximately 2.5 billion and reach 60% of the total population (Ramaswami et al. 2016). Rapid urbanization can alter the surface radiation, humidity, and thermal environment in urban areas (Tran et al. 2017;Yao et al. 2022) and cause the urban heat island (UHI) effects (Zhou, Wang, and Cadenasso 2017). UHIs adversely affect the net primary production of ecosystems (Mohammad and Goswami 2021) and can aggravate respiratory diseases and even cause death (Ren et al. 2018). For example, a study indicated that for every degree of temperature increase, the mortality rate of urban residents increased by 4% (Peng et al. 2012). Quantifying the spatiotemporal variations of UHIs and identifying their drivers can aid decision-makers in developing and executing rational policies for urban sustainability .
SUHI results from a modification of the energy balance at the urban surface (Weng 2009). Compared with in-situ air temperature measurements, LST data from the Moderate Resolution Imaging Spectroradiometer (MODIS) sensor had wider sources, broader coverage, and steadier periodicity for monitoring SUHI (Voogt and Oke 2003;Stewart 2011). Previous studies related to LST data have been conducted in different cities revealing, that the daytime SUHI intensity (SUHII) in major cities generally reaches 1.5°C (Li, Zha, and Wang 2020b). Previous studies on SUHIs have primarily focused on specific regions, such as mega cities or major cities with well-developed economies and rapid urbanization (Imhoff et al. 2010). However, small and medium cities have been neglected, which can occupy a larger proportion of total cities (Fei and Zhao 2019). As these studies only selected large cities and did not examine the full coverage of city size, they concluded that city size had no significant effect (Schwarz, Lautenbach, and Seppelt 2011). However, it has been suggested that SUHIs expand with the growth of urban size. It is even unclear about the relationship between city size and SUHIs and whether there is a threshold range for this relationship between SUHIs and urban area size. (Liu and Morawska 2020). Thus, some conclusions about large and mega cities may mislead decisionmakers about the sustainable development of medium and small cities (Elmes et al. 2020). As urbanization continues, the distances between urbanized areas have shrunk or even disappeared, and thus, SUHIs are no longer restricted to single cities (Arshad et al. 2022). Therefore, using urbanized areas with full sample coverage as the research object can decrease the inaccuracy of the results owing to sampling surveys in the study area and more accurately estimate the SUHIs of cities of different sizes (Imhoff et al. 2010).
In addition, past studies have explored the relationships between the spatiotemporal variations of SUHIs and associated factors. SUHIs were strongly linked to anthropogenic heat release and urban form (Kim and Brown 2021). Vegetation, temperature, and precipitation also influence the SUHIs, while wind speed contributes less (Peng et al. 2019b). Urban landscapes are characterized by complex spatial heterogeneity with different Land-cover/land-use (LCLU) types (Talukdar et al. 2022). Therefore, the urban landscape structure also could influence the variations of heat islands, for example, LCLU and the shape of impervious surface patches. Past studies have analyzed the detailed effects of these factors on SUHII through statistical methods such as correlation analysis, stepwise regression analysis, and multivariate statistical analysis . These methods measure the interrelationship between variables, suitable for time series, without considering the spatial association between regions (Naikoo et al. 2022). However, previous studies have indicated that LST tends to correlate spatially with LSTs observed nearby owing to atmospheric flows (Song et al. 2014). Statistical analysis ignores the spatial correlation of the SUHI distributions, leading to evident limitations. Therefore, some researchers have selected spatial autocorrelation models to eliminate the interactions between geographic units. In addition, these models can only focus on the impact of the elements on SUHI alone, while urban development does have feedback effects on natural conditions (e.g. background climate). This indirect relationship of drivers to SUHI also works on the SUHI change. There is still a lack of investigation of the interaction of these internal factors on SUHII.
Apparently, in the context of rapid urbanization, the spatiotemporal patterns and potential mechanisms of SUHII in different-sized cities need to be comprehensively evaluated. In this study, we examined urbanized areas in mainland China and analyzed the temporal change tendency and spatial patterns of SUHII, and the associated factors. The objectives of the study were to (1) illustrate the urbanization process with the characteristics of impervious surfaces in China from 2003 to 2015; (2) characterize the SUHII evolution at diurnal, seasonal, and annual scales and analyze SUHI spatial patterns in urbanized areas of different sizes; (3) determine the relationships between SUHII and variables, particularly for urbanized area size, using spatial autoregressive models. Structural equation modeling (SEM) was used further to test the causal relationship between drivers and SUHII.

Study area
China has a vast territory, covering an area of 9.6 million km 2 . The study area is the urbanized area of China in this paper ( Figure 1). The terrain of China gradually decreased from west to east. Affected by the monsoons, the climate in China differs greatly between seasons and contains five climate zones, ranging from tropical to cold temperate zones. The average temperature in China decreases from south to north, while precipitation decreases from the southeast coast to the northwest inland. The socioeconomic development in China also presents significant diversity. China has urbanized rapidly owing to accelerated industrialization. By 2018, the proportion of the urbanized population increased from 17.92% to 59.58%. According to the empirical fact of global urbanization (Kuang et al. 2016), China is still in a rapid development stage. However, the urbanization scale varies greatly within China because of different natural and economic conditions. The remarkable regional differences in urbanization and natural ecosystems make China an ideal study area for exploring the urbanization and thermal environmental effects it generates. To investigate the geographical variations in SUHI, China was divided into seven regions according to its geographical location and urbanization: South China, Southwest China, Central China, East China, Northeast China, North China, and Northwest China (Figure 1). The South China, Southwest China, Central China, and East China have typical humid-hot climate. Northeast China, North China, and Northwest China have typical humidcold, subhumid/semiarid-temperate, and arid climates.

Defining urbanization
The impervious surface data used in our study were derived from high-resolution the annual artificial impervious area dynamics at a 30 m resolution (GAIA), http://data.ess.tsinghua.edu.cn.).The data were interpreted from Landsat imagery (TM/ ETM+/OLI) available on the Google Earth Engine with an overall accuracy of over 90% (Peng et al. 2019a).
Urbanization has profoundly altered climate change. Urbanization has been widely simulated using impervious surface area, population density, and building density indicators. As SUHIs occur primarily because of the impervious surfaces absorbing radiation and emitting heat to the surroundings, our study selected the impervious surface area to quantify the urbanization (Dong et al. 2020).
In our study, the city is an abstract concept relative to the countryside, rather than just the concept of a city in administrative planning. Therefore, this study considered continuous impervious surface pixels as urbanized areas (Li et al. 2020a). We selected impervious surface data, performed the priority search algorithm, and defined urban areas as impervious surface patches around county administrative centers of each city for the years 2003, 2005, 2010, and 2015. According to the impervious surface area, cities were classified into five categories: mega urbanized areas (area > 1000 km 2 ), super urbanized areas (400 < area ≤ 1000 km 2 ), large urbanized areas (200 < area ≤ 400 km 2 ), medium urbanized areas (50 < area ≤ 200 km 2 ), petty urbanized areas (area ≤ 50 km 2 ) (Notice on Adjusting the Criteria for City Size Classification, Notice of the State Council, 2015).

Quantifying the surface urban heat island effect
MODIS LST data were selected for this study (Bo et al. 2013). Numerous SUHI studies have documented the spatial patterns of single intra-city SUHII that utilize the resolution and long time series of Landsat TM/ ETM+ (Shen et al. 2016). The variability and characteristics of SUHII from the national to the global scale have also been studied using the MODIS LST product (Clinton and Gong 2013). MODIS has stable transit moments, providing LST data of spatiotemporal continuity and integrity at typical moments. MODIS LST data are appropriate for large-scale SUHI studies.
We used the MODIS Aqua/Terra 8-day composite (MYD11A2 and MOD11A2) LST products from 2003 to 2015 (http://ladsweb.nascom.nasa.gov/data/search. html). Terra satellite overpass time is at 10:30 and 22:30, and Aqua satellite overpass time is at 13:30 and 1:30 (local solar time). Using only one transit time LST data may increase the error of LST, so daytime LSTs were obtained by averaging the LST at 10:30 and 13:30, while nighttime LSTs were obtained by averaging the LST at 22:30 and 1:30 (Hu et al. 2014;Peng et al. 2014). The errors of the MODIS LST products are usually within ± 2 K. Specifically, the MODIS LST products had already been geometric correction and atmospheric correction. However, MODIS LST is a clear-sky instantaneous LST. Hence, we averaged the LST difference between the urban area and the surrounding area on a monthly basis to obtain the SUHI for each urbanized area. Furthermore, the annual mean values were calculated for spring, summer, autumn, and winter, respectively (Hu and Brunsell 2013;Sharma and Joshi 2014). In this study, we defined spring as March-May, summer as June-August, autumn as September-October, and winter as December-February. However, the Aqua satellite was launched in May 2002, and our study period included the winter months of December-February and the spring months of March-May. Due to data limitations, the assessment period was set from 2003 to 2015.
The magnitude of the SUHI can be quantitatively expressed by the SUHII. This study defined SUHII as the LST difference between the urban and surrounding rural areas. We used a buffer method with an area equal in size to its corresponding urban area, i.e. always keeping the area of each urban area the same as its buffer area, which was defined as the rural area surrounding these urban areas. In the determination of the buffers, we estimated the impervious surface patches as circles to obtain the buffer radius. This approach to defining buffer zones is more consistent than previous studies, which selected regions within cyclic annular distances to the urban extents (Li, Chen, and Gao 2022). Urban pixels in the buffer zone were removed before calculating SUHII. Then the SUHII was calculated as the difference in the average LST in the urban and rural area, i.e. surrounding buffer area. The following formula can calculate SUHII: The seasonally averaged SUHIIs for each urban area were obtained in 2003, 2005, 2010, and 2015. The present study calculated the SUHII of impervious surface area at the patch-scale rather than a city-scale. And we select urban patches with SUHII > 0 by extracting the boundary of impermeable surfaces to determine whether urbanization has a heating effect.

Drivers of surface urban heat island effect particularly urban size
Influencing drivers of multiple dimensions were selected, covering urban development, underlying surface characteristics, climatic backgrounds, geographic position, and community economy. Considering the magnitude of this study and availability of data, the following datasets were assembled to examine their relationship between drivers and SUHII: elevation, air temperature, precipitation, NDVI (Normalized Difference Vegetation Index), longitude and latitude, compactness index (CI), urban size, GDP (Gross Domestic Product), population density and night lights (NL). Natural factors include elevation, air temperature, precipitation, NDVI, longitude, and latitude. The compactness index, population density, urban size, GDP and night lights represent anthropogenic factors (Feng et al. 2021). The specific descriptions of the data are presented in Table 1.

Spatial autoregressive model (SAR) methods
Our study employs a spatial autoregressive model (SAR) to evaluate the relationship between SUHII and all independent variables. The spatial autoregressive model can be expressed by the following equation (1): Where y is the dependent variable and represents SUHII, x is the independent variable and represents the relevant factors impacting the SUHII, β is the regression coefficient of the independent variables, α is the intercept, μ is a random error term, ε is a random error obeying a mean of 0 and a variance of δ 2 , W 1 is the spatial weight matrix of dependent variables, W 1 is the spatial weights matrix of random error, ρ regression coefficient of the spatial lag term, and λ is the regression coefficient of the spatial error term. Spatial observations of the dependent and independent variables of the ordinary least squares model (OLS) were not affected by spatial variability. In the spatial error model, there was no spatial dependence of the spatially dependent variables. In the spatial lag model, the observations of the dependent variables are related to both the corresponding independent variables and the dependent variables in the neighboring regions.
We constructed an OLS model using SUHII as the dependent variable. The observations and regression errors are independent of the OLS model. Therefore, in the second step, we identified whether the spatial dependence of the data affected the OLS model results using the global univariate Moran's I-statistic to test the presence of spatial autocorrelation of the residuals. Finally, if there was spatial dependence in the OLS model residuals, the Lagrange multiplier (LM), robust Lagrange multiplier (RLM), log likelihood (LogL), Akaike information criterion (AIC), and Schwartz criterion (SC) were used to determine whether a spatial error model or spatial lag model should be selected; the higher the LogL, LM, and RLM values and the lower the AIC and SC values, the better the model. In our study, the CI of urban impervious surface, population density (POP), nighttime lighting (NL), elevation, urban size, GDP, longitude (Long.) and latitude (Lati.), NDVI, precipitation (PREC), and background temperature (TEMP) were selected as the influencing drivers. However, a more detailed collinearity diagnostics indicated collinearity between latitude and background temperature, and then we removed latitude from the variables (Peng et al. 2019b). Because each influence driver has  Table  S1. The goodness-of-fit of the spatial error model was better than that of the spatial lag model. Thus, our study used the spatial error model to explore the concurrent effects of the ten independent variables on SUHI. The ten independent variables were divided into two main categories.

Structural equation model (SEM) methods
The SAR model can only provide information about the coefficient of influence of each driver on the SUHI individually. Furthermore, we used the structural equation model (SEM) to construct an integrated framework of the associations between anthropogenic and natural drivers and SUHII. Unlike conventional analytical methods, SEM allows quantifying the direct causal relationships between these drivers and SUHI, as well as the driver-mediated indirect relationships. These provided illustrations and insights into the mechanisms of how drivers affect SUHI. First, we proposed a hypothetical model defining the expected relationship between SUHII and the ten drivers according to scientific experience and background knowledge (Wu, Li, and Li 2021). The correlation between the drivers by Pearson's correlation coefficient and the results of the spatial statistical model were considered when constructing the hypothetical model. In our study, we assumed that longitude, elevation, urban size, night lights, population density, precipitation, air temperature, GDP, CI, and NDVI would have effects on SUHII and they would interaction with each other. We then tested whether important pathways were overlooked (modification indices ≥ 4) and whether the existing pathways were significant (p ≤ 0.05). We tested if important paths were missing, as well as removed insignificant pathways, and revised the a prior model. To estimate the fit of the SEM model, we used the values of Tucker-Lewis index (TLI), normed fit index (NFI), and comparative fit index (CFI). In our study, these indexes should be greater than 0.9 (Bentler 1980). SEM was conducted in AMOS 26 (IBM Corporation, Meadville, PA, USA) and then the results were drawn artificially by the authors in the software draw.io.

Dynamics of impervious surface area in China from 2003-2015
Overall, the impervious surface area exhibited a rapidly increasing trend from 2003 to 2015, with a growth rate of 2.15 times (Figure 2a). In 2003, the impervious surface area in China was dominated by petty and medium-sized urbanized areas, with proportions of 54.83% and 21.11%, respectively. However, by 2015, the proportion of mega urbanized areas increased from 9.69% in 2003 to 18.32% in 2015 (Figure 2c). Urban agglomerations in the mega urbanized areas are adjacent in space, and urbanization proceeded with a high degree of spatial symbiosis and consistency, particularly in Beijing-Tianjin-Hebei, the Yangtze River Delta, and the Pearl River Delta urban agglomerations. Meanwhile, petty urbanized areas decreased from 54.83% in 2003 to 40.14% in 2015 as they transformed to medium urbanized areas and others. Moreover, the major expansion pattern of mega urbanized areas, super urbanized areas, and large urbanized areas extended outward from the former center; however, the pattern of petty urbanized areas and medium urbanized areas was in-fill expansion (Figure 2c).
Meanwhile, the dynamics of urban expansion across cities with different urbanized area sizes over time were delineated (Figure 2b). The sharpest increase in the impervious surface area occurred in the period 2005-2010, and the growth rate began to  (Figure 3). In our study, petty urbanized areas had the largest proportion of area for all regions other than South China; the proportion of mega urbanized areas in South China was the highest owing to the rapid development of urban agglomeration in the Pearl River Delta. Second, we concluded that the proportion of petty urbanized areas declined, whereas large and medium urbanized areas remained invariant, and that of super and mega urbanized areas increased. This characteristic was evident in East China. Petty urbanized areas in East China expanded and transformed into mega and super urbanized areas owing to the expansion of the Yangtze River Delta urban agglomeration. Table 3 showed the average, maximum and minimum values of seasonal and diurnal SUHI across China in the period 2003-2015. In this study, the annual SUHII  was stronger in the daytime (1.12 ± 0.60°C) than the in nighttime (0.91 ± 0.53°C), particularly in spring, summer, and autumn. The maximum value of SUHII reached up to 7.52°C in the summer daytime. However, the SUHII at night was stronger during daytime in winter. The negative values of SUHII (SUHII < 0°C) were more found during the daytime in winter. The SUHII also has significant spatiotemporal variations among the seven geographical regions (Figures 4 and 5). It demonstrated a larger spatial variation in the day (0.21 to 1.69°C) than at night (0.41 to 1.40°C). The median SUHII was higher in South and East China in the daytime, where the climate is warmer and more humid. Higher median SUHII at night was observed primarily in the colder zones in Northeast China, where lower temperatures and precipitation characterize the background climate. In particular, the greatest SUHII was observed in Northeast China on summer days (1.69°C) and winter nights (1.41°C). SUHII also differed greatly by season, with the largest SUHII in summer and the lowest in winter ( Figure 6). Furthermore, the diurnal difference in SUHII was greater in spring and summer than in autumn and winter. In summary, the average daytime SUHII was stronger than at night, and SUHII was stronger in summer than in other seasons from 2003-2015. Seasonal trends in large and asymmetric SUHIIs were observed during the day and night. The difference between daytime and nighttime SUHII is greatest in summer (0.74°C) and least in autumn (0.12°C).

Spatiotemporal patterns of the SUHII from 2003 to 2015
The SUHII for urbanized areas of different sizes from 2003 to 2015 is shown in Figure 7. We observed that SUHI existed in every type of urbanized area, and the petty urbanized areas had the lowest SUHI intensity. During the day, the magnitude of SUHII differed considerably in urbanized areas of different sizes,  ranging from 0.36°C to 3.22°C from petty urbanized areas to mega urbanized areas. However, at night, the distinction between the different sizes of urbanized areas was minor, within 1°C. The average amplitude of SUHIs was remarkably asymmetric across urban sizes. Our analysis indicated that the strongest SUHIs occurred during the day in the largest urbanized areas, the mega urbanized areas. However, the nighttime SUHIs did not precisely follow the urban scale.

Drivers for the spatial distribution of SUHII
The coefficients of the spatial regression model indicated that among the anthropogenic drivers, urban size, night lights, and population density had a strong positive correlation with SUHI. In contrast, CI negatively correlated with SUHII at day and night (Figure 8). The larger the impervious surface area and the higher the NL index, the more intense the SUHI; however, the more compact the impervious surface boundaries, the lower the heat island intensity in the urbanized areas. A nonsignificant negative correlation was observed between SUHII and GDP during the day, while a negative correlation was observed during the night. Meanwhile, the negative impact of GDP on SUHII is more significant in winter. Among the natural factors, NDVI had the most apparent negative correlation with SUHII. Intra-annually, anthropogenic factors contributed more to SUHII  variations than land surface descriptors and background climate. However, the above results only explained the impact of individual drivers on SUHI. We selected the summer with the strongest SUHII for structural equation modeling (SEM). The results of structural equation modeling showed that interactions existed within the anthropogenic and natural drivers (Figure 9). In general, both anthropogenic and natural drivers had direct or indirect impacts on SUHII at night. However, there was no direct influence of GDP on the daytime SUHII. In the daytime SUHII driving model, urban size had the largest coefficient of direct contribution to the heat island, reaching 0.63. TEMP and NDVI had a direct negative contribution to the daytime SUHII. Also, they were influenced by NL and thus had an indirect contribution to the daytime SUHII. CI, the inner-city landscape indicator, only contributed indirectly to the daytime SUHII by affecting NL. POP, NL, and PREC positively directly contributed to the daytime SUHII. In the model driving the nighttime SUHII, the direct contribution of urban size to the SUHII was lower than the daytime, being 0.38. PREC and NDVI had a negative direct contribution to the nighttime SUHII. In contrast, TEMP did not play a significant role in the nighttime SUHII. POP, GDP,   and NL all directly affected the nighttime SUHII, while CI had a direct negative effect on it. Urban size had the strongest direct contribution to SUHII. We proposed that the total area of the urbanized areas is the most critical driver in explaining the variations in SUHII. As the urbanized area size increased, the SUHII for both day and night increased. However, it is unclear whether there is a threshold range for this relationship between SUHIs and urban area size. To explore the potential impact of urbanized area size on SUHIs, we further analyzed the comprehensive relationships between urban size and SUHII. The relationships between SUHII and urban size were significantly and positively loglinear ( Figure 10). Consequently, the larger the urbanized area, the stronger the SUHI during daytime and nighttime. However, the drastic increase in SUHIs due to the expansion of urban size was not infinite. As predicted, there was a gradient effect in the response of SUHIs to urbanization, and there was a certain threshold range. First, SUHII increased rapidly when the urban size grew from 0 to 200 km 2 , from petty to medium urbanized areas. Afterward, the growth rate of SUHII started decelerating when the urban size increased from 200 to 400 km 2 , which is from medium to the largest urbanized area. Ultimately, the heating effect of urban expansion flattened with the increasing urban area when the size was larger than 400 km 2 . Additionally, we found a stronger positive relationship between  SUHII and urbanized area size in spring and summer than in autumn and winter.

Discussion
This research indicates that China's total impervious surface area in China doubled during 2003-2015, consistent with other studies (Kuang et al. 2016). Our observations confirmed that the growth rate of impervious areas first increased and then decreased in China from 2003 to 2015. With rapid urbanization, the most drastic increase in the impervious surface area occurred during 2005-2010 and began to decline after 2010. In this study, the proportion of petty urbanized areas decreased, whereas large and medium urbanized areas remained constant, with a significant increase in super and mega urbanized areas. Super and mega urbanized areas are the first rapidly developing regions in a country (Dou and Kuang 2020;Kuang 2020). They are often considered as "hubs" of national development and are most vulnerable to policy preferences and resource concentration . Hence socioeconomic developments and policy orientations have an inequality effect on urban expansion (Kuang 2020).
Our results showed that urban expansion was unequal across various regions in China. Petty urbanized areas had the largest proportion of all regions except South China, which exhibited a downward trend. Notably, some regions are inherently well urbanized, particularly in South China, where the proportion of mega urbanized areas is over 50% and increases annually. We consider that as land urbanization proceeds, the petty urbanized areas were consolidated into large, super, or mega urbanized areas in various regions, implying less fragmented urbanization. We compared our results with those of other studies (Weng et al. 2019) and found that small cities experienced greater expansion than large cities in China and the United States over the past ten years. Particularly, petty and medium urbanized areas in western China should further develop compact and efficient urbanized patterns. Previous studies have also proved that conventional urban-rural boundaries are obsolete in modern densely populated regions (Kuang et al. 2016). With the development of urbanization, the distance between urban areas tends to decrease and even disappear as urban land expansion accelerates. The expansion of the impervious surface area can explain the variations in SUHIs during rapid urbanization. Larger urban areas tend to have large areas of impervious surfaces that can absorb more solar radiation and store anthropogenic heat release, resulting in higher LSTs (Stewart and Oke 2012). Impervious surface expansion caused by SUHIs cannot be estimated by city boundaries of administrative districts but should be addressed in the patch scale of the impermeable surface.
In our study, spatial and temporal trends of the SUHII across urbanized areas were systematically analyzed at the national scale. Solar radiation is strongest in summer owing to direct sun irradiation on the Tropic of Cancer (Zhao et al. 2020). It causes more sensible but less latent heat fluxes and leads to summer experiencing the strongest SUHII (Imhoff et al. 2010), which is similar to the results of other SUHI studies (Sun et al. 2019).And the weakest SUHI was present during the daytime in winter with a minimum value of −1.99°C.At the patch scale, the distribution of SUHII was not uniform because the standard deviation values were reasonably high. Due to many patches, the SUHII was lower than the studies mentioned above, which focused on large cities or major cities. In this study, we evaluated SUHII from the impervious surface area in all cities (Kalnay and Cai 2003). We also observed significant differences in SUHII among geographical regions ( Figure 5). The SUHII was stronger in higher latitude regions, such as Northeast and Northwest China. This phenomenon was due to the typical arid, semi-arid, and semi-humid temperate climate of Northeast and Northwest China (Si et al. 2022), where the seasonal variations in vegetation activity are the largest. It is consistent with previous findings that the UHI is stronger in cold-or high-latitude climates (Feng et al. 2021).
Our results confirmed the different underlying mechanisms between day and night (Figure 9). Specifically, the daytime SUHII is related closely to the background precipitation or air temperature. Inconsistent with the results of (Zhou et al. 2016), the temperature in this study showed a negative direct contribution to SUHII, which corresponded well with previous findings (Yao et al. 2018). We believed that the effect of precipitation and temperature on SUHI is mainly dependent on soil moisture content. Soils in rural areas have a higher water holding capacity than impervious surface areas. However, low precipitation and high temperatures during the summer daytime will reduce soil moisture content in rural areas (Kuang et al. 2022). LSTs rise rapidly in rural areas during the daytime, thus reducing the SUHII. The soil moisture is generally higher in southern regions than in the north due to the humid climate, which causes lower SUHII at night (Figure 5b). Comparatively, precipitation amplifies the difference in soil moisture, thus increasing daytime SUHII and decreasing nighttime SUHII. Vegetation can increase the latent heat fluxes through transpiration and have a cooling effect on surface temperature, thereby mitigating daytime SUHII (Fang, Wu, and He 2021). At the same time, the increase in vegetation cover amplifies the shading effect and reduces the heat stored during the day, thus reducing nighttime SUHII.
The urban size contributes positively and directly to daytime and nighttime SUHII. Night light represents anthropogenic heat emissions. This direct anthropogenic heat flux also helps maintain SUHII (Peng et al. 2012). Areas with higher population densities normally have more indirect energy consumption, such as buildings and pavements; hence, more energy is released . In particularly, we found that nighttime lights and population density had a stronger effect on SUHII during the daytime in summer, suggesting that heat release caused by air conditioning usage can increase SUHII. CI can measure the complexity of urban landscape shape, and the structural equation modeling indicated that CI contributes indirectly to the daytime SUHII by influencing energy consumption. Boundary complexity contributes greatly to the solar energy absorption and emission ability of urbanized areas. Anthropogenic heat release is an important factor affecting SUHI (Voogt and Oke 2003). Human activities directly affect SUHII variations, and they also act as the indirect driver of SUHII through feedback on the urban climate. Overall, urban areas with a larger scale, complicated boundaries, and high intensity of socioeconomic development have higher SUHII. Meanwhile, some studies demonstrated that threedimensional urban morphology also has an impact on SUHI formation, particularly urban configuration and morphological dimensions Fu and Weng 2018).
Previous studies also confirmed the regulation of urban size on SUHI . The impacts of urban size on SUHII are complicated and involve population and building densities. It is well known that urbanization has significantly enhanced SUHIs; however, research on the thermal environmental impact of urban size is insufficient (Peng et al. 2012). Our study identified that spatial and temporal variations in SUHIs are extremely significant at different urban sizes from the patch scale. SEM revealed that urban area size makes a direct contribution to the SUHI, more than other drivers in both the daytime and nighttime. We found that the relationship between SUHII and urban size is significantly and positively log-linear. Land urbanization increased LSTs on an annual mean scale, regardless of the urbanized area size. However, the increase in SUHII caused by urban expansion is not unlimited. Petty urbanized areas, which are smaller than 50 km 2 , do present SUHIs, although they are weaker than in other types of regions. Petty and medium urbanized areas, which are the smallest in size, had the lowest SUHII, which confirms the enhanced heating effect of urban expansion (Kuang et al. 2022). A gradient effect exists in the response of SUHIs to urbanization. With the increasing size of urbanized areas, SUHII exhibits a logarithmic increase and reaches stability after a certain size. When the urbanized area is larger than 400 km 2 , this heating effect may not be sufficiently significant or fade away, particularly at nighttime. Conversely, the logarithmic growth relationships of SUHI's response to urban size were more significant during the daytime (Figure 10). There are several explanations for this: first, the increase in the density of urban infrastructure and human activities with the increase in urban size is not infinite (Georgescu et al. 2014). Second, the cooling effect is reduced because small cities lose more urban green space than large cities during urban expansion (Yao et al. 2017;Dong et al. 2022). Third, the increased SUHIs due to urban expansion could alter the local wind system, thereby mitigating the heating effect of urban expansion (Brazel et al. 2007). To mitigate the SUHIs increase, we advocate that the extent of the impervious surface area should be greater than 400 km 2 when cities are developing.
Our study had some limitations that should be acknowledged. Although we processed the cloud cover of remote sensing images in the data preprocessing, the LST values of some areas were still disturbed by cloud cover. We will try to employ the all-weather seamless LST product in a subsequent study. To match the temporal scale of all data, our study only investigated change patterns in the urban impervious surface area and SUHIs for four-time points : 2003, 2005, 2010, and 2015. Therefore, the time scale portrayal was not sufficient, and in the future, longer LST time series should be combined to study the variations in SUHIs per year. However, it was not easy to exploring the actual relationships for each variable, city, and period. In addition, while correlation analyses between SUHII and associated drivers can be conducted across space and time, we only analyzed the relationships across space. The temporal differences in the influencing factors should be investigated in future studies to develop policies to cope with SUHIs under different urbanization processes. Additionally, only macro factors were considered in this study, and future studies could select factors such as the shape, number, and structure of impervious surface patches at the landscape scale.

Conclusion
This study systematically quantified the dynamics of impervious surface areas and the spatiotemporal patterns of SUHIs in China from 2003 to 2015 using multisource remote sensing data. It also constructed a framework for assessing the underlying contributions of drivers to SUHIs using the spatial autoregressive model and structural equation model, and further investigated the relationship between urban size and SUHIs.
Our results showed that impervious areas had an increasing trend over the period 2003-2015, with accelerated growth from 2003 to 2010. And the expansion of petty urbanized areas was the dominant contributor to urban expansion. SUHII was more intense in summer than in other seasons and higher during the day than at night. Furthermore, SUHIs existed in low-ranking urbanized areas (i.e. petty urbanized areas). SUHI ranged from 0.36°C to 3.22°C for different sizes of the urbanized areas. The results also showed that urban size was the most crucial driver in explaining SUHII variations. In addition, the dominant drivers of daytime SUHII were background temperature and NDVI, while the dominant drivers of nighttime SUHII were population density and precipitation. It should be noted that natural factors are also influenced by anthropogenic factors as mediators and thus contribute indirectly to SUHII change. There are interactions between the anthropogenic and natural drivers, which complicates the mechanism of SUHI variations. It is undeniable that vegetation has an evident cooling effect. Urban green spaces should be increased and anthropogenic heat emissions should be reduced to mitigate SUHIs. Finally, we observed a statistically significant and positive loglinear relationship between SUHIs and urbanized area size. The heating effect of urban expansion flattened out with increasing urban areas when the urbanized area was larger than 400 km 2 . Hence, we consider that the impervious surface area of the city should exceed 400 km 2 during urban development.
These results may have policy implications for city managers. It is necessary to reduce the fragmented distribution of small towns to confront SUHI effects. Our results illustrate the spatial correlation between the surface temperature and urban development situation, as well as the cooling effects of the green space, which have important implications for urban planning to mitigate SUHI effects.