Urban Expansion in China: Spatiotemporal Dynamics and Determinants

Abstract

China’s urbanization has attracted many scholars’ attention due to its significant impact on socioeconomic sustainability. Many studies have explored the spatial pattern and effects of the factors influencing urban expansion. However, the spatiotemporal dynamics integrating spatial and temporal dimensions and the spatial scales of the influencing factors are always ignored. This study applied the framework of exploratory space–time data analysis (ESTDA) to investigate the spatiotemporal dynamics of urban expansion across 342 cities in China from 1990 to 2017 and, further, used multiscale geographically weighted regression (MGWR) to estimate the effects of influencing factors on urban expansion. We found that urban expansion had an obvious south–north division, and yet the effects of influencing factors usually showed an east–west division. We also found that the dynamic local spatial dependency of urban expansion was accompanied by a volatile coevolution process and inclined to transfer from heterogeneity to homogeneity, and homogeneity tended to be stable. The coevolution of urban expansion between cities and other neighboring ones became stronger with increases in time and regional integration. These findings support the use of customized urban planning for specific regions in different spatial dependence to improve land-use efficiency and coordinate regional development.

1. Introduction

China has been experiencing dramatic urbanization since the intensifying economic reforms of the early 1990s. Urbanization has greatly changed human settlement and the aggregation of resources and population, promoting economic development and the expansion of urban areas. Therefore, recent studies have considered urbanization as the fourth transition in China [1,2]. Meanwhile, urban expansion, fulfilling the demands of production, and living significantly impacts future sustainability [3]. Urban expansion has occupied large natural areas [4,5,6], and the loss of natural areas threatens terrestrial biodiversity [7] and ecosystem services [8,9,10]. Urban sprawl also worsens people’s living environments via air pollution [11,12], water pollution [13,14], and urban heat islands [15,16], which endanger people’s health. Many aspects of urbanization, however, remain unexplored [17]. Clarifying the spatiotemporal dynamics of urban expansion and its influencing factors are fundamental to improving land-use efficiency and coordinating regional development.

Many studies have explored the spatial patterns of urban expansion, and they have mostly used exploratory spatial data analysis (ESDA) methods, such as Getis’s Gi* [18], Moran’s I [19], the local indicator of spatial association (LISA) [20], and the standard deviational ellipse (SDE) [21]. These methods have revealed the spatial aggregation and heterogeneity of urban expansion. For example, Wu et al. [22] discovered that China’s urban expansion showed a dramatic club convergence and that the east–west division was obvious. Xu et al. [23] applied the SDE to investigate the principal orientation and direction of urban expansion in Guangzhou (China) and found that the directions of urban expansion varied with scales. As stated above, the ESDA methods provide effective tools with which to explore the spatial pattern of urban expansion, but struggle to convey the properties of spatial correlations and the time dimension [24]. Exploratory space–time data analysis (ESTDA), however, integrates the time dimension into spatial analyses and visualizes spatiotemporal interaction via graph theory [25,26]. In addition, ESTDA methods have been used in many fields, such as poverty reduction [27], carbon footprint [28], and air pollution [29].

Multiple factors influence urban expansion, and the identification of influencing factors can reveal urban expansion mechanisms and provide essential information with which to coordinate regional urbanization. Many researchers have investigated the determinants of urban expansion based on natural and economic factors. Natural factors, such as elevation, topography, and slope, determine the starting point of economic development, and the long-term invariant promotes or restricts urban expansion [30,31]. Social and economic development are the internal forces of urban expansion, and the aggregation of population and resources requires a sufficient urban area. Many studies have considered population and gross domestic product (GDP) as the basic socioeconomic factors, because populations are the dominant participants in society and the economy [22,32] and GDP is associated with many aspects of socioeconomic development [33].

The effects of influencing factors on urban expansion, however, seem to vary across studies. Some of the literature has attributed the varying effects to spatial heterogeneity and compared the differences between regions, urban sizes, and administrative levels [32], but the influencing factors’ spatial scale is always omitted. The heterogeneity of urban expansion is the result of the superposition of different factors with different influencing scales. Geographically weighted regression (GWR) considers the spatially varying effects of influencing factors [22] but assumes that all factors operate at the same scale. Multiscale geographically weighted regression (MGWR) relaxes the assumption, allows each factor to work at a specific scale, and improves GWR to perform a more reliable estimation [34].

This study aims to explore the spatiotemporal dynamics of urban expansion in China and examine the effects of the influencing factors. This paper is divided into five sections. The Section 2 describes the data and methods. The Section 3 presents the major findings on the spatiotemporal dynamics of urban expansion and the spatial heterogeneity of the influencing factors. The Section 4 compares the major findings with other studies, proposes policy implications, and discusses their limitations. The Section 5 offers a short conclusion.

2. Data and Methods

This study’s framework is shown in Figure 1. Based on urban expansion datasets containing 342 cities from 1990 to 2017, we first applied the ESTDA framework to investigate the dynamics of the spatial pattern of urban expansion, reveal the transitions of local spatial structures, and visualize spatial coevolution between cities and their neighbors. Then, we applied MGWR to estimate the effects of various influencing factors on urban expansion.

2.1. Data Sources

2.1.1. Data of Urban Areas

The urban areas in existence between 1990 and 2017 were selected from a long-term (1978–2017) data package created by Gong et al. [35]. This data package applies reliable impervious surface mapping algorithms to Landsat imagery and night-time light (NTL) data based on the Google Earth Engine (GEE) platform. The data resolution was 60 m in 1978 and 30 m during the 1985–2017 period, and the overall accuracy exceeds 90%.

We extracted the urban area of 342 cities at the prefecture level and above and separated the entire study period into three subperiods (1990–2000, 2000–2010, and 2010–2017). The urban expansion rate (UER) eliminates the effect of city size and was used in ESTDA. The following equation was used to calculate UER:

$$UER=((Are{a}_{ed}/Are{a}_{st})-1)/(ed-st)$$

(1)

where $Are{a}_{st}$ and $Are{a}_{ed}$ refer to the urban area at the start and end times, respectively.

2.1.2. Natural and Socioeconomic Data

The influencing factors include natural factors and socioeconomic factors. For natural factors, elevation, topography, and slope were used. The topography was represented by the standard error of elevation, and the mean values of the three indicators reflect the natural characteristics of each city. For socioeconomic factors, population and GDP were used. The two indicators denote the living demand and agglomeration capability of each city. Both natural and socioeconomic factors were derived from the raster dataset, and

Table 1. Influencing factors of urban expansion.

Category	Variable	Description	Sources
Physical factors	Elevation	The dataset has a 1 km resolution, and the unit is the meter.	The dataset is provided by China Resources and Environmental Science Data Center (http://www.resdc.cn/, accessed on 24 February 2022).
	Topography	The standard error of elevation.
	Slope	The dataset has a 1.8 km resolution, and the unit is the degree.	The dataset is provided by National Tibetan Plateau Data Center (http://data.tpdc.ac.cn, accessed on 24 February 2022).
Socioeconomic factors	Population	The dataset has a 1 km resolution at 1995, 2005, and 2015, and the unit is 10,000 people/km2.	The dataset is provided by China Resources and Environmental Science Data Center (http://www.resdc.cn/, accessed on 24 February 2022).
	GDP	The dataset has a 1 km resolution at 1995, 2005, and 2015, and the unit is RMB 10,000/km2.	The dataset is provided by China Resources and Environmental Science Data Center (http://www.resdc.cn/, accessed on 24 February 2022).

provides detailed descriptions of these variables with corresponding sources. Population and GDP from 1995, 2005, and 2015 represent the three periods (1990–2000, 2000–2010, and 2010–2017), respectively.

2.2. Methodology

2.2.1. Exploratory Space–Time Data Analysis (ESTDA)

The ESTDA is a collection of spatial data analysis methods and mainly includes the Global Moran’s I index, LISA, LISA time path, space–time transition, and correlation network methods [25,26]. The Global Moran’s I index and the LISA index measure spatial autocorrelation on the global and local scales, respectively. The geometric features of the LISA time path depict local spatial dependence and the spatial volatility of urban expansion. The space–time transition describes the transfer characteristics of urban expansion. In addition, the correlation network is used to visualize the local spatial coevolution of urban expansion. (1) Moran’s I and LISA

Tobler’s first law of geography indicates that “everything is related to everything else, but near things are more related to each other” [36]. This is spatial autocorrelation. We applied Moran’s I [19] and LISA [20] to explore the global and local spatial autocorrelation of urban expansion, respectively.

The Moran’s I is given by

$$\mathrm{Moran}’\mathrm{s}\mathrm{I}=\frac{N}{{\sum }_{ij}{w}_{ij}}\frac{{\sum }_i{\sum }_j{w}_{ij}\left(x_i-\overline{x}\right)\left(x_j-\overline{x}\right)}{{\sum }_i{\left(x_i-\overline{x}\right)}^2}$$

(2)

where $x_i$ and $x_j$ represent the urban expansion rate of city i and city j; $\overline{x}$ is the mean value of urban expansion of all cities; $w_{i,j}$ is a weight matrix, which considers the spatial effects (Queen adjacent was used in this study); and n denotes the city numbers. When the global Moran’s I has a positive value, the cities and their neighbors have a positive spatial association; otherwise, there is a negative spatial association.

LISA is given by

$$I_i=z_i\sum _j{w}_{ij}{z}_j$$

(3)

where $z_i$ and $z_j$ represent the standardized values of the urban expansion of city i and city j, respectively, and $w_{i,j}$ is the row-standardized contiguity matrix. Four clusters on a Moran scatter plot can be divided based on $I_i$.

(2) LISA time path

The LISA time path [25] adds the temporal dimension to the LISA analysis to reveal the migration of spatial elements in the Moran scatter plot. The pairwise movement of the variable and its spatial lag suggest the dynamics of local spatial dependence. The path of the city i over time can be represented as [$(y_{i,1},y{l}_{i,1}),(y_{i,2},y{l}_{i,2}),\dots ,(y_{i,t},y{l}_{i,t})$]. $y_{i,t}$ is the urban expansion of city i at time t, and $y{l}_{i,t}$ is city i’s spatial lag at time t. The geometric features of the LISA time path include relative length and tortuosity.

The relative length can be represented as follows:

$${\mathsf{\Gamma }}_i=\frac{N\ast {\sum }_{t=1}^{T-1}d\left(L_{i,t},L_{i,t+1}\right)}{{\sum }_{i=1}^N{\sum }_{t=1}^{T-1}d\left(L_{i,t},L_{i,t+1}\right)}$$

(4)

where $L_{i,t}$ represents the coordinate ($y_{i,t},y{l}_{i,t}$) on the Moran scatter plot, implying the location of the city i at time t. $d(L_{i,t},L_{i,t+1})$ is the distance between location $L_{i,t}$ and $L_{i,t+1}$. N is the number of cities. When a city’s moving distance is longer than average, the relative length is larger than 1.

The tortuosity can be represented as follows:

$${\mathsf{\Delta }}_i=\frac{{\sum }_{t=1}^{T-1}d\left(L_{i,t},L_{i,t+1}\right)}{d\left(L_{i,1},L_{i,T}\right)}$$

(5)

where ${\mathsf{\Delta }}_i$ refers to the city i’s tortuosity on the Moran scatter plot over time. A larger ${\mathsf{\Delta }}_i$ indicates a more curved path and cities and their neighbors with more incongruity in terms of urban expansion.

(3) Space–time transition

The space–time transition explores the conversion of the local spatial dependence forms of urban expansion over time on the Moran scatter plot [25]. The Moran scatter plot has four classes: High–High (first quadrant), Low–High (second quadrant), Low–Low (third quadrant), and High–Low (fourth quadrant). There are 16 transitions among the four classes, which can be classified into four types (

Table 2. Space–time transition types.

Types	Description	Transitions
I	Self-transition with neighborhood stabilization	$H{H}_t\to L{H}_{t+1},L{H}_t\to H{H}_{t+1},H{L}_t\to L{L}_{t+1},L{L}_t\to H{L}_{t+1}$
II	Self-stabilization with neighborhood transition	$H{H}_t\to H{L}_{t+1},L{H}_t\to L{L}_{t+1},H{L}_t\to H{H}_{t+1},L{L}_t\to L{H}_{t+1}$
III	Self-transition with neighborhood transition	$H{H}_t\to L{L}_{t+1},L{H}_t\to H{L}_{t+1},H{L}_t\to L{H}_{t+1},L{L}_t\to H{H}_{t+1}$
IV	Self-stabilization with neighborhood stabilization	$H{H}_t\to H{H}_{t+1},L{H}_t\to L{H}_{t+1},H{L}_t\to H{L}_{t+1},L{L}_t\to L{L}_{t+1}$

).

We further defined a spatial homogeneous tendency indicator (SHTI) to measure the tendency to shift from heterogeneity (HL and LH) to homogeneity (HH and LL) during transition. The SHTI can be represented as follows:

$$SHTI=\left(Hm-Ht\right)/2$$

(6)

where Hm refers to the probability of transition from heterogeneity (HL and LH) to homogeneity (HH and LL), and Ht expresses the opposite transition probability. The SHTI value ranges between −1 and 1. When SHTI is greater than 0, this indicates that transitions have a homogenous tendency. When SHTI is less than 0, this suggests that transitions have a heterogeneous tendency.

(4) Correlation Network

Spatial coevolution can be represented by graph theory [25]. The correlation coefficient between two neighboring cities’ urban expansion was adopted to measure coevolution, linking the centroids of cities to construct the spatial coevolution network of urban expansion. The correlation links are mapped by a specified criterion (

Table 3. Coevolution strength.

ID	Correlation Coefficient	Coevolution Strength
1	−1∼−0.5	Strong decoupling
2	−0.5∼0	Low decoupling
3	0.1∼0.5	Low coevolution
4	0.5∼0.8	Moderate coevolution
5	0.8∼1	Strong coevolution

), so regionally shared temporal dynamics can be identified.

2.2.2. Multiscale Geographically Weighted Regression (MGWR)

Multiscale geographically weighted regression (MGWR) relaxes the assumption that all variables should be at the same spatial scale, as in classical GWR, and allows for different variables at different spatial scales by deriving an optimal bandwidth vector [34]. This model can improve the credibility of estimation and provide valuable information on the scales of different variables. This study, therefore, applied MGWR to estimate the spatiotemporal heterogeneity of the factors influencing urban expansion, and the function is expressed below:

$$y_i=\sum _1^k{\beta }_{bwj}({\mu }_i,v_i)x_{ij}+{ϵ}_i$$

(7)

where $y_i$ indicates urban expansion in city i, $x_{ij}$ indicates the influencing factor j in city i, and ${\beta }_{bwj}({\mu }_i,v_i)$ denotes the coefficient of the city i for location $({\mu }_i,v_i)$, in which bwj is the optimal bandwidth of influencing factor j.

The local parameter ${\beta }_{bwj}({\mu }_i,v_i)$ is estimated by specifying spatial weight matrices, and the estimate can be described as follows:

$$\widehat{{\beta }_{bwj}}({\mu }_i,v_i)={\left(X_j^{\prime }{W}_{bwj}({\mu }_i,v_i)X_j\right)}^{-1}{X}_j^{\prime }{W}_{bwj}({\mu }_i,v_i)y$$

(8)

where $x_j$ represents the influencing factor j, $W_{bwj}({\mu }_i,v_i)$ denotes the spatial weight of j at city i, and y is the observation of urban expansion.

In MGWR models, each influencing factor has a unique spatial scale, and the optimal bandwidth is detected via the corrected Akaike information criterion (AICc). We used MGWR software to calculate the results in this study (https://sgsup.asu.edu/sparc/mgwr, accessed on 24 February 2022).

3. Results

3.1. Spatiotemporal Pattern of Urban Expansion across China

On a national scale, urban areas increased 4.5 times, from 32,263 km² in 1990 to 146,073 km² in 2017 (Figure 2). We divided the 28 years into three periods and observed that the urban expansion rate accelerated (

Table 4. Urban expansion rate in the three periods.

	1990–2000	2000–2010	2010–2017
Average rate (%)	4.5	5.6	8.8
Area (km2)/year	1599	3346	7427

). The first period (1990–2000) had a fluctuating urban expansion rate, with an average yearly expansion of 4.5%, namely, 1599 km². The second period (2000–2010) witnessed the average annual acceleration of the expansion rate, at 5.6%, and enlarged urban areas of 3346 km² each year, more than double that of the first period. The third period (2010–2017) saw the highest expansion rate, at 8.8%, with 7427 km² of annual expansion area. Furthermore, the center of gravity of all urban areas continued to shift toward the southeast, with a strong southerly movement between 1990 and 2001, before slowing between 2001 and 2010, and remaining steady between 2011 and 2017.

We employed the Gini coefficient and Moran’s I to investigate the differences in spatial autocorrelation of urban expansion between cities. As shown in Figure 3, the Gini coefficient shows a downward trend over time, whereas the trend for Moran’s I is upward. In addition, all values of Moran’s I are positive. The two indices indicate that the gap in urban expansion is narrowing and the correlation between cities and their neighbors is increasing.

We mapped the urban expansion in three periods (1990–2000, 2000–2010, and 2010–2017) and explored the local spatial dependence by using LISA. As shown in Figure 4a, the south–north division in terms of urban expansion remains stable throughout the three periods, and aggregation patterns become more and more obvious. The southern region has a high urban expansion rate, and the Yangtze River Delta and the Sichuan Basin region lead the way in this regard. The northern region experiences low urban expansion, and the northwest region provided a typical example of such a low urban expansion. The LISA cluster maps (Figure 4b) identify the cold spots and hot spots, and urban expansion shows club convergence in China. The Low–Low cluster was located in the northern region, and the High–High cluster was transferred from the east coast to the southwest. These results suggest that urban expansion remained regionally unbalanced, although the Gini coefficient shows a decreasing trend (Figure 3).

3.2. Dynamics of the Local Spatial Dependence of Urban Expansion

The relative length (RL) of the LISA time path reveals the dynamic local spatial dependence, and the longer the relative length, the stronger the dynamics. As shown in Figure 5, the relative length was classified into four levels: low (RL < 0.6), lower middle (0.6 < RL < 1), higher middle (1 < RL < 1.4), and high (RL > 1.4), and the heterogeneity between the eastern and western regions is obvious and stable. The cities with high relative lengths were mainly aggregated in the western region, whereas cities with low relative lengths were mainly distributed in the east. The northeast grew from a low level to a high level over the three periods. In addition, the number of cities with a relative length greater than 1 decreased over time, with 126, 110, and 106 cities for each period. Overall, these results indicate that the local spatial dependence in the western region has more dynamics than that in the eastern region, and that the northeastern region has moved from stability to dynamics.

The tortuosity (To) of the LISA time path reflects the volatility of the local spatial dependence of urban expansion by detecting the change in dynamics between cities and their neighbors. The tortuosity was also classified into four levels: low (To < 4), lower middle (4 < To < 8), higher middle (8 < To < 16), and high (To > 16). The number of cities with high or higher middle tortuosity was 116, 130, and 106 for each period, respectively. Figure 5 shows that the cities’ tortuosity changed dramatically in the western region throughout the three periods and the high tortuosity expanded from the western border region in the first period to the western region in the second period, and then transformed into low tortuosity in the third period. For the eastern region, the pattern of mixing high and low tortuosity persisted over the three periods. In summary, the second period had the strongest volatility in the dynamics of local spatial structures, and the third period had the weakest volatility. From a spatial perspective, the western region experienced more volatility than the eastern region.

3.3. Transition of Local Spatial Dependence of Urban Expansion

The LISA transition matrix reveals the transfer characteristics of local spatial dependence among the four states (HH, HL, LH, and LL) (

Table 5. Spatiotemporal transition matrix, types and SHTI.

Period	t/t+1	HH	HL	LH	LL	Type	Type	SHTI
2000/2010	HH	0.61	0.05	0.13	0.2	I	0.2	0.32
	HL	0.17	0.48	0.06	0.29	II	0.14
	LH	0.32	0.06	0.38	0.23	III	0.09
	LL	0.03	0.08	0.11	0.78	IV	0.56
2010/2017	HH	0.36	0.08	0.21	0.35	I	0.28	0.37
	HL	0.13	0.38	0.02	0.47	II	0.14
	LH	0.36	0.02	0.36	0.26	III	0.14
	LL	0.15	0.09	0.1	0.66	IV	0.44

). In the first transition, from the first period to the second period, the transition from LL to LL had the maximum possibility of 78%, and the transition from HH to HH had the second-maximum possibility of 61%. Regarding transition types, the cities with a type IV transition (self-stabilization with neighborhood stabilization) had the highest proportion of 56%, indicating that the spatial dependence of urban expansion between cities and their neighbors was inclined to remain stable. As the spatial homogeneous tendency indicator (SHTI) was 0.32, which is larger than 0, the local heterogeneous spatial structures (HL and LH) tended to be homogeneous (HH and LL).

In the second transition, from the second period to the third period, the transfer characteristics were similar to those in the first transition. The highest probability of transition still occurred in the transition from LL to LL, with 66%, and the transition of type IV remained the most probable transition, at 44%. In addition, the SHTI was larger than 0, at 0.37. These results indicate that the spatial aggregation of urban expansion depends on the local stable spatial dependence, especially the Low–Low spatial structure, and that the local spatial dependence was inclined to transfer from heterogeneity to homogeneity.

3.4. Coevolution of Urban Expansion among Cities

The spatiotemporal correlation network illustrates the coevolution of urban expansion in neighboring cities (Figure 6). The network reflects the fact that the coevolution has become strong, and the proportion of cities with a strong pairwise correlation was 0.16, 0.19, and 0.27 in the three periods, respectively. In the first period, strong coevolution was mainly distributed on the east coast, and the difference between the east coast and the inland was obvious. In the second period, the coevolution strengthened in the whole of China, especially in central China. This indicates that cities and their neighbors tend to be integrated. The third period saw strong coevolution, aggregated in the Yangtze River Delta, the Pearl River Delta, and the Chengdu-Chongqing region, and this pattern supports the national plan for developing major functional zones [37].

3.5. Factors Influencing Urban Expansion

The MGWR model was used to estimate the influencing factors of urban expansion.

Table 6. The model fit metrics in different periods.

Model Index	1995	2005	2015
R-squared	0.55	0.62	0.7
AICc	715.33	656.98	573.68
RSS	155.5	131.11	102.76

reports the model fit metrics, and the R-square values in 1995, 2005, and 2015 are greater than 0.5. Figure 7 depicts the spatial heterogeneity of each factor’s coefficients, with most patterns showing an east–west division.

The GDP coefficients indicate that the correlations between economies and urban expansion were positive in the three years and exhibit an obvious east–west division (Figure 7). Economic development had a greater impact in the west than in the east. The disparities in GDP impacts between the east and west could result from China’s regional development policies. The Great Western Development, beginning in 1999, attracted a large amount of investment and improved infrastructure in the western region. These measures encouraged economic development and urban expansion in the west. The coefficients were lower in northeast China compared with the other regions because the urban system there is mature and stable, and the Revitalization of Northeast China, proposed in 2004, aimed to reverse the development path dependency on heavy industry and restructure the industries in that area. The bandwidth of GDP used in the MGWR model was 340 in 1995 and 44 in 2005 and 2015, demonstrating that the stable scales have become small and economic effects tend to include spatial aggregation.

The estimated population maps suggest that stratified heterogeneity significantly increased over time and the coefficient pattern of heterogeneity being high in the east and low in the west became clear in 2015. The significant coefficients were concentrated in the central region in 1995, and the population was aggregated in this region. The significant coefficient expanded to the east and center of China in 2005 because the income gap between the eastern and western regions drove large numbers of the population to migrate from the west to the east. The Rise of Middle China strategy also strengthened the aggregation capacity of the central region. Some western cities exhibited a negative relationship between population and urban expansion in 1995 and 2005, but the coefficient was insignificant. This, however, reflects the fact that urban expansion in these cities was driven by economic development rather than population aggregation in this period. The bandwidth of the population increased in 2015, indicating that the influencing scale expanded. This change was accelerated by inland infrastructure improvements and population migration.

Many cities’ coefficients for natural factors were insignificant. For elevation, cities in the western regions had significant positive coefficients, but cities in the eastern regions had insignificant negative coefficients. The distinction between the east and the west was obvious and stable. The great western development strategy promoted urban expansion in high-altitude areas, and this contributed to positive relations in western regions. For topography, the coefficients were usually negative and insignificant. As indicated by the bandwidth, the influencing scale of topography remained stable. Regarding slope, the coefficients were positive in all three years but only significant in 2005. The slope does not limit urban expansion, especially when flat land is protected for cultivation or fully exploited. A noteworthy illustration is that urban areas severely expanded in the hilly southeast area.

4. Discussion

This study investigated the spatiotemporal dynamics of urban expansion across 342 cities from 1990 to 2017 in China and further applied MGWR to characterize the spatial pattern of multiple natural and socioeconomic factors. We discovered that urban expansion shows a south–north division, and regions with a high urban expansion are clustered in the Yangtze River Delta and the Sichuan Basin regions. Nevertheless, the spatial patterns of influencing factors usually show an east–west division. In addition, the local spatial dependence of urban expansion between cities and their neighbors tends to be stable and is inclined toward homogeneity. The number of cities with pairwise coevolution in terms of urban expansion increased throughout the three periods, and local integration became obvious.

4.1. Contributions of This Study

Regarding the spatial pattern of urban expansion, this study is consistent with other studies [22,30]. For example, regions with severe urban expansion were identified in the Yangtze River Delta and the Sichuan Basin. However, this study further explores the spatiotemporal dynamics of urban expansion. Although some studies compared the spatial patterns between different periods [22,31], the dynamics and transfer characteristics were ignored. This study focuses on the spatiotemporal dynamics of urban expansion to compensate for the lack of other studies on this topic.

However, the effects of the influencing factors estimated in our study are inconsistent with those of previous studies. For example, we found that economic development always encourages urban expansion, but Wu et al. [22] report that, in some cities, there is a negative relationship between economic development and urban expansion. These differences can be attributed to the model and dataset used. Wu’s study applied GTWR with a fixed bandwidth and only included 286 cities, leaving out numerous western cities.

4.2. Mechanisms of Urban Expansion

Natural and socioeconomic factors influence the spatiotemporal dynamics of urban expansion, and the bandwidths of influencing factors reflect the influence scale change across periods. Natural factors are the substrate of city development, and their influence scales are stable despite economic development (Figure 7). The heterogeneity of natural factors introduces the original hierarchy across cities [38,39,40]. The original hierarchy induces the flow of resources and population, enhancing cities’ interaction [41,42] and leading the regions to convergence according to neoclassical economics [43,44]. This is why the local spatial dependence of urban expansion is inclined to transfer from heterogeneity to homogeneity. In this process, as indicated by the change in bandwidths, the regional economic center is formed, and population migration increases.

The spatial aggregation of coevolutionary cities in the network map is similar to that of the major function-oriented zones in China [37]. The spatial distribution reveals the relationship between urbanization, economic growth, and government planning. Urban planning usually focuses on individual cities [37,45]. The coevolution within spatial aggregations suggests that urban planning should go beyond individual cities [46]. China, however, lacks uniform standards to coordinate urban expansion across cities. Planning on a regional scale lags behind the development of urban aggregation [37]. Urban planning on a regional scale can lead to urbanization sustainability, which can improve resource efficiency and coordinate regional development [47].

4.3. Policy Implications

These findings serve as a guide for Chinese policymakers. Some cities saw a dramatic fluctuation in urban expansion as a result of disruptions in land-use planning and implementation. Strictly implementing relevant land-use planning may ensure orderly urban land expansion and improve land-use efficiency. Regional integration has become a trend. The customized urban expansion pattern requires targeted policies, and these customized policies improve upon the measures enacted by a unified policy [48]. Policies focused on individual cities may not achieve optimal results, and a comprehensive plan based on a regional scale must be considered.

5. Conclusions

China has been experiencing rapid urbanization since the deepening economic reforms of the early 1990s, and the urban expansion presents spatiotemporal heterogeneity. In order to investigate the spatiotemporal dynamics of urban expansion in China and estimate the effects of natural and socioeconomic factors on urban expansion, we first explored the spatiotemperal pattern of urban expansion across 342 cities from 1990 to 2017, then investigated the dynamical features with ESTDA. We further applied MGWR to characterize the spatial pattern of the relationship between multiple natural and socioeconomic factors and urban expansion.

This study found that the urban expansion rate accelerated and that the south–north division of urban expansion remained stable. The total urban areas increased from 32,263 km² in 1990 to 146,073 km² in 2017, and the urban expansion rate grew from 4.5% during 1990–2000 and 5.6% during 2000–2010, to 8.8% during 2010–2017. In terms of spatial pattern, the southern region had a high urban expansion rate, whereas the northern region experienced low urban expansion. Furthermore, the movement of the center of gravity of all urban areas continued to shift toward the southeast, particularly with a southerly movement.

Regarding the dynamic features of urban expansion, the spatial pattern of relative length and tortuosity of the LISA time path indicated that the western region had more dynamics and volatility than that in the eastern region. The transition matrix reveals that the cities with type IV transitions (self-stabilization with neighborhood stabilization) had the highest proportion (56% in the first transition and 44% in the second transition), which indicates that the spatial dependence of urban expansion is inclined to remain stable. The SHTIs were larger than zero in the two transitions, and local heterogeneous spatial structures tended to be homogeneous. In the correlation network, the proportion of cities with a strong pairwise association grew, and regional integration became apparent.

Different natural and socioeconomic factors have different effects on urban expansion, and most of them showed an east–west division. Specifically, the correlations between GDP and urban expansion were positive, despite the spatiotemporal variations. With population migration driven by the income gap between the east and the west, the influence of the population on urban expansion increased in the east. Natural factors were insignificant in many cities, and the limiting effects of natural factors on urban expansion were not significant.

Our empirical study clarified the spatiotemporal dynamics of urban expansion and validated the spatiotemporal heterogeneous effects of natural and socioeconomic factors in China. These findings could help with customized urban planning for specific regions in different spatial dependences and the coevolution of urban expansion. However, this research has certain drawbacks. First, urbanization is a complicated process and is influenced by many factors. In addition to the natural and socioeconomic factors, other factors, such as historical and political factors, should also be comprehensively considered. In addition, spatiotemporal dynamics at different urban sizes and development levels should be investigated in future work.