Evolution Process and Land Use/Land Cover Response of Urban–Rural Space in Wuhan under Polycentric Structure

Abstract

Polycentric development facilitates urban–rural spatial reshaping and land use/land cover (LULC) protection. Previous studies have predominantly focused on urban areas, with spatial delineation methods biased towards the macro-level, lacking a holistic perspective that situates them within the urban–rural spatial framework. This study proposes a spatial delineation framework that is applicable to the polycentric structure, taking into account the social, economic, and natural characteristics of urbanization. It employs semivariance analysis and spatial continuous wavelet transform (SCWT) to analyze the effects of polycentric development on the urban–rural space of Wuhan from 2012 to 2021 and applies a land use transition matrix, landscape indices, and bivariate spatial autocorrelation to quantify the responses and differences of LULC within urban–rural space. The results indicate that $600\mathrm{m}\times 600\mathrm{m}$ is the best scale for exhibiting the multidimensional characterization of urbanization. The polycentric structure alleviates the compact development of the central city, and it drives rapid expansion at the urban–rural fringe, exacerbating the spatial heterogeneity in LULC change pattern, spatial configuration, and urbanization response within urban–rural spaces. The overall effects of urbanization on LULC are relatively weak along the urban–rural gradient, experiencing a transition from positive to negative and back to positive. This study employs a novel spatial delineation framework to depict the polycentric transformation of metropolitan areas and provides valuable insights for regional planning and ecological conservation in the urban–rural fringe.

1. Introduction

Rapid urbanization has exacerbated urban–rural spatial conflicts and posed significant challenges to environmental conservation. Since the 21st century, China has experienced unprecedented urbanization, with the urbanization rate reaching 64.72% in 2021 [1], and is projected to further increase to 80% by 2050 [2]. Urban expansion has encroached on a substantial amount of arable and ecological land, resulting in chaotic urban–rural spatial patterns and difficulties in effectively controlling development disparities. As of 2021, the urban–rural per capita income ratio remains as high as 2.5. Urban expansion has also caused drastic changes in land use/land cover (LULC), triggering issues such as the heat island effect, land degradation, and vegetation loss. Therefore, there is an urgent need to improve urban development strategies to foster urban–rural spatial reshaping and ecological spatial governance.

Polycentricity is a significant trend in urbanization, connecting previously adjacent but independent urban residential areas into a comprehensive urban system. This integration fosters a more balanced spatial structure and reduces agglomeration diseconomies, contributing to the resolution of the aforementioned issues [3]. As a spatial strategy for sustainable development, polycentricity has been widely applied in the planning and design of city agglomerations and town clusters, prompting an increasing number of cities to shift towards a polycentric structure. This transformation is especially obvious in metropolitan areas such as Wuhan City, which is manifested as the development and nurturing of edge cities and satellite towns, along with their functional linkage and spatial integration [4]. Wuhan is a core city in the Yangtze River Economic Belt development strategy. At the beginning of the 20th century, its urban area was primarily confined to the “two rivers and three towns”—Hankou, Hanyang, and Wuchang districts at the confluence of the Han and Yangtze Rivers. In the process of regional integration over the past two decades, it has gradually evolved into a polycentric and clustered metropolitan area, with its economic output rising to the eighth largest among Chinese cities. Nevertheless, it still faces challenges such as mixed urban land use, imbalanced urban–rural development, and ecological fragmentation. Undoubtedly, polycentric development enhances regional competitiveness but introduces uncertainties to the urbanization process. The government necessitates assessing the implementation effects of polycentric strategies, including their impact on the urban–rural spatial morphology and the functionality of LULC [5]. Therefore, dynamically delineating the urban–rural space of Wuhan and analyzing the interactions and differences of LULC within these spaces contribute to a comprehensive understanding of polycentric practices and urban–rural transformations in metropolitan areas.

With the rapid advancement of GIS and RS technologies, numerous delineation methods have been proposed for different city spatial structures. In the monocentric model, urban–rural fringe has received more attention due to its dynamic complexity. Relevant studies have shifted from qualitative judgment to quantitative analysis, primarily employing methods such as comprehensive analysis [6], threshold setting [7], mutation detection [8], and spatial clustering [9], combined with land use, nighttime lights, and point of interest data [10,11]. These methods have been proven to have good identification performance, but they remain to be validated in complex polycentric city scenarios. In the polycentric model, more attention is paid to the location and structural identification of urban centers [12], which primarily utilize remote sensing data and employ methods such as minimum cut-off point (MCOP) [13], parametric and non-parametric estimation [14], and spatial statistics [15]. These studies are predominantly conducted at the meso–macro-scale, where the setting of urban–rural space is based on the percentage of impervious surfaces within grid cells [16] and buffer distances [17] (as urban–rural space is typically correlated with urban development intensity). Moreover, in the case analysis of specific regions, urban–rural space is often artificially defined [18]. Overall, there is a lack of a more comprehensive and detailed approach to delineate urban–rural space in polycentric structures. Additionally, another key issue is the “modifiable areal unit problem”, since these studies are dependent on the grid scale and constrained by the coarse resolution of remote sensing data [17]. The spatial representation of urban development is affected by scale effects, and, prior to conducting spatial delineation, it is essential to quantify the dependency on grid cells. Additionally, it is necessary to fully leverage the instrumental role of geographic big data for spatial structure research at fine scales, adopting higher precision and more multidimensional indicator data.

Numerous studies have quantified the urbanization impact on LULC by developing a suite of indicators that characterize urbanization, land use, and landscape patterns [19,20,21], combined with statistical methods such as linear regression [22] and spatial autocorrelation [23]. The research focus is expanding from urban to suburban areas [24], as they exhibit more direct, intense, and relevant responses to urbanization. Nonetheless, few studies have analyzed the effects of urbanization across different urban–rural spaces and their development contexts. The urbanization process involves the entire range of urban–rural spaces, even rural areas, whether the impacts are direct or indirect, linear or non-linear [25,26]. It is necessary to explore the effects and differences of LULC driven by urbanization within a broader spatial framework.

In summary, this study aims to analyze the shaping process of urban–rural space by polycentric development and quantify the responses and variations of LULC within this space. The main contributions are as follows: (1) Proposing a spatial delineation framework suitable for a polycentric structure to reveal the expansion process and spatiotemporal patterns of urban–rural space; (2) Analyzing the spatiotemporal dynamics of LULC in different urban–rural spaces and quantifying their urbanization responses. This study takes Wuhan as the research object. Firstly, considering the social, economic, and natural characteristics of urbanization, it constructed an identification framework based on multi-source data, used semi-variogram analysis to determine the optimal grid scale, and combined principal component analysis and spatial continuous wavelet transform (SCWT) to determine the urban–rural spatial pattern and its evolution process from 2012 to 2021. Then, a land use transfer matrix and landscape index are used to analyze the change models of LULC within the urban–rural space, and ecological risk index and bivariate spatial autocorrelation are employed to quantify the effects and variations of urbanization on LULC. This study evaluates the performance of other delineation methods in the polycentric structure, elucidates the positive contributions of planning and construction to polycentric growth, and encapsulates the spatiotemporal heterogeneity of LULC along the urban–rural gradient. It provides valuable insights into the polycentric transformation of metropolitan areas and the ecological preservation of peripheral areas.

2. Materials and Methods

2.1. Study Area and Data

Wuhan (29°58’∼31°22’ N, 113°41’∼115°05’ E) is the geographic center of the Yangtze River Economic Belt and the capital of Hubei Province. Wuhan covers an area of 8569.15 km², with jurisdiction over 7 main-urban districts and 6 distant-urban districts (Figure 1). In 2021, it hosted a permanent population of 13.65 million, with an urbanization rate of 84.56%, and its GDP reached CNY 1.77 trillion. As the development core of the Yangtze River Midstream Urban Agglomeration, Wuhan has experienced continuous adjustments in land use and spatial planning policies. With a strengthened agglomeration effect of population and economic factors, as well as a noticeable trend towards a multi-centered and clustered spatial layout, it has gradually developed into an internationalized metropolitan area, leading to significant changes in its spatial form, geographical landscape, and urban–rural structure. Therefore, understanding the development and evolution of Wuhan’s urban–rural space is of great significance for the territorial spatial optimization and rural development in the Wuhan metropolitan area.

This study mainly employs multi-source data, land use data, and statistical data, collected in 2012, 2016, and 2021. Among them, the Point of Interest (POI) data are obtained through the Amap API (https://www.amap.com) on 1 June of the corresponding year, which mainly involves the tertiary industry (catering, shopping, educational institutions, leisure and entertainment, life services, automobile services, hotels, finance) and the secondary industry (factories and mines, industrial parks, agriculture, forestry, fishing- and animal husbandry-related industries, furniture manufacturing, high-tech industrial parks) [27]; the road network data are from OpenStreetMap (http://www.openstreetmap.org), which are mainly extracted from primary, secondary, and tertiary roads; and the population data are obtained from the WorldPop program (http://www.worldpop.org), with a resolution of $100\mathrm{m}\times 100\mathrm{m}$. Land use data are derived from the Landsat-derived annual China land cover dataset (CLCD) [28], with a resolution of $30\mathrm{m}\times 30\mathrm{m}$. Among the statistical data, vector data such as administrative divisions and water systems are obtained from the National Geomatics Center of China (https://www.webmap.cn), and statistical yearbooks and bulletins are obtained from the Wuhan Bureau of Statistics (http://tjj.wuhan.gov.cn). For the pre-processing of multi-source data, the required industrial POI data are obtained by merging, filtering, and removing duplicate and ambiguous point data. The road network data are obtained by removing duplicated roads, converting double lines to single lines, and conducting topology checking and connectivity checking. Population data are assigned to study grids in different scales based on administrative divisions using the “Zoning Statistics” tool in ArcGIS (v 10.8). All data are projected to the WGS_1984_UTM_Zone_49 N coordinate system.

2.2. Study Methods

This study mainly encompasses three sections (Figure 2): urban–rural spatial division under the polycentric structure, a spatiotemporal evolution analysis of urban–rural space, and a comprehensive analysis of LULC.

2.2.1. Construction of Index System

Urbanization is the main force driving the remodeling of urban–rural space, which is subject to the interaction of social, economic, natural, and other factors. This study primarily characterizes urbanization under the polycentric mode from five dimensions: industry, economy, population, transportation, and land use (

Table 1. Index system for urbanization.

Dimension	Index	Description	Weight
Population	Population density	Population in the study grids	0.1893
Economy	NTL value	NTL value in the study grids	0.1858
Industry	POI density ratio	POI density ratio between tertiary and secondary industries in the study grids	0.1736
Land use	Land use intensity	Total use intensity of LULCs in the study grids	0.2096
Transportation	Road network density	Total length of roads in the study grids, including highways, as well as primary, secondary, and tertiary roads	0.2416

).

(1) Industry: there are spatial characteristics in the distribution of the secondary and tertiary industries, which are guided by policies such as urban planning and influenced by factors such as rents and pollution [27]. The service industry is concentrated in the urban polycentric areas, while the manufacturing industry is mostly scattered in its periphery, and, when urban expansion occurs, the latter moves outward before the former. Therefore, the POI density ratio between secondary and tertiary industries is used to reflect the industrial gradient differences in urban–rural spaces;

(2) Economy and Population: the nighttime light (NTL) value and population density are employed, which have been widely proven to exhibit strong correlations with regional productivity and human activities [29];

(3) Transportation: the improvement of transportation and other public infrastructure is a prerequisite for urban outward expansion, particularly in supporting polycenter development and the conduction of urban functions. Hence, the road network density is used to reflect the infrastructure construction and its resulting resource location effect;

(4) Land Use: the most direct manifestation of urbanization is the rapid land use/land cover change (LUCC), so the land use intensity (LUI) is adopted to reflect the differences in land development in urban–rural space [20], and its calculation formula is given in Equation (1). In this study, principal component analysis (PCA) is applied to determine the indicator weights, then linear weighting is used to obtain the urbanization attribute values. It reduces the multicollinearity between indicators, and the statistical results of PCA are shown in the Supplementary Materials.

$$LUI=\sum _{i=1}^4{\displaystyle \frac{L{A}_i}{{\sum }_{i=1}^{n}L{A}_i}}\times D_i$$

(1)

where $L{A}_{(i,t)}$ represents the area of land use type i, n denotes the number of land use types in the study area, and $D_i$ indicates the weighted value for each land use type. Unused land is assigned a weighted value of 1; construction land 4; forestland, grassland, water areas, and wetlands 2; while cultivated land is assigned as 3.

2.2.2. Determination of Grid Scale

The spatial representation of urbanization is affected by scale effects, and this study employs a combined quantitative and qualitative approach to determine the optimal grid scale under the MAUP. Among them, the quantitative method is semivariance analysis (

Table 2. Semi-variogram parameters at different grid scales.

Size (m)	${\mathit{C}}_0$	${\mathit{C}}_0+\mathit{C}$	${\mathit{C}}_0/({\mathit{C}}_0+\mathit{C})$	A	${\mathit{R}}^2$	$\mathit{RSS}(\times {10}^{-8})$
$200\times 200$	0.00020	0.00041	0.50	17590	0.924	3.348
$400\times 400$	0.00030	0.00093	0.32	18460	0.959	8.302
$600\times 600$	0.00037	0.00165	0.22	19490	0.974	14.51
$800\times 800$	0.00047	0.00350	0.13	19460	0.965	32.03

Note: There are three basic parameters in the semi-variogram model [30]. $C_0$ represents the random error introduced by grid division, while $C_0+C$ reflects the degree of spatial heterogeneity for urban attributes. $C_0/(C_0+C)$ quantifies the nugget effect, which means spatial correlation, and a larger $C_0/(C_0+C)$ represents a weaker spatial correlation. The range A represents the maximum distance of the correlation.

), which reflects the spatial heterogeneity under different grid scales through the parameters $C_0/(C_0+C)$, $R^2$, and $RSS$. Taking 2021 as an example, we compare its representation and fitting effects in the following grid scales separately (Figure 3). At $200\mathrm{m}\times 200\mathrm{m}$ and $400\mathrm{m}\times 400\mathrm{m}$ scales, urbanization attributes appear fragmented and scattered, with significant spatial random errors ($C_0/(C_0+C)>25\%$). At an $800\mathrm{m}\times 800\mathrm{m}$ scale, the overall characteristics of urbanization can be displayed, but local details are ignored, so it is difficult to precisely delineate the spatial extent. Additionally, its error also increases accordingly, nearly 2 to 10 times that of other scales. In contrast, the $600\mathrm{m}\times 600\mathrm{m}$ scale reflects the spatial heterogeneity of urbanization features, with a lower nugget coefficient and the highest $R^2$ (0.974). Therefore, the grid scale is determined to be $600\mathrm{m}\times 600\mathrm{m}$ finally, and the semi-variogram model is as follows:

$$\gamma \left(d\right)={\displaystyle \frac{1}{2N\left(d\right)}}\sum _{i=1}^{N\left(d\right)}{[Z\left(x_i\right)-Z(x_i+d)]}^2$$

(2)

where $\gamma \left(d\right)$ is the semi-variogram value; d is the sample point spacing; $N\left(d\right)$ is the paired numbers of sample points with d distance; and $Z\left(x_i\right)$ and $Z(x_i+d)$ are the regionalization variable $Z\left(x_i\right)$ at the spatial locations $x_i$ and $x_i+d$, respectively, which represents the integrated results of the index system.

2.2.3. Division of Urban–Rural Space

Urbanization attributes generally exhibit gradient characteristics at the numerical level and fluctuation characteristics in the numerical structure [31], which is crucial for delineating urban–rural spaces. Mutation detection is widely used in the monocentric model, and we verify its applicability in this study (see the Supplementary Materials), and utilize the spatial continuous wavelet transform (SCWT) to identify the mutation locations of urbanization attributes in geographic space. The main steps are as follows:

(1) Given that the main center of Wuhan is merged from the main-urban districts with the municipal government as the center, while the sub-centers are mostly distributed in the distant-urban districts. Therefore, the locations of the municipal and distant-urban district governments are designated as the central points for the main center and sub-centers, respectively (a total of 8 centers);

(2) A total of 360 cross-sectional lines covering the study area are drawn from the central points, with an angular interval of 1°. Then, these lines are intersected with the urbanization attribute layers of other grids, resulting in a total of 8 sets of 360 spatial data sequences [32]. These data sequences contain rich urban–rural gradient information, and we use SCWT to detect the mutation points in the sequences;

(3) The mutation points are corrected in different ways according to the Supplementary Materials;

(4) The corrected mutation point groups (Figure 4a) form two concentric rings around the urban center, with the inner and outer rings corresponding to urban and rural/natural areas, respectively, while the intermediate layer corresponds to the urban–rural fringe. In addition, there are also individual rings that mainly correspond to the urban–rural fringe. We connect these rings accordingly and ultimately determine the urban–rural spatial extent of Wuhan in 2012, 2016, and 2021. The formula for SCWT is as follows:

$$SCWT=S\left(x\right)\phi \left(x\right)={\displaystyle \frac{1}{\sqrt{a}}}{\int }_0^{lx}S\left(x\right)\phi \left({\displaystyle \frac{x-b}{a}}\right)dx$$

(3)

where $S\left(x\right)$ and $\phi \left(x\right)$ represent the spatial data sequence and $Daubechies$ wavelet function, respectively [33]. a and b represent the scale factor and the shift factor, which reflect the frequency range and the translation distance of the wavelet, respectively. Notably, SCWT requires balancing the effect by the scale factor a on the wavelet coefficients, as too small of a value of a contains more noise and too large of a value of a lacks valid information. Referring to the application of the variance method [34], we compare the coefficient variance at different scales (Figure 4b) to determine the optimal scale for a as 3.

2.2.4. Comprehensive Analysis of LULC

In terms of LULC change, this study utilizes the land use transfer matrix (LUTM) to reflect the change pattern of LULC in urban–rural space. Due to the continuous movements of urban–rural spatial, the transfer analysis at each stage is based on the end range [35]. Then, it employs landscape index to quantify the spatial characteristics of LULC along the urban–rural gradient, and selected indices include the Area Percentage of Landscape (PLAND), Patch Density (PD), Landscape Shape Index (LSI), and Aggregation Index (AI). These indices, respectively, describe the richness, fragmentation, shape complexity, and connectivity of LULC, enabling the sensitive detection of external disturbances and their spatiotemporal heterogeneity.

In terms of LULC response, it is affected by land encroachment, habitat destruction due to urban expansion, as well as ecological conservation efforts associated with urban development. Therefore, this study employs the ecological risk index (ERI) to quantify the positive and negative impacts on LULC. This index consists of the landscape interference index and the landscape vulnerability index, effectively reflecting both external pressures and the resilience of LULC itself [36]. Subsequently, the bivariate local spatial autocorrelation is utilized to assess the extent of LULC response to urbanization [37]. Compared to traditional linear regression models, it demonstrates the spatial correlation between variables through spatial clustering, which is particularly suitable for identifying the local correlation patterns of complex urban forms under the polycentric model. Their formulas are as follows:

$$\begin{array}{c}\hfill ERI=IV\end{array}$$

(4)

$$\begin{array}{c}\hfill I=aF+bS+cD\end{array}$$

(5)

$$\begin{array}{c}\hfill V=\sum _{k=1}^n{P}_k{V}_k\end{array}$$

(6)

where $ERI$ is the ecological risk index; I is the landscape disturbance index; F, S, and D are the normalized landscape fragmentation index (PD), landscape separation index (SPLIT), and landscape dominance index (DIVISION) in the study grid, respectively; a, b, and c are the weights assigned to 0.5, 0.3, and 0.2, respectively; V represents the landscape vulnerability index; and $P_k$ and $V_k$ denote the area of land use type k within the study grid and the corresponding vulnerability coefficient (see the Supplementary Materials), respectively.

$$I_i=Z_{xi}\sum _{j=1,j\ne i}^N{w}_{ij}{Z}_{yj}$$

(7)

where $I_i$ is the local Moran’s I (LISA), $w_{ij}$ is the spatial weight matrix, and $Z_{xi}$, $Z_{yj}$ represent the standardized Z-scores values of the variables x (urbanization attribute value) and y (eco-risk index) for the study grids i and j, respectively.

3. Results

3.1. Division Results and Evolution Process of Urban–Rural Space

Wuhan exhibits a “one core and multiple centers” distribution pattern, with urban expansion dominated by edge and outlying (Figure 5). In 2012, in the vast rural/natural background, the urban area was mainly concentrated along the two rivers (Yangtze and Han), while the urban–rural fringe distributed around them in “wheel-shaped clusters”. In 2016, the urban area had relatively limited expansion, but the urban–rural fringe extended continuously from the original boundaries and formed a “dispersed cluster” around the government seats of the distant-urban districts. In 2021, substantial rural areas were developed as construction land, while the main city center was interconnected with surrounding sub-centers through the urban–rural fringe and axially extended to the Tianhe Airport in the northwest, gradually forming a closely connected “belt-shaped city”.

The main characteristics of urban–rural space include compact urban development, the extensive expansion of urban–rural fringe, and consequent encroachment on agricultural/natural land. This study utilizes urban–rural fringe to reflect the evolution process of urban–rural space, as it is more dynamic than urban area. It is noteworthy that there are differences in urban expansion between the main and distant urban districts (Figure 6a). The main-urban districts follow a “circle development and clustered layout” pattern, with a large urban area and high proportion, but gradually approaching saturation. For instance, Jianghan District has completely transformed into urban area during the study period. The urban expansion of the main-urban districts primarily occurs in Hongshan and Hanyang Districts, with the urban–rural fringe expanding by 65.99 km² and 12.92 km², respectively. Meanwhile, the distant-urban districts construct an “axial development” township space based on transportation arteries. Although they start late, the overall improvement is fast. The urban expansion of the distant-urban districts mainly occurs in Huangpi and Jiangxia District, with the urban–rural fringe expanding by 173.32 km² and 125.10 km², respectively, nearly three times that of other main-urban districts. The directional analysis (Figure 6b) demonstrates that urbanization is accelerating, with a rapid fringe expansion. In 2012–2016, the urban area and urban–rural fringe expanded mainly to the southeast, with rates of 1.85 km²/a and 16.10 km²/a, respectively. In 2016–2021, due to previous development, the urban area continued to extend in the southeast, with the rate increasing to 5.06 km²/a. Meanwhile, the urban–rural fringe shifts towards the north and south, maintaining high growth rates of 18.14 km²/a and 17.92 km²/a, respectively.

3.2. LULC Changes in Urban–Rural Space

The transfer of LULC exhibits a drastic and frequent trend (Figure 7a). In 2012–2016, the most drastic transfer occurred in the rural area, mainly involving the transfer out of forest and construction land, and the transfer in of waters, with areas of 153.46 km², 89.09 km², and 107.60 km², respectively. In contrast, the urban area and rural–urban fringes show less change. In 2016–2021, the rural area continued in the previous stage, and forest, construction land, and waters still dominated, with transfer areas of 156.62 km², 78.32 km², and 129.846 km², respectively. It is noteworthy that, for the first time in the rural area, there has been a shift from construction land to cropland, with an area of 61.63 km². This mainly occurred in the fragmented edges of the rural periphery, possibly due to the expropriation of collective construction land and unused land by Wuhan in 2019. The transfer rate in the urban–rural fringe has significantly accelerated, where the area of cropland converted to construction land increased from 39.76 km² to 81.23 km². Overall, large-scale land transfer primarily occurs in the rural area, dominated by cultivated land and cropland; small-scale but high-proportion transfer mainly occurs in the fringe, dominated by water and forest.

The spatial configuration of LULC displays significant variations within urban–rural space (Figure 7b). The PLAND of cropland remains high at around 26.49 in the rural area, while the PLAND of construction land gradually increases to 9.17 in the fringe. The PD of cropland is relatively stable in the rural area but rapidly increases to 1.34 in the fringe. The LSI and AI reveal increasing landscape heterogeneity of LULC along the urban–rural gradient, with the highest LSI and AI observed for cropland in the rural area and for construction land in the urban area, reaching 123.33 and 96.55 in 2021, respectively.

3.3. LULC Responses in Urban–Rural Space

The LISA cluster (Figure 8) shows significant spatial heterogeneity between LULC and urbanization, although the overall correlation is weak. In the main urban center, there is a negative correlation and an enhanced trend of “High-low” clustering, with Moran’s I reaching −0.366 in 2021. This suggests that high urbanization reduces ecological risks in this region and stimulates positive responses in LULC. On its periphery—specifically, the urban–rural fringe and adjacent rural area—there are more complex spatial correlations, characterized by network-like “High-High” clustering and corridor-like “Low-High” clustering. According to clustering patterns and previous analyses, these are caused by the inefficient use of construction land, direct encroachment on ecological land, and fragmentation of cropland. These indicate that, regardless of the correlation, urban expansion exacerbates ecological risks in these regions. In the urban sub-centers, the correlation is not significant, maintaining a delicate balance between the positive effects in the main center and the negative impacts in the peripheral areas. In the vast rural areas, LULC is virtually unaffected by urbanization (Moran’s I = 0.014), which explains the insignificant correlation within urban–rural space. Nevertheless, there are large-scale “Low-Low” and “Low-High” clusterings, specifically in the northern and central parts of Wuhan, which contribute to resisting human disturbances and regulating the ecosystem.

4. Discussions

4.1. Performance of Other Delineation Methods in the Polycentric Structure

This study compares the classical methods (spatial clustering and threshold determination) based on overall distribution and spatial details, and pays special attention to their practical performance in the complex multi-center scenario of Wuhan. The clustering method is able to tap into the intrinsic structure and regularity of urban–rural space, divide it into multiple highly similar clusters, and have wide applicability at the meso–macro-scale. The threshold method determines the spatial extent by repeatedly comparing the delineation effects at different thresholds, which is one of the most commonly used methods in the early monocentric model.

The results indicate that the clustering method may produce inaccurate delineation results when dealing with spatial data disturbed by water systems. Additionally, it leads to coarse estimations near urban sub-centers. For example (Figure 9(e-1,e-2)), in the Huangpi District, influenced by the Fu River and its wetlands, the clustering method divides half of Tianhe Airport into rural areas. In the government seat of Jiangxia District, where there are rural settlements in the vicinity of the urban sub-center, the clustering method spans over the dense barriers of forests and farmland and identifies the entire periphery, including villages, as urban–rural fringe. Moreover, due to the nature of unsupervised learning algorithms, it is relatively difficult to establish effective external correction mechanisms to correct these errors, and additional criteria and bases reduce its applicability. Composite indicators increase the noise of spatial data, thereby increasing the complexity of the clustering method, which is more pronounced in the thresholding method. Therefore, a single metric (land use intensity) is used instead to assess the performance of the thresholding method [38]. Figure 9b shows that the extracted spatial extents are relatively scattered and, due to the dependence on construction land, it tends to overestimate when delineating urban and rural–urban fringe. Most importantly, a single threshold setting makes it difficult to identify heterogeneous areas within the urban–rural space because it has almost no ability to recognize polycenters.

In contrast, this study disaggregates urban centers into multiple strictly monocentric regions and obtains a comprehensive urban–rural spatial extent by employing mutation detection and considering the multidimensional characteristics of urbanization. Although the correction of mutation points remains complex and intricate, the semi-quantitative approach helps to adopt unified criteria to address the complex environment in Wuhan, thereby providing relatively accurate results for subsequent LULC analysis.

4.2. Impact of Planning and Construction on Polycentric Expansion

Urban planning and construction are important factors driving urbanization, directly and indirectly influencing the spatiotemporal patterns of urban–rural space and the evolution of its internal LULC. Wuhan has implemented various measures to promote polycentric growth, so this study primarily discusses the positive contributions of urban planning, industrial development, and transportation (Figure 10).

In terms of urban planning and industrial layout, Wuhan focuses on the development of the central city and its surrounding areas (such as the construction of the Greater Wuhan Metropolitan Area and new town clusters). In 2013, Dongxihu District established the Airport Economic and Technological Development Zone, forming a national foodstuff industry zone dominated by food processing, modern logistics, and light industries. Subsequently, the East Lake High-tech Development Zone expanded eastward along Gaoxin Road and Wuhuang Road, gradually forming a high-tech industrial cluster dominated by optoelectronics, biomedicine, and mechatronics. The increase in the proportion of secondary and tertiary industries stimulates the demand for industrial land and population agglomeration, and accelerates urban construction and living space expansion. As a result, in 2012–2016, the urban–rural fringe mainly expanded towards the Dongxihu and Hongshan District, while the changes in the sub-centers were relatively small.

In terms of transportation location, Wuhan is strengthening the radiation capacity of the main center and the connection with the sub-centers through the construction of international hubs and suburban railways. In 2017, the Airport New City in the Huangpi District was positioned as the suburban core and accelerated the construction of the “aviation-subway hub” modern comprehensive transportation system relying on Tianhe Airport, which has accelerated the urbanization process in the northern region. In 2018, the East Lake High-tech Development Zone opened several transportation arteries, such as the Zhifang Line, driving the development of high-tech industries in the southern region, including medicine and electromechanics. Road transportation, serving as the framework and bridge for economic development, supports the urban expansion towards the south and corridor-like connectivity with the sub-center in Jiangxia District via the fringe.

4.3. Spatiotemporal Heterogeneity of LULC along Urban–Rural Gradient

From urban area to urban–rural fringe and further to rural/natural areas, the interaction of LULCs along the gradient creates changing and distinct landscape types and also serves as a significant characteristic of the urban–rural gradient. LULC exhibits notable differences in urban–rural spaces, necessitating its placement in a gradient framework (Figure 11). Previous studies have shown that fragmented and complex-shaped LULCs increase from urban area to urban–rural fringe, and then rapidly decline in densely rural agricultural regions [39]. In this study, this trend is primarily driven by cropland, as it demonstrates a rising and then declining trend, while other LULCs exhibit a monotonic increase. In Wuhan, cropland is a crucial agricultural resource, covering 64.64% of the total area and serving as a spatial continuum that almost encompasses all LULCs, playing a significant role in guaranteeing food security and enhancing ecological resilience. Nevertheless, cropland in the urban–rural fringe is facing the highest PD and AI, and this trend is still increasing. This suggests that urban sprawl in the fringe is still only partially controlled, resulting in the significant destruction of agricultural land and ecological systems, with the negative impacts persisting. Overall, the urbanization response of LULC along the urban–rural gradient experiences a shift from positive to negative and then back to positive.

From the perspective of urban planning and ecological conservation, efforts should be made to establish a negative correlation between LULC changes and urbanization in urban areas and urban–rural fringe by implementing effective spatial governance to address the “human-land” conflict. In the rural area, a positive correlation between the two should be promoted, by avoiding the radiating “destructive force”. This viewpoint aligns with the study on Mediterranean cities [40], which indicates that moving towards more competitive and sustainable polycentric development requires specialized and differentiated planning solutions. Although there are relatively complete systems for cropland regulation and ecological protection in the master planning of Wuhan, the more pressured peripheral landscapes have not received the same attention as the urban landscapes when implementing in the district and county. Therefore, ecological conservation efforts should enhance risk identification and project supervision, considering regional differences in scale and resources as well as trade-offs between multiple interests [41], to establish a transitional landscape that integrates urban development and regional characteristics under an ecological protection framework [42]. For villages or regions significantly impacted by urbanization, it should involve cropland fragmentation in the planning assessment indicators, and carry out the comprehensive land consolidation [43].

5. Conclusions

This study proposes a spatial delineation framework suitable for the polycentric structure, which accounts for the social, economic, and natural characteristics of urbanization to meticulously depict the polycentric transformation and to discuss the effects and gradient differences of polycentric development on LULC from a holistic urban–rural perspective. The key findings are as follows: (1) Under the complex urban scenario of polycentricity, the grid scale of $600\mathrm{m}\times 600\mathrm{m}$ effectively captures the multidimensional characteristics of urbanization, and the mutation detection method can be applied to the spatial delineation of polycentric structure; (2) The urbanization process of Wuhan is accelerating, with the development of the central city approaching saturation. The polycentric structure propels the rapid expansion of the urban–rural fringe, which mainly expands along the north and south directions, at the rates of 15.81 km²/a and 13.69 km²/a, respectively; (3) LULC changes are more drastic and frequent, which differ significantly within the urban–rural space. Large-scale land conversions mainly occur in rural areas, dominated by cropland and construction land; small-scale, high-proportion land conversions mainly occur in the urban–rural fringe, dominated by waters and forest. The overall effects of urbanization on LULC are relatively weak, along the urban–rural gradient, experiencing a transition from positive to negative and back to positive.

Compared to other studies, these findings reveal that utilizing multi-source data and comprehensive indicators facilitates the spatial delineation of the polycentric structure. More importantly, they highlight that the focus should be shifted beyond urban areas—polycentric development is not always beneficial. The government should consider the affordability of urban periphery when introducing measures to promote polycentric development; otherwise, the gap between the center and periphery may widen. This study expands our understanding of polycentric development, provides crucial insights into the interactions between urbanization and land use in complex urban scenarios, and offers a theoretical basis for polycentric practices and ecological space governance in metropolitan areas. Nevertheless, there are several limitations. Firstly, in terms of spatial delineation, although mutation detection provides relatively accurate results for subsequent LULC analysis, the locations of urban primary/secondary centers are predetermined, and the individual delineation of each center may not be applicable at the meso–macro-scale. Secondly, we only analyze the relationship between urbanization and LULC in the context of the urban–rural spatial framework. Natural resources, urban morphology, industrial development, topography, and elevation are factors that should be further considered to explain the driving forces and patterns behind spatial variations, thereby enabling differentiated regulation and conservation strategies.