Identification and Analysis of Production–Living–Ecological Space Based on Multi-Source Geospatial Data: A Case Study of Xuzhou City

Abstract

Effective production, living, and ecological space allocation is essential for improving and optimizing urban space development. In this study, we proposed a production–living–ecological space (PLES) identification method based on Point of Interest (POI) data and China Land Cover Dataset (CLCD) to identify PLESs in Xuzhou City for the years 2012, 2018, and 2022, with an average recognition accuracy of 89.81%. Moreover, the land-use transfer matrix, center of gravity migration, and Geo-detector were used to reveal the spatiotemporal pattern evolution of PLESs. The results showed that: (1) The distribution of PLESs presented significant differentiation between Urban Built-Up Area (UBUA) and Non-Urban Built-Up Area (NUBUA). UBUA was mainly composed of living spaces, while NUBUA was primarily characterized by production–ecological spaces. (2) The intensive utilization of urban land led to an increase in the area of multifunctional spaces, while the complexity of urban space increased. (3) During 2012 to 2022, the center of gravity of PLESs remained relatively stable. The moving distances were all less than 1 km (except for ecological space from 2012 to 2018). (4) The evolution of PLESs was closely linked with socio-economic factors, and the interactions between the factors also had a significant driving effect on PLESs.

1. Introduction

In recent years, the rapid urbanization and economic development in China have significantly transformed land use and allocation patterns, leading to profound changes in the structure and functions of regional landscapes [1]. As a country with a large population and limited land resources, China faces pressing challenges in balancing urban expansion with territorial spatial pattern preservation [2]. Under this circumstance, optimizing the allocation and coordinated development of production, living, and ecological spaces has emerged as a crucial strategy for addressing resource and environmental challenges and achieving sustainable regional development [3]. The balanced development of the territorial spatial pattern and its function is an important prerequisite for the realization of regional high-quality development [4]. In 2012, the 18th National Congress of the Communist Party of China proposed the objective of “building the production space that is intensive and efficient, the living space that is livable and appropriate, and an ecological space that is beautiful and unspoiled”, which introduced the concept of the development significance of production, living, and ecological spaces [5,6]. These three kinds of spaces are abbreviated as Production–Living–Ecological Space (PLES), which serves as a vital indicator in land planning and development [7]. Production space supports economic activities such as agriculture, industry, and commerce [8,9]; living space accommodates human habitation and essential services [10]; and ecological space contributes to environmental regulation and resource conservation [11]. In recent years, studies on urban spaces have increasingly focused on balancing ecological protection with sustainable development. Yuan et al. [12] developed natural ecological networks to evaluate the ecological system characteristics, providing a scientific basis for balancing ecological protection with regional development in Changzhi City. Similarly, Bai et al. [13] demonstrated that the coordinated allocation of urban space and resources is essential for achieving sustainable urban agglomeration development from the contributions of infrastructure, industrial layout, public services, and environmental services. Their study further validated the importance of balanced urban spatial development by incorporating infrastructure, industrial layout, public services, and environmental services. Production, living, and ecological spaces are both independent and interdependent, and their relationships shape the functional patterns and development potential of a region [14]. However, under rapid urbanization, these spaces often face imbalances. For instance, large areas of agricultural land have been converted to urban construction land [15], and human activities have severely impacted the stability of urban ecosystems [16], resulting in the expansion of production and living spaces at the expense of ecological spaces, which weakens ecological functions and sustainability [17]. The interactions among production, living, and ecological spaces are often imbalanced, which could lead to ecological degradation and inefficient resource utilization. Therefore, monitoring and analyzing the spatiotemporal changes in these spaces is essential for providing scientific support for land management and urban planning, helping to alleviate the environmental pressures caused by land resource scarcity. Huang et al. [18] selected land-use data and constructed a regional spatial model to reveal the ecological environment quality and the effects of land-use transformation in the PLESs of Baishan City. Aruban and Liu [19] discussed the evolution of the spatiotemporal pattern of PLESs to offer a basis for the effective use and rational development of rural land in Duolun County. The identification and analysis of PLESs from diverse perspectives provide deeper insights into urban growth patterns, land-use optimization, and the sustainable development of urban areas. However, most studies [20,21] have focused on the identification and optimization of single-function spaces, often neglecting the synergy, interrelationships, and functional balance among different types of spaces. In reality, the various functional spaces involved in urban development are often interwoven and interdependent. Focusing solely on the optimization of a single space is insufficient to resolve the conflicts and trade-offs among different spaces effectively. Therefore, systematic and comprehensive identification of PLESs and analyzing their spatiotemporal pattern evolution are critical for achieving efficient land management and sustainable development.

The previous research on PLESs has focused on identification and classification [22,23], the spatiotemporal pattern evolution [24], evolution characteristics [25], driving force exploration [26,27], and spatial simulation and optimization [28]. Hu et al. [29] selected land-use/cover data and developed a classifying assessment index system for analyzing the spatial distribution of PLESs in Qinghai. From the perspective of PLESs, Liang et al. [30] explored the evolutionary trend of land-use change and landscape ecological risk in an ecologically fragile area the based on the land cover data. These studies also provided targeted policy recommendations to support regional sustainable development, ecological environment protection, and land and space optimization based on their identification and analysis of PLESs. Under these circumstances, effective PLES identification is a necessary condition for achieving a fair allocation of PLESs. For the identification and classification of PLESs, it mostly relies on land-use classification data, which offer limited representation of land functions and have difficulty fully reflecting the diversity and dynamic change characteristics of spatial functions in previous research. To address these limitations, multi-source geospatial data have been widely applied across various fields due to their high accuracy, comprehensiveness, multidimensionality, and timeliness. Yin et al. [31] analyzed the economy–society–environment sustainability in the Yangtze River delta urban agglomeration based on the multi-source geospatial data. Qin et al. also used the data to assess subsidence risk in urban areas [32]. The rapid development and widespread application of multi-source geospatial data have provided new technical support for the detailed study of PLESs [33]. The rise of geographic big data (such as Point of Interest (POI) data) provides new opportunities and perspectives for the realization of the simulation and inference of geographic systems and the exploration of the development rules and trends of complex urban space [34]. POI data not only offer precise information about urban functional layouts, but also reflect the intensity and distribution of human activities [35]. When integrated with remote sensing imagery, socioeconomic statistics, land cover data, and other sources, it becomes possible to identify and monitor PLESs across multiple scales dynamically [36]. Therefore, conducting refined identification and spatial pattern evolution analysis of PLESs based on multi-source data holds significant theoretical value and practical implications.

Different from the previous studies that identify and analyze PLESs in urban cities, this study focused on Xuzhou City and divided it into Urban Built-Up Area (UBUA) and Non-Urban Built-Up Area (NUBUA) using the Densi-graph method and POI data. In UBUA, a PLES identification method based on POI data and CLCD was proposed to identify PLESs in Xuzhou City for the years 2012, 2018, and 2022. This method achieved a more detailed and nuanced classification of PLESs and successfully captured its diverse functional characteristics. However, POI data are more densely distributed in UBUAs and sparsely distributed in NUBUAs [37]. Due to the sparse distribution of POI data in NUBUAs, the CLCD was applied to identify PLESs, which could improve the reliability and accuracy of the identification results. Furthermore, based on the PLES identification results, the spatiotemporal pattern evolution and driving mechanism of PLES distribution in Xuzhou City were discussed. The findings can provide a scientific reference for the optimization of regional spatial resource allocation in Xuzhou City.

Descriptions of the materials (including the study area and data sources) and methods used in this paper are provided in Section 2. The distributions of the PLESs classified by our method are described in Section 3. Moreover, the spatiotemporal pattern evolution of PLESs and the significance and limitations of this study are discussed in Section 4. Finally, the main conclusions are drawn in Section 5.

2. Materials and Methods

2.1. Overview of the Study Area

The city of Xuzhou (116°22′ E~118°40′ E, 33°43′ N~34°58′ N) is situated in northwestern Jiangsu Province, covering an area of 11,765 km² (Figure 1). It consists of the central urban districts (Tongshan, Quanshan, Yunlong, and Gulou), Jiawang District, Feng County, Pei County, Suining County, Xinyi City, and Pizhou City. As one of the national old industrial bases and the first resource-depleted cities, Xuzhou has steadily advanced industrial, ecological, urban, and social transformations in revitalization. According to the China Resource-Based Old Industrial City Transformation Development Index Report (2022), Xuzhou City ranked first among the top 20 cities in China for its comprehensive transformation performance as a resource-based old industrial city. However, Xuzhou City is under the challenges of enhancing autonomous innovation, alleviating resource and environmental constraints. Consequently, it is very important to detect the evolution and interrelationships of the PLESs in Xuzhou City, which will help to scientifically make regional spatial planning and guide the coordinated development of regional environment, economy, and society.

2.2. Data

2.2.1. Data Sources

The primary data used in this study include POI data, CLCD, and Open Street Map (OSM) data. The POI data were obtained from the Automatic Navigation Map (AMAP) for Xuzhou City in 2012, 2018, and 2022. A total of 654,328 elements were obtained, with 82,289 collected in 2012, 317,119 in 2018, and 254,920 in 2022. Furthermore, a variety of socio-economic and natural geographic data were incorporated in this study. More details of these data are provided in

Table 1. Data description.

Data Type	Data Name	Data Source
Point data	POI	AMAP inside
Land cover data	CLCD	Wuhan University
Administrative division spatial data	Basic geographic data	National geographic Information Resource Catalog Service System
Map data	Open Street Map (OSM)	Open-source website: OSM
Statistical data	Area of built districts	Xuzhou statistical yearbook
Precipitation data	1 km monthly mean temperature dataset for China	Resource and environment science data platform
Temperature data	1 km monthly precipitation dataset for China	National earth system science data center
Terrain data	DEM	NASA
Socio-economic data	GDP	Xuzhou statistical yearbook
	Demographic data	Xuzhou statistical yearbook
Net Primary Productivity (NPP) data	Data product of MOD17A3H.006	PIE-Engine

.

2.2.2. Data Preprocessing

According to the AMAP classification system, the POI data were divided into three levels, including base type, subtype, and category. Each POI dataset includes 9 attributes: name, province, city, address, longitude, latitude, and the three-level categories. Due to the extensive number of categories and the inherent redundancy and overlap within the data, duplicate POIs were eliminated based on identical names and geographic proximity. The POI data, initially comprising nine attributes, were filtered to retain only six key attributes: name, longitude, latitude, and the three-level categories, while excluding the province, city, and address fields. POI elements with misclassified categories were relabeled based on their category attributes. Additionally, referencing the “Current land use classification (GB/T 21010-2017)” [38] and the “Code for classification of urban land use and planning standards of development land (GB 50137-2011)” [39], the POI data were reclassified based on the living, production, and ecological functions.

The POI data captured both the spatial location and attribute information of geographic entities. However, significant differences in area were observed among various types of entities, such as supermarkets, schools, gymnasiums, and zoos. To address this, relevance (P1) and influence (P2) were adopted [34], based on the classification of living, production, and ecological functions, to calculate a comprehensive weight (P). Relevance (P1) was determined for each POI element using the analytic hierarchy process (AHP) [40], while influence (P2) was assigned a score ranging from 0 to 100 based on the public perception of the area of geographic entities, with scoring criteria referenced from previous studies [41,42]. Finally, the comprehensive weight (P) for each POI element was obtained by multiplying relevance (P1) and influence (P2). The results of the reclassification and comprehensive weight calculations are presented in

Table 2. Reclassification and comprehensive weight calculation results of POI data.

Spatial Type	Sub-Categories	POI Type	P1	P2	P
Living space	Service-oriented living space	Food and Beverage Services	0.0280	10	0.028
		Shopping Services	0.0242	20	0.484
		Hotel Accommodation	0.0356	10	0.356
		Lifestyle Services	0.0954	15	1.431
		Automotive Services	0.0108	10	0.108
		Public Facilities	0.0122	15	0.183
		Science, Education and Culture	0.1882	25	4.705
		Sports and Leisure	0.0216	10	0.216
		Healthcare Services	0.1721	20	3.442
	Residential living space	Business Residences	0.0501	50	2.505
Production space	Service-oriented production space	Business Districts	0.0463	30	1.389
		Financial Institutions	0.0439	20	0.878
	Enterprise-oriented production space	Industrial Parks	0.0297	70	2.079
		Corporate Enterprises	0.0632	50	3.160
	Transportation-oriented production space	Transportation Facilities	0.0874	15	1.311
		Roadside Facilities	0.0140	10	0.140
Ecological space	Service-oriented ecological space	Parks and Green Spaces	0.0374	90	3.366
		Scenic Spots and Landmarks	0.0274	90	2.466
		Cultural Landscapes	0.0126	30	0.378

.

2.3. Methods

In this study, the methodology is primarily divided into four parts. Firstly, POI data were applied to extract the urban built-up area by the Densi-graph method. Xuzhou City was divided into Urban Built-Up Areas (UBUAs) and Non-Urban Built-Up Areas (NUBUAs). Secondly, the identification units were delineated by using the OpenStreetMap (OSM) data for UBUAs. Meanwhile, a PLES identification method that combined POI data and CLCD was proposed within UBUAs. In contrast, the classification and identification of PLESs in NUBUAs primarily relied on the CLCD, which provided comprehensive land cover information. Finally, the spatial distribution and temporal evolution of PLESs in Xuzhou City were discussed and analyzed. The detailed flowchart of the research process is shown in Figure 2.

2.3.1. Extraction of Urban Built-Up Area

POI data are closely related to human activities, and it reflects the spatial concentration of these activities. The POI data is densely distributed in UBUA. While at the boundaries of UBUA, the number of POIs significantly decreases. Based on this characteristic, the Densi-graph method proposed by Xu et al. [43], which have proven effective for extracting the boundaries of urban built-up areas.

The main idea of the Densi-graph method is as follows: First, the Kernel Density Estimation (KDE) [44,45] was applied with a search radius set to 5 m to generate the kernel density map for POI data, which provides the contour lines of the POIs. And then, the kernel density values (d) and the theoretical radius ($\surd S_d$) of the area (S_d) enclosed by the corresponding closed curves are respectively used as the X and Y values to plot the Densi-graph curve. By analyzing the changes in the curve, the threshold at which the overall decline in POI density leads to an irreversible upward trend is identified. The POI density value corresponding to this threshold was taken as the critical value (r) for the urban built-up area boundary. The area enclosed by the contour line at this value (r) was defined as the urban built-up area. The main calculation formula of this method is shown in Formula (1):

$$D=lim{\displaystyle \frac{\mathrm{d}\left(\sqrt{{∆S}_d}\right)}{\mathrm{d}d}}$$

(1)

where d represents the density value, $S_d$ is the area enclosed by the corresponding closed curve, $\sqrt{{∆S}_d}$is the theoretical radius of the area defined by the closed curve, and D is the derivative value of the theoretical radius increment. When D = 0, the density curve is uniform, when D > 0, the density curve is expanding, and when D < 0, the density curve is contracting.

2.3.2. Unit Division for PLES Identification in UBUAs

The identification units (hereinafter referred to as units) represent the smallest spatial divisions used to categorize different areas. Defining these units is a crucial step in identifying PLESs within UBUAs. In recent research, two primary methods are commonly used to define units: grid-based division [46,47,48] and road network-based division [49,50,51]. In reality, urban block forms are typically represented as irregular polygons, and grid-based division may lead to discrepancies between the shapes and sizes of the grid units and the actual functional zones. In contrast, road network-based division considers the actual layout of urban roads, resulting in units with irregular geometric shapes that more accurately reflect the true boundaries and structures of functional zones. Therefore, this study employed the road network-based division method to define the units.

2.3.3. PLES Identification in UBUAs

The Spatial Join tool in ArcMap 10.2 was used to spatially connected the POI data with the units in UBUAs. The attributes of production, living, and ecological functions were assigned to each corresponding unit. In addition, the number of POI datapoints that belong to each unit was counted. Based on the specific types and quantities of living, production, and ecological proportions within each unit, the proportions of different functions for each unit were calculated using Formulas (2)–(4):

$$N_i=C_i\times P_i , i=1, 2, 3,\dots \dots ,n$$

(2)

$$S=\sum _i^n{N}_i ,i=1, 2, 3,\dots \dots ,n$$

(3)

$$R_i={\displaystyle \frac{S_i}{\sum _i^n{S}_i}} , i=1, 2, 3,\dots \dots n$$

(4)

where $N_i$ represents the comprehensive count of the i-th type of POI within the unit, $C_i$ is the number of POI data points of the i-th type of POI data in the unit, $P_i$ indicates the comprehensive weight of the i-th type of POI data, $S$ represents the combined value of living, production, and ecological functions in each unit, $S_i$ denotes the proportion value of the i-th type functions within the unit, and $R_i$ is the value of the proportion of the i-th type functional element in the unit.

Relying on the calculated proportion of POI data assigned to different PLES functions, the spatial type can be determined within each unit. First, the proportions of living, production, and ecological functions within each unit were calculated. Second, the spatial identification of the units was conducted based on a specified classification standard. It is noteworthy that null areas may appear in the identification results, requiring further supplementary identification using CLCD. The spatial types were identified and classified into 10 categories of PLES types in UBUAs. The specific spatial categories are listed in

Table 3. The description of different spatial types of PLESs in this study.

Function Type	Spatial Type	Description
Single-function space	Living space (LS)	One function of PLES dominates or only one function exists in the unit
	Production space (PS)
	Ecological space (ES)
Dual-function space	Sub-living space (Sub-LS)	Two functions of PLES are predominant or only two functions are present in similar proportion in the unit
	Sub-production space (Sub-PS)
	Sub-ecological space (Sub-ES)
	Living–production space (LPS)	Two functions of PLES are predominant or only two functions are present in significantly different proportion in the unit
	Living–ecological space (LES)
	Production–ecological space (PES)
Mixed-function space	Mixed space (MS)	Three functions of PLES with similar proportions or one of the functions is approximately the sum of the other two similar proportions in the unit

.

The process for identifying PLESs, as detailed in Figure 3, is presented in the following steps:

• When a unit contains no POI data, it is classified as “No data”, and its final spatial type is determined by supplementary identification using the CLCD.
• If POI data exist within the unit, the proportions of each function are calculated and ranked as primary, secondary, and tertiary proportions, represented by $R_i$, $R_j$, and $R_k$, respectively.

  (1)

  When $R_i\ge 70\%$, the unit is determined as a single-function space (LS, PS, or ES) dominated by the function of $R_i$.

  (2)

  When$50\%\le R_i<70\%$, $R_j$ is compared to 30%, resulting in four possible identification scenarios:

• If $R_j<30\%$ and $R_k\ge 20\%$, the unit may be defined as MS.
• If $R_j<30\%$, $R_k<20\%$ and the relative difference between $R_j$ and $R_k$ is $\ge 10\%$, the unit is possibly classified as SLS, SPS, or SES based on the function of $R_i$. If the difference is less than 10%, the unit may be determined to be a single-function space (LS, PS, or ES) dominated by the function of $R_i$.
• If $R_j\ge 30\%$ and $R_k\ge 10\%$, the unit is possibly classified as SLS, SPS, or SES according to the $R_i$.
• If $R_j\ge 30\%$, $R_k<10\%$, and the relative difference between $R_i$ and $R_j$ is $\ge 20\%$, the unit is possibly classified as SLS, SPS, or SES based on the function of $R_i$. If the difference is less than 20%, the unit is possibly defined as the LPS, LES, or PES depending on the functions of $R_i$ and $R_j$.

  (3)

  When $R_i<50\%$, $R_j$ is compared to 30%, resulting in three possible identification scenarios:

• If $R_j\ge 30\%$ and $R_k\ge 20\%$, the unit is classified as MS.
• If $R_j\ge 30\%$ and $R_k<20\%$, the unit is likely classified as LPS, LES, or PES, depending on the functions of $R_i$ and $R_j$.
• If $R_j<30\%$, the unit is classified as MS.

2.3.4. PLES Identification in NUBUAs

In NUBUAs, the drastic decrease in the number of POIs compromises the accuracy and reliability of PLES identification results. To enhance the precision of these results and avoid errors caused by data sparsity, the CLCD was utilized for identifying PLESs within NUBUAs. Subsequently, a classification framework for NUBUAs was developed by reclassifying and reorganizing the CLCD from the perspective of PLESs. The classification criteria for this system are detailed in

Table 4. Classification framework of PLESs in NUBUAs.

Spatial Type	Dominant Function	CLCD Type
Living–production space	Rural or urban living	Impervious
Production–ecological space	Crops	Cropland
Ecological space	Aquatic or greenery ecosystem	Forest, Shrub, Grassland, Water, Barren, Wetland

. The CLCD offered valuable insights into land cover types, which were essential for identifying and classifying areas that support ecological and productive functions.

3. Results

3.1. Extraction Results and Verification of UBUAs

3.1.1. Extraction Results

Based on the methods described in Section 2.3.1, the urban built-up areas of Xuzhou in 2012, 2018, and 2022 were extracted using POI data. The extraction results are presented in Figure 4, which clearly illustrates a discernible trend of urban built-up area expansion in Xuzhou City.

3.1.2. Verification of Extraction Results

After obtaining the built-up area extraction results, accuracy validation was necessary to verify their reliability and to determine whether the results can be used to divide Xuzhou City into UBUAs and NUBUAs for subsequent PLES identification. Firstly, the areas of built districts from the Xuzhou Statistical Yearbook were used as reference data. Secondly, the relative errors between the extracted area and the statistical data were calculated. The pertinent data presented in

Table 5. Accuracy validation results for urban built-up area extraction.

Year	Extracted Area (km2)	Statistical Area (km2)	Error (%)
2012	421.16	403.1	4.48
2018	462.07	454.7	1.62
2022	533.04	497.1	7.23

showed that the relative errors between the extracted area and the statistical area were consistently below 10%, indicating that the extraction results were reliable.

3.2. Identification Results and Verification of PLESs in Xuzhou City

3.2.1. Identification Results

After UBUAs were extracted, Xuzhou City was divided into UBUAs and NUBUAs for the identification of PLESs. By applying the methods described in Section 2.3.3 and Section 2.3.4, the PLESs were identified for the years 2012, 2018, and 2022. The results are shown in Figure 5. From 2012 to 2022, the PLESs in Xuzhou City changed significantly. It is obvious that the living–production space in NUBUA was expanding, especially around the boundary of UBUAs. This phenomenon revealed the gradual urbanization of the spatial layout in NUBUAs with the advancement of the transformation of industries and urbanization processes in Xuzhou City.

3.2.2. Verification of Identification Results of PLESs in UBUAs

The objective is to evaluate the accuracy and applicability of the proposed method for identifying PLESs in UBUAs. Based on the identification results, 50% of the units for each year were randomly selected as validation units and were compared with the PLES classifications derived from visual interpretation of remote sensing imagery in UBUAs using a confusion matrix [52]. In 2012, the identification result achieved an accuracy of 91.58%, with a kappa coefficient of 0.88. For 2018, the accuracy was slightly lower at 88.42%, with a kappa coefficient of 0.84. In 2022, the accuracy improved to 89.42%, with a kappa coefficient of 0.86. Overall, the average accuracy across the three years was 89.81%, with kappa coefficients consistently above 0.86. The accuracy validation results confirmed the reliability and effectiveness of our method for identifying PLESs in UBUAs over time. The validation results in 2022 are shown in

Table 6. The confusion matrix of PLES verification results for UBUAs in 2022.

PLES 1	LS	PS	ES	Sub-LS	Sub-PS	Sub-ES	LPS	LES	PES	MS	N 2	PA (%)
LS	427	5	2	10	1	3	35	0	1	2	486	87.86
PS	2	119	0	0	0	0	0	0	2	0	123	96.74
ES	0	0	79	0	0	0	0	0	0	0	79	100
Sub-LS	0	0	0	72	0	0	0	0	0	0	72	100
Sub-PS	0	1	0	0	48	0	1	0	0	0	50	96.00
Sub-ES	0	0	6	0	0	30	0	0	1	0	37	81.08
LPS	31	4	4	5	6	1	330	6	1	1	389	84.83
LES	0	0	4	1	0	1	0	24	0	0	30	80.00
PES	0	0	0	0	0	1	0	0	24	0	25	96.00
MS	0	0	0	0	0	0	0	0	0	13	13	100
N2	460	129	95	88	55	36	366	30	29	16	1304	-
UA (%)	92.83	92.25	83.16	81.82	87.28	83.33	90.16	80.00	82.76	81.25	-	-
OA			89.42%			Kappa			0.86

¹ The horizontal table header is the comparison data obtained through visual interpretation of remote sensing images; the vertical table header is the PLES results obtained based on POI data identification. ² The number of validation units.

.

According to the validation results in Table 6, we further analyzed the misclassified units and found that 138 units out of the 1304 selected units exhibited discrepancies between the types identified using POI data and the types visually interpreted from remote sensing imagery. Among these misclassified units, 28 units (counted 20.29% for our validation units) were located near the boundaries of built-up areas, where the sparse distribution of POIs led to inaccuracies in the identification results of PLESs. The remaining 110 misclassified units were primarily concentrated in areas with no data. When supplemented with CLCD for identification, these units were often classified as similar functional spaces. Similarly, these phenomena were also present in the identification results of PLESs in 2012 and 2018. However, the validation results indicate an average accuracy of 89.81% and an average kappa coefficient of 0.86. Conclusively, the PLES identification results were reliable enough to provide effective support for the identification of PLESs in UBUAs.

3.3. Identification Results in UBUAs

In this section, the detailed presentation of the identification results of the PLESs in UBUAs (Figure 6) are provided to enhance comprehension of the distributions and interactions of the PLESs. The area of each type of PLES is calculated and summarized in Table 6. From the figure, it could be observed that the distribution of PLESs in UBUAs presented diversified characteristics. The living spaces, sub-living spaces, and living–production spaces represented the main components of the PLESs in Xuzhou City. Living and sub-living spaces primarily supported the daily activities of urban residents, serving as an essential foundation for residential and service functions. Living–production spaces combined residential and production functions to fulfill the accommodation requirements of the population. Meanwhile, they also provided support for economic activities. These spaces constitute the primary spatial distribution within UBUAs, which reflected the coordination and diversity in the utilization of space in Xuzhou City.

According to the

Table 7. The areas of PLESs in UBUAs (km²).

PLES Type	Xuzhou City			The Central Urban Districts
	2012	2018	2022	2012	2018	2022
LS	263.23	223.72	261.39	70.99	66.00	70.48
PS	19.99	18.97	17.45	10.28	13.48	9.11
ES	5.11	11.76	10.67	3.37	6.14	10.44
Sub-LS	51.21	65.11	66.01	22.91	30.31	28.00
Sub-PS	9.60	16.14	18.38	4.59	7.09	10.77
Sub-ES	1.50	4.25	5.20	0.91	1.17	4.56
LPS	53.55	97.69	117.07	20.97	32.56	37.91
LES	12.39	15.87	17.36	6.64	10.87	9.16
PES	1.88	4.23	5.41	0.75	2.09	3.60
MS	2.71	4.33	14.10	1.86	2.29	8.61

, the area of living space and production space dominated by a single function showed a decreasing tendency in UBUAs over the past decades. In contrast, the area of multifunctional space and mixed-function space increased significantly. Notably, the growth in living–production spaces and mixed spaces were particularly remarkable, with increases of 63.52 km² and 11.39 km². Through the active promotion of emerging industries, an increase in the proportion of production functions was achieved in Xuzhou City. It marked a phased success in its industrial transformation.

3.4. Identification Results in NUBUAs

The results for identification of PLESs in NUBUAs are presented in

Table 8. The areas of PLESs in NUBUAs (km²).

PLES Type	PES			LPS			ES
Year	2012	2018	2022	2012	2018	2022	2012	2018	2022
The central urban district	1554.97	1500.12	1473.79	581.42	618.15	637.51	102.41	76.42	66.74
Jiawang	408.78	394.24	390.56	166.94	183.11	191.13	29.68	28.92	26.05
Pizhou	1623.50	1592.04	1576.17	376.76	411.12	425.83	49.79	46.88	43.55
Peixian	939.08	920.52	907.39	259.48	273.92	275.87	13.36	14.48	14.33
Xinyi	1166.35	1139.12	1125.35	269.45	296.98	309.88	122.37	121.86	118.41
Fengxian	1134.43	1108.67	1096.91	268.36	295.12	303.76	11.22	12.71	12.36
Suining	1405.58	1369.18	1348.75	300.71	334.19	345.21	22.02	24.29	23.83

.

(1)

The scale of ecological space decreased, and the protection of ecological land was under pressure.

The area of ecological space decreased significantly from 2012 to 2022, especially in districts with rapid urban expansion, such as the central urban districts and Xinyi City. The ecological space in the central urban districts decreased from 102.41 km² in 2012 to 66.74 km² in 2022, indicating a reduction rate of 34.8%. Similarly, Xinyi City experienced a reduction of 3.2%. The ecological space was increasingly compressed by urban development and production demands, which reflected weakened ecological protection in areas with rapid development pressures. Nevertheless, Pei County and Feng County maintained relatively stable ecological space proportions, suggesting that localized conservation measures may have been effective.

(2)

The scale of production–ecological space decreased due to urban expansion and land-use transformation.

From 2012 to 2022, the production–ecological space in Suining County decreased by 569.8 km², which accounted for an annual average reduction of 56.98 km², while Pizhou City decreased by 469.2 km², with an annual reduction of 46.92 km². Feng County, Xinyi City, and Pei County also showed decreases of 375.2 km², 409.8 km², and 316.9 km² over the same period, respectively. The most notable annual reduction rate occurred in Suining County, where its production–ecological space decreased by approximately 4.05% over the decade. The consistency of total scale change in production–living space and change in its internal structure showed that Xuzhou City, part of the agricultural land in Xuzhou City, was converted into construction land, which has led to the weakening of regional agricultural production and ecological functions, and has increased the pressure on land resource allocation.

(3)

The scale of living–production space increased significantly due to urbanization and rural development.

The living–production space in all districts of Xuzhou presented an expansion trend. In the central urban districts, the living–production space reached 637.51 km² over the decade, marking an increase of 56.09 km² from 2012, with a growth rate of 9.65%. In Pizhou City, the living–production space increased by 425.83 km², representing a growth rate of 13.02%. Feng County and Xinyi City also increased by 35.4 km² and 40.43 km², representing increases of 13.2% and 15.01%, respectively. The expansion of the total area of living space revealed that rapid development and urbanization occurred in Xuzhou. It also reflected the increasing emphasis on improving residential and production conditions in NUBUAs.

4. Discussion

The spatial evolution of PLESs in Xuzhou City has been influenced by various factors, including urbanization, economic restructuring, and environmental changes. It is crucial to investigate how PLES shave changed over time and to identify the key driving forces behind these transformations in order to support the development of effective urban planning strategies. In particular, investigating the migration of the center of gravity and the driving factors affecting spatial transitions can provide valuable insights into the complex relationship between urban sprawl and land use. This study has discussed the spatiotemporal evolution patterns of PLESs in Xuzhou and examined the driving factors behind these changes by focusing on quantitative spatial transitions, analyzing migration in the centers of gravity, and investigating driving factors. These findings are expected to provide insights that contribute to the optimization of land-use strategies and the promotion of balanced urban–rural integration, which are essential for sustainable regional development.

4.1. Quantitative Description of Spatial Transitions in PLESs

The land-use transfer matrix [53,54,55] is primarily used to evaluate the rates and directions of transitions between various land-use categories, while also uncovering the interconnections and dynamic trends among these categories [56]. Similarly, the principle of the land-use transfer matrix is applicable to exploring the internal correlations and change trends in PLESs, which can provide a quantitative description of the transfer of various spatial types. The land-use transfer matrix formula is shown in Formula (5).

$${\mathrm{S}}_{\mathrm{i}\mathrm{j}}=\left[\begin{array}{cccc}{\mathrm{S}}_{11}& {\mathrm{S}}_{12}& \cdots & {\mathrm{S}}_{1\mathrm{n}}\\ {\mathrm{S}}_{21}& {\mathrm{S}}_{22}& \cdots & {\mathrm{S}}_{2\mathrm{n}}\\ \cdots & \cdots & \cdots & \cdots \\ {\mathrm{S}}_{\mathrm{n}1}& {\mathrm{S}}_{\mathrm{n}2}& \cdots & {\mathrm{S}}_{\mathrm{n}\mathrm{n}}\end{array}\right]$$

(5)

where S represents the area transferred by the type of PLES; n represents the number of PLES classifications in the study area; i and j indicate the types of PLES before and after the transfer, respectively; and S_ij represents the area of PLES converted from Class i to Class j.

Based on the land-use transfer matrix, the changes in various types of PLESs within Xuzhou were determined.

Table 9. The transfer matrix of PLESs in UBUAs from 2012 to 2022 (km²).

Year	PLES Type	LS	PS	ES	Sub-LS	Sub-PS	Sub-ES	LPS	LES	PES	MS
2012–2018	LS	-	5.21	6.08	32.95	0.14	0	60.07	0	15.18	0.35
	PS	10.01	-	0.12	2	0	0	4.64	0.12	3.78	0
	ES	0.2	0.01	-	0	0	0	0.1	4.33	0	0
	Sub-LS	28.69	3.75	1.45	-	0	0.21	4.94	6.67	4.32	0.06
	Sub-PS	1.83	0.4	0.09	0.44	-	0	4.02	0	2.52	0
	Sub-ES	0.04	0.39	0.11	0.57	0	-	0.18	0.22	0	0
	LPS	27.61	5.51	0.64	8.75	0	0.04	-	3.27	0.26	2.66
	LES	0.15	0	0	0	0	0	2.24	-	0	0
	PES	0.51	0.16	0.07	0.07	0	0	0.61	0	-	0
	MS	0	0	2.57	0	0	0	0	0	0	-
2018–2022	LS	-	3.36	0.08	18.04	6.58	0.06	17.67	16.97	0.34	6.84
	PS	16.15	-	0	2.88	4.03	0.01	1.57	0	0	0
	ES	5.57	0	-	2.09	0.1	0.05	1	0.01	0.03	2.67
	Sub-LS	30.75	0.61	0.01	-	0.35	3.78	10.06	0.03	0.01	1.39
	Sub-PS	0	0	0.14	0	-	0	0	0	0	0
	Sub-ES	0	0.04	0	0	0	-	0	0	0	0
	LPS	54.21	2.32	0.34	15.89	2.18	0.06	-	0.04	0.12	0.07
	LES	3.29	0.12	0	0	0	0	7.87	-	4.33	0
	PES	15.56	2.09	0.01	0.74	1.62	0	6.46	0	-	0
	MS	0.01	0.24	0	0	0	0	2.74	0	0	-

illustrates the specific patterns of spatial transitions among various PLES types within UBUAs. Meanwhile, the corresponding transitions in NUBUAs are presented in

Table 10. The transfer matrix of PLESs in NUBUAs from 2012 to 2022 (km²).

Year	PLES Type	PES	LPS	ES
2012–2018	PES	-	1684.07	308.72
	LPS	598.05	-	97.37
	ES	43.54	12.74	-
2018–2022	PES	-	501.63	129.40
	LPS	346.09	-	21.20
	ES	77.19	27.99	-

.

From 2012 to 2022, the spatial transitions of PLESs in UBUAs showed variable tendencies. Based on the transition matrix in Table 9, living spaces expanded significantly over the period. From 2012 to 2018, approximately 27.61 km² of living spaces were converted into multifunctional spaces such as living–production space, which reflected the growing urban demand for integrated spaces. This trend continued from 2018 to 2022, though at a reduced intensity, indicating a slowing yet persistent shift of living spaces into multifunctional spaces. Based on the findings from these two periods, the conversion of PLES types in UBUAs of Xuzhou City primarily occurred among living space, sub-living space, and living–production space. This phenomenon was largely driven by differences in the proportion of production functions within the identified units, which resulted in similar spatial transformations.

As the most significant area changes, the living space primarily transitioned from production–ecological and ecological spaces. On the basis of requirements on urbanization and industrialization processes in Xuzhou City, the area of living–production space expanded substantially, increasing from 2547.46 km² in 2012 to 3827.48 km² in 2022. From 2012 to 2018, production–ecological space accounted for the largest transfer-out area at 1992.79 km², with 84.51% (1684.07 km²) transitioning into living–production space. In comparison, ecological space had the smallest transfer-out area, totaling just 56.28 km². From 2018 to 2022, the transfer-out area of production–ecological space declined significantly to 631.03 km², while living–production space absorbed the majority of this transfer, including 501.63 km² from production–ecological space and 27.99 km² from ecological space.

4.2. Center of Gravity Migration of PLESs in NUBUAs

The center of gravity migration model in PLESs is the regional weighted average center, which could effectively illustrate the overall patterns of spatial changes in regional land use. Many scholars have used the center of gravity migration models to calculate changes in vegetation [57], land use [58], soil erosion [59], and other factors, providing insights into spatial patterns and change dynamics across various domains. The mathematical expressions of the center of gravity migration model are presented in Formulas (6) and (7).

$$\overline{X}={\displaystyle \frac{\sum _{i=1}^n{M}_i{X}_i}{\sum _{i=1}^n{M}_i}}$$

(6)

$$\overline{Y}={\displaystyle \frac{\sum _{i=1}^n{M}_i{Y}_i}{\sum _{i=1}^n{M}_i}}$$

(7)

where $\overline{X}$ is the latitude of the center of gravity of one type of PLES cover; $\overline{Y}$ is the longitude of the center of gravity of one type of PLES cover; $X_i$ is the latitude of the i-th data unit; $Y_i$ is the longitude of the i-th data unit; and $M_i$ represents the area of the PLES resource in the i-th data unit.

Figure 7 illustrates that during the study period, the centers of gravity of living–production space and production–ecological space were concentrated in the Jiawang District. Meanwhile, the ecological space was focused in the Tongshan District. Specifically, the center of gravity of living–production space migrated to the northwest from 2012 to 2018, and the migration distance was approximately 640.63 m; from 2018 to 2022, it migrated to the northeast with a migration distance of approximately 145.18 m. The production–ecological space mainly migrated to the east, and the migration distances in the two periods from 2012 to 2018 and 2018 to 2022 were approximately 184.89 m and 53.98 m, respectively. The ecological space migrated approximately 1702.23 m to the southwest during 2012–2018, and approximately 391.27 m to the southwest from 2018 to 2022.

4.3. Driving Forces of the Spatiotemporal Evolution in PLESs

4.3.1. Selection of Driving Factors

Based on prior research findings [58,60], eight factors were selected to explore the driving forces behind the trade-offs and synergies in the spatial distribution of PLESs. These factors, categorized into topography, socio-economic conditions, and natural climate, are illustrated in Figure 8. The topographic factors include Digital Elevation Model (DEM, X1), slope (X2), and aspect (X3). The DEM is a representation of the terrain’s elevation, which directly influences water drainage, soil erosion, and vegetation distribution [61]. Slope and aspect significantly affect land suitability for various activities. Steep slopes and unfavorable aspects often limit agricultural development and urban expansion [62]. Social and economic factors, including GDP (X4) and population (X5), reflect the influence of human activities. GDP serves as an indicator of regional economic development, correlating with urbanization intensity, infrastructure growth, and industrial activities [63]. Similarly, population serves as a fundamental determinant of spatial demand [64]. The natural climate factors, including precipitation (X6), temperature (X7), and NPP (X8), are essential environmental indicators that influence regional land-use patterns, resource allocation, and the spatial distribution of economic activities. Precipitation impacts the availability of water resources, which are vital for agricultural production, ecosystem health, and urban water supply [65,66]. Temperature influences ecosystem growth [67] and agricultural productivity, while also influencing human living conditions [68]. In urban areas, elevated temperatures can exacerbate the urban heat island effect, further affecting the demand for urban spaces [69]. NPP represents fundamental natural ecosystem functions, and serves as a crucial indicator for evaluating ecosystem responses, offering an effective measure of the ecological impacts of urbanization [70]. To quantify these driving factors, the natural breaks (Jenks) method, provided by the Reclassify tool of ArcMap 10.2, was employed to individually discretize each factor into nine levels, with values assigned on a scale from 1 to 9.

4.3.2. Driving Force Exploration with Geo-Detector

Geo-detector is a tool proposed by Wang et al. for detecting spatially stratified heterogeneity in 2016 [71]. It consists of four components: factor detection, ecological detection, interaction detection, and risk detection. Geo-detector [72,73] is utilized to analyze spatial disparities in geographic phenomena and to reveal the fundamental factors driving these variations, and it is widely utilized in the field of ecological research. In this study, the discretized values of each driving factor were used as independent variables, and the PLES utilization index calculated by the calculation formula of the land-use intensity served as the dependent variable. These values were provided as input data in a version of Excel Geo-detector, with further details available at http://geodetector.cn/. Finally, the results of factor detection and interaction detection were primarily selected to investigate the driving forces of PLESs. The calculation formulas of Geo-detector are presented in Formula (8).

$$q=1-{\displaystyle \frac{SSW}{SST}}=1-{\displaystyle \frac{\sum _{h=1}^L{N}_h{\sigma }_h^2}{N{\sigma }^2}}$$

(8)

where q [0–1] value presents the contributions of the driving factors to PLESs. A higher q corresponds to the greater influence of the driving factor X. $L$ represents the number of categories for the selected driving factor. $N_h$ and $N$ represent the unit numbers of stratum h and the whole region of interest, respectively. ${\sigma }_h^2$ and ${\sigma }^2$ represent the variance of PLESs for stratum h and the whole area, respectively.

Based on the results of the factor detector, the average q values of various factors, which were calculated from the q values of 2012, 2018, and 2022, are shown in Figure 9. The primary driving forces behind the evolution of PLES types in Xuzhou were socio-economic factors, with population and GDP showing mean q-values of 0.601 and 0.492, respectively. Belonging to the natural climate factors, the annual average precipitation, temperature, and NPP demonstrated a balanced influence on the spatial distribution of PLESs in Xuzhou City. In contrast, the influence of topographic factors was relatively weak. Among the topographic factors, the effect of the DEM was 0.111, while the slope and aspect decreased by 0.013 and 0.001, respectively. The influence of different factors showed significant differentiation in the evolution of PLESs. Exploring the reasons for this phenomenon, it is obvious that rapid economic and social development affected the spatial distribution of urban functional spaces. On the other hand, approximately 90% of Xuzhou City’s area is covered by flat plains. Therefore, the influence of topographic factors on PLES changes was relatively weak compared to socio-economic and natural climate factors.

The interaction detection results (Figure 10) indicated that the distribution and changes of PLESs in Xuzhou City were also influenced by multiple factors together. Nonlinear enhancements and two-factor enhancement effects were observed for each factor in the spatiotemporal evolution of PLESs. Specifically, the two-factor enhancement effect in 2012 was weaker compared to that in 2018 and 2022, with the maximum interaction effect observed at 0.7 (between population and temperature). In contrast, the maximum interaction effects in 2018 and 2022 increased to 0.96 and 0.97, respectively, driven by the interaction between GDP and population. In the development of the PLES pattern, the driving forces of socio-economic and natural climate factors has become increasingly significant.

4.4. Limitation and Prospective

One of the main aspects of our study is the identification of PLESs in Xuzhou City based on combining the POI data and CLCD. POI data that provide abundant information for urban functional zoning and spatial structure analysis are widely available across the country. Meanwhile, the CLCD offers broad coverage and updates regularly, which could provide accurate representations of land-use types across various regions. With the advancement of geographic information sharing and open data policies in China, the availability of these datasets is steadily increasing. Therefore, the method for identifying PLESs based on POI data and CLCD is applicable to most comparable regions. However, due to the heavy reliance on the availability and accuracy of data sources, which could potentially affect the comprehensiveness of the results in areas with limited or inconsistent data coverage, the classification accuracy of our method is influenced by the resolution and reliability of the POI and CLCD data. Under the circumstances, integrating additional high-precision datasets (such as satellite imagery or remote sensing data) to improve classification accuracy is a meaningful research direction.

Furthermore, the detailed analysis of the accuracy validation results for our proposed method in UBUAs revealed instances of misclassification near the boundaries of built-up areas, primarily due to sparse POI distribution causing inaccuracies in PLES identification, as well as misclassifications in areas with no data, where units were frequently identified as similar functional spaces. In the future, incorporating additional data for the correction and verification of misclassifications to enhance the reliability of our study will remain a critical challenge that needs to be resolved.

In addition, the comprehensiveness of selected driving factors determines the accuracy of driving force analysis. The selection of driving factors in this study focused on topography, natural climate, and socio-economic aspects, while neglecting institutional factors such as government policies. These elements often play a crucial role in shaping land-use patterns and spatial configurations. The scope of driving factors should be broadened by incorporating policy orientations, planning strategies, and other institutional influences. This would provide a more comprehensive understanding of the mechanisms driving PLES distribution and evolution, and offer stronger support for evidence-based spatial planning and governance.

5. Conclusions

With the rapid advancement of industrialization and urbanization, the rational distribution and development of urban spaces have been profoundly impacted and shocked. Exploring the distribution patterns and interconnections of PLESs is highly significant for spatial planning and zoning management. In this paper, considering the spatial characteristics of our study area, Xuzhou City was separated into UBUAs and NUBUAs using the Densi-graph method and POI data. Within UBUAs, a PLES identification method based on POI data and CLCD was proposed to identify PLESs in Xuzhou City for the years 2012, 2018, and 2022, with an average recognition accuracy of 89.81%. It proves the feasibility and potential of combining POI data and CLCD in PLES spatial identification, which also provides a theoretical and methodological basis for the exploration of urban spatial planning. However, due to the sparse distribution of POI data in NUBUAs, the CLCD was applied to identify PLESs, which improved the reliability and accuracy of the identification results. Particularly, based on the PLES identification results, the spatiotemporal evolution and driving mechanism of PLES distribution in Xuzhou City were further discussed. The following conclusions were drawn from this study:

(1)

Based on the PLES classification results, UBUAs were mainly composed of living spaces, while NUBUAs were primarily characterized by production–ecological spaces. Additionally, the area of multifunctional space increased from 10,876.98 km² in 2012 to 10,898.55 km² in 2018, and further expanded to 10,941.15 km² by 2022. This continuous growth in multifunctional space indicates the increasingly integrated relationship be-tween living, production, and ecological activities, as well as the growing interaction of spatial functions within Xuzhou City.

(2)

According to the results of the spatial transfer matrix in PLESs, the most significant scale of conversion in UBUAs was the transformation of living space into living–production space. In NUBUAs, the area of ecological space transferred shows a significant increase, growing from 56.28 km² to 105.18 km². However, compared to the areas transferred out from living–production and production–ecological spaces, the transfer-out area of ecological space remains the smallest. The effective ecological restoration and protection efforts of Xuzhou City have reduced the loss of ecological space and have contributed to enhancing environmental sustainability.

(3)

The centers of gravity for living–production space and production–ecological space in the years 2012, 2018, and 2022 consistently remained within Jiawang District, while the centers of gravity for ecological space were located in Tongshan District. From 2012 to 2018, the moving distances for living–production space, production–ecological space, and ecological space were 640.63 m, 184.89 m, and 1702.23 m, respectively. However, the moving distances were reduced to 145.18 m, 53.98 m, and 391.27 m from 2018 to 2022, respectively. During 2012 to 2022, the center of gravity of PLESs remained relatively stable.

(4)

Driving force analysis revealed that socio-economic and natural climate factors were the main drivers of the evolution of PLESs. The primary social factors influencing the spatial variations of PLESs were GDP and population; the primary natural factor influencing the geographical variations of PLESs is precipitation. The interaction between GDP and population had the greatest effect on the spatial changes of PLESs.

In conclusion, this study presents a detailed, scientific, and practical method for identifying PLESs, contributing valuable insights for the sustainable development of urban spatial planning. The spatiotemporal evolution patterns and driving forces of PLESs were analyzed, providing a comprehensive understanding of the factors influencing PLES distribution. This understanding is crucial for guiding urban planning and land-use strategies in rapidly growing cities. By identifying the trade-offs and synergies among different spatial functions, this study offers key insights to optimize land allocation, ensuring a balanced and sustainable approach to development and ecological conservation. Additionally, the findings provide a theoretical foundation and important guidance for managing the spatial balance of Xuzhou City, while also offering a scientific framework for managing geospatial resources in comparable regions.